Propagator.hpp

Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file Propagator.hpp

// This file is part of the Acts project.
//
// Copyright (C) 2016-2023 CERN for the benefit of the Acts project
//
// This Source Code Form is subject to the terms of the Mozilla Public
// License, v. 2.0. If a copy of the MPL was not distributed with this
// file, You can obtain one at http://mozilla.org/MPL/2.0/.

#pragma once

// clang-format off
// Workaround for building on clang+libstdc++. Must be the first include.
#include "Acts/Utilities/detail/ReferenceWrapperAnyCompat.hpp"
// clang-format on

#include "Acts/Definitions/Algebra.hpp"
#include "Acts/Definitions/PdgParticle.hpp"
#include "Acts/Definitions/Units.hpp"
#include "Acts/EventData/TrackParametersConcept.hpp"
#include "Acts/Geometry/GeometryContext.hpp"
#include "Acts/MagneticField/MagneticFieldContext.hpp"
#include "Acts/Propagator/AbortList.hpp"
#include "Acts/Propagator/ActionList.hpp"
#include "Acts/Propagator/StandardAborters.hpp"
#include "Acts/Propagator/StepperConcept.hpp"
#include "Acts/Propagator/detail/ParameterTraits.hpp"
#include "Acts/Propagator/detail/VoidPropagatorComponents.hpp"
#include "Acts/Utilities/Logger.hpp"
#include "Acts/Utilities/Result.hpp"

#include <optional>

namespace Acts {

template <typename parameters_t, typename... result_list>
struct PropagatorResult : private detail::Extendable<result_list...> {
  using detail::Extendable<result_list...>::get;

  std::optional<parameters_t> endParameters = std::nullopt;

  std::optional<BoundMatrix> transportJacobian = std::nullopt;

  unsigned int steps = 0;

  double pathLength = 0.;
};

struct PropagatorPlainOptions {
  Direction direction = Direction::Forward;

  unsigned int maxSteps = 1000;

  unsigned int maxRungeKuttaStepTrials = 10000;

  double maxStepSize = std::numeric_limits<double>::max();

  double pathLimit = std::numeric_limits<double>::max();

  double targetTolerance = s_onSurfaceTolerance;

  bool loopProtection = true;
  double loopFraction = 0.5;  

  // Configurations for Stepper
  double tolerance = 1e-4;

  double stepSizeCutOff = 0.;
};

template <typename action_list_t = ActionList<>,
          typename aborter_list_t = AbortList<>>
struct PropagatorOptions : public PropagatorPlainOptions {
  using action_list_type = action_list_t;
  using aborter_list_type = aborter_list_t;

  PropagatorOptions() = delete;

  PropagatorOptions(
      const PropagatorOptions<action_list_t, aborter_list_t>& po) = default;

  PropagatorOptions(const GeometryContext& gctx,
                    const MagneticFieldContext& mctx)
      : geoContext(gctx), magFieldContext(mctx) {}

  template <typename extended_aborter_list_t>
  PropagatorOptions<action_list_t, extended_aborter_list_t> extend(
      extended_aborter_list_t aborters) const {
    PropagatorOptions<action_list_t, extended_aborter_list_t> eoptions(
        geoContext, magFieldContext);
    // Copy the options over
    eoptions.direction = direction;
    eoptions.maxSteps = maxSteps;
    eoptions.maxRungeKuttaStepTrials = maxRungeKuttaStepTrials;
    eoptions.maxStepSize = maxStepSize;
    eoptions.targetTolerance = targetTolerance;
    eoptions.pathLimit = direction * std::abs(pathLimit);
    eoptions.loopProtection = loopProtection;
    eoptions.loopFraction = loopFraction;

    // Stepper options
    eoptions.tolerance = tolerance;
    eoptions.stepSizeCutOff = stepSizeCutOff;
    // Action / abort list
    eoptions.actionList = std::move(actionList);
    eoptions.abortList = std::move(aborters);
    // And return the options
    return eoptions;
  }

  void setPlainOptions(const PropagatorPlainOptions& pOptions) {
    // Copy the options over
    direction = pOptions.direction;
    maxSteps = pOptions.maxSteps;
    maxRungeKuttaStepTrials = pOptions.maxRungeKuttaStepTrials;
    maxStepSize = pOptions.maxStepSize;
    targetTolerance = pOptions.targetTolerance;
    pathLimit = direction * std::abs(pOptions.pathLimit);
    loopProtection = pOptions.loopProtection;
    loopFraction = pOptions.loopFraction;
    tolerance = pOptions.tolerance;
    stepSizeCutOff = pOptions.stepSizeCutOff;
  }

  action_list_t actionList;

  aborter_list_t abortList;

  std::reference_wrapper<const GeometryContext> geoContext;

  std::reference_wrapper<const MagneticFieldContext> magFieldContext;
};

template <typename stepper_t, typename navigator_t = detail::VoidNavigator>
class Propagator final {
  using StepperBoundTrackParameters =
      detail::stepper_bound_parameters_type_t<stepper_t>;
  static_assert(
      Concepts::BoundTrackParametersConcept<StepperBoundTrackParameters>,
      "Stepper bound track parameters do not fulfill bound "
      "parameters concept.");

  using StepperCurvilinearTrackParameters =
      detail::stepper_curvilinear_parameters_type_t<stepper_t>;
  static_assert(
      Concepts::BoundTrackParametersConcept<StepperCurvilinearTrackParameters>,
      "Stepper bound track parameters do not fulfill bound "
      "parameters concept.");

  using Jacobian = BoundMatrix;
  using BoundState = std::tuple<StepperBoundTrackParameters, Jacobian, double>;
  using CurvilinearState =
      std::tuple<StepperCurvilinearTrackParameters, Jacobian, double>;

  static_assert(StepperStateConcept<typename stepper_t::State>,
                "Stepper does not fulfill stepper concept.");
  static_assert(StepperConcept<stepper_t>,
                "Stepper does not fulfill stepper concept.");

 public:
  using Stepper = stepper_t;

  using Navigator = navigator_t;

  using StepperState = typename Stepper::State;

  using NavigatorState = typename navigator_t::State;

  explicit Propagator(stepper_t stepper, navigator_t navigator = navigator_t(),
                      std::shared_ptr<const Logger> _logger =
                          getDefaultLogger("Propagator", Acts::Logging::INFO))
      : m_stepper(std::move(stepper)),
        m_navigator(std::move(navigator)),
        m_logger{std::move(_logger)} {}

  template <typename propagator_options_t>
  struct State {
    State(const propagator_options_t& topts, StepperState steppingIn,
          NavigatorState navigationIn)
        : options(topts),
          stepping{std::move(steppingIn)},
          navigation{std::move(navigationIn)},
          geoContext(topts.geoContext) {}

    propagator_options_t options;

    StepperState stepping;

    NavigatorState navigation;

    std::reference_wrapper<const GeometryContext> geoContext;
  };

 private:
  template <typename parameters_t, typename action_list_t>
  struct result_type_helper {
    template <typename... args>
    using this_result_type = PropagatorResult<parameters_t, args...>;

    using type = typename action_list_t::template result_type<this_result_type>;
  };

 public:
  template <typename parameters_t, typename action_list_t>
  using action_list_t_result_t =
      typename result_type_helper<parameters_t, action_list_t>::type;

 private:
  template <typename result_t, typename propagator_state_t>
  Result<void> propagate_impl(propagator_state_t& state,
                              result_t& result) const;

 public:
  template <typename parameters_t, typename propagator_options_t,
            typename path_aborter_t = PathLimitReached>
  Result<
      action_list_t_result_t<StepperCurvilinearTrackParameters,
                             typename propagator_options_t::action_list_type>>
  propagate(const parameters_t& start, const propagator_options_t& options,
            bool makeCurvilinear = true) const;

  template <typename parameters_t, typename propagator_options_t,
            typename path_aborter_t = PathLimitReached>
  Result<
      action_list_t_result_t<StepperCurvilinearTrackParameters,
                             typename propagator_options_t::action_list_type>>
  propagate(
      const parameters_t& start, const propagator_options_t& options,
      bool makeCurvilinear,
      action_list_t_result_t<StepperCurvilinearTrackParameters,
                             typename propagator_options_t::action_list_type>&&
          inputResult) const;

  template <typename parameters_t, typename propagator_options_t,
            typename target_aborter_t = SurfaceReached,
            typename path_aborter_t = PathLimitReached>
  Result<
      action_list_t_result_t<StepperBoundTrackParameters,
                             typename propagator_options_t::action_list_type>>
  propagate(const parameters_t& start, const Surface& target,
            const propagator_options_t& options) const;

  template <typename parameters_t, typename propagator_options_t,
            typename target_aborter_t = SurfaceReached,
            typename path_aborter_t = PathLimitReached>
  Result<
      action_list_t_result_t<StepperBoundTrackParameters,
                             typename propagator_options_t::action_list_type>>
  propagate(
      const parameters_t& start, const Surface& target,
      const propagator_options_t& options,
      action_list_t_result_t<StepperBoundTrackParameters,
                             typename propagator_options_t::action_list_type>
          inputResult) const;

 private:
  const Logger& logger() const { return *m_logger; }

  stepper_t m_stepper;

  navigator_t m_navigator;

  std::shared_ptr<const Logger> m_logger;
};

}  // namespace Acts

#include "Acts/Propagator/Propagator.ipp"