Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-08-28 08:11:33

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/Definitions/Algebra.hpp"
0012 #include "Acts/Definitions/Direction.hpp"
0013 #include "Acts/Definitions/TrackParametrization.hpp"
0014 #include "Acts/EventData/MultiComponentTrackParameters.hpp"
0015 #include "Acts/EventData/TrackParameters.hpp"
0016 #include "Acts/EventData/detail/CorrectedTransformationFreeToBound.hpp"
0017 #include "Acts/MagneticField/MagneticFieldProvider.hpp"
0018 #include "Acts/Material/IVolumeMaterial.hpp"
0019 #include "Acts/Propagator/ConstrainedStep.hpp"
0020 #include "Acts/Propagator/StepperConcept.hpp"
0021 #include "Acts/Propagator/StepperOptions.hpp"
0022 #include "Acts/Propagator/StepperStatistics.hpp"
0023 #include "Acts/Propagator/detail/LoopStepperUtils.hpp"
0024 #include "Acts/Propagator/detail/SteppingHelper.hpp"
0025 #include "Acts/Surfaces/Surface.hpp"
0026 #include "Acts/Utilities/Intersection.hpp"
0027 #include "Acts/Utilities/Logger.hpp"
0028 #include "Acts/Utilities/Result.hpp"
0029 
0030 #include <algorithm>
0031 #include <cmath>
0032 #include <cstddef>
0033 #include <limits>
0034 #include <sstream>
0035 #include <vector>
0036 
0037 #include <boost/container/small_vector.hpp>
0038 
0039 namespace Acts {
0040 
0041 using namespace Acts::UnitLiterals;
0042 
0043 namespace detail {
0044 
0045 struct MaxMomentumComponent {
0046   template <typename component_range_t>
0047   auto operator()(const component_range_t& cmps) const {
0048     return std::max_element(cmps.begin(), cmps.end(),
0049                             [&](const auto& a, const auto& b) {
0050                               return std::abs(a.state.pars[eFreeQOverP]) >
0051                                      std::abs(b.state.pars[eFreeQOverP]);
0052                             });
0053   }
0054 };
0055 
0056 struct MaxWeightComponent {
0057   template <typename component_range_t>
0058   auto operator()(const component_range_t& cmps) {
0059     return std::max_element(
0060         cmps.begin(), cmps.end(),
0061         [&](const auto& a, const auto& b) { return a.weight < b.weight; });
0062   }
0063 };
0064 
0065 template <typename component_chooser_t>
0066 struct SingleComponentReducer {
0067   template <typename stepper_state_t>
0068   static Vector3 position(const stepper_state_t& s) {
0069     return component_chooser_t{}(s.components)
0070         ->state.pars.template segment<3>(eFreePos0);
0071   }
0072 
0073   template <typename stepper_state_t>
0074   static Vector3 direction(const stepper_state_t& s) {
0075     return component_chooser_t{}(s.components)
0076         ->state.pars.template segment<3>(eFreeDir0);
0077   }
0078 
0079   template <typename stepper_state_t>
0080   static double qOverP(const stepper_state_t& s) {
0081     const auto cmp = component_chooser_t{}(s.components);
0082     return cmp->state.pars[eFreeQOverP];
0083   }
0084 
0085   template <typename stepper_state_t>
0086   static double absoluteMomentum(const stepper_state_t& s) {
0087     const auto cmp = component_chooser_t{}(s.components);
0088     return s.particleHypothesis.extractMomentum(cmp->state.pars[eFreeQOverP]);
0089   }
0090 
0091   template <typename stepper_state_t>
0092   static Vector3 momentum(const stepper_state_t& s) {
0093     const auto cmp = component_chooser_t{}(s.components);
0094     return s.particleHypothesis.extractMomentum(cmp->state.pars[eFreeQOverP]) *
0095            cmp->state.pars.template segment<3>(eFreeDir0);
0096   }
0097 
0098   template <typename stepper_state_t>
0099   static double charge(const stepper_state_t& s) {
0100     const auto cmp = component_chooser_t{}(s.components);
0101     return s.particleHypothesis.extractCharge(cmp->state.pars[eFreeQOverP]);
0102   }
0103 
0104   template <typename stepper_state_t>
0105   static double time(const stepper_state_t& s) {
0106     return component_chooser_t{}(s.components)->state.pars[eFreeTime];
0107   }
0108 
0109   template <typename stepper_state_t>
0110   static FreeVector pars(const stepper_state_t& s) {
0111     return component_chooser_t{}(s.components)->state.pars;
0112   }
0113 
0114   template <typename stepper_state_t>
0115   static FreeVector cov(const stepper_state_t& s) {
0116     return component_chooser_t{}(s.components)->state.cov;
0117   }
0118 };
0119 
0120 }  // namespace detail
0121 
0122 using MaxMomentumReducerLoop =
0123     detail::SingleComponentReducer<detail::MaxMomentumComponent>;
0124 using MaxWeightReducerLoop =
0125     detail::SingleComponentReducer<detail::MaxWeightComponent>;
0126 
0127 /// @brief Stepper based on a single-component stepper, but can handle
0128 /// Multi-Component Tracks (e.g., for the GSF). Internally, this only
0129 /// manages a vector of states of the single stepper. This simplifies
0130 /// implementation, but has several drawbacks:
0131 /// * There are certain redundancies between the global State and the
0132 /// component states
0133 /// * The components do not share a single magnetic-field-cache
0134 /// @tparam sstepper_t The single-component stepper type to use
0135 /// @tparam component_reducer_t How to map the multi-component state to a single
0136 /// component
0137 template <Concepts::SingleStepper single_stepper_t,
0138           typename component_reducer_t = MaxWeightReducerLoop>
0139 class MultiStepperLoop : public single_stepper_t {
0140   /// Limits the number of steps after at least one component reached the
0141   /// surface
0142   std::size_t m_stepLimitAfterFirstComponentOnSurface = 50;
0143 
0144   /// The logger (used if no logger is provided by caller of methods)
0145   std::unique_ptr<const Acts::Logger> m_logger;
0146 
0147   /// Small vector type for speeding up some computations where we need to
0148   /// accumulate stuff of components. We think 16 is a reasonable amount here.
0149   template <typename T>
0150   using SmallVector = boost::container::small_vector<T, 16>;
0151 
0152  public:
0153   /// @brief Typedef to the Single-Component Eigen Stepper
0154   using SingleStepper = single_stepper_t;
0155 
0156   /// @brief Typedef to the Single-Component Stepper Options
0157   using SingleOptions = typename SingleStepper::Options;
0158 
0159   /// @brief Typedef to the State of the single component Stepper
0160   using SingleState = typename SingleStepper::State;
0161 
0162   /// @brief Typedef to the Config of the single component Stepper
0163   using SingleConfig = typename SingleStepper::Config;
0164 
0165   /// @brief Use the definitions from the Single-stepper
0166   using typename SingleStepper::Covariance;
0167   using typename SingleStepper::Jacobian;
0168 
0169   /// @brief Define an own bound state
0170   using BoundState =
0171       std::tuple<MultiComponentBoundTrackParameters, Jacobian, double>;
0172 
0173   /// @brief The reducer type
0174   using Reducer = component_reducer_t;
0175 
0176   /// @brief How many components can this stepper manage?
0177   static constexpr int maxComponents = std::numeric_limits<int>::max();
0178 
0179   struct Config : public SingleStepper::Config {
0180     /// Limits the number of steps after at least one component reached the
0181     /// surface
0182     std::size_t stepLimitAfterFirstComponentOnSurface = 50;
0183   };
0184 
0185   struct Options : public SingleOptions {
0186     Options(const GeometryContext& gctx, const MagneticFieldContext& mctx)
0187         : SingleOptions(gctx, mctx) {}
0188 
0189     void setPlainOptions(const StepperPlainOptions& options) {
0190       static_cast<StepperPlainOptions&>(*this) = options;
0191     }
0192   };
0193 
0194   struct State {
0195     /// The struct that stores the individual components
0196     struct Component {
0197       SingleState state;
0198       double weight;
0199       IntersectionStatus status;
0200 
0201       Component(SingleState state_, double weight_, IntersectionStatus status_)
0202           : state(std::move(state_)), weight(weight_), status(status_) {}
0203     };
0204 
0205     Options options;
0206 
0207     /// Particle hypothesis
0208     ParticleHypothesis particleHypothesis = ParticleHypothesis::pion();
0209 
0210     /// The components of which the state consists
0211     SmallVector<Component> components;
0212 
0213     bool covTransport = false;
0214     double pathAccumulated = 0.;
0215     std::size_t steps = 0;
0216 
0217     /// Step-limit counter which limits the number of steps when one component
0218     /// reached a surface
0219     std::optional<std::size_t> stepCounterAfterFirstComponentOnSurface;
0220 
0221     /// The stepper statistics
0222     StepperStatistics statistics;
0223 
0224     /// Constructor from the initial bound track parameters
0225     ///
0226     /// @param [in] optionsIn is the options object for the stepper
0227     ///
0228     /// @note the covariance matrix is copied when needed
0229     explicit State(const Options& optionsIn) : options(optionsIn) {}
0230   };
0231 
0232   /// Constructor from a magnetic field and a optionally provided Logger
0233   /// TODO this requires that every stepper can be constructed like this...
0234   explicit MultiStepperLoop(std::shared_ptr<const MagneticFieldProvider> bField,
0235                             std::unique_ptr<const Logger> logger =
0236                                 getDefaultLogger("GSF", Logging::INFO))
0237       : SingleStepper(std::move(bField)), m_logger(std::move(logger)) {}
0238 
0239   /// Constructor from a configuration and optionally provided Logger
0240   explicit MultiStepperLoop(const Config& config,
0241                             std::unique_ptr<const Logger> logger =
0242                                 getDefaultLogger("MultiStepperLoop",
0243                                                  Logging::INFO))
0244       : SingleStepper(config),
0245         m_stepLimitAfterFirstComponentOnSurface(
0246             config.stepLimitAfterFirstComponentOnSurface),
0247         m_logger(std::move(logger)) {}
0248 
0249   State makeState(const Options& options) const {
0250     State state(options);
0251     return state;
0252   }
0253 
0254   void initialize(State& state,
0255                   const MultiComponentBoundTrackParameters& par) const {
0256     if (par.components().empty()) {
0257       throw std::invalid_argument(
0258           "Cannot construct MultiEigenStepperLoop::State with empty "
0259           "multi-component parameters");
0260     }
0261 
0262     state.particleHypothesis = par.particleHypothesis();
0263 
0264     const auto surface = par.referenceSurface().getSharedPtr();
0265 
0266     for (auto i = 0ul; i < par.components().size(); ++i) {
0267       const auto& [weight, singlePars] = par[i];
0268       auto& cmp =
0269           state.components.emplace_back(SingleStepper::makeState(state.options),
0270                                         weight, IntersectionStatus::onSurface);
0271       SingleStepper::initialize(cmp.state, singlePars);
0272     }
0273 
0274     if (std::get<2>(par.components().front())) {
0275       state.covTransport = true;
0276     }
0277   }
0278 
0279   /// A proxy struct which allows access to a single component of the
0280   /// multi-component state. It has the semantics of a const reference, i.e.
0281   /// it requires a const reference of the single-component state it
0282   /// represents
0283   using ConstComponentProxy =
0284       detail::LoopComponentProxyBase<const typename State::Component,
0285                                      MultiStepperLoop>;
0286 
0287   /// A proxy struct which allows access to a single component of the
0288   /// multi-component state. It has the semantics of a mutable reference, i.e.
0289   /// it requires a mutable reference of the single-component state it
0290   /// represents
0291   using ComponentProxy =
0292       detail::LoopComponentProxy<typename State::Component, MultiStepperLoop>;
0293 
0294   /// Creates an iterable which can be plugged into a range-based for-loop to
0295   /// iterate over components
0296   /// @note Use a for-loop with by-value semantics, since the Iterable returns a
0297   /// proxy internally holding a reference
0298   auto componentIterable(State& state) const {
0299     struct Iterator {
0300       using difference_type [[maybe_unused]] = std::ptrdiff_t;
0301       using value_type [[maybe_unused]] = ComponentProxy;
0302       using reference [[maybe_unused]] = ComponentProxy;
0303       using pointer [[maybe_unused]] = void;
0304       using iterator_category [[maybe_unused]] = std::forward_iterator_tag;
0305 
0306       typename decltype(state.components)::iterator it;
0307       const State& s;
0308 
0309       // clang-format off
0310       auto& operator++() { ++it; return *this; }
0311       auto operator==(const Iterator& other) const { return it == other.it; }
0312       auto operator*() const { return ComponentProxy(*it, s); }
0313       // clang-format on
0314     };
0315 
0316     struct Iterable {
0317       State& s;
0318 
0319       // clang-format off
0320       auto begin() { return Iterator{s.components.begin(), s}; }
0321       auto end() { return Iterator{s.components.end(), s}; }
0322       // clang-format on
0323     };
0324 
0325     return Iterable{state};
0326   }
0327 
0328   /// Creates an constant iterable which can be plugged into a range-based
0329   /// for-loop to iterate over components
0330   /// @note Use a for-loop with by-value semantics, since the Iterable returns a
0331   /// proxy internally holding a reference
0332   auto constComponentIterable(const State& state) const {
0333     struct ConstIterator {
0334       using difference_type [[maybe_unused]] = std::ptrdiff_t;
0335       using value_type [[maybe_unused]] = ConstComponentProxy;
0336       using reference [[maybe_unused]] = ConstComponentProxy;
0337       using pointer [[maybe_unused]] = void;
0338       using iterator_category [[maybe_unused]] = std::forward_iterator_tag;
0339 
0340       typename decltype(state.components)::const_iterator it;
0341       const State& s;
0342 
0343       // clang-format off
0344       auto& operator++() { ++it; return *this; }
0345       auto operator==(const ConstIterator& other) const { return it == other.it; }
0346       auto operator*() const { return ConstComponentProxy{*it}; }
0347       // clang-format on
0348     };
0349 
0350     struct Iterable {
0351       const State& s;
0352 
0353       // clang-format off
0354       auto begin() const { return ConstIterator{s.components.cbegin(), s}; }
0355       auto end() const { return ConstIterator{s.components.cend(), s}; }
0356       // clang-format on
0357     };
0358 
0359     return Iterable{state};
0360   }
0361 
0362   /// Get the number of components
0363   ///
0364   /// @param state [in,out] The stepping state (thread-local cache)
0365   std::size_t numberComponents(const State& state) const {
0366     return state.components.size();
0367   }
0368 
0369   /// Remove missed components from the component state
0370   ///
0371   /// @param state [in,out] The stepping state (thread-local cache)
0372   void removeMissedComponents(State& state) const {
0373     auto new_end = std::remove_if(
0374         state.components.begin(), state.components.end(), [](const auto& cmp) {
0375           return cmp.status == IntersectionStatus::unreachable;
0376         });
0377 
0378     state.components.erase(new_end, state.components.end());
0379   }
0380 
0381   /// Reweight the components
0382   ///
0383   /// @param [in,out] state The stepping state (thread-local cache)
0384   void reweightComponents(State& state) const {
0385     double sumOfWeights = 0.0;
0386     for (const auto& cmp : state.components) {
0387       sumOfWeights += cmp.weight;
0388     }
0389     for (auto& cmp : state.components) {
0390       cmp.weight /= sumOfWeights;
0391     }
0392   }
0393 
0394   /// Reset the number of components
0395   ///
0396   /// @param [in,out] state  The stepping state (thread-local cache)
0397   void clearComponents(State& state) const { state.components.clear(); }
0398 
0399   /// Add a component to the Multistepper
0400   ///
0401   /// @param [in,out] state  The stepping state (thread-local cache)
0402   /// @param [in] pars Parameters of the component to add
0403   /// @param [in] weight Weight of the component to add
0404   ///
0405   /// @note: It is not ensured that the weights are normalized afterwards
0406   /// @note This function makes no garantuees about how new components are
0407   /// initialized, it is up to the caller to ensure that all components are
0408   /// valid in the end.
0409   /// @note The returned component-proxy is only garantueed to be valid until
0410   /// the component number is again modified
0411   Result<ComponentProxy> addComponent(State& state,
0412                                       const BoundTrackParameters& pars,
0413                                       double weight) const {
0414     auto& cmp =
0415         state.components.emplace_back(SingleStepper::makeState(state.options),
0416                                       weight, IntersectionStatus::onSurface);
0417     SingleStepper::initialize(cmp.state, pars);
0418 
0419     return ComponentProxy{state.components.back(), state};
0420   }
0421 
0422   /// Get the field for the stepping, it checks first if the access is still
0423   /// within the Cell, and updates the cell if necessary.
0424   ///
0425   /// @param [in,out] state is the propagation state associated with the track
0426   ///                 the magnetic field cell is used (and potentially updated)
0427   /// @param [in] pos is the field position
0428   ///
0429   /// @note This uses the cache of the first component stored in the state
0430   Result<Vector3> getField(State& state, const Vector3& pos) const {
0431     return SingleStepper::getField(state.components.front().state, pos);
0432   }
0433 
0434   /// Global particle position accessor
0435   ///
0436   /// @param state [in] The stepping state (thread-local cache)
0437   Vector3 position(const State& state) const {
0438     return Reducer::position(state);
0439   }
0440 
0441   /// Momentum direction accessor
0442   ///
0443   /// @param state [in] The stepping state (thread-local cache)
0444   Vector3 direction(const State& state) const {
0445     return Reducer::direction(state);
0446   }
0447 
0448   /// QoP access
0449   ///
0450   /// @param state [in] The stepping state (thread-local cache)
0451   double qOverP(const State& state) const { return Reducer::qOverP(state); }
0452 
0453   /// Absolute momentum accessor
0454   ///
0455   /// @param state [in] The stepping state (thread-local cache)
0456   double absoluteMomentum(const State& state) const {
0457     return Reducer::absoluteMomentum(state);
0458   }
0459 
0460   /// Momentum accessor
0461   ///
0462   /// @param state [in] The stepping state (thread-local cache)
0463   Vector3 momentum(const State& state) const {
0464     return Reducer::momentum(state);
0465   }
0466 
0467   /// Charge access
0468   ///
0469   /// @param state [in] The stepping state (thread-local cache)
0470   double charge(const State& state) const { return Reducer::charge(state); }
0471 
0472   /// Particle hypothesis
0473   ///
0474   /// @param state [in] The stepping state (thread-local cache)
0475   ParticleHypothesis particleHypothesis(const State& state) const {
0476     return state.particleHypothesis;
0477   }
0478 
0479   /// Time access
0480   ///
0481   /// @param state [in] The stepping state (thread-local cache)
0482   double time(const State& state) const { return Reducer::time(state); }
0483 
0484   /// Update surface status
0485   ///
0486   /// It checks the status to the reference surface & updates
0487   /// the step size accordingly
0488   ///
0489   /// @param [in,out] state The stepping state (thread-local cache)
0490   /// @param [in] surface The surface provided
0491   /// @param [in] index The surface intersection index
0492   /// @param [in] navDir The navigation direction
0493   /// @param [in] boundaryTolerance The boundary check for this status update
0494   /// @param [in] surfaceTolerance Surface tolerance used for intersection
0495   /// @param [in] stype The step size type to be set
0496   /// @param [in] logger A @c Logger instance
0497   IntersectionStatus updateSurfaceStatus(
0498       State& state, const Surface& surface, std::uint8_t index,
0499       Direction navDir, const BoundaryTolerance& boundaryTolerance,
0500       double surfaceTolerance, ConstrainedStep::Type stype,
0501       const Logger& logger = getDummyLogger()) const {
0502     using Status = IntersectionStatus;
0503 
0504     std::array<int, 3> counts = {0, 0, 0};
0505 
0506     for (auto& component : state.components) {
0507       component.status = detail::updateSingleSurfaceStatus<SingleStepper>(
0508           *this, component.state, surface, index, navDir, boundaryTolerance,
0509           surfaceTolerance, stype, logger);
0510       ++counts[static_cast<std::size_t>(component.status)];
0511     }
0512 
0513     // If at least one component is on a surface, we can remove all missed
0514     // components before the step. If not, we must keep them for the case that
0515     // all components miss and we need to retarget
0516     if (counts[static_cast<std::size_t>(Status::onSurface)] > 0) {
0517       removeMissedComponents(state);
0518       reweightComponents(state);
0519     }
0520 
0521     ACTS_VERBOSE("Component status wrt "
0522                  << surface.geometryId() << " at {"
0523                  << surface.center(state.options.geoContext).transpose()
0524                  << "}:\t" << [&]() {
0525                       std::stringstream ss;
0526                       for (auto& component : state.components) {
0527                         ss << component.status << "\t";
0528                       }
0529                       return ss.str();
0530                     }());
0531 
0532     // Switch on stepCounter if one or more components reached a surface, but
0533     // some are still in progress of reaching the surface
0534     if (!state.stepCounterAfterFirstComponentOnSurface &&
0535         counts[static_cast<std::size_t>(Status::onSurface)] > 0 &&
0536         counts[static_cast<std::size_t>(Status::reachable)] > 0) {
0537       state.stepCounterAfterFirstComponentOnSurface = 0;
0538       ACTS_VERBOSE("started stepCounterAfterFirstComponentOnSurface");
0539     }
0540 
0541     // If there are no components onSurface, but the counter is switched on
0542     // (e.g., if the navigator changes the target surface), we need to switch it
0543     // off again
0544     if (state.stepCounterAfterFirstComponentOnSurface &&
0545         counts[static_cast<std::size_t>(Status::onSurface)] == 0) {
0546       state.stepCounterAfterFirstComponentOnSurface.reset();
0547       ACTS_VERBOSE("switch off stepCounterAfterFirstComponentOnSurface");
0548     }
0549 
0550     // This is a 'any_of' criterium. As long as any of the components has a
0551     // certain state, this determines the total state (in the order of a
0552     // somewhat importance)
0553     if (counts[static_cast<std::size_t>(Status::reachable)] > 0) {
0554       return Status::reachable;
0555     } else if (counts[static_cast<std::size_t>(Status::onSurface)] > 0) {
0556       state.stepCounterAfterFirstComponentOnSurface.reset();
0557       return Status::onSurface;
0558     } else {
0559       return Status::unreachable;
0560     }
0561   }
0562 
0563   /// Update step size
0564   ///
0565   /// This method intersects the provided surface and update the navigation
0566   /// step estimation accordingly (hence it changes the state). It also
0567   /// returns the status of the intersection to trigger onSurface in case
0568   /// the surface is reached.
0569   ///
0570   /// @param state [in,out] The stepping state (thread-local cache)
0571   /// @param oIntersection [in] The ObjectIntersection to layer, boundary, etc
0572   /// @param direction [in] The propagation direction
0573   /// @param stype [in] The step size type to be set
0574   template <typename object_intersection_t>
0575   void updateStepSize(State& state, const object_intersection_t& oIntersection,
0576                       Direction direction, ConstrainedStep::Type stype) const {
0577     const Surface& surface = *oIntersection.object();
0578 
0579     for (auto& component : state.components) {
0580       auto intersection = surface.intersect(
0581           component.state.options.geoContext,
0582           SingleStepper::position(component.state),
0583           direction * SingleStepper::direction(component.state),
0584           BoundaryTolerance::None())[oIntersection.index()];
0585 
0586       SingleStepper::updateStepSize(component.state, intersection, direction,
0587                                     stype);
0588     }
0589   }
0590 
0591   /// Update step size - explicitly with a double
0592   ///
0593   /// @param state [in,out] The stepping state (thread-local cache)
0594   /// @param stepSize [in] The step size value
0595   /// @param stype [in] The step size type to be set
0596   void updateStepSize(State& state, double stepSize,
0597                       ConstrainedStep::Type stype) const {
0598     for (auto& component : state.components) {
0599       SingleStepper::updateStepSize(component.state, stepSize, stype);
0600     }
0601   }
0602 
0603   /// Get the step size
0604   ///
0605   /// @param state [in] The stepping state (thread-local cache)
0606   /// @param stype [in] The step size type to be returned
0607   /// @note This returns the smallest step size of all components. It uses
0608   /// std::abs for comparison to handle backward propagation and negative
0609   /// step sizes correctly.
0610   double getStepSize(const State& state, ConstrainedStep::Type stype) const {
0611     return std::min_element(state.components.begin(), state.components.end(),
0612                             [=](const auto& a, const auto& b) {
0613                               return std::abs(a.state.stepSize.value(stype)) <
0614                                      std::abs(b.state.stepSize.value(stype));
0615                             })
0616         ->state.stepSize.value(stype);
0617   }
0618 
0619   /// Release the step-size for all components
0620   ///
0621   /// @param state [in,out] The stepping state (thread-local cache)
0622   /// @param [in] stype The step size type to be released
0623   void releaseStepSize(State& state, ConstrainedStep::Type stype) const {
0624     for (auto& component : state.components) {
0625       SingleStepper::releaseStepSize(component.state, stype);
0626     }
0627   }
0628 
0629   /// Output the Step Size of all components into one std::string
0630   ///
0631   /// @param state [in,out] The stepping state (thread-local cache)
0632   std::string outputStepSize(const State& state) const {
0633     std::stringstream ss;
0634     for (const auto& component : state.components) {
0635       ss << component.state.stepSize.toString() << " || ";
0636     }
0637 
0638     return ss.str();
0639   }
0640 
0641   /// Create and return the bound state at the current position
0642   ///
0643   /// @brief This transports (if necessary) the covariance
0644   /// to the surface and creates a bound state. It does not check
0645   /// if the transported state is at the surface, this needs to
0646   /// be guaranteed by the propagator.
0647   /// @note This is done by combining the gaussian mixture on the specified
0648   /// surface. If the conversion to bound states of some components
0649   /// fails, these components are ignored unless all components fail. In this
0650   /// case an error code is returned.
0651   ///
0652   /// @param [in] state State that will be presented as @c BoundState
0653   /// @param [in] surface The surface to which we bind the state
0654   /// @param [in] transportCov Flag steering covariance transport
0655   /// @param [in] freeToBoundCorrection Flag steering non-linear correction during global to local correction
0656   ///
0657   /// @return A bound state:
0658   ///   - the parameters at the surface
0659   ///   - the stepwise jacobian towards it (from last bound)
0660   ///   - and the path length (from start - for ordering)
0661   Result<BoundState> boundState(
0662       State& state, const Surface& surface, bool transportCov = true,
0663       const FreeToBoundCorrection& freeToBoundCorrection =
0664           FreeToBoundCorrection(false)) const;
0665 
0666   /// @brief If necessary fill additional members needed for curvilinearState
0667   ///
0668   /// Compute path length derivatives in case they have not been computed
0669   /// yet, which is the case if no step has been executed yet.
0670   ///
0671   /// @param [in, out] state The stepping state (thread-local cache)
0672   /// @return true if nothing is missing after this call, false otherwise.
0673   bool prepareCurvilinearState(State& state) const {
0674     (void)state;
0675     return true;
0676   }
0677 
0678   /// Create and return a curvilinear state at the current position
0679   ///
0680   /// @brief This transports (if necessary) the covariance
0681   /// to the current position and creates a curvilinear state.
0682   /// @note This is done as a simple average over the free representation
0683   /// and covariance of the components.
0684   ///
0685   /// @param [in] state State that will be presented as @c CurvilinearState
0686   /// @param [in] transportCov Flag steering covariance transport
0687   ///
0688   /// @return A curvilinear state:
0689   ///   - the curvilinear parameters at given position
0690   ///   - the stepweise jacobian towards it (from last bound)
0691   ///   - and the path length (from start - for ordering)
0692   BoundState curvilinearState(State& state, bool transportCov = true) const;
0693 
0694   /// Method for on-demand transport of the covariance
0695   /// to a new curvilinear frame at current  position,
0696   /// or direction of the state
0697   ///
0698   /// @param [in,out] state State of the stepper
0699   void transportCovarianceToCurvilinear(State& state) const {
0700     for (auto& component : state.components) {
0701       SingleStepper::transportCovarianceToCurvilinear(component.state);
0702     }
0703   }
0704 
0705   /// Method for on-demand transport of the covariance
0706   /// to a new curvilinear frame at current position,
0707   /// or direction of the state
0708   ///
0709   /// @tparam surface_t the Surface type
0710   ///
0711   /// @param [in,out] state State of the stepper
0712   /// @param [in] surface is the surface to which the covariance is forwarded
0713   /// @param [in] freeToBoundCorrection Flag steering non-linear correction during global to local correction
0714   /// to
0715   /// @note no check is done if the position is actually on the surface
0716   void transportCovarianceToBound(
0717       State& state, const Surface& surface,
0718       const FreeToBoundCorrection& freeToBoundCorrection =
0719           FreeToBoundCorrection(false)) const {
0720     for (auto& component : state.components) {
0721       SingleStepper::transportCovarianceToBound(component.state, surface,
0722                                                 freeToBoundCorrection);
0723     }
0724   }
0725 
0726   /// Perform a Runge-Kutta track parameter propagation step
0727   ///
0728   /// @param [in,out] state The state of the stepper
0729   /// @param propDir is the direction of propagation
0730   /// @param material is the material properties
0731   /// @return the result of the step
0732   ///
0733   /// The state contains the desired step size. It can be negative during
0734   /// backwards track propagation, and since we're using an adaptive
0735   /// algorithm, it can be modified by the stepper class during propagation.
0736   Result<double> step(State& state, Direction propDir,
0737                       const IVolumeMaterial* material) const;
0738 };
0739 
0740 }  // namespace Acts
0741 
0742 #include "MultiStepperLoop.ipp"