Acts/Propagator/Propagator.ipp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #include "Acts/EventData/TrackParametersConcept.hpp"
0010 #include "Acts/Propagator/ActorList.hpp"
0011 #include "Acts/Propagator/ConstrainedStep.hpp"
0012 #include "Acts/Propagator/NavigationTarget.hpp"
0013 #include "Acts/Propagator/PropagatorError.hpp"
0014 #include "Acts/Propagator/StandardAborters.hpp"
0015 #include "Acts/Propagator/detail/LoopProtection.hpp"
0016 #include "Acts/Surfaces/BoundaryTolerance.hpp"
0017 #include "Acts/Utilities/Intersection.hpp"
0018
0019 #include <concepts>
0020
0021 template <typename S, typename N>
0022 template <typename propagator_state_t>
0023 Acts::Result<void> Acts::Propagator<S, N>::propagate(
0024     propagator_state_t& state) const {
0025   ACTS_VERBOSE("Entering propagation.");
0026
0027   state.stage = PropagatorStage::prePropagation;
0028
0029   // Pre-Propagation: call to the actor list, abort condition check
0030   state.options.actorList.act(state, m_stepper, m_navigator, logger());
0031
0032   if (state.options.actorList.checkAbort(state, m_stepper, m_navigator,
0033                                          logger())) {
0034     ACTS_VERBOSE("Propagation terminated without going into stepping loop.");
0035
0036     state.stage = PropagatorStage::postPropagation;
0037
0038     state.options.actorList.act(state, m_stepper, m_navigator, logger());
0039
0040     return Result<void>::success();
0041   }
0042
0043   auto getNextTarget = [&]() -> Result<NavigationTarget> {
0044     for (unsigned int i = 0; i < state.options.maxTargetSkipping; ++i) {
0045       NavigationTarget nextTarget = m_navigator.nextTarget(
0046           state.navigation, state.position, state.direction);
0047       if (nextTarget.isNone()) {
0048         return NavigationTarget::None();
0049       }
0050       IntersectionStatus preStepSurfaceStatus = m_stepper.updateSurfaceStatus(
0051           state.stepping, *nextTarget.surface,
0052           nextTarget.surfaceIntersectionIndex, state.options.direction,
0053           nextTarget.boundaryTolerance, state.options.surfaceTolerance,
0054           ConstrainedStep::Type::Navigator, logger());
0055       if (preStepSurfaceStatus == IntersectionStatus::reachable ||
0056           preStepSurfaceStatus == IntersectionStatus::onSurface) {
0057         return nextTarget;
0058       }
0059     }
0060
0061     ACTS_ERROR("getNextTarget failed to find a valid target surface after "
0062                << state.options.maxTargetSkipping << " attempts.");
0063     return Result<NavigationTarget>::failure(
0064         PropagatorError::NextTargetLimitReached);
0065   };
0066
0067   // priming error condition
0068   bool terminatedNormally = false;
0069
0070   // Pre-Stepping: target setting
0071   state.stage = PropagatorStage::preStep;
0072
0073   Result<NavigationTarget> nextTargetResult = getNextTarget();
0074   if (!nextTargetResult.ok()) {
0075     return nextTargetResult.error();
0076   }
0077   NavigationTarget nextTarget = *nextTargetResult;
0078
0079   ACTS_VERBOSE("Starting stepping loop.");
0080
0081   // Stepping loop
0082   for (; state.steps < state.options.maxSteps; ++state.steps) {
0083     // Perform a step
0084     Result<double> res =
0085         m_stepper.step(state.stepping, state.options.direction,
0086                        m_navigator.currentVolumeMaterial(state.navigation));
0087     if (!res.ok()) {
0088       ACTS_ERROR("Step failed with " << res.error() << ": "
0089                                      << res.error().message());
0090       // pass error to caller
0091       return res.error();
0092     }
0093     // Accumulate the path length
0094     state.pathLength += *res;
0095     // Update the position and direction
0096     state.position = m_stepper.position(state.stepping);
0097     state.direction =
0098         state.options.direction * m_stepper.direction(state.stepping);
0099
0100     ACTS_VERBOSE("Step with size " << *res << " performed. We are now at "
0101                                    << state.position.transpose()
0102                                    << " with direction "
0103                                    << state.direction.transpose());
0104
0105     // release actor and aborter constrains after step was performed
0106     m_stepper.releaseStepSize(state.stepping, ConstrainedStep::Type::Navigator);
0107     m_stepper.releaseStepSize(state.stepping, ConstrainedStep::Type::Actor);
0108
0109     // Post-stepping: check target status, call actors, check abort conditions
0110     state.stage = PropagatorStage::postStep;
0111
0112     if (!nextTarget.isNone()) {
0113       IntersectionStatus postStepSurfaceStatus = m_stepper.updateSurfaceStatus(
0114           state.stepping, *nextTarget.surface,
0115           nextTarget.surfaceIntersectionIndex, state.options.direction,
0116           nextTarget.boundaryTolerance, state.options.surfaceTolerance,
0117           ConstrainedStep::Type::Navigator, logger());
0118       if (postStepSurfaceStatus == IntersectionStatus::onSurface) {
0119         m_navigator.handleSurfaceReached(state.navigation, state.position,
0120                                          state.direction, *nextTarget.surface);
0121       }
0122       if (postStepSurfaceStatus != IntersectionStatus::reachable) {
0123         nextTarget = NavigationTarget::None();
0124       }
0125     }
0126
0127     state.options.actorList.act(state, m_stepper, m_navigator, logger());
0128
0129     if (state.options.actorList.checkAbort(state, m_stepper, m_navigator,
0130                                            logger())) {
0131       terminatedNormally = true;
0132       break;
0133     }
0134
0135     // Update the position and direction because actors might have changed it
0136     state.position = m_stepper.position(state.stepping);
0137     state.direction =
0138         state.options.direction * m_stepper.direction(state.stepping);
0139
0140     // Pre-Stepping: target setting
0141     state.stage = PropagatorStage::preStep;
0142
0143     if (!nextTarget.isNone() &&
0144         !m_navigator.checkTargetValid(state.navigation, state.position,
0145                                       state.direction)) {
0146       ACTS_VERBOSE("Target is not valid anymore.");
0147       nextTarget = NavigationTarget::None();
0148     }
0149
0150     if (nextTarget.isNone()) {
0151       // navigator step constraint is not valid anymore
0152       m_stepper.releaseStepSize(state.stepping,
0153                                 ConstrainedStep::Type::Navigator);
0154
0155       nextTargetResult = getNextTarget();
0156       if (!nextTargetResult.ok()) {
0157         return nextTargetResult.error();
0158       }
0159       nextTarget = *nextTargetResult;
0160     }
0161   }  // end of stepping loop
0162
0163   // check if we didn't terminate normally via aborters
0164   if (!terminatedNormally) {
0165     ACTS_ERROR("Propagation reached the step count limit of "
0166                << state.options.maxSteps << " (did " << state.steps
0167                << " steps)");
0168     return PropagatorError::StepCountLimitReached;
0169   }
0170
0171   ACTS_VERBOSE("Stepping loop done.");
0172
0173   state.stage = PropagatorStage::postPropagation;
0174
0175   // Post-stepping call to the actor list
0176   state.options.actorList.act(state, m_stepper, m_navigator, logger());
0177
0178   return Result<void>::success();
0179 }
0180
0181 template <typename S, typename N>
0182 template <typename parameters_t, typename propagator_options_t,
0183           typename path_aborter_t>
0184 auto Acts::Propagator<S, N>::propagate(const parameters_t& start,
0185                                        const propagator_options_t& options,
0186                                        bool makeCurvilinear) const
0187     -> Result<
0188         actor_list_t_result_t<StepperCurvilinearTrackParameters,
0189                               typename propagator_options_t::actor_list_type>> {
0190   static_assert(std::copy_constructible<StepperCurvilinearTrackParameters>,
0191                 "return track parameter type must be copy-constructible");
0192
0193   auto state = makeState<propagator_options_t, path_aborter_t>(options);
0194
0195   auto initRes =
0196       initialize<decltype(state), parameters_t, path_aborter_t>(state, start);
0197   if (!initRes.ok()) {
0198     return initRes.error();
0199   }
0200
0201   // Perform the actual propagation
0202   auto propagationResult = propagate(state);
0203
0204   return makeResult(std::move(state), propagationResult, options,
0205                     makeCurvilinear);
0206 }
0207
0208 template <typename S, typename N>
0209 template <typename parameters_t, typename propagator_options_t,
0210           typename target_aborter_t, typename path_aborter_t>
0211 auto Acts::Propagator<S, N>::propagate(
0212     const parameters_t& start, const Surface& target,
0213     const propagator_options_t& options) const
0214     -> Result<
0215         actor_list_t_result_t<StepperBoundTrackParameters,
0216                               typename propagator_options_t::actor_list_type>> {
0217   static_assert(BoundTrackParametersConcept<parameters_t>,
0218                 "Parameters do not fulfill bound parameters concept.");
0219
0220   auto state =
0221       makeState<propagator_options_t, target_aborter_t, path_aborter_t>(
0222           target, options);
0223
0224   auto initRes =
0225       initialize<decltype(state), parameters_t, path_aborter_t>(state, start);
0226   if (!initRes.ok()) {
0227     return initRes.error();
0228   }
0229
0230   // Perform the actual propagation
0231   auto propagationResult = propagate(state);
0232
0233   return makeResult(std::move(state), propagationResult, target, options);
0234 }
0235
0236 template <typename S, typename N>
0237 template <typename propagator_options_t, typename path_aborter_t>
0238 auto Acts::Propagator<S, N>::makeState(
0239     const propagator_options_t& options) const {
0240   // Type of track parameters produced by the propagation
0241   using ReturnParameterType = StepperCurvilinearTrackParameters;
0242
0243   static_assert(std::copy_constructible<ReturnParameterType>,
0244                 "return track parameter type must be copy-constructible");
0245
0246   // Expand the actor list with a path aborter
0247   path_aborter_t pathAborter;
0248   pathAborter.internalLimit = options.pathLimit;
0249
0250   auto actorList = options.actorList.append(pathAborter);
0251
0252   // Create the extended options and declare their type
0253   auto eOptions = options.extend(actorList);
0254
0255   using OptionsType = decltype(eOptions);
0256   using StateType =
0257       actor_list_t_state_t<OptionsType,
0258                            typename propagator_options_t::actor_list_type>;
0259
0260   StateType state{eOptions, m_stepper.makeState(eOptions.stepping),
0261                   m_navigator.makeState(eOptions.navigation)};
0262
0263   return state;
0264 }
0265
0266 template <typename S, typename N>
0267 template <typename propagator_options_t, typename target_aborter_t,
0268           typename path_aborter_t>
0269 auto Acts::Propagator<S, N>::makeState(
0270     const Surface& target, const propagator_options_t& options) const {
0271   // Expand the actor list with a target and path aborter
0272   target_aborter_t targetAborter;
0273   targetAborter.surface = &target;
0274   path_aborter_t pathAborter;
0275   pathAborter.internalLimit = options.pathLimit;
0276
0277   auto actorList = options.actorList.append(targetAborter, pathAborter);
0278
0279   // Create the extended options and declare their type
0280   auto eOptions = options.extend(actorList);
0281   eOptions.navigation.targetSurface = &target;
0282
0283   using OptionsType = decltype(eOptions);
0284   using StateType =
0285       actor_list_t_state_t<OptionsType,
0286                            typename propagator_options_t::actor_list_type>;
0287
0288   StateType state{eOptions, m_stepper.makeState(eOptions.stepping),
0289                   m_navigator.makeState(eOptions.navigation)};
0290
0291   return state;
0292 }
0293
0294 template <typename S, typename N>
0295 template <typename propagator_state_t, typename parameters_t,
0296           typename path_aborter_t>
0297 Acts::Result<void> Acts::Propagator<S, N>::initialize(
0298     propagator_state_t& state, const parameters_t& start) const {
0299   static_assert(BoundTrackParametersConcept<parameters_t>,
0300                 "Parameters do not fulfill bound parameters concept.");
0301
0302   m_stepper.initialize(state.stepping, start);
0303
0304   state.position = m_stepper.position(state.stepping);
0305   state.direction =
0306       state.options.direction * m_stepper.direction(state.stepping);
0307
0308   state.navigation.options.startSurface = &start.referenceSurface();
0309
0310   // Navigator initialize state call
0311   auto navInitRes =
0312       m_navigator.initialize(state.navigation, state.position, state.direction,
0313                              state.options.direction);
0314   if (!navInitRes.ok()) {
0315     return navInitRes.error();
0316   }
0317
0318   // Apply the loop protection - it resets the internal path limit
0319   detail::setupLoopProtection(
0320       state, m_stepper, state.options.actorList.template get<path_aborter_t>(),
0321       false, logger());
0322
0323   return Result<void>::success();
0324 }
0325
0326 template <typename S, typename N>
0327 template <typename propagator_state_t, typename propagator_options_t>
0328 auto Acts::Propagator<S, N>::makeResult(propagator_state_t state,
0329                                         Result<void> propagationResult,
0330                                         const propagator_options_t& /*options*/,
0331                                         bool makeCurvilinear) const
0332     -> Result<
0333         actor_list_t_result_t<StepperCurvilinearTrackParameters,
0334                               typename propagator_options_t::actor_list_type>> {
0335   // Type of track parameters produced by the propagation
0336   using ReturnParameterType = StepperCurvilinearTrackParameters;
0337
0338   static_assert(std::copy_constructible<ReturnParameterType>,
0339                 "return track parameter type must be copy-constructible");
0340
0341   // Type of the full propagation result, including output from actors
0342   using ResultType =
0343       actor_list_t_result_t<ReturnParameterType,
0344                             typename propagator_options_t::actor_list_type>;
0345
0346   if (!propagationResult.ok()) {
0347     return propagationResult.error();
0348   }
0349
0350   ResultType result{};
0351   moveStateToResult(state, result);
0352
0353   if (makeCurvilinear) {
0354     if (!m_stepper.prepareCurvilinearState(state.stepping)) {
0355       // information to compute curvilinearState is incomplete.
0356       return propagationResult.error();
0357     }
0358     /// Convert into return type and fill the result object
0359     auto curvState = m_stepper.curvilinearState(state.stepping);
0360     // Fill the end parameters
0361     result.endParameters =
0362         std::get<StepperCurvilinearTrackParameters>(curvState);
0363     // Only fill the transport jacobian when covariance transport was done
0364     if (state.stepping.covTransport) {
0365       result.transportJacobian = std::get<Jacobian>(curvState);
0366     }
0367   }
0368
0369   return Result<ResultType>::success(std::move(result));
0370 }
0371
0372 template <typename S, typename N>
0373 template <typename propagator_state_t, typename propagator_options_t>
0374 auto Acts::Propagator<S, N>::makeResult(
0375     propagator_state_t state, Result<void> propagationResult,
0376     const Surface& target, const propagator_options_t& /*options*/) const
0377     -> Result<
0378         actor_list_t_result_t<StepperBoundTrackParameters,
0379                               typename propagator_options_t::actor_list_type>> {
0380   // Type of track parameters produced at the end of the propagation
0381   using ReturnParameterType = StepperBoundTrackParameters;
0382
0383   static_assert(std::copy_constructible<ReturnParameterType>,
0384                 "return track parameter type must be copy-constructible");
0385
0386   // Type of the full propagation result, including output from actors
0387   using ResultType =
0388       actor_list_t_result_t<ReturnParameterType,
0389                             typename propagator_options_t::actor_list_type>;
0390
0391   if (!propagationResult.ok()) {
0392     return propagationResult.error();
0393   }
0394
0395   ResultType result{};
0396   moveStateToResult(state, result);
0397
0398   // Compute the final results and mark the propagation as successful
0399   auto bsRes = m_stepper.boundState(state.stepping, target);
0400   if (!bsRes.ok()) {
0401     return bsRes.error();
0402   }
0403   const auto& bs = *bsRes;
0404
0405   // Fill the end parameters
0406   result.endParameters = std::get<StepperBoundTrackParameters>(bs);
0407   // Only fill the transport jacobian when covariance transport was done
0408   if (state.stepping.covTransport) {
0409     result.transportJacobian = std::get<Jacobian>(bs);
0410   }
0411   return Result<ResultType>::success(std::move(result));
0412 }
0413
0414 template <typename S, typename N>
0415 template <typename propagator_state_t, typename result_t>
0416 void Acts::Propagator<S, N>::moveStateToResult(propagator_state_t& state,
0417                                                result_t& result) const {
0418   result.tuple() = std::move(state.tuple());
0419
0420   result.steps = state.steps;
0421   result.pathLength = state.pathLength;
0422
0423   result.statistics.stepping = state.stepping.statistics;
0424   result.statistics.navigation = state.navigation.statistics;
0425 }
0426
0427 template <typename derived_t>
0428 Acts::Result<Acts::BoundTrackParameters>
0429 Acts::detail::BasePropagatorHelper<derived_t>::propagateToSurface(
0430     const BoundTrackParameters& start, const Surface& target,
0431     const Options& options) const {
0432   using ResultType = Result<typename derived_t::template actor_list_t_result_t<
0433       BoundTrackParameters, ActorList<>>>;
0434   using DerivedOptions = typename derived_t::template Options<>;
0435
0436   DerivedOptions derivedOptions(options);
0437
0438   // dummy initialization
0439   ResultType res = ResultType::failure(PropagatorError::Failure);
0440
0441   // Due to the geometry of the perigee surface the overstepping tolerance
0442   // is sometimes not met.
0443   if (target.type() == Surface::SurfaceType::Perigee) {
0444     res = static_cast<const derived_t*>(this)
0445               ->template propagate<BoundTrackParameters, DerivedOptions,
0446                                    ForcedSurfaceReached, PathLimitReached>(
0447                   start, target, derivedOptions);
0448   } else {
0449     res = static_cast<const derived_t*>(this)
0450               ->template propagate<BoundTrackParameters, DerivedOptions,
0451                                    SurfaceReached, PathLimitReached>(
0452                   start, target, derivedOptions);
0453   }
0454
0455   if (!res.ok()) {
0456     return res.error();
0457   }
0458
0459   // Without errors we can expect a valid endParameters when propagating to a
0460   // target surface
0461   assert((*res).endParameters);
0462   return std::move((*res).endParameters.value());
0463 }