Acts/Propagator/RiddersPropagator.hpp

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/ActorList.hpp"
0012 #include "Acts/Propagator/PropagatorOptions.hpp"
0013 #include "Acts/Propagator/VoidNavigator.hpp"
0014 #include "Acts/Surfaces/Surface.hpp"
0015
0016 namespace Acts {
0017
0018 /// Options for Ridders-based covariance propagation.
0019 template <typename propagator_t, typename actor_list_t = ActorList<>>
0020 struct RiddersPropagatorOptions
0021     : public propagator_t::template Options<actor_list_t> {
0022   /// Base type
0023   using base_type = propagator_t::template Options<actor_list_t>;
0024
0025   /// Stepper options type
0026   using stepper_options_type = typename base_type::stepper_options_type;
0027   /// Navigator options type
0028   using navigator_options_type = typename base_type::navigator_options_type;
0029   /// Actor list type
0030   using actor_list_type = actor_list_t;
0031
0032   /// PropagatorOptions with context
0033   /// @param gctx Geometry context
0034   /// @param mctx Magnetic field context
0035   RiddersPropagatorOptions(const GeometryContext& gctx,
0036                            const MagneticFieldContext& mctx)
0037       : base_type(gctx, mctx) {}
0038
0039   /// PropagatorOptions with context and plain options
0040   /// @param pOptions Plain propagator options
0041   explicit RiddersPropagatorOptions(const PropagatorPlainOptions& pOptions)
0042       : base_type(pOptions) {}
0043
0044   using base_type::operator PropagatorPlainOptions;
0045
0046   /// @brief Expand the options with extended actors
0047   ///
0048   /// @tparam extended_actor_list_t Type of the new actor list
0049   ///
0050   /// @param extendedActorList The new actor list to be used (internally)
0051   /// @return New options with extended actor list
0052   template <typename extended_actor_list_t>
0053   RiddersPropagatorOptions<propagator_t, extended_actor_list_t> extend(
0054       extended_actor_list_t extendedActorList) const {
0055     RiddersPropagatorOptions<propagator_t, extended_actor_list_t> eoptions(
0056         base_type::geoContext, base_type::magFieldContext);
0057
0058     static_cast<decltype(eoptions)::base_type&>(eoptions) =
0059         base_type::extend(std::move(extendedActorList));
0060
0061     return eoptions;
0062   }
0063
0064   using base_type::setPlainOptions;
0065
0066   /// Initial scale for the deviation of the individual bound track parameters
0067   BoundVector deviationScale = {1e-4, 1e-4, 1e-4, 1e-4, 1e-4, 1e-4};
0068
0069   /// Different factors applied to the initial scale to create the
0070   /// deviations of the individual bound track parameters. The resulting
0071   /// function value deviations are then fitted to a line to determine the
0072   /// first order derivatives of the final parameters wrt. the initial
0073   /// parameters.
0074   std::vector<double> deviationFactors = {-2, -1, 1, 2};
0075 };
0076
0077 /// @brief This class performs the Ridders algorithm to estimate the propagation
0078 /// of the covariance to a certain point in space.
0079 ///
0080 /// The algorithm is based on the small deviations of the start parameters based
0081 /// on their uncertainty at the beginning of the propagation. This deviation is
0082 /// represented here by a vector of relative deviations of these parameters and
0083 /// fix for all parameters. So, a common choice has to be found that is able to
0084 /// actually fit into the order of magnitude of the uncertainty of each
0085 /// parameter. Using these deviations, the propagation is repeated multiple
0086 /// times and the final covariance matrix at a given target surface is
0087 /// afterwards evaluated by first order derivatives of the final state
0088 /// parameters wrt. the initial parameters. Therefore this evaluation represents
0089 /// a first order approximation of the transport jacobian. Since performing
0090 /// multiple propagations and a numerical evaluation of the covariance requires
0091 /// more time than a single propagation towards a target + a common propagation
0092 /// of the covariance, this class just serves to verify the results of the
0093 /// latter classes.
0094 template <typename propagator_t>
0095 class RiddersPropagator {
0096  public:
0097   /// Type of the propagator used internally
0098   using Propagator = propagator_t;
0099
0100   /// Type of the stepper in use
0101   using Stepper = Propagator::Stepper;
0102
0103   /// Type of the navigator in use
0104   using Navigator = Propagator::Navigator;
0105
0106   /// Type of the stepper options
0107   using StepperOptions = Propagator::StepperOptions;
0108
0109   /// Type of the navigator options
0110   using NavigatorOptions = Propagator::NavigatorOptions;
0111
0112   /// Type of the propagator options with actor list
0113   template <typename actor_list_t = ActorList<>>
0114   using Options = RiddersPropagatorOptions<propagator_t, actor_list_t>;
0115
0116   /// Type of the stepper state
0117   using StepperState = Propagator::StepperState;
0118
0119   /// Type of the navigator state
0120   using NavigatorState = Propagator::NavigatorState;
0121
0122   /// Type of the propagation state derived from the propagator options
0123   /// @tparam propagator_options_t Type of the propagator options
0124   template <typename propagator_options_t>
0125   using State = Propagator::template State<propagator_options_t>;
0126
0127   /// Type of the propagation result derived from the propagator options
0128   /// @tparam propagator_options_t Type of the propagator options
0129   template <typename propagator_options_t>
0130   using ResultType =
0131       typename Propagator::template ResultType<propagator_options_t>;
0132
0133   /// @brief Constructor using a propagator
0134   ///
0135   /// @param [in] propagator The propagator to use
0136   explicit RiddersPropagator(propagator_t propagator)
0137       : m_propagator(std::move(propagator)) {}
0138
0139   /// @brief Propagation method targeting curvilinear parameters
0140   ///
0141   /// @tparam parameters_t Type of the start parameters
0142   /// @tparam propagator_options_t Type of the propagator options
0143   ///
0144   /// @param [in] start Start parameters
0145   /// @param [in] options Options of the propagations
0146   ///
0147   /// @return Result of the propagation
0148   template <typename parameters_t, typename propagator_options_t>
0149   Result<ResultType<propagator_options_t>> propagate(
0150       const parameters_t& start, const propagator_options_t& options) const;
0151
0152   /// @brief Propagation method targeting bound parameters
0153   ///
0154   /// @tparam parameters_t Type of the start parameters
0155   /// @tparam propagator_options_t Type of the propagator options
0156   ///
0157   /// @param [in] start Start parameters
0158   /// @param [in] target The target surface
0159   /// @param [in] options Options of the propagations
0160   ///
0161   /// @return Result of the propagation
0162   /// @note If the target surface is a disc, the resulting covariance may be
0163   /// inconsistent. In this case a zero matrix is returned.
0164   template <typename parameters_t, typename propagator_options_t>
0165   Result<ResultType<propagator_options_t>> propagate(
0166       const parameters_t& start, const Surface& target,
0167       const propagator_options_t& options) const;
0168
0169  private:
0170   /// Does the actual ridders propagation by wiggling the parameters and
0171   /// propagating again. This function is called from the different
0172   /// propagation overloads in order to deduplicate code.
0173   ///
0174   /// @param [in] options Options of the propagations
0175   /// @param [in] start Start parameters
0176   /// @param [in] nominalResult The result of the nominal propagation
0177   template <typename parameters_t, typename propagator_options_t>
0178   BoundMatrix wiggleAndCalculateJacobian(
0179       const parameters_t& start, const propagator_options_t& options,
0180       const ResultType<propagator_options_t>& nominalResult) const;
0181
0182   /// @brief This function wiggles one dimension of the starting parameters,
0183   /// performs the propagation to a surface and collects for each change of the
0184   /// start parameters the slope
0185   ///
0186   /// @tparam options_t PropagatorOptions object
0187   /// @tparam parameters_t Type of the parameters to start the propagation with
0188   ///
0189   /// @param [in] options Options do define how to wiggle
0190   /// @param [in] start Start parameters which will be modified
0191   /// @param [in] param Index to get the parameter that will be modified
0192   /// @param [in] target Target surface
0193   /// @param [in] nominal Nominal end parameters
0194   /// @param [in] deviations Vector of deviations
0195   ///
0196   /// @return Vector containing each slope
0197   template <typename propagator_options_t, typename parameters_t>
0198   std::vector<BoundVector> wiggleParameter(
0199       const propagator_options_t& options, const parameters_t& start,
0200       unsigned int param, const Surface& target, const BoundVector& nominal,
0201       const std::vector<double>& deviations) const;
0202
0203   propagator_t m_propagator;
0204 };
0205
0206 }  // namespace Acts
0207
0208 #include "Acts/Propagator/RiddersPropagator.ipp"