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