Core/Propagator/AtlasStepperTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2020 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 http://mozilla.org/MPL/2.0/.
0008
0009 #include <boost/test/test_tools.hpp>
0010 #include <boost/test/unit_test.hpp>
0011
0012 #include "Acts/Definitions/Algebra.hpp"
0013 #include "Acts/Definitions/Common.hpp"
0014 #include "Acts/Definitions/Direction.hpp"
0015 #include "Acts/Definitions/TrackParametrization.hpp"
0016 #include "Acts/Definitions/Units.hpp"
0017 #include "Acts/EventData/Charge.hpp"
0018 #include "Acts/EventData/GenericBoundTrackParameters.hpp"
0019 #include "Acts/EventData/GenericCurvilinearTrackParameters.hpp"
0020 #include "Acts/EventData/ParticleHypothesis.hpp"
0021 #include "Acts/EventData/TrackParameters.hpp"
0022 #include "Acts/EventData/TransformationHelpers.hpp"
0023 #include "Acts/Geometry/GeometryContext.hpp"
0024 #include "Acts/MagneticField/ConstantBField.hpp"
0025 #include "Acts/MagneticField/MagneticFieldContext.hpp"
0026 #include "Acts/MagneticField/MagneticFieldProvider.hpp"
0027 #include "Acts/Propagator/AtlasStepper.hpp"
0028 #include "Acts/Propagator/ConstrainedStep.hpp"
0029 #include "Acts/Surfaces/BoundaryCheck.hpp"
0030 #include "Acts/Surfaces/DiscSurface.hpp"
0031 #include "Acts/Surfaces/PerigeeSurface.hpp"
0032 #include "Acts/Surfaces/PlaneSurface.hpp"
0033 #include "Acts/Surfaces/StrawSurface.hpp"
0034 #include "Acts/Surfaces/Surface.hpp"
0035 #include "Acts/Tests/CommonHelpers/Assertions.hpp"
0036 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0037 #include "Acts/Utilities/Helpers.hpp"
0038 #include "Acts/Utilities/Result.hpp"
0039
0040 #include <algorithm>
0041 #include <array>
0042 #include <functional>
0043 #include <iterator>
0044 #include <limits>
0045 #include <memory>
0046 #include <optional>
0047 #include <tuple>
0048 #include <type_traits>
0049 #include <utility>
0050
0051 namespace Acts::Test {
0052
0053 using namespace Acts::UnitLiterals;
0054 using Acts::VectorHelpers::makeVector4;
0055 using Covariance = BoundSquareMatrix;
0056 using Jacobian = BoundMatrix;
0057 using Stepper = Acts::AtlasStepper;
0058
0059 /// Simplified propagator state.
0060 struct MockPropagatorState {
0061   MockPropagatorState(Stepper::State stepperState)
0062       : stepping(std::move(stepperState)) {}
0063
0064   /// Stepper state.
0065   Stepper::State stepping;
0066   /// Propagator options with only the relevant components.
0067   struct {
0068     double stepTolerance = 10_um;
0069     Direction direction = Direction::Backward;
0070   } options;
0071 };
0072
0073 struct MockNavigator {};
0074
0075 static constexpr MockNavigator navigator;
0076
0077 // epsilon for floating point comparisons
0078 static constexpr auto eps = 1024 * std::numeric_limits<double>::epsilon();
0079
0080 // propagation settings
0081 static constexpr auto stepSize = 10_mm;
0082 static constexpr auto tolerance = 10_um;
0083 static constexpr Direction navDir = Direction::Backward;
0084 static auto magneticField =
0085     std::make_shared<ConstantBField>(Vector3(0.1_T, -0.2_T, 2_T));
0086
0087 // initial parameter state
0088 static const Vector4 pos4(1_mm, -1_mm, 2_mm, 2_ns);
0089 static const Vector3 pos = pos4.segment<3>(ePos0);
0090 static const auto time = pos4[eTime];
0091 static const Vector3 unitDir = Vector3(-2, 2, 1).normalized();
0092 static constexpr auto absMom = 1_GeV;
0093 static constexpr auto charge = -1_e;
0094 static const auto particleHypothesis = ParticleHypothesis::pion();
0095 static const Covariance cov = Covariance::Identity();
0096
0097 // context objects
0098 static const GeometryContext geoCtx;
0099 static const MagneticFieldContext magCtx;
0100
0101 BOOST_AUTO_TEST_SUITE(AtlasStepper)
0102
0103 // test state construction from parameters w/o covariance
0104 BOOST_AUTO_TEST_CASE(ConstructState) {
0105   Stepper::State state(
0106       geoCtx, magneticField->makeCache(magCtx),
0107       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, std::nullopt,
0108                                  particleHypothesis),
0109       stepSize, tolerance);
0110
0111   BOOST_CHECK(!state.covTransport);
0112   BOOST_CHECK_EQUAL(state.covariance, nullptr);
0113   BOOST_CHECK_EQUAL(state.pVector[0], pos.x());
0114   BOOST_CHECK_EQUAL(state.pVector[1], pos.y());
0115   BOOST_CHECK_EQUAL(state.pVector[2], pos.z());
0116   BOOST_CHECK_EQUAL(state.pVector[3], time);
0117   CHECK_CLOSE_ABS(state.pVector[4], unitDir.x(), eps);
0118   CHECK_CLOSE_ABS(state.pVector[5], unitDir.y(), eps);
0119   CHECK_CLOSE_ABS(state.pVector[6], unitDir.z(), eps);
0120   BOOST_CHECK_EQUAL(state.pVector[7], charge / absMom);
0121   BOOST_CHECK_EQUAL(state.pathAccumulated, 0.);
0122   BOOST_CHECK_EQUAL(state.stepSize.value(), stepSize);
0123   BOOST_CHECK_EQUAL(state.previousStepSize, 0.);
0124   BOOST_CHECK_EQUAL(state.tolerance, tolerance);
0125 }
0126
0127 // test state construction from parameters w/ covariance
0128 BOOST_AUTO_TEST_CASE(ConstructStateWithCovariance) {
0129   Stepper::State state(
0130       geoCtx, magneticField->makeCache(magCtx),
0131       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0132                                  particleHypothesis),
0133       stepSize, tolerance);
0134
0135   BOOST_CHECK(state.covTransport);
0136   BOOST_CHECK_EQUAL(*state.covariance, cov);
0137   BOOST_CHECK_EQUAL(state.pVector[0], pos.x());
0138   BOOST_CHECK_EQUAL(state.pVector[1], pos.y());
0139   BOOST_CHECK_EQUAL(state.pVector[2], pos.z());
0140   BOOST_CHECK_EQUAL(state.pVector[3], time);
0141   CHECK_CLOSE_ABS(state.pVector[4], unitDir.x(), eps);
0142   CHECK_CLOSE_ABS(state.pVector[5], unitDir.y(), eps);
0143   CHECK_CLOSE_ABS(state.pVector[6], unitDir.z(), eps);
0144   BOOST_CHECK_EQUAL(state.pVector[7], charge / absMom);
0145   BOOST_CHECK_EQUAL(state.pathAccumulated, 0.);
0146   BOOST_CHECK_EQUAL(state.stepSize.value(), stepSize);
0147   BOOST_CHECK_EQUAL(state.previousStepSize, 0.);
0148   BOOST_CHECK_EQUAL(state.tolerance, tolerance);
0149 }
0150
0151 // test stepper getters for particle state
0152 BOOST_AUTO_TEST_CASE(Getters) {
0153   Stepper stepper(magneticField);
0154   Stepper::State state(
0155       geoCtx, magneticField->makeCache(magCtx),
0156       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0157                                  particleHypothesis),
0158       stepSize, tolerance);
0159
0160   CHECK_CLOSE_ABS(stepper.position(state), pos, eps);
0161   CHECK_CLOSE_ABS(stepper.time(state), time, eps);
0162   CHECK_CLOSE_ABS(stepper.direction(state), unitDir, eps);
0163   CHECK_CLOSE_ABS(stepper.absoluteMomentum(state), absMom, eps);
0164   BOOST_CHECK_EQUAL(stepper.charge(state), charge);
0165 }
0166
0167 // test stepper update methods with bound state as input
0168 BOOST_AUTO_TEST_CASE(UpdateFromBound) {
0169   Stepper stepper(magneticField);
0170   Stepper::State state(
0171       geoCtx, magneticField->makeCache(magCtx),
0172       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0173                                  particleHypothesis),
0174       stepSize, tolerance);
0175
0176   auto newPos4 = (pos4 + Vector4(1_mm, 2_mm, 3_mm, 20_ns)).eval();
0177   auto newPos = newPos4.segment<3>(ePos0);
0178   auto newTime = newPos4[eTime];
0179   auto newUnitDir = (unitDir + Vector3(1, -1, -1)).normalized();
0180   auto newAbsMom = 0.9 * absMom;
0181
0182   // example surface and bound parameters at the updated position
0183   auto plane = Surface::makeShared<PlaneSurface>(newPos, newUnitDir);
0184   auto params =
0185       BoundTrackParameters::create(plane, geoCtx, newPos4, newUnitDir,
0186                                    charge / absMom, cov, particleHypothesis)
0187           .value();
0188   FreeVector freeParams;
0189   freeParams[eFreePos0] = newPos4[ePos0];
0190   freeParams[eFreePos1] = newPos4[ePos1];
0191   freeParams[eFreePos2] = newPos4[ePos2];
0192   freeParams[eFreeTime] = newPos4[eTime];
0193   freeParams[eFreeDir0] = newUnitDir[eMom0];
0194   freeParams[eFreeDir1] = newUnitDir[eMom1];
0195   freeParams[eFreeDir2] = newUnitDir[eMom2];
0196   freeParams[eFreeQOverP] = charge / newAbsMom;
0197
0198   // WARNING for some reason there seems to be an additional flag that makes
0199   //         the update method not do anything when it is set. Why?
0200   state.state_ready = false;
0201   BOOST_CHECK(params.covariance().has_value());
0202   stepper.update(state, freeParams, params.parameters(), *params.covariance(),
0203                  *plane);
0204   CHECK_CLOSE_ABS(stepper.position(state), newPos, eps);
0205   CHECK_CLOSE_ABS(stepper.time(state), newTime, eps);
0206   CHECK_CLOSE_ABS(stepper.direction(state), newUnitDir, eps);
0207   CHECK_CLOSE_ABS(stepper.absoluteMomentum(state), newAbsMom, eps);
0208   BOOST_CHECK_EQUAL(stepper.charge(state), charge);
0209 }
0210
0211 // test stepper update methods with individual components as input
0212 BOOST_AUTO_TEST_CASE(UpdateFromComponents) {
0213   Stepper stepper(magneticField);
0214   Stepper::State state(
0215       geoCtx, magneticField->makeCache(magCtx),
0216       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0217                                  particleHypothesis),
0218       stepSize, tolerance);
0219
0220   auto newPos = (pos + Vector3(1_mm, 2_mm, 3_mm)).eval();
0221   auto newTime = time + 20_ns;
0222   auto newUnitDir = (unitDir + Vector3(1, -1, -1)).normalized();
0223   auto newAbsMom = 0.9 * absMom;
0224
0225   stepper.update(state, newPos, newUnitDir, charge / newAbsMom, newTime);
0226   CHECK_CLOSE_ABS(stepper.position(state), newPos, eps);
0227   CHECK_CLOSE_ABS(stepper.time(state), newTime, eps);
0228   CHECK_CLOSE_ABS(stepper.direction(state), newUnitDir, eps);
0229   CHECK_CLOSE_ABS(stepper.absoluteMomentum(state), newAbsMom, eps);
0230   BOOST_CHECK_EQUAL(stepper.charge(state), charge);
0231 }
0232
0233 // test building a bound state object from the stepper state
0234 BOOST_AUTO_TEST_CASE(BuildBound) {
0235   Stepper stepper(magneticField);
0236   Stepper::State state(
0237       geoCtx, magneticField->makeCache(magCtx),
0238       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0239                                  particleHypothesis),
0240       stepSize, tolerance);
0241   // example surface at the current state position
0242   auto plane = Surface::makeShared<PlaneSurface>(pos, unitDir);
0243
0244   auto&& [pars, jac, pathLength] = stepper.boundState(state, *plane).value();
0245   // check parameters
0246   CHECK_CLOSE_ABS(pars.position(geoCtx), pos, eps);
0247   CHECK_CLOSE_ABS(pars.time(), time, eps);
0248   CHECK_CLOSE_ABS(pars.momentum(), absMom * unitDir, eps);
0249   BOOST_CHECK_EQUAL(pars.charge(), charge);
0250   BOOST_CHECK(pars.covariance().has_value());
0251   BOOST_CHECK_NE(*pars.covariance(), cov);
0252   // check Jacobian. should be identity since there was no propagation yet
0253   CHECK_CLOSE_ABS(jac, Jacobian(Jacobian::Identity()), eps);
0254   // check propagation length
0255   CHECK_CLOSE_ABS(pathLength, 0., eps);
0256 }
0257
0258 // test building a curvilinear state object from the stepper state
0259 BOOST_AUTO_TEST_CASE(BuildCurvilinear) {
0260   Stepper stepper(magneticField);
0261   Stepper::State state(
0262       geoCtx, magneticField->makeCache(magCtx),
0263       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0264                                  particleHypothesis),
0265       stepSize, tolerance);
0266
0267   auto&& [pars, jac, pathLength] = stepper.curvilinearState(state);
0268   // check parameters
0269   CHECK_CLOSE_ABS(pars.position(geoCtx), pos, eps);
0270   CHECK_CLOSE_ABS(pars.time(), time, eps);
0271   CHECK_CLOSE_ABS(pars.momentum(), absMom * unitDir, eps);
0272   BOOST_CHECK_EQUAL(pars.charge(), charge);
0273   BOOST_CHECK(pars.covariance().has_value());
0274   BOOST_CHECK_NE(*pars.covariance(), cov);
0275   // check Jacobian. should be identity since there was no propagation yet
0276   CHECK_CLOSE_ABS(jac, Jacobian(Jacobian::Identity()), eps);
0277   // check propagation length
0278   CHECK_CLOSE_ABS(pathLength, 0., eps);
0279 }
0280
0281 // test step method without covariance transport
0282 BOOST_AUTO_TEST_CASE(Step) {
0283   Stepper stepper(magneticField);
0284   MockPropagatorState state(
0285       Stepper::State(geoCtx, magneticField->makeCache(magCtx),
0286                      CurvilinearTrackParameters(pos4, unitDir, charge / absMom,
0287                                                 cov, particleHypothesis),
0288                      stepSize, tolerance));
0289   state.stepping.covTransport = false;
0290
0291   // ensure step does not result in an error
0292   auto res = stepper.step(state, navigator);
0293   BOOST_CHECK(res.ok());
0294
0295   // extract the actual step size
0296   auto h = res.value();
0297   BOOST_CHECK_EQUAL(state.stepping.stepSize.value(), stepSize);
0298   BOOST_CHECK_EQUAL(state.stepping.stepSize.value(), h * navDir);
0299
0300   // check that the position has moved
0301   auto deltaPos = (stepper.position(state.stepping) - pos).eval();
0302   BOOST_CHECK_LT(0, deltaPos.norm());
0303   // check that time has changed
0304   auto deltaTime = stepper.time(state.stepping) - time;
0305   BOOST_CHECK_LT(0, std::abs(deltaTime));
0306   // check that the direction was rotated
0307   auto projDir = stepper.direction(state.stepping).dot(unitDir);
0308   BOOST_CHECK_LT(projDir, 1);
0309
0310   // momentum and charge should not change
0311   CHECK_CLOSE_ABS(stepper.absoluteMomentum(state.stepping), absMom, eps);
0312   BOOST_CHECK_EQUAL(stepper.charge(state.stepping), charge);
0313 }
0314
0315 // test step method with covariance transport
0316 BOOST_AUTO_TEST_CASE(StepWithCovariance) {
0317   Stepper stepper(magneticField);
0318   MockPropagatorState state(
0319       Stepper::State(geoCtx, magneticField->makeCache(magCtx),
0320                      CurvilinearTrackParameters(pos4, unitDir, charge / absMom,
0321                                                 cov, particleHypothesis),
0322                      stepSize, tolerance));
0323   state.stepping.covTransport = true;
0324
0325   // ensure step does not result in an error
0326   auto res = stepper.step(state, navigator);
0327   BOOST_CHECK(res.ok());
0328
0329   // extract the actual step size
0330   auto h = res.value();
0331   BOOST_CHECK_EQUAL(state.stepping.stepSize.value(), stepSize);
0332   BOOST_CHECK_EQUAL(state.stepping.stepSize.value(), h * navDir);
0333
0334   // check that the position has moved
0335   auto deltaPos = (stepper.position(state.stepping) - pos).eval();
0336   BOOST_CHECK_LT(0, deltaPos.norm());
0337   // check that time has changed
0338   auto deltaTime = stepper.time(state.stepping) - time;
0339   BOOST_CHECK_LT(0, std::abs(deltaTime));
0340   // check that the direction was rotated
0341   auto projDir = stepper.direction(state.stepping).dot(unitDir);
0342   BOOST_CHECK_LT(projDir, 1);
0343
0344   // momentum and charge should not change
0345   CHECK_CLOSE_ABS(stepper.absoluteMomentum(state.stepping), absMom, eps);
0346   BOOST_CHECK_EQUAL(stepper.charge(state.stepping), charge);
0347
0348   stepper.transportCovarianceToCurvilinear(state.stepping);
0349   BOOST_CHECK_NE(state.stepping.cov, cov);
0350 }
0351
0352 // test state reset method
0353 BOOST_AUTO_TEST_CASE(Reset) {
0354   Stepper stepper(magneticField);
0355   MockPropagatorState state(
0356       Stepper::State(geoCtx, magneticField->makeCache(magCtx),
0357                      CurvilinearTrackParameters(pos4, unitDir, charge / absMom,
0358                                                 cov, particleHypothesis),
0359                      stepSize, tolerance));
0360   state.stepping.covTransport = true;
0361
0362   // ensure step does not result in an error
0363   stepper.step(state, navigator);
0364
0365   // Construct the parameters
0366   Vector3 newPos(1.5, -2.5, 3.5);
0367   auto newAbsMom = 4.2 * absMom;
0368   double newTime = 7.5;
0369   double newCharge = 1.;
0370   BoundSquareMatrix newCov = 8.5 * Covariance::Identity();
0371   CurvilinearTrackParameters cp(makeVector4(newPos, newTime), unitDir,
0372                                 newCharge / newAbsMom, newCov,
0373                                 particleHypothesis);
0374   FreeVector freeParams = transformBoundToFreeParameters(
0375       cp.referenceSurface(), geoCtx, cp.parameters());
0376   Direction navDir = Direction::Forward;
0377   double stepSize = -256.;
0378
0379   auto copyState = [&](auto& field, const auto& other) {
0380     using field_t = std::decay_t<decltype(field)>;
0381     std::decay_t<decltype(other)> copy(geoCtx, field.makeCache(magCtx), cp,
0382                                        stepSize, tolerance);
0383
0384     copy.state_ready = other.state_ready;
0385     copy.useJacobian = other.useJacobian;
0386     copy.step = other.step;
0387     copy.maxPathLength = other.maxPathLength;
0388     copy.mcondition = other.mcondition;
0389     copy.needgradient = other.needgradient;
0390     copy.newfield = other.newfield;
0391     copy.field = other.field;
0392     copy.pVector = other.pVector;
0393     std::copy(std::begin(other.parameters), std::end(other.parameters),
0394               std::begin(copy.parameters));
0395     copy.covariance = other.covariance;
0396     copy.covTransport = other.covTransport;
0397     std::copy(std::begin(other.jacobian), std::end(other.jacobian),
0398               std::begin(copy.jacobian));
0399     copy.pathAccumulated = other.pathAccumulated;
0400     copy.stepSize = other.stepSize;
0401     copy.previousStepSize = other.previousStepSize;
0402     copy.tolerance = other.tolerance;
0403
0404     copy.fieldCache = MagneticFieldProvider::Cache(
0405         std::in_place_type<typename field_t::Cache>,
0406         other.fieldCache.template as<typename field_t::Cache>());
0407
0408     copy.geoContext = other.geoContext;
0409     copy.debug = other.debug;
0410     copy.debugString = other.debugString;
0411     copy.debugPfxWidth = other.debugPfxWidth;
0412     copy.debugMsgWidth = other.debugMsgWidth;
0413
0414     return copy;
0415   };
0416
0417   // Reset all possible parameters
0418   Stepper::State stateCopy(copyState(*magneticField, state.stepping));
0419   BOOST_CHECK(cp.covariance().has_value());
0420   stepper.resetState(stateCopy, cp.parameters(), *cp.covariance(),
0421                      cp.referenceSurface(), stepSize);
0422   // Test all components
0423   BOOST_CHECK(stateCopy.covTransport);
0424   BOOST_CHECK_EQUAL(*stateCopy.covariance, newCov);
0425   BOOST_CHECK_EQUAL(stepper.position(stateCopy),
0426                     freeParams.template segment<3>(eFreePos0));
0427   BOOST_CHECK_EQUAL(stepper.direction(stateCopy),
0428                     freeParams.template segment<3>(eFreeDir0).normalized());
0429   BOOST_CHECK_EQUAL(stepper.absoluteMomentum(stateCopy),
0430                     std::abs(1. / freeParams[eFreeQOverP]));
0431   BOOST_CHECK_EQUAL(stepper.charge(stateCopy), stepper.charge(state.stepping));
0432   BOOST_CHECK_EQUAL(stepper.time(stateCopy), freeParams[eFreeTime]);
0433   BOOST_CHECK_EQUAL(stateCopy.pathAccumulated, 0.);
0434   BOOST_CHECK_EQUAL(stateCopy.stepSize.value(), navDir * stepSize);
0435   BOOST_CHECK_EQUAL(stateCopy.previousStepSize,
0436                     state.stepping.previousStepSize);
0437   BOOST_CHECK_EQUAL(stateCopy.tolerance, state.stepping.tolerance);
0438
0439   // Reset all possible parameters except the step size
0440   stateCopy = copyState(*magneticField, state.stepping);
0441   stepper.resetState(stateCopy, cp.parameters(), *cp.covariance(),
0442                      cp.referenceSurface());
0443   // Test all components
0444   BOOST_CHECK(stateCopy.covTransport);
0445   BOOST_CHECK_EQUAL(*stateCopy.covariance, newCov);
0446   BOOST_CHECK_EQUAL(stepper.position(stateCopy),
0447                     freeParams.template segment<3>(eFreePos0));
0448   BOOST_CHECK_EQUAL(stepper.direction(stateCopy),
0449                     freeParams.template segment<3>(eFreeDir0).normalized());
0450   BOOST_CHECK_EQUAL(stepper.absoluteMomentum(stateCopy),
0451                     std::abs(1. / freeParams[eFreeQOverP]));
0452   BOOST_CHECK_EQUAL(stepper.charge(stateCopy), stepper.charge(state.stepping));
0453   BOOST_CHECK_EQUAL(stepper.time(stateCopy), freeParams[eFreeTime]);
0454   BOOST_CHECK_EQUAL(stateCopy.pathAccumulated, 0.);
0455   BOOST_CHECK_EQUAL(stateCopy.stepSize.value(),
0456                     std::numeric_limits<double>::max());
0457   BOOST_CHECK_EQUAL(stateCopy.previousStepSize,
0458                     state.stepping.previousStepSize);
0459   BOOST_CHECK_EQUAL(stateCopy.tolerance, state.stepping.tolerance);
0460
0461   // Reset the least amount of parameters
0462   stateCopy = copyState(*magneticField, state.stepping);
0463   stepper.resetState(stateCopy, cp.parameters(), *cp.covariance(),
0464                      cp.referenceSurface());
0465   // Test all components
0466   BOOST_CHECK(stateCopy.covTransport);
0467   BOOST_CHECK_EQUAL(*stateCopy.covariance, newCov);
0468   BOOST_CHECK_EQUAL(stepper.position(stateCopy),
0469                     freeParams.template segment<3>(eFreePos0));
0470   BOOST_CHECK_EQUAL(stepper.direction(stateCopy),
0471                     freeParams.template segment<3>(eFreeDir0).normalized());
0472   BOOST_CHECK_EQUAL(stepper.absoluteMomentum(stateCopy),
0473                     std::abs(1. / freeParams[eFreeQOverP]));
0474   BOOST_CHECK_EQUAL(stepper.charge(stateCopy), stepper.charge(state.stepping));
0475   BOOST_CHECK_EQUAL(stepper.time(stateCopy), freeParams[eFreeTime]);
0476   BOOST_CHECK_EQUAL(stateCopy.pathAccumulated, 0.);
0477   BOOST_CHECK_EQUAL(stateCopy.stepSize.value(),
0478                     std::numeric_limits<double>::max());
0479   BOOST_CHECK_EQUAL(stateCopy.previousStepSize,
0480                     state.stepping.previousStepSize);
0481   BOOST_CHECK_EQUAL(stateCopy.tolerance, state.stepping.tolerance);
0482
0483   // Reset using different surface shapes
0484   // 1) Disc surface
0485   // Setting some parameters
0486   newPos << 0.5, -1.5, 0.;
0487   newAbsMom *= 1.23;
0488   newTime = 8.4;
0489   newCharge = -1.;
0490   newCov = 10.9 * Covariance::Identity();
0491   Transform3 trafo = Transform3::Identity();
0492   auto disc = Surface::makeShared<DiscSurface>(trafo);
0493   auto boundDisc = BoundTrackParameters::create(
0494                        disc, geoCtx, makeVector4(newPos, newTime), unitDir,
0495                        newCharge / newAbsMom, newCov, particleHypothesis)
0496                        .value();
0497
0498   // Reset the state and test
0499   Stepper::State stateDisc = copyState(*magneticField, state.stepping);
0500   BOOST_CHECK(boundDisc.covariance().has_value());
0501   stepper.resetState(stateDisc, boundDisc.parameters(), *boundDisc.covariance(),
0502                      boundDisc.referenceSurface());
0503
0504   CHECK_NE_COLLECTIONS(stateDisc.pVector, stateCopy.pVector);
0505   CHECK_NE_COLLECTIONS(stateDisc.pVector, state.stepping.pVector);
0506
0507   // 2) Perigee surface
0508   // Setting some parameters
0509   newPos << -2.06155, -2.06155, 3.5;
0510   newAbsMom *= 0.45;
0511   newTime = 2.3;
0512   newCharge = 1.;
0513   newCov = 8.7 * Covariance::Identity();
0514   auto perigee = Surface::makeShared<PerigeeSurface>(trafo);
0515   auto boundPerigee =
0516       BoundTrackParameters::create(
0517           perigee, geoCtx, makeVector4(newPos, newTime), unitDir,
0518           newCharge / newAbsMom, newCov, particleHypothesis)
0519           .value();
0520
0521   // Reset the state and test
0522   Stepper::State statePerigee = copyState(*magneticField, state.stepping);
0523   BOOST_CHECK(boundPerigee.covariance().has_value());
0524   stepper.resetState(statePerigee, boundPerigee.parameters(),
0525                      *boundPerigee.covariance(),
0526                      boundPerigee.referenceSurface());
0527   CHECK_NE_COLLECTIONS(statePerigee.pVector, stateCopy.pVector);
0528   CHECK_NE_COLLECTIONS(statePerigee.pVector, state.stepping.pVector);
0529   CHECK_NE_COLLECTIONS(statePerigee.pVector, stateDisc.pVector);
0530
0531   // 3) Straw surface
0532   // Use the same parameters as for previous Perigee surface
0533   auto straw = Surface::makeShared<StrawSurface>(trafo);
0534   auto boundStraw = BoundTrackParameters::create(
0535                         straw, geoCtx, makeVector4(newPos, newTime), unitDir,
0536                         newCharge / newAbsMom, newCov, particleHypothesis)
0537                         .value();
0538
0539   // Reset the state and test
0540   Stepper::State stateStraw = copyState(*magneticField, state.stepping);
0541   BOOST_CHECK(boundStraw.covariance().has_value());
0542   stepper.resetState(stateStraw, boundStraw.parameters(),
0543                      *boundStraw.covariance(), boundStraw.referenceSurface());
0544   CHECK_NE_COLLECTIONS(stateStraw.pVector, stateCopy.pVector);
0545   CHECK_NE_COLLECTIONS(stateStraw.pVector, state.stepping.pVector);
0546   CHECK_NE_COLLECTIONS(stateStraw.pVector, stateDisc.pVector);
0547   BOOST_CHECK_EQUAL_COLLECTIONS(
0548       stateStraw.pVector.begin(), stateStraw.pVector.end(),
0549       statePerigee.pVector.begin(), statePerigee.pVector.end());
0550 }
0551
0552 BOOST_AUTO_TEST_CASE(StepSize) {
0553   Stepper stepper(magneticField);
0554   Stepper::State state(
0555       geoCtx, magneticField->makeCache(magCtx),
0556       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0557                                  particleHypothesis),
0558       stepSize, tolerance);
0559
0560   // TODO figure out why this fails and what it should be
0561   // BOOST_CHECK_EQUAL(stepper.overstepLimit(state), tolerance);
0562
0563   stepper.updateStepSize(state, -5_cm, ConstrainedStep::actor);
0564   BOOST_CHECK_EQUAL(state.previousStepSize, stepSize);
0565   BOOST_CHECK_EQUAL(state.stepSize.value(), -5_cm);
0566
0567   stepper.releaseStepSize(state, ConstrainedStep::actor);
0568   BOOST_CHECK_EQUAL(state.stepSize.value(), stepSize);
0569 }
0570
0571 // test step size modification with target surfaces
0572 BOOST_AUTO_TEST_CASE(StepSizeSurface) {
0573   Stepper stepper(magneticField);
0574   Stepper::State state(
0575       geoCtx, magneticField->makeCache(magCtx),
0576       CurvilinearTrackParameters(pos4, unitDir, charge / absMom, cov,
0577                                  particleHypothesis),
0578       stepSize, tolerance);
0579
0580   auto distance = 10_mm;
0581   auto target = Surface::makeShared<PlaneSurface>(
0582       pos + navDir * distance * unitDir, unitDir);
0583
0584   stepper.updateSurfaceStatus(state, *target, 0, navDir, BoundaryCheck(false));
0585   BOOST_CHECK_EQUAL(state.stepSize.value(ConstrainedStep::actor), distance);
0586
0587   // test the step size modification in the context of a surface
0588   stepper.updateStepSize(
0589       state,
0590       target
0591           ->intersect(state.geoContext, stepper.position(state),
0592                       navDir * stepper.direction(state), BoundaryCheck(false))
0593           .closest(),
0594       navDir, false);
0595   BOOST_CHECK_EQUAL(state.stepSize.value(), distance);
0596
0597   // start with a different step size
0598   state.stepSize.setUser(navDir * stepSize);
0599   stepper.updateStepSize(
0600       state,
0601       target
0602           ->intersect(state.geoContext, stepper.position(state),
0603                       navDir * stepper.direction(state), BoundaryCheck(false))
0604           .closest(),
0605       navDir, true);
0606   BOOST_CHECK_EQUAL(state.stepSize.value(), navDir * stepSize);
0607 }
0608
0609 BOOST_AUTO_TEST_SUITE_END()
0610
0611 }  // namespace Acts::Test