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