ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
StraightLineStepper.hpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file StraightLineStepper.hpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2016-2019 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #pragma once
10 
11 #include <cmath>
12 #include <functional>
13 
24 
25 namespace Acts {
26 
33  public:
36  using BoundState = std::tuple<BoundParameters, Jacobian, double>;
37  using CurvilinearState = std::tuple<CurvilinearParameters, Jacobian, double>;
38  using BField = NullBField;
39 
42  struct State {
44  State() = delete;
45 
56  template <typename parameters_t>
57  explicit State(std::reference_wrapper<const GeometryContext> gctx,
58  std::reference_wrapper<const MagneticFieldContext> /*mctx*/,
59  const parameters_t& par, NavigationDirection ndir = forward,
60  double ssize = std::numeric_limits<double>::max(),
61  double stolerance = s_onSurfaceTolerance)
62  : pos(par.position()),
63  dir(par.momentum().normalized()),
64  p(par.momentum().norm()),
65  q(par.charge()),
66  t(par.time()),
67  navDir(ndir),
68  stepSize(ndir * std::abs(ssize)),
69  tolerance(stolerance),
70  geoContext(gctx) {
71  if (par.covariance()) {
72  // Get the reference surface for navigation
73  const auto& surface = par.referenceSurface();
74  // set the covariance transport flag to true and copy
75  covTransport = true;
76  cov = BoundSymMatrix(*par.covariance());
77  surface.initJacobianToGlobal(gctx, jacToGlobal, pos, dir,
78  par.parameters());
79  }
80  }
81 
83  BoundToFreeMatrix jacToGlobal = BoundToFreeMatrix::Zero();
84 
86  FreeMatrix jacTransport = FreeMatrix::Identity();
87 
89  Jacobian jacobian = Jacobian::Identity();
90 
92  FreeVector derivative = FreeVector::Zero();
93 
95  bool covTransport = false;
96  Covariance cov = Covariance::Zero();
97 
99  Vector3D pos = Vector3D(0., 0., 0.);
100 
102  Vector3D dir = Vector3D(1., 0., 0.);
103 
105  double p = 0.;
106 
108  double q = 0.;
109 
111  double t = 0.;
112 
115 
117  double pathAccumulated = 0.;
118 
121 
122  // Previous step size for overstep estimation (ignored for SL stepper)
123  double previousStepSize = 0.;
124 
127 
128  // Cache the geometry context of this propagation
129  std::reference_wrapper<const GeometryContext> geoContext;
130  };
131 
134  using state_type = State;
135 
137  StraightLineStepper() = default;
138 
144  Vector3D getField(State& /*state*/, const Vector3D& /*pos*/) const {
145  // get the field from the cell
146  return Vector3D(0., 0., 0.);
147  }
148 
152  Vector3D position(const State& state) const { return state.pos; }
153 
157  Vector3D direction(const State& state) const { return state.dir; }
158 
162  double momentum(const State& state) const { return state.p; }
163 
167  double charge(const State& state) const { return state.q; }
168 
172  double time(const State& state) const { return state.t; }
173 
177  double overstepLimit(const State& /*state*/) const {
178  return s_onSurfaceTolerance;
179  }
180 
192  const BoundaryCheck& bcheck) const {
193  return detail::updateSingleSurfaceStatus<StraightLineStepper>(
194  *this, state, surface, bcheck);
195  }
196 
205  template <typename object_intersection_t>
206  void updateStepSize(State& state, const object_intersection_t& oIntersection,
207  bool release = true) const {
208  detail::updateSingleStepSize<StraightLineStepper>(state, oIntersection,
209  release);
210  }
211 
217  void setStepSize(State& state, double stepSize,
219  state.previousStepSize = state.stepSize;
220  state.stepSize.update(stepSize, stype, true);
221  }
222 
226  void releaseStepSize(State& state) const {
228  }
229 
233  std::string outputStepSize(const State& state) const {
234  return state.stepSize.toString();
235  }
236 
249  BoundState boundState(State& state, const Surface& surface) const;
250 
262 
267  void update(State& state, const BoundParameters& pars) const;
268 
276  void update(State& state, const Vector3D& uposition,
277  const Vector3D& udirection, double up, double time) const;
278 
284  void covarianceTransport(State& state) const;
285 
297  void covarianceTransport(State& state, const Surface& surface) const;
298 
309  template <typename propagator_state_t>
310  Result<double> step(propagator_state_t& state) const {
311  // use the adjusted step size
312  const auto h = state.stepping.stepSize;
313  // time propagates along distance as 1/b = sqrt(1 + m²/p²)
314  const auto dtds = std::hypot(1., state.options.mass / state.stepping.p);
315  // Update the track parameters according to the equations of motion
316  state.stepping.pos += h * state.stepping.dir;
317  state.stepping.t += h * dtds;
318  // Propagate the jacobian
319  if (state.stepping.covTransport) {
320  // The step transport matrix in global coordinates
321  FreeMatrix D = FreeMatrix::Identity();
322  D.block<3, 3>(0, 4) = ActsSymMatrixD<3>::Identity() * h;
323  // Extend the calculation by the time propagation
324  // Evaluate dt/dlambda
325  D(3, 7) = h * state.options.mass * state.options.mass *
326  (state.stepping.q == 0. ? 1. : state.stepping.q) /
327  (state.stepping.p * dtds);
328  // Set the derivative factor the time
329  state.stepping.derivative(3) = dtds;
330  // Update jacobian and derivative
331  state.stepping.jacTransport = D * state.stepping.jacTransport;
332  state.stepping.derivative.template head<3>() = state.stepping.dir;
333  }
334  // state the path length
335  state.stepping.pathAccumulated += h;
336 
337  // return h
338  return h;
339  }
340 };
341 
342 } // namespace Acts