EIC 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/EventData/NeutralTrackParameters.hpp"
12 #include "Acts/EventData/TrackParameters.hpp"
13 #include "Acts/EventData/detail/TransformationBoundToFree.hpp"
14 #include "Acts/Propagator/StraightLineStepper.hpp"
15 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
16 #include "Acts/Utilities/Helpers.hpp"
17 
18 #include <limits>
19 
20 namespace tt = boost::test_tools;
21 using Acts::VectorHelpers::makeVector4;
22 
23 namespace Acts {
24 namespace Test {
25 
26 using Covariance = BoundSymMatrix;
27 using Jacobian = BoundMatrix;
28 
30 struct PropState {
32  PropState(StraightLineStepper::State sState) : stepping(sState) {}
34  StraightLineStepper::State stepping;
36  struct {
37  double mass = 42.;
38  } options;
39 };
40 
41 static constexpr auto eps = 2 * std::numeric_limits<double>::epsilon();
42 
44 BOOST_AUTO_TEST_CASE(straight_line_stepper_state_test) {
45  // Set up some variables
46  GeometryContext tgContext = GeometryContext();
47  MagneticFieldContext mfContext = MagneticFieldContext();
48  NavigationDirection ndir = backward;
49  double stepSize = 123.;
50  double tolerance = 234.;
51 
52  Vector3D pos(1., 2., 3.);
53  Vector3D dir(4., 5., 6.);
54  double time = 7.;
55  double absMom = 8.;
56  double charge = -1.;
57 
58  // Test charged parameters without covariance matrix
59  CurvilinearTrackParameters cp(makeVector4(pos, time), dir, absMom, charge);
60  StraightLineStepper::State slsState(tgContext, mfContext, cp, ndir, stepSize,
61  tolerance);
62 
63  // Test the result & compare with the input/test for reasonable members
64  BOOST_CHECK_EQUAL(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
65  BOOST_CHECK_EQUAL(slsState.jacTransport, FreeMatrix::Identity());
66  BOOST_CHECK_EQUAL(slsState.derivative, FreeVector::Zero());
67  BOOST_CHECK(!slsState.covTransport);
68  BOOST_CHECK_EQUAL(slsState.cov, Covariance::Zero());
69  CHECK_CLOSE_OR_SMALL(slsState.pos, pos, eps, eps);
70  CHECK_CLOSE_OR_SMALL(slsState.dir, dir.normalized(), eps, eps);
71  CHECK_CLOSE_REL(slsState.p, absMom, eps);
72  BOOST_CHECK_EQUAL(slsState.q, charge);
73  CHECK_CLOSE_OR_SMALL(slsState.t, time, eps, eps);
74  BOOST_CHECK_EQUAL(slsState.navDir, ndir);
75  BOOST_CHECK_EQUAL(slsState.pathAccumulated, 0.);
76  BOOST_CHECK_EQUAL(slsState.stepSize, ndir * stepSize);
77  BOOST_CHECK_EQUAL(slsState.previousStepSize, 0.);
78  BOOST_CHECK_EQUAL(slsState.tolerance, tolerance);
79 
80  // Test without charge and covariance matrix
81  NeutralCurvilinearTrackParameters ncp(makeVector4(pos, time), dir,
82  1 / absMom);
83  slsState = StraightLineStepper::State(tgContext, mfContext, ncp, ndir,
84  stepSize, tolerance);
85  BOOST_CHECK_EQUAL(slsState.q, 0.);
86 
87  // Test with covariance matrix
88  Covariance cov = 8. * Covariance::Identity();
89  ncp = NeutralCurvilinearTrackParameters(makeVector4(pos, time), dir,
90  1 / absMom, cov);
91  slsState = StraightLineStepper::State(tgContext, mfContext, ncp, ndir,
92  stepSize, tolerance);
93  BOOST_CHECK_NE(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
94  BOOST_CHECK(slsState.covTransport);
95  BOOST_CHECK_EQUAL(slsState.cov, cov);
96 }
97 
100 BOOST_AUTO_TEST_CASE(straight_line_stepper_test) {
101  // Set up some variables for the state
102  GeometryContext tgContext = GeometryContext();
103  MagneticFieldContext mfContext = MagneticFieldContext();
104  NavigationDirection ndir = backward;
105  double stepSize = 123.;
106  double tolerance = 234.;
107 
108  // Construct the parameters
109  Vector3D pos(1., 2., 3.);
110  Vector3D dir = Vector3D(4., 5., 6.).normalized();
111  double time = 7.;
112  double absMom = 8.;
113  double charge = -1.;
114  Covariance cov = 8. * Covariance::Identity();
115  CurvilinearTrackParameters cp(makeVector4(pos, time), dir, charge / absMom,
116  cov);
117 
118  // Build the state and the stepper
119  StraightLineStepper::State slsState(tgContext, mfContext, cp, ndir, stepSize,
120  tolerance);
121  StraightLineStepper sls;
122 
123  // Test the getters
124  BOOST_CHECK_EQUAL(sls.position(slsState), slsState.pos);
125  BOOST_CHECK_EQUAL(sls.direction(slsState), slsState.dir);
126  BOOST_CHECK_EQUAL(sls.momentum(slsState), slsState.p);
127  BOOST_CHECK_EQUAL(sls.charge(slsState), slsState.q);
128  BOOST_CHECK_EQUAL(sls.time(slsState), slsState.t);
129 
130  //~ BOOST_CHECK_EQUAL(sls.overstepLimit(slsState), tolerance);
131 
132  // Step size modifies
133  const std::string originalStepSize = slsState.stepSize.toString();
134 
135  sls.setStepSize(slsState, 1337.);
136  BOOST_CHECK_EQUAL(slsState.previousStepSize, ndir * stepSize);
137  BOOST_CHECK_EQUAL(slsState.stepSize, 1337.);
138 
139  sls.releaseStepSize(slsState);
140  BOOST_CHECK_EQUAL(slsState.stepSize, -123.);
141  BOOST_CHECK_EQUAL(sls.outputStepSize(slsState), originalStepSize);
142 
143  // Test the curvilinear state construction
144  auto curvState = sls.curvilinearState(slsState);
145  auto curvPars = std::get<0>(curvState);
146  CHECK_CLOSE_ABS(curvPars.position(tgContext), cp.position(tgContext), 1e-6);
147  CHECK_CLOSE_ABS(curvPars.momentum(), cp.momentum(), 1e-6);
148  CHECK_CLOSE_ABS(curvPars.charge(), cp.charge(), 1e-6);
149  CHECK_CLOSE_ABS(curvPars.time(), cp.time(), 1e-6);
150  BOOST_CHECK(curvPars.covariance().has_value());
151  BOOST_CHECK_NE(*curvPars.covariance(), cov);
152  CHECK_CLOSE_COVARIANCE(std::get<1>(curvState),
153  BoundMatrix(BoundMatrix::Identity()), 1e-6);
154  CHECK_CLOSE_ABS(std::get<2>(curvState), 0., 1e-6);
155 
156  // Test the update method
157  Vector3D newPos(2., 4., 8.);
158  Vector3D newMom(3., 9., 27.);
159  double newTime(321.);
160  sls.update(slsState, newPos, newMom.normalized(), newMom.norm(), newTime);
161  BOOST_CHECK_EQUAL(slsState.pos, newPos);
162  BOOST_CHECK_EQUAL(slsState.dir, newMom.normalized());
163  BOOST_CHECK_EQUAL(slsState.p, newMom.norm());
164  BOOST_CHECK_EQUAL(slsState.q, charge);
165  BOOST_CHECK_EQUAL(slsState.t, newTime);
166 
167  // The covariance transport
168  slsState.cov = cov;
169  sls.covarianceTransport(slsState);
170  BOOST_CHECK_NE(slsState.cov, cov);
171  BOOST_CHECK_NE(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
172  BOOST_CHECK_EQUAL(slsState.jacTransport, FreeMatrix::Identity());
173  BOOST_CHECK_EQUAL(slsState.derivative, FreeVector::Zero());
174 
175  // Perform a step without and with covariance transport
176  slsState.cov = cov;
177  PropState ps(slsState);
178 
179  ps.stepping.covTransport = false;
180  double h = sls.step(ps).value();
181  BOOST_CHECK_EQUAL(ps.stepping.stepSize, ndir * stepSize);
182  BOOST_CHECK_EQUAL(ps.stepping.stepSize, h);
183  CHECK_CLOSE_COVARIANCE(ps.stepping.cov, cov, 1e-6);
184  BOOST_CHECK_GT(ps.stepping.pos.norm(), newPos.norm());
185  BOOST_CHECK_EQUAL(ps.stepping.dir, newMom.normalized());
186  BOOST_CHECK_EQUAL(ps.stepping.p, newMom.norm());
187  BOOST_CHECK_EQUAL(ps.stepping.q, charge);
188  BOOST_CHECK_LT(ps.stepping.t, newTime);
189  BOOST_CHECK_EQUAL(ps.stepping.derivative, FreeVector::Zero());
190  BOOST_CHECK_EQUAL(ps.stepping.jacTransport, FreeMatrix::Identity());
191 
192  ps.stepping.covTransport = true;
193  double h2 = sls.step(ps).value();
194  BOOST_CHECK_EQUAL(ps.stepping.stepSize, ndir * stepSize);
195  BOOST_CHECK_EQUAL(h2, h);
196  CHECK_CLOSE_COVARIANCE(ps.stepping.cov, cov, 1e-6);
197  BOOST_CHECK_GT(ps.stepping.pos.norm(), newPos.norm());
198  BOOST_CHECK_EQUAL(ps.stepping.dir, newMom.normalized());
199  BOOST_CHECK_EQUAL(ps.stepping.p, newMom.norm());
200  BOOST_CHECK_EQUAL(ps.stepping.q, charge);
201  BOOST_CHECK_LT(ps.stepping.t, newTime);
202  BOOST_CHECK_NE(ps.stepping.derivative, FreeVector::Zero());
203  BOOST_CHECK_NE(ps.stepping.jacTransport, FreeMatrix::Identity());
204 
206  // Construct the parameters
207  Vector3D pos2(1.5, -2.5, 3.5);
208  Vector3D dir2 = Vector3D(4.5, -5.5, 6.5).normalized();
209  double time2 = 7.5;
210  double absMom2 = 8.5;
211  double charge2 = 1.;
212  BoundSymMatrix cov2 = 8.5 * Covariance::Identity();
213  CurvilinearTrackParameters cp2(makeVector4(pos2, time2), dir2, absMom2,
214  charge2, cov2);
215  FreeVector freeParams = detail::transformBoundToFreeParameters(
216  cp2.referenceSurface(), tgContext, cp2.parameters());
217  ndir = forward;
218  double stepSize2 = -2. * stepSize;
219 
220  // Reset all possible parameters
221  StraightLineStepper::State slsStateCopy(ps.stepping);
222  sls.resetState(slsStateCopy, cp2.parameters(), *cp2.covariance(),
223  cp2.referenceSurface(), ndir, stepSize2);
224  // Test all components
225  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, BoundToFreeMatrix::Zero());
226  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, ps.stepping.jacToGlobal);
227  BOOST_CHECK_EQUAL(slsStateCopy.jacTransport, FreeMatrix::Identity());
228  BOOST_CHECK_EQUAL(slsStateCopy.derivative, FreeVector::Zero());
229  BOOST_CHECK(slsStateCopy.covTransport);
230  BOOST_CHECK_EQUAL(slsStateCopy.cov, cov2);
231  BOOST_CHECK_EQUAL(slsStateCopy.pos,
232  freeParams.template segment<3>(eFreePos0));
233  BOOST_CHECK_EQUAL(slsStateCopy.dir,
234  freeParams.template segment<3>(eFreeDir0).normalized());
235  BOOST_CHECK_EQUAL(slsStateCopy.p, std::abs(1. / freeParams[eFreeQOverP]));
236  BOOST_CHECK_EQUAL(slsStateCopy.q, ps.stepping.q);
237  BOOST_CHECK_EQUAL(slsStateCopy.t, freeParams[eFreeTime]);
238  BOOST_CHECK_EQUAL(slsStateCopy.navDir, ndir);
239  BOOST_CHECK_EQUAL(slsStateCopy.pathAccumulated, 0.);
240  BOOST_CHECK_EQUAL(slsStateCopy.stepSize, ndir * stepSize2);
241  BOOST_CHECK_EQUAL(slsStateCopy.previousStepSize,
242  ps.stepping.previousStepSize);
243  BOOST_CHECK_EQUAL(slsStateCopy.tolerance, ps.stepping.tolerance);
244 
245  // Reset all possible parameters except the step size
246  slsStateCopy = ps.stepping;
247  sls.resetState(slsStateCopy, cp2.parameters(), *cp2.covariance(),
248  cp2.referenceSurface(), ndir);
249  // Test all components
250  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, BoundToFreeMatrix::Zero());
251  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, ps.stepping.jacToGlobal);
252  BOOST_CHECK_EQUAL(slsStateCopy.jacTransport, FreeMatrix::Identity());
253  BOOST_CHECK_EQUAL(slsStateCopy.derivative, FreeVector::Zero());
254  BOOST_CHECK(slsStateCopy.covTransport);
255  BOOST_CHECK_EQUAL(slsStateCopy.cov, cov2);
256  BOOST_CHECK_EQUAL(slsStateCopy.pos,
257  freeParams.template segment<3>(eFreePos0));
258  BOOST_CHECK_EQUAL(slsStateCopy.dir,
259  freeParams.template segment<3>(eFreeDir0));
260  BOOST_CHECK_EQUAL(slsStateCopy.p, std::abs(1. / freeParams[eFreeQOverP]));
261  BOOST_CHECK_EQUAL(slsStateCopy.q, ps.stepping.q);
262  BOOST_CHECK_EQUAL(slsStateCopy.t, freeParams[eFreeTime]);
263  BOOST_CHECK_EQUAL(slsStateCopy.navDir, ndir);
264  BOOST_CHECK_EQUAL(slsStateCopy.pathAccumulated, 0.);
265  BOOST_CHECK_EQUAL(slsStateCopy.stepSize,
266  ndir * std::numeric_limits<double>::max());
267  BOOST_CHECK_EQUAL(slsStateCopy.previousStepSize,
268  ps.stepping.previousStepSize);
269  BOOST_CHECK_EQUAL(slsStateCopy.tolerance, ps.stepping.tolerance);
270 
271  // Reset the least amount of parameters
272  slsStateCopy = ps.stepping;
273  sls.resetState(slsStateCopy, cp2.parameters(), *cp2.covariance(),
274  cp2.referenceSurface());
275  // Test all components
276  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, BoundToFreeMatrix::Zero());
277  BOOST_CHECK_NE(slsStateCopy.jacToGlobal, ps.stepping.jacToGlobal);
278  BOOST_CHECK_EQUAL(slsStateCopy.jacTransport, FreeMatrix::Identity());
279  BOOST_CHECK_EQUAL(slsStateCopy.derivative, FreeVector::Zero());
280  BOOST_CHECK(slsStateCopy.covTransport);
281  BOOST_CHECK_EQUAL(slsStateCopy.cov, cov2);
282  BOOST_CHECK_EQUAL(slsStateCopy.pos,
283  freeParams.template segment<3>(eFreePos0));
284  BOOST_CHECK_EQUAL(slsStateCopy.dir,
285  freeParams.template segment<3>(eFreeDir0).normalized());
286  BOOST_CHECK_EQUAL(slsStateCopy.p, std::abs(1. / freeParams[eFreeQOverP]));
287  BOOST_CHECK_EQUAL(slsStateCopy.q, ps.stepping.q);
288  BOOST_CHECK_EQUAL(slsStateCopy.t, freeParams[eFreeTime]);
289  BOOST_CHECK_EQUAL(slsStateCopy.navDir, forward);
290  BOOST_CHECK_EQUAL(slsStateCopy.pathAccumulated, 0.);
291  BOOST_CHECK_EQUAL(slsStateCopy.stepSize, std::numeric_limits<double>::max());
292  BOOST_CHECK_EQUAL(slsStateCopy.previousStepSize,
293  ps.stepping.previousStepSize);
294  BOOST_CHECK_EQUAL(slsStateCopy.tolerance, ps.stepping.tolerance);
295 
297  auto plane = Surface::makeShared<PlaneSurface>(pos, dir);
298  BoundTrackParameters bp(plane, tgContext, makeVector4(pos, time), dir,
299  charge / absMom, cov);
300  slsState = StraightLineStepper::State(tgContext, mfContext, cp, ndir,
301  stepSize, tolerance);
302 
303  // Test the intersection in the context of a surface
304  auto targetSurface =
305  Surface::makeShared<PlaneSurface>(pos + ndir * 2. * dir, dir);
306  sls.updateSurfaceStatus(slsState, *targetSurface, BoundaryCheck(false));
307  CHECK_CLOSE_ABS(slsState.stepSize.value(ConstrainedStep::actor), ndir * 2.,
308  1e-6);
309 
310  // Test the step size modification in the context of a surface
311  sls.updateStepSize(
312  slsState,
313  targetSurface->intersect(slsState.geoContext, slsState.pos,
314  slsState.navDir * slsState.dir, false),
315  false);
316  CHECK_CLOSE_ABS(slsState.stepSize, 2, 1e-6);
317  slsState.stepSize = ndir * stepSize;
318  sls.updateStepSize(
319  slsState,
320  targetSurface->intersect(slsState.geoContext, slsState.pos,
321  slsState.navDir * slsState.dir, false),
322  true);
323  CHECK_CLOSE_ABS(slsState.stepSize, 2, 1e-6);
324 
325  // Test the bound state construction
326  auto boundState = sls.boundState(slsState, *plane);
327  auto boundPars = std::get<0>(boundState);
328  CHECK_CLOSE_ABS(boundPars.position(tgContext), bp.position(tgContext), 1e-6);
329  CHECK_CLOSE_ABS(boundPars.momentum(), bp.momentum(), 1e-6);
330  CHECK_CLOSE_ABS(boundPars.charge(), bp.charge(), 1e-6);
331  CHECK_CLOSE_ABS(boundPars.time(), bp.time(), 1e-6);
332  BOOST_CHECK(boundPars.covariance().has_value());
333  BOOST_CHECK_NE(*boundPars.covariance(), cov);
334  CHECK_CLOSE_COVARIANCE(std::get<1>(boundState),
335  BoundMatrix(BoundMatrix::Identity()), 1e-6);
336  CHECK_CLOSE_ABS(std::get<2>(boundState), 0., 1e-6);
337 
338  // Transport the covariance in the context of a surface
339  sls.covarianceTransport(slsState, *plane);
340  BOOST_CHECK_NE(slsState.cov, cov);
341  BOOST_CHECK_NE(slsState.jacToGlobal, BoundToFreeMatrix::Zero());
342  BOOST_CHECK_EQUAL(slsState.jacTransport, FreeMatrix::Identity());
343  BOOST_CHECK_EQUAL(slsState.derivative, FreeVector::Zero());
344 
345  // Update in context of a surface
346  freeParams = detail::transformBoundToFreeParameters(
347  bp.referenceSurface(), tgContext, bp.parameters());
348  freeParams.segment<3>(eFreePos0) *= 2;
349  freeParams[eFreeTime] *= 2;
350  freeParams.segment<3>(eFreeDir0) *= 2;
351  freeParams[eFreeQOverP] *= -0.5;
352 
353  sls.update(slsState, freeParams, 2 * (*bp.covariance()));
354  CHECK_CLOSE_OR_SMALL(slsState.pos, 2. * pos, eps, eps);
355  CHECK_CLOSE_OR_SMALL(slsState.dir, dir, eps, eps);
356  CHECK_CLOSE_REL(slsState.p, 2. * absMom, eps);
357  BOOST_CHECK_EQUAL(slsState.q, 1. * charge);
358  CHECK_CLOSE_OR_SMALL(slsState.t, 2. * time, eps, eps);
359  CHECK_CLOSE_COVARIANCE(slsState.cov, Covariance(2. * cov), 1e-6);
360 }
361 } // namespace Test
362 } // namespace Acts