Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
StraightLineStepperTests.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file StraightLineStepperTests.cpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2020 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 #include <boost/test/unit_test.hpp>
10 
11 #include "Acts/Definitions/Algebra.hpp"
12 #include "Acts/Definitions/Direction.hpp"
13 #include "Acts/Definitions/TrackParametrization.hpp"
14 #include "Acts/EventData/Charge.hpp"
15 #include "Acts/EventData/GenericBoundTrackParameters.hpp"
16 #include "Acts/EventData/GenericCurvilinearTrackParameters.hpp"
17 #include "Acts/EventData/TrackParameters.hpp"
18 #include "Acts/EventData/detail/CorrectedTransformationFreeToBound.hpp"
19 #include "Acts/EventData/detail/TransformationBoundToFree.hpp"
20 #include "Acts/Geometry/GeometryContext.hpp"
21 #include "Acts/MagneticField/MagneticFieldContext.hpp"
22 #include "Acts/Propagator/ConstrainedStep.hpp"
23 #include "Acts/Propagator/StraightLineStepper.hpp"
24 #include "Acts/Surfaces/BoundaryCheck.hpp"
25 #include "Acts/Surfaces/PlaneSurface.hpp"
26 #include "Acts/Surfaces/Surface.hpp"
27 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
28 #include "Acts/Utilities/Helpers.hpp"
29 #include "Acts/Utilities/Result.hpp"
30 
31 #include <functional>
32 #include <limits>
33 #include <memory>
34 #include <optional>
35 #include <string>
36 #include <tuple>
37 #include <utility>
38 
39 namespace tt = boost::test_tools;
40 using Acts::VectorHelpers::makeVector4;
41 
42 namespace Acts {
43 namespace Test {
44 
45 using Covariance = BoundSquareMatrix;
46 using Jacobian = BoundMatrix;
47 
49 struct PropState {
51  PropState(Direction direction, StraightLineStepper::State sState)
52  : stepping(std::move(sState)) {
53  options.direction = direction;
54  }
56  StraightLineStepper::State stepping;
58  struct {
59  Direction direction = Direction::Forward;
60  } options;
61 };
62 
63 struct MockNavigator {};
64 
65 static constexpr MockNavigator mockNavigator;
66 
67 static constexpr auto eps = 2 * std::numeric_limits<double>::epsilon();
68 
70 BOOST_AUTO_TEST_CASE(straight_line_stepper_state_test) {
71  // Set up some variables
72  GeometryContext tgContext = GeometryContext();
73  MagneticFieldContext mfContext = MagneticFieldContext();
74  double stepSize = 123.;
75  double tolerance = 234.;
76 
77  Vector3 pos(1., 2., 3.);
78  Vector3 dir(4., 5., 6.);
79  double time = 7.;
80  double absMom = 8.;
81  double charge = -1.;
82 
83  // Test charged parameters without covariance matrix
84  CurvilinearTrackParameters cp(makeVector4(pos, time), dir, charge / absMom,
85  std::nullopt, ParticleHypothesis::pion());
86  StraightLineStepper::State slsState(tgContext, mfContext, cp, stepSize,
87  tolerance);
88 
89  StraightLineStepper sls;
90 
91  // Test the result & compare with the input/test for reasonable members
92  BOOST_CHECK_EQUAL(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
93  BOOST_CHECK_EQUAL(slsState.jacTransport, FreeMatrix::Identity());
94  BOOST_CHECK_EQUAL(slsState.derivative, FreeVector::Zero());
95  BOOST_CHECK(!slsState.covTransport);
96  BOOST_CHECK_EQUAL(slsState.cov, Covariance::Zero());
97  CHECK_CLOSE_OR_SMALL(sls.position(slsState), pos, eps, eps);
98  CHECK_CLOSE_OR_SMALL(sls.direction(slsState), dir.normalized(), eps, eps);
99  CHECK_CLOSE_REL(sls.absoluteMomentum(slsState), absMom, eps);
100  BOOST_CHECK_EQUAL(sls.charge(slsState), charge);
101  CHECK_CLOSE_OR_SMALL(sls.time(slsState), time, eps, eps);
102  BOOST_CHECK_EQUAL(slsState.pathAccumulated, 0.);
103  BOOST_CHECK_EQUAL(slsState.stepSize.value(), stepSize);
104  BOOST_CHECK_EQUAL(slsState.previousStepSize, 0.);
105  BOOST_CHECK_EQUAL(slsState.tolerance, tolerance);
106 
107  // Test without charge and covariance matrix
108  CurvilinearTrackParameters ncp(makeVector4(pos, time), dir, 1 / absMom,
109  std::nullopt, ParticleHypothesis::pion0());
110  slsState = StraightLineStepper::State(tgContext, mfContext, ncp, stepSize,
111  tolerance);
112 
113  // Test with covariance matrix
114  Covariance cov = 8. * Covariance::Identity();
115  ncp = CurvilinearTrackParameters(makeVector4(pos, time), dir, 1 / absMom, cov,
116  ParticleHypothesis::pion0());
117  slsState = StraightLineStepper::State(tgContext, mfContext, ncp, stepSize,
118  tolerance);
119  BOOST_CHECK_NE(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
120  BOOST_CHECK(slsState.covTransport);
121  BOOST_CHECK_EQUAL(slsState.cov, cov);
122 }
123 
126 BOOST_AUTO_TEST_CASE(straight_line_stepper_test) {
127  // Set up some variables for the state
128  GeometryContext tgContext = GeometryContext();
129  MagneticFieldContext mfContext = MagneticFieldContext();
130  Direction navDir = Direction::Backward;
131  double stepSize = 123.;
132  double tolerance = 234.;
133 
134  // Construct the parameters
135  Vector3 pos(1., 2., 3.);
136  Vector3 dir = Vector3(4., 5., 6.).normalized();
137  double time = 7.;
138  double absMom = 8.;
139  double charge = -1.;
140  Covariance cov = 8. * Covariance::Identity();
141  CurvilinearTrackParameters cp(makeVector4(pos, time), dir, charge / absMom,
142  cov, ParticleHypothesis::pion());
143 
144  // Build the state and the stepper
145  StraightLineStepper::State slsState(tgContext, mfContext, cp, stepSize,
146  tolerance);
147  StraightLineStepper sls;
148 
149  // Test the getters
150  CHECK_CLOSE_ABS(sls.position(slsState), pos, 1e-6);
151  CHECK_CLOSE_ABS(sls.direction(slsState), dir, 1e-6);
152  CHECK_CLOSE_ABS(sls.absoluteMomentum(slsState), absMom, 1e-6);
153  BOOST_CHECK_EQUAL(sls.charge(slsState), charge);
154  BOOST_CHECK_EQUAL(sls.time(slsState), time);
155 
156  //~ BOOST_CHECK_EQUAL(sls.overstepLimit(slsState), tolerance);
157 
158  // Step size modifies
159  const std::string originalStepSize = slsState.stepSize.toString();
160 
161  sls.setStepSize(slsState, -1337.);
162  BOOST_CHECK_EQUAL(slsState.previousStepSize, stepSize);
163  BOOST_CHECK_EQUAL(slsState.stepSize.value(), -1337.);
164 
165  sls.releaseStepSize(slsState);
166  BOOST_CHECK_EQUAL(slsState.stepSize.value(), stepSize);
167  BOOST_CHECK_EQUAL(sls.outputStepSize(slsState), originalStepSize);
168 
169  // Test the curvilinear state construction
170  auto curvState = sls.curvilinearState(slsState);
171  auto curvPars = std::get<0>(curvState);
172  CHECK_CLOSE_ABS(curvPars.position(tgContext), cp.position(tgContext), 1e-6);
173  CHECK_CLOSE_ABS(curvPars.absoluteMomentum(), cp.absoluteMomentum(), 1e-6);
174  CHECK_CLOSE_ABS(curvPars.charge(), cp.charge(), 1e-6);
175  CHECK_CLOSE_ABS(curvPars.time(), cp.time(), 1e-6);
176  BOOST_CHECK(curvPars.covariance().has_value());
177  BOOST_CHECK_NE(*curvPars.covariance(), cov);
178  CHECK_CLOSE_COVARIANCE(std::get<1>(curvState),
179  BoundMatrix(BoundMatrix::Identity()), 1e-6);
180  CHECK_CLOSE_ABS(std::get<2>(curvState), 0., 1e-6);
181 
182  // Test the update method
183  Vector3 newPos(2., 4., 8.);
184  Vector3 newMom(3., 9., 27.);
185  double newTime(321.);
186  sls.update(slsState, newPos, newMom.normalized(), charge / newMom.norm(),
187  newTime);
188  CHECK_CLOSE_ABS(sls.position(slsState), newPos, 1e-6);
189  CHECK_CLOSE_ABS(sls.direction(slsState), newMom.normalized(), 1e-6);
190  CHECK_CLOSE_ABS(sls.absoluteMomentum(slsState), newMom.norm(), 1e-6);
191  BOOST_CHECK_EQUAL(sls.charge(slsState), charge);
192  BOOST_CHECK_EQUAL(sls.time(slsState), newTime);
193 
194  // The covariance transport
195  slsState.cov = cov;
196  sls.transportCovarianceToCurvilinear(slsState);
197  BOOST_CHECK_NE(slsState.cov, cov);
198  BOOST_CHECK_NE(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
199  BOOST_CHECK_EQUAL(slsState.jacTransport, FreeMatrix::Identity());
200  BOOST_CHECK_EQUAL(slsState.derivative, FreeVector::Zero());
201 
202  // Perform a step without and with covariance transport
203  slsState.cov = cov;
204  PropState ps(navDir, slsState);
205 
206  ps.stepping.covTransport = false;
207  double h = sls.step(ps, mockNavigator).value();
208  BOOST_CHECK_EQUAL(ps.stepping.stepSize.value(), stepSize);
209  BOOST_CHECK_EQUAL(ps.stepping.stepSize.value(), h * navDir);
210  CHECK_CLOSE_COVARIANCE(ps.stepping.cov, cov, 1e-6);
211  BOOST_CHECK_GT(sls.position(ps.stepping).norm(), newPos.norm());
212  CHECK_CLOSE_ABS(sls.direction(ps.stepping), newMom.normalized(), 1e-6);
213  CHECK_CLOSE_ABS(sls.absoluteMomentum(ps.stepping), newMom.norm(), 1e-6);
214  CHECK_CLOSE_ABS(sls.charge(ps.stepping), charge, 1e-6);
215  BOOST_CHECK_LT(sls.time(ps.stepping), newTime);
216  BOOST_CHECK_EQUAL(ps.stepping.derivative, FreeVector::Zero());
217  BOOST_CHECK_EQUAL(ps.stepping.jacTransport, FreeMatrix::Identity());
218 
219  ps.stepping.covTransport = true;
220  double h2 = sls.step(ps, mockNavigator).value();
221  BOOST_CHECK_EQUAL(ps.stepping.stepSize.value(), stepSize);
222  BOOST_CHECK_EQUAL(h2, h);
223  CHECK_CLOSE_COVARIANCE(ps.stepping.cov, cov, 1e-6);
224  BOOST_CHECK_GT(sls.position(ps.stepping).norm(), newPos.norm());
225  CHECK_CLOSE_ABS(sls.direction(ps.stepping), newMom.normalized(), 1e-6);
226  CHECK_CLOSE_ABS(sls.absoluteMomentum(ps.stepping), newMom.norm(), 1e-6);
227  CHECK_CLOSE_ABS(sls.charge(ps.stepping), charge, 1e-6);
228  BOOST_CHECK_LT(sls.time(ps.stepping), newTime);
229  BOOST_CHECK_NE(ps.stepping.derivative, FreeVector::Zero());
230  BOOST_CHECK_NE(ps.stepping.jacTransport, FreeMatrix::Identity());
231 
233  // Construct the parameters
234  Vector3 pos2(1.5, -2.5, 3.5);
235  Vector3 dir2 = Vector3(4.5, -5.5, 6.5).normalized();
236  double time2 = 7.5;
237  double absMom2 = 8.5;
238  double charge2 = 1.;
239  BoundSquareMatrix cov2 = 8.5 * Covariance::Identity();
240  CurvilinearTrackParameters cp2(makeVector4(pos2, time2), dir2,
241  charge2 / absMom2, cov2,
242  ParticleHypothesis::pion());
243  BOOST_CHECK(cp2.covariance().has_value());
244  FreeVector freeParams = detail::transformBoundToFreeParameters(
245  cp2.referenceSurface(), tgContext, cp2.parameters());
246  navDir = Direction::Forward;
247  double stepSize2 = -2. * stepSize;
248 
249  // Reset all possible parameters
250  StraightLineStepper::State slsStateCopy(ps.stepping);
251  sls.resetState(slsStateCopy, cp2.parameters(), *cp2.covariance(),
252  cp2.referenceSurface(), stepSize2);
253  // Test all components
254  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, BoundToFreeMatrix::Zero());
255  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, ps.stepping.jacToGlobal);
256  BOOST_CHECK_EQUAL(slsStateCopy.jacTransport, FreeMatrix::Identity());
257  BOOST_CHECK_EQUAL(slsStateCopy.derivative, FreeVector::Zero());
258  BOOST_CHECK(slsStateCopy.covTransport);
259  BOOST_CHECK_EQUAL(slsStateCopy.cov, cov2);
260  CHECK_CLOSE_ABS(sls.position(slsStateCopy),
261  freeParams.template segment<3>(eFreePos0), 1e-6);
262  CHECK_CLOSE_ABS(sls.direction(slsStateCopy),
263  freeParams.template segment<3>(eFreeDir0).normalized(), 1e-6);
264  CHECK_CLOSE_ABS(sls.absoluteMomentum(slsStateCopy),
265  std::abs(1. / freeParams[eFreeQOverP]), 1e-6);
266  CHECK_CLOSE_ABS(sls.charge(slsStateCopy), -sls.charge(ps.stepping), 1e-6);
267  CHECK_CLOSE_ABS(sls.time(slsStateCopy), freeParams[eFreeTime], 1e-6);
268  BOOST_CHECK_EQUAL(slsStateCopy.pathAccumulated, 0.);
269  BOOST_CHECK_EQUAL(slsStateCopy.stepSize.value(), navDir * stepSize2);
270  BOOST_CHECK_EQUAL(slsStateCopy.previousStepSize,
271  ps.stepping.previousStepSize);
272  BOOST_CHECK_EQUAL(slsStateCopy.tolerance, ps.stepping.tolerance);
273 
274  // Reset all possible parameters except the step size
275  slsStateCopy = ps.stepping;
276  sls.resetState(slsStateCopy, cp2.parameters(), *cp2.covariance(),
277  cp2.referenceSurface());
278  // Test all components
279  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, BoundToFreeMatrix::Zero());
280  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, ps.stepping.jacToGlobal);
281  BOOST_CHECK_EQUAL(slsStateCopy.jacTransport, FreeMatrix::Identity());
282  BOOST_CHECK_EQUAL(slsStateCopy.derivative, FreeVector::Zero());
283  BOOST_CHECK(slsStateCopy.covTransport);
284  BOOST_CHECK_EQUAL(slsStateCopy.cov, cov2);
285  CHECK_CLOSE_ABS(sls.position(slsStateCopy),
286  freeParams.template segment<3>(eFreePos0), 1e-6);
287  CHECK_CLOSE_ABS(sls.direction(slsStateCopy),
288  freeParams.template segment<3>(eFreeDir0), 1e-6);
289  CHECK_CLOSE_ABS(sls.absoluteMomentum(slsStateCopy),
290  std::abs(1. / freeParams[eFreeQOverP]), 1e-6);
291  CHECK_CLOSE_ABS(sls.charge(slsStateCopy), -sls.charge(ps.stepping), 1e-6);
292  CHECK_CLOSE_ABS(sls.time(slsStateCopy), freeParams[eFreeTime], 1e-6);
293  BOOST_CHECK_EQUAL(slsStateCopy.pathAccumulated, 0.);
294  BOOST_CHECK_EQUAL(slsStateCopy.stepSize.value(),
295  std::numeric_limits<double>::max());
296  BOOST_CHECK_EQUAL(slsStateCopy.previousStepSize,
297  ps.stepping.previousStepSize);
298  BOOST_CHECK_EQUAL(slsStateCopy.tolerance, ps.stepping.tolerance);
299 
300  // Reset the least amount of parameters
301  slsStateCopy = ps.stepping;
302  sls.resetState(slsStateCopy, cp2.parameters(), *cp2.covariance(),
303  cp2.referenceSurface());
304  // Test all components
305  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, BoundToFreeMatrix::Zero());
306  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, ps.stepping.jacToGlobal);
307  BOOST_CHECK_EQUAL(slsStateCopy.jacTransport, FreeMatrix::Identity());
308  BOOST_CHECK_EQUAL(slsStateCopy.derivative, FreeVector::Zero());
309  BOOST_CHECK(slsStateCopy.covTransport);
310  BOOST_CHECK_EQUAL(slsStateCopy.cov, cov2);
311  CHECK_CLOSE_ABS(sls.position(slsStateCopy),
312  freeParams.template segment<3>(eFreePos0), 1e-6);
313  CHECK_CLOSE_ABS(sls.direction(slsStateCopy),
314  freeParams.template segment<3>(eFreeDir0).normalized(), 1e-6);
315  CHECK_CLOSE_ABS(sls.absoluteMomentum(slsStateCopy),
316  std::abs(1. / freeParams[eFreeQOverP]), 1e-6);
317  CHECK_CLOSE_ABS(sls.charge(slsStateCopy), -sls.charge(ps.stepping), 1e-6);
318  CHECK_CLOSE_ABS(sls.time(slsStateCopy), freeParams[eFreeTime], 1e-6);
319  BOOST_CHECK_EQUAL(slsStateCopy.pathAccumulated, 0.);
320  BOOST_CHECK_EQUAL(slsStateCopy.stepSize.value(),
321  std::numeric_limits<double>::max());
322  BOOST_CHECK_EQUAL(slsStateCopy.previousStepSize,
323  ps.stepping.previousStepSize);
324  BOOST_CHECK_EQUAL(slsStateCopy.tolerance, ps.stepping.tolerance);
325 
327  auto plane = Surface::makeShared<PlaneSurface>(pos, dir);
328  auto bp = BoundTrackParameters::create(
329  plane, tgContext, makeVector4(pos, time), dir, charge / absMom,
330  cov, ParticleHypothesis::pion())
331  .value();
332  slsState =
333  StraightLineStepper::State(tgContext, mfContext, cp, stepSize, tolerance);
334 
335  // Test the intersection in the context of a surface
336  auto targetSurface =
337  Surface::makeShared<PlaneSurface>(pos + navDir * 2. * dir, dir);
338  sls.updateSurfaceStatus(slsState, *targetSurface, navDir,
339  BoundaryCheck(false));
340  CHECK_CLOSE_ABS(slsState.stepSize.value(ConstrainedStep::actor), navDir * 2.,
341  1e-6);
342 
343  // Test the step size modification in the context of a surface
344  sls.updateStepSize(
345  slsState,
346  targetSurface
347  ->intersect(slsState.geoContext, sls.position(slsState),
348  navDir * sls.direction(slsState), false)
349  .closest(),
350  navDir, false);
351  CHECK_CLOSE_ABS(slsState.stepSize.value(), 2, 1e-6);
352  slsState.stepSize.setUser(navDir * stepSize);
353  sls.updateStepSize(
354  slsState,
355  targetSurface
356  ->intersect(slsState.geoContext, sls.position(slsState),
357  navDir * sls.direction(slsState), false)
358  .closest(),
359  navDir, true);
360  CHECK_CLOSE_ABS(slsState.stepSize.value(), 2, 1e-6);
361 
362  // Test the bound state construction
363  FreeToBoundCorrection freeToBoundCorrection(false);
364  auto boundState =
365  sls.boundState(slsState, *plane, true, freeToBoundCorrection).value();
366  auto boundPars = std::get<0>(boundState);
367  CHECK_CLOSE_ABS(boundPars.position(tgContext), bp.position(tgContext), 1e-6);
368  CHECK_CLOSE_ABS(boundPars.momentum(), bp.momentum(), 1e-6);
369  CHECK_CLOSE_ABS(boundPars.charge(), bp.charge(), 1e-6);
370  CHECK_CLOSE_ABS(boundPars.time(), bp.time(), 1e-6);
371  BOOST_CHECK(boundPars.covariance().has_value());
372  BOOST_CHECK_NE(*boundPars.covariance(), cov);
373  CHECK_CLOSE_COVARIANCE(std::get<1>(boundState),
374  BoundMatrix(BoundMatrix::Identity()), 1e-6);
375  CHECK_CLOSE_ABS(std::get<2>(boundState), 0., 1e-6);
376 
377  // Transport the covariance in the context of a surface
378  sls.transportCovarianceToBound(slsState, *plane, freeToBoundCorrection);
379  BOOST_CHECK_NE(slsState.cov, cov);
380  BOOST_CHECK_NE(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
381  BOOST_CHECK_EQUAL(slsState.jacTransport, FreeMatrix::Identity());
382  BOOST_CHECK_EQUAL(slsState.derivative, FreeVector::Zero());
383 
384  // Update in context of a surface
385  freeParams = detail::transformBoundToFreeParameters(
386  bp.referenceSurface(), tgContext, bp.parameters());
387  freeParams.segment<3>(eFreePos0) *= 2;
388  freeParams[eFreeTime] *= 2;
389 
390  BOOST_CHECK(bp.covariance().has_value());
391  sls.update(slsState, freeParams, bp.parameters(), 2 * (*bp.covariance()),
392  *plane);
393  CHECK_CLOSE_OR_SMALL(sls.position(slsState), 2. * pos, eps, eps);
394  BOOST_CHECK_EQUAL(sls.charge(slsState), 1. * charge);
395  CHECK_CLOSE_OR_SMALL(sls.time(slsState), 2. * time, eps, eps);
396  CHECK_CLOSE_COVARIANCE(slsState.cov, Covariance(2. * cov), 1e-6);
397 }
398 } // namespace Test
399 } // namespace Acts