Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-11-12 09:43:43

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/.
0008 
0009 #include <boost/test/unit_test.hpp>
0010 
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Definitions/Units.hpp"
0013 #include "Acts/EventData/Types.hpp"
0014 #include "Acts/EventData/VectorMultiTrajectory.hpp"
0015 #include "Acts/EventData/VectorTrackContainer.hpp"
0016 #include "Acts/EventData/detail/TestSourceLink.hpp"
0017 #include "Acts/Geometry/CuboidVolumeBuilder.hpp"
0018 #include "Acts/Geometry/TrackingGeometry.hpp"
0019 #include "Acts/Geometry/TrackingGeometryBuilder.hpp"
0020 #include "Acts/Material/HomogeneousSurfaceMaterial.hpp"
0021 #include "Acts/Material/MaterialSlab.hpp"
0022 #include "Acts/Propagator/EigenStepper.hpp"
0023 #include "Acts/Propagator/Navigator.hpp"
0024 #include "Acts/Propagator/Propagator.hpp"
0025 #include "Acts/Propagator/StraightLineStepper.hpp"
0026 #include "Acts/Surfaces/PlaneSurface.hpp"
0027 #include "Acts/Surfaces/RectangleBounds.hpp"
0028 #include "Acts/TrackFitting/GlobalChiSquareFitter.hpp"
0029 #include "Acts/Utilities/Logger.hpp"
0030 #include "Acts/Visualization/EventDataView3D.hpp"
0031 #include "Acts/Visualization/GeometryView3D.hpp"
0032 #include "Acts/Visualization/ObjVisualization3D.hpp"
0033 #include "ActsTests/CommonHelpers/DetectorElementStub.hpp"
0034 #include "ActsTests/CommonHelpers/MeasurementsCreator.hpp"
0035 #include "ActsTests/CommonHelpers/PredefinedMaterials.hpp"
0036 
0037 #include <numbers>
0038 #include <vector>
0039 
0040 #include "FitterTestsCommon.hpp"
0041 
0042 using namespace Acts;
0043 using namespace Acts::Experimental;
0044 using namespace Acts::UnitLiterals;
0045 using namespace Acts::detail::Test;
0046 
0047 Logging::Level logLevel = Logging::VERBOSE;
0048 const auto gx2fLogger = getDefaultLogger("Gx2f", logLevel);
0049 
0050 namespace ActsTests {
0051 
0052 /// @brief Helper function to visualise measurements in a 3D environment.
0053 ///
0054 /// This function iterates through the provided measurements and visualises each
0055 /// one using the specified 3D visualisation helper. The visualisation takes
0056 /// into account the surface transformations and localisation errors.
0057 ///
0058 /// @param helper The 3D visualisation helper used to draw the measurements.
0059 /// @param measurements A collection of measurements to be visualised, containing source links with parameters and covariance information.
0060 /// @param geometry A shared pointer to the constant tracking geometry used to find surfaces associated with measurements.
0061 /// @param geoCtx The geometry context used for transformations and accessing geometry-related information.
0062 /// @param locErrorScale Scaling factor for localisation errors. Default value is 1.0.
0063 /// @param viewConfig Configuration settings for the visualisation view. Default value is s_viewMeasurement.
0064 static void drawMeasurements(
0065     IVisualization3D& helper, const Measurements& measurements,
0066     const std::shared_ptr<const TrackingGeometry>& geometry,
0067     const GeometryContext& geoCtx, double locErrorScale = 1.,
0068     const ViewConfig& viewConfig = s_viewMeasurement) {
0069   std::cout << "\n*** Draw measurements ***\n" << std::endl;
0070 
0071   for (auto& singleMeasurement : measurements.sourceLinks) {
0072     auto cov = singleMeasurement.covariance;
0073     auto lposition = singleMeasurement.parameters;
0074 
0075     auto surf = geometry->findSurface(singleMeasurement.m_geometryId);
0076     auto transf = surf->localToGlobalTransform(geoCtx);
0077 
0078     EventDataView3D::drawMeasurement(helper, lposition, cov, transf,
0079                                      locErrorScale, viewConfig);
0080   }
0081 }
0082 
0083 //// Construct initial track parameters.
0084 BoundTrackParameters makeParameters(
0085     const double x = 0.0_m, const double y = 0.0_m, const double z = 0.0_m,
0086     const double w = 42_ns, const double phi = 0_degree,
0087     const double theta = 90_degree, const double p = 2_GeV,
0088     const double q = 1_e,
0089     const ParticleHypothesis& particleHypothesis = ParticleHypothesis::pion()) {
0090   // create covariance matrix from reasonable standard deviations
0091   BoundVector stddev;
0092   stddev[eBoundLoc0] = 100_um;
0093   stddev[eBoundLoc1] = 100_um;
0094   stddev[eBoundTime] = 25_ns;
0095   stddev[eBoundPhi] = 2_degree;
0096   stddev[eBoundTheta] = 2_degree;
0097   stddev[eBoundQOverP] = 1 / 100_GeV;
0098   const BoundMatrix cov = stddev.cwiseProduct(stddev).asDiagonal();
0099   // define a track in the transverse plane along x
0100   const Vector4 mPos4(x, y, z, w);
0101   return BoundTrackParameters::createCurvilinear(mPos4, phi, theta, q / p, cov,
0102                                                  particleHypothesis);
0103 }
0104 
0105 static std::vector<SourceLink> prepareSourceLinks(
0106     const std::vector<TestSourceLink>& sourceLinks) {
0107   std::vector<SourceLink> result;
0108   std::transform(sourceLinks.begin(), sourceLinks.end(),
0109                  std::back_inserter(result),
0110                  [](const auto& sl) { return SourceLink{sl}; });
0111   return result;
0112 }
0113 
0114 /// @brief Create a simple telescope detector
0115 ///
0116 /// @param geoCtx
0117 /// @param nSurfaces Number of surfaces
0118 /// @param surfaceIndexWithMaterial A set of index of the material surfaces
0119 std::shared_ptr<const TrackingGeometry> makeToyDetector(
0120     const GeometryContext& geoCtx, const std::size_t nSurfaces = 5,
0121     const std::set<std::size_t>& surfaceIndexWithMaterial = {}) {
0122   if (nSurfaces < 1) {
0123     throw std::invalid_argument("At least 1 surfaces needs to be created.");
0124   }
0125 
0126   // Define the dimensions of the square surfaces
0127   const double halfSizeSurface = 1_m;
0128 
0129   // Rotation of the surfaces around the y-axis
0130   const double rotationAngle = std::numbers::pi / 2.;
0131   const Vector3 xPos(std::cos(rotationAngle), 0., std::sin(rotationAngle));
0132   const Vector3 yPos(0., 1., 0.);
0133   const Vector3 zPos(-std::sin(rotationAngle), 0., std::cos(rotationAngle));
0134 
0135   // Construct builder
0136   CuboidVolumeBuilder cvb;
0137 
0138   // Create configurations for surfaces
0139   std::vector<CuboidVolumeBuilder::SurfaceConfig> surfaceConfig;
0140   for (std::size_t surfPos = 1; surfPos <= nSurfaces; surfPos++) {
0141     // Position of the surfaces
0142     CuboidVolumeBuilder::SurfaceConfig cfg;
0143     cfg.position = {surfPos * UnitConstants::m, 0., 0.};
0144 
0145     // Rotation of the surfaces
0146     cfg.rotation.col(0) = xPos;
0147     cfg.rotation.col(1) = yPos;
0148     cfg.rotation.col(2) = zPos;
0149 
0150     // Boundaries of the surfaces (shape)
0151     cfg.rBounds = std::make_shared<const RectangleBounds>(
0152         RectangleBounds(halfSizeSurface, halfSizeSurface));
0153 
0154     // Add material only for selected surfaces
0155     if (surfaceIndexWithMaterial.contains(surfPos)) {
0156       // Material of the surfaces
0157       MaterialSlab matProp(makeSilicon(), 5_mm);
0158       cfg.surMat = std::make_shared<HomogeneousSurfaceMaterial>(matProp);
0159     }
0160 
0161     // Thickness of the detector element
0162     cfg.thickness = 1_um;
0163 
0164     cfg.detElementConstructor =
0165         [](const Transform3& trans,
0166            const std::shared_ptr<const RectangleBounds>& bounds,
0167            double thickness) {
0168           return new DetectorElementStub(trans, bounds, thickness);
0169         };
0170     surfaceConfig.push_back(cfg);
0171   }
0172 
0173   // Build layer configurations
0174   std::vector<CuboidVolumeBuilder::LayerConfig> layerConfig;
0175   for (auto& sCfg : surfaceConfig) {
0176     CuboidVolumeBuilder::LayerConfig cfg;
0177     cfg.surfaceCfg = {sCfg};
0178     cfg.active = true;
0179     cfg.envelopeX = {-0.1_mm, 0.1_mm};
0180     cfg.envelopeY = {-0.1_mm, 0.1_mm};
0181     cfg.envelopeZ = {-0.1_mm, 0.1_mm};
0182     layerConfig.push_back(cfg);
0183   }
0184 
0185   // Inner Volume - Build volume configuration
0186   CuboidVolumeBuilder::VolumeConfig volumeConfig;
0187   volumeConfig.length = {(nSurfaces + 1) * 1_m, 2 * halfSizeSurface,
0188                          2 * halfSizeSurface};
0189   volumeConfig.position = {volumeConfig.length.x() / 2, 0., 0.};
0190   volumeConfig.layerCfg = layerConfig;
0191   volumeConfig.name = "TestVolume";
0192 
0193   // Outer volume - Build TrackingGeometry configuration
0194   CuboidVolumeBuilder::Config config;
0195   config.length = {(nSurfaces + 1) * 1_m, 2 * halfSizeSurface,
0196                    2 * halfSizeSurface};
0197   config.position = {volumeConfig.length.x() / 2, 0., 0.};
0198   config.volumeCfg = {volumeConfig};
0199 
0200   cvb.setConfig(config);
0201 
0202   TrackingGeometryBuilder::Config tgbCfg;
0203 
0204   tgbCfg.trackingVolumeBuilders.push_back(
0205       [=](const auto& context, const auto& inner, const auto&) {
0206         return cvb.trackingVolume(context, inner, nullptr);
0207       });
0208 
0209   TrackingGeometryBuilder tgb(tgbCfg);
0210 
0211   std::unique_ptr<const TrackingGeometry> detector =
0212       tgb.trackingGeometry(geoCtx);
0213   return detector;
0214 }
0215 
0216 struct Detector {
0217   // geometry
0218   std::shared_ptr<const TrackingGeometry> geometry;
0219 };
0220 
0221 BOOST_AUTO_TEST_SUITE(TrackFittingSuite)
0222 
0223 ACTS_LOCAL_LOGGER(getDefaultLogger("Gx2fTests", logLevel))
0224 
0225 // Context objects
0226 const auto geoCtx = GeometryContext::dangerouslyDefaultConstruct();
0227 const MagneticFieldContext magCtx;
0228 const CalibrationContext calCtx;
0229 
0230 // Measurement resolutions
0231 const MeasurementResolution resPixel = {MeasurementType::eLoc01,
0232                                         {25_um, 50_um}};
0233 const MeasurementResolution resStrip0 = {MeasurementType::eLoc0, {25_um}};
0234 const MeasurementResolution resStrip1 = {MeasurementType::eLoc1, {50_um}};
0235 const MeasurementResolutionMap resMapAllPixel = {
0236     {GeometryIdentifier().withVolume(0), resPixel}};
0237 
0238 // This test checks if the call to the fitter works and no errors occur in the
0239 // framework, without fitting and updating any parameters
0240 BOOST_AUTO_TEST_CASE(NoFit) {
0241   ACTS_INFO("*** Test: NoFit -- Start");
0242 
0243   std::default_random_engine rng(42);
0244 
0245   ACTS_DEBUG("Create the detector");
0246   const std::size_t nSurfaces = 5;
0247   Detector detector;
0248   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0249 
0250   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0251   const auto parametersMeasurements = makeParameters();
0252   const auto startParametersFit = makeParameters(
0253       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0254 
0255   ACTS_DEBUG("Create the measurements");
0256   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0257   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0258   const auto measurements =
0259       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0260                          resMapAllPixel, rng);
0261   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0262   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0263 
0264   ACTS_DEBUG("Set up the fitter");
0265   // Reuse the SimPropagator, since we will not actually use it
0266   using Gx2Fitter = Gx2Fitter<SimPropagator, VectorMultiTrajectory>;
0267   const Gx2Fitter fitter(simPropagator, gx2fLogger->clone());
0268 
0269   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0270 
0271   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0272   extensions.calibrator
0273       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0274   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0275   extensions.surfaceAccessor
0276       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0277 
0278   Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0279                                PropagatorPlainOptions(geoCtx, magCtx), rSurface,
0280                                false, false, FreeToBoundCorrection(false), 0,
0281                                0);
0282 
0283   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0284 
0285   ACTS_DEBUG("Fit the track");
0286   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0287   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0288   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0289                               startParametersFit, gx2fOptions, tracks);
0290 
0291   BOOST_REQUIRE(res.ok());
0292 
0293   const auto& track = *res;
0294   BOOST_CHECK_EQUAL(track.tipIndex(), kTrackIndexInvalid);
0295   BOOST_CHECK(track.hasReferenceSurface());
0296 
0297   // Track quantities
0298   BOOST_CHECK_EQUAL(track.chi2(), 0.);
0299   BOOST_CHECK_EQUAL(track.nDoF(), 0u);
0300   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
0301   BOOST_CHECK_EQUAL(track.nMeasurements(), 0u);
0302   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
0303   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
0304 
0305   // Parameters
0306   BOOST_CHECK_EQUAL(track.parameters(), startParametersFit.parameters());
0307   BOOST_CHECK_EQUAL(track.covariance(), BoundMatrix::Identity());
0308 
0309   // Convergence
0310   BOOST_CHECK_EQUAL(
0311       (track.template component<
0312           std::uint32_t, hashString(Gx2fConstants::gx2fnUpdateColumn)>()),
0313       0);
0314 
0315   ACTS_INFO("*** Test: NoFit -- Finish");
0316 }
0317 
0318 BOOST_AUTO_TEST_CASE(Fit5Iterations) {
0319   ACTS_INFO("*** Test: Fit5Iterations -- Start");
0320 
0321   std::default_random_engine rng(42);
0322 
0323   ACTS_DEBUG("Create the detector");
0324   const std::size_t nSurfaces = 5;
0325   Detector detector;
0326   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0327 
0328   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0329   const auto parametersMeasurements = makeParameters();
0330   // Non-default hypothesis to check that it is propagated to the output track
0331   const auto startParametersFit =
0332       makeParameters(7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV,
0333                      1_e, ParticleHypothesis::muon());
0334 
0335   ACTS_DEBUG("Create the measurements");
0336   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0337   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0338   const auto measurements =
0339       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0340                          resMapAllPixel, rng);
0341   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0342   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0343 
0344   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0345 
0346   ACTS_DEBUG("Set up the fitter");
0347   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0348 
0349   using RecoStepper = EigenStepper<>;
0350   const auto recoPropagator =
0351       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0352 
0353   using RecoPropagator = decltype(recoPropagator);
0354   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
0355   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
0356 
0357   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0358   extensions.calibrator
0359       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0360   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0361   extensions.surfaceAccessor
0362       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0363 
0364   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0365                                      PropagatorPlainOptions(geoCtx, magCtx),
0366                                      rSurface, false, false,
0367                                      FreeToBoundCorrection(false), 5, 0);
0368 
0369   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0370 
0371   ACTS_DEBUG("Fit the track");
0372   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0373   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0374   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0375                               startParametersFit, gx2fOptions, tracks);
0376 
0377   BOOST_REQUIRE(res.ok());
0378 
0379   const auto& track = *res;
0380 
0381   BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
0382   BOOST_CHECK(track.hasReferenceSurface());
0383 
0384   // Track quantities
0385   CHECK_CLOSE_ABS(track.chi2(), 8., 2.);
0386   BOOST_CHECK_EQUAL(track.nDoF(), nSurfaces * 2);
0387   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
0388   BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
0389   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
0390   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
0391   BOOST_CHECK_EQUAL(track.particleHypothesis(),
0392                     startParametersFit.particleHypothesis());
0393 
0394   // Parameters
0395   // We need quite coarse checks here, since on different builds
0396   // the created measurements differ in the randomness
0397   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 7e0);
0398   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 6e0);
0399   BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-5, 1e3);
0400   BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], std::numbers::pi / 2,
0401                     1e-3);
0402   BOOST_CHECK_EQUAL(track.parameters()[eBoundQOverP], 1);
0403   BOOST_CHECK_CLOSE(track.parameters()[eBoundTime],
0404                     startParametersFit.parameters()[eBoundTime], 1e-6);
0405   BOOST_CHECK_CLOSE(track.covariance().determinant(), 1e-27, 4e0);
0406 
0407   // Convergence
0408   BOOST_CHECK_EQUAL(
0409       (track.template component<
0410           std::uint32_t, hashString(Gx2fConstants::gx2fnUpdateColumn)>()),
0411       4);
0412 
0413   ACTS_INFO("*** Test: Fit5Iterations -- Finish");
0414 }
0415 
0416 // This test checks that the fitter reaches a measurement on a free surface,
0417 // which is not part of the tracking geometry
0418 BOOST_AUTO_TEST_CASE(FreeMeasurementSurface) {
0419   ACTS_INFO("*** Test: FreeMeasurementSurface -- Start");
0420 
0421   std::default_random_engine rng(42);
0422 
0423   ACTS_DEBUG("Create the detector");
0424   const std::size_t nSurfaces = 5;
0425   Detector detector;
0426   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0427 
0428   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0429   const auto parametersMeasurements = makeParameters();
0430   const auto startParametersFit = makeParameters(
0431       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0432 
0433   ACTS_DEBUG("Create the measurements");
0434   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0435   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0436   const auto measurements =
0437       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0438                          resMapAllPixel, rng);
0439   auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0440   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0441 
0442   ACTS_DEBUG("Add a measurement on a free surface between two layers");
0443   const Surface* firstSurface = detector.geometry->findSurface(
0444       measurements.sourceLinks.front().m_geometryId);
0445   BOOST_REQUIRE(firstSurface != nullptr);
0446   Transform3 freeTransform = Transform3::Identity();
0447   freeTransform.linear() =
0448       firstSurface->localToGlobalTransform(geoCtx).linear();
0449   freeTransform.translation() = Vector3(2.5_m, 0., 0.);
0450   const auto freeSurface = Surface::makeShared<PlaneSurface>(
0451       freeTransform, std::make_shared<RectangleBounds>(1_m, 1_m));
0452   BOOST_REQUIRE_EQUAL(freeSurface->geometryId(), GeometryIdentifier{});
0453   sourceLinks.emplace_back(
0454       TestSourceLink(eBoundLoc0, eBoundLoc1, Vector2::Zero(),
0455                      Vector2(25_um * 25_um, 50_um * 50_um).asDiagonal()));
0456 
0457   // Resolve the source link without a geometry id to the free surface
0458   struct FreeSurfaceAccessor {
0459     const TrackingGeometry* geometry{};
0460     const Surface* freeSurface{};
0461 
0462     const Surface* operator()(const SourceLink& sourceLink) const {
0463       const GeometryIdentifier geoId =
0464           sourceLink.get<TestSourceLink>().m_geometryId;
0465       return geoId == GeometryIdentifier{} ? freeSurface
0466                                            : geometry->findSurface(geoId);
0467     }
0468   };
0469   FreeSurfaceAccessor surfaceAccessor{detector.geometry.get(),
0470                                       freeSurface.get()};
0471 
0472   ACTS_DEBUG("Set up the fitter");
0473   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0474 
0475   using RecoStepper = EigenStepper<>;
0476   const auto recoPropagator =
0477       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0478 
0479   using RecoPropagator = decltype(recoPropagator);
0480   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
0481   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
0482 
0483   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0484   extensions.calibrator
0485       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0486   extensions.surfaceAccessor.connect<&FreeSurfaceAccessor::operator()>(
0487       &surfaceAccessor);
0488 
0489   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0490                                      PropagatorPlainOptions(geoCtx, magCtx),
0491                                      rSurface, false, false,
0492                                      FreeToBoundCorrection(false), 5, 0);
0493 
0494   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0495 
0496   ACTS_DEBUG("Fit the track");
0497   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0498                               startParametersFit, gx2fOptions, tracks);
0499 
0500   BOOST_REQUIRE(res.ok());
0501 
0502   const auto& track = *res;
0503   BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces + 1);
0504   BOOST_CHECK_EQUAL(track.nDoF(), (nSurfaces + 1) * 2);
0505   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
0506   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
0507 
0508   // The free surface carries a measurement state
0509   bool freeSurfaceMeasured = false;
0510   for (const auto trackState : track.trackStatesReversed()) {
0511     if (&trackState.referenceSurface() == freeSurface.get()) {
0512       freeSurfaceMeasured = trackState.typeFlags().isMeasurement();
0513     }
0514   }
0515   BOOST_CHECK(freeSurfaceMeasured);
0516 
0517   ACTS_INFO("*** Test: FreeMeasurementSurface -- Finish");
0518 }
0519 
0520 BOOST_AUTO_TEST_CASE(MixedDetector) {
0521   ACTS_INFO("*** Test: MixedDetector -- Start");
0522 
0523   std::default_random_engine rng(42);
0524 
0525   ACTS_DEBUG("Create the detector");
0526   const std::size_t nSurfaces = 7;
0527   Detector detector;
0528   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0529 
0530   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0531   const auto parametersMeasurements = makeParameters();
0532   const auto startParametersFit = makeParameters(
0533       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0534 
0535   ACTS_DEBUG("Create the measurements");
0536   const MeasurementResolutionMap resMap = {
0537       {GeometryIdentifier().withVolume(2).withLayer(2), resPixel},
0538       {GeometryIdentifier().withVolume(2).withLayer(4), resStrip0},
0539       {GeometryIdentifier().withVolume(2).withLayer(6), resStrip1},
0540       {GeometryIdentifier().withVolume(2).withLayer(8), resPixel},
0541       {GeometryIdentifier().withVolume(2).withLayer(10), resStrip0},
0542       {GeometryIdentifier().withVolume(2).withLayer(12), resStrip1},
0543       {GeometryIdentifier().withVolume(2).withLayer(14), resPixel},
0544   };
0545 
0546   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0547   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0548   const auto measurements = createMeasurements(
0549       simPropagator, geoCtx, magCtx, parametersMeasurements, resMap, rng);
0550   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0551   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0552 
0553   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0554 
0555   ACTS_DEBUG("Set up the fitter");
0556   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0557 
0558   using RecoStepper = EigenStepper<>;
0559   const auto recoPropagator =
0560       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0561 
0562   using RecoPropagator = decltype(recoPropagator);
0563   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
0564   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
0565 
0566   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0567   extensions.calibrator
0568       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0569   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0570   extensions.surfaceAccessor
0571       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0572 
0573   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0574                                      PropagatorPlainOptions(geoCtx, magCtx),
0575                                      rSurface, false, false,
0576                                      FreeToBoundCorrection(false), 5, 0);
0577 
0578   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0579 
0580   ACTS_DEBUG("Fit the track");
0581   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0582   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0583   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0584                               startParametersFit, gx2fOptions, tracks);
0585 
0586   BOOST_REQUIRE(res.ok());
0587 
0588   const auto& track = *res;
0589   BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
0590   BOOST_CHECK(track.hasReferenceSurface());
0591 
0592   // Track quantities
0593   CHECK_CLOSE_ABS(track.chi2(), 8.5, 4.);
0594   BOOST_CHECK_EQUAL(track.nDoF(), 10u);
0595   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
0596   BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
0597   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
0598   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
0599 
0600   // Parameters
0601   // We need quite coarse checks here, since on different builds
0602   // the created measurements differ in the randomness
0603   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 7e0);
0604   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 6e0);
0605   BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-5, 1e3);
0606   BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], std::numbers::pi / 2,
0607                     1e-3);
0608   BOOST_CHECK_EQUAL(track.parameters()[eBoundQOverP], 1);
0609   BOOST_CHECK_CLOSE(track.parameters()[eBoundTime],
0610                     startParametersFit.parameters()[eBoundTime], 1e-6);
0611   BOOST_CHECK_CLOSE(track.covariance().determinant(), 2e-28, 1e0);
0612 
0613   // Convergence
0614   BOOST_CHECK_EQUAL(
0615       (track.template component<
0616           std::uint32_t, hashString(Gx2fConstants::gx2fnUpdateColumn)>()),
0617       4);
0618 
0619   ACTS_INFO("*** Test: MixedDetector -- Finish");
0620 }
0621 
0622 // This test checks if we can fit QOverP, when a magnetic field is introduced
0623 BOOST_AUTO_TEST_CASE(FitWithBfield) {
0624   ACTS_INFO("*** Test: FitWithBfield -- Start");
0625 
0626   std::default_random_engine rng(42);
0627 
0628   ACTS_DEBUG("Create the detector");
0629   const std::size_t nSurfaces = 5;
0630   Detector detector;
0631   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0632 
0633   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0634   const auto parametersMeasurements = makeParameters();
0635   const auto startParametersFit = makeParameters(
0636       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0637 
0638   ACTS_DEBUG("Create the measurements");
0639   using SimStepper = EigenStepper<>;
0640   const auto simPropagator =
0641       makeConstantFieldPropagator<SimStepper>(detector.geometry, 0.3_T);
0642 
0643   const auto measurements =
0644       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0645                          resMapAllPixel, rng);
0646 
0647   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0648   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0649 
0650   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0651 
0652   ACTS_DEBUG("Set up the fitter");
0653   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0654 
0655   // Reuse the SimPropagator, since it already uses the EigenStepper<>
0656   using SimPropagator = decltype(simPropagator);
0657   using Gx2Fitter = Gx2Fitter<SimPropagator, VectorMultiTrajectory>;
0658   const Gx2Fitter fitter(simPropagator, gx2fLogger->clone());
0659 
0660   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0661   extensions.calibrator
0662       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0663   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0664   extensions.surfaceAccessor
0665       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0666 
0667   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0668                                      PropagatorPlainOptions(geoCtx, magCtx),
0669                                      rSurface, false, false,
0670                                      FreeToBoundCorrection(false), 5, 0);
0671 
0672   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0673 
0674   ACTS_DEBUG("Fit the track");
0675   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0676   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0677   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0678                               startParametersFit, gx2fOptions, tracks);
0679 
0680   BOOST_REQUIRE(res.ok());
0681 
0682   const auto& track = *res;
0683   BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
0684   BOOST_CHECK(track.hasReferenceSurface());
0685 
0686   // Track quantities
0687   CHECK_CLOSE_ABS(track.chi2(), 7.5, 1.5);
0688   BOOST_CHECK_EQUAL(track.nDoF(), nSurfaces * 2);
0689   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
0690   BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
0691   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
0692   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
0693 
0694   // Parameters
0695   // We need quite coarse checks here, since on different builds
0696   // the created measurements differ in the randomness
0697   // TODO investigate further the reference values for eBoundPhi and det(cov)
0698   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 8e0);
0699   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 6e0);
0700   BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-4, 1e3);
0701   BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], std::numbers::pi / 2,
0702                     1e-3);
0703   BOOST_CHECK_CLOSE(track.parameters()[eBoundQOverP], 0.5, 2e-1);
0704   BOOST_CHECK_CLOSE(track.parameters()[eBoundTime],
0705                     startParametersFit.parameters()[eBoundTime], 1e-6);
0706   BOOST_CHECK_CLOSE(track.covariance().determinant(), 8e-35, 4e0);
0707 
0708   // Convergence
0709   BOOST_CHECK_EQUAL(
0710       (track.template component<
0711           std::uint32_t, hashString(Gx2fConstants::gx2fnUpdateColumn)>()),
0712       4);
0713 
0714   ACTS_INFO("*** Test: FitWithBfield -- Finish");
0715 }
0716 
0717 BOOST_AUTO_TEST_CASE(relChi2changeCutOff) {
0718   ACTS_INFO("*** Test: relChi2changeCutOff -- Start");
0719 
0720   std::default_random_engine rng(42);
0721 
0722   ACTS_DEBUG("Create the detector");
0723   const std::size_t nSurfaces = 5;
0724   Detector detector;
0725   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0726 
0727   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0728   const auto parametersMeasurements = makeParameters();
0729   const auto startParametersFit = makeParameters(
0730       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0731 
0732   ACTS_DEBUG("Create the measurements");
0733   // simulation propagator
0734   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0735   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0736   const auto measurements =
0737       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0738                          resMapAllPixel, rng);
0739   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0740   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0741 
0742   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0743 
0744   ACTS_DEBUG("Set up the fitter");
0745   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0746 
0747   using RecoStepper = EigenStepper<>;
0748   const auto recoPropagator =
0749       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0750 
0751   using RecoPropagator = decltype(recoPropagator);
0752   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
0753   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
0754 
0755   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0756   extensions.calibrator
0757       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0758   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0759   extensions.surfaceAccessor
0760       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0761 
0762   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0763                                      PropagatorPlainOptions(geoCtx, magCtx),
0764                                      rSurface, false, false,
0765                                      FreeToBoundCorrection(false), 500, 1e-5);
0766 
0767   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0768 
0769   ACTS_DEBUG("Fit the track");
0770   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0771   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0772   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0773                               startParametersFit, gx2fOptions, tracks);
0774 
0775   BOOST_REQUIRE(res.ok());
0776 
0777   const auto& track = *res;
0778   BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
0779   BOOST_CHECK(track.hasReferenceSurface());
0780 
0781   // Track quantities
0782   CHECK_CLOSE_ABS(track.chi2(), 8., 2.);
0783   BOOST_CHECK_EQUAL(track.nDoF(), nSurfaces * 2);
0784   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
0785   BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
0786   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
0787   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
0788 
0789   // Parameters
0790   // We need quite coarse checks here, since on different builds
0791   // the created measurements differ in the randomness
0792   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 7e0);
0793   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 6e0);
0794   BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-5, 1e3);
0795   BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], std::numbers::pi / 2,
0796                     1e-3);
0797   BOOST_CHECK_EQUAL(track.parameters()[eBoundQOverP], 1);
0798   BOOST_CHECK_CLOSE(track.parameters()[eBoundTime],
0799                     startParametersFit.parameters()[eBoundTime], 1e-6);
0800   BOOST_CHECK_CLOSE(track.covariance().determinant(), 1e-27, 4e0);
0801 
0802   // Convergence
0803   BOOST_CHECK_LT(
0804       (track.template component<
0805           std::uint32_t, hashString(Gx2fConstants::gx2fnUpdateColumn)>()),
0806       10);
0807 
0808   ACTS_INFO("*** Test: relChi2changeCutOff -- Finish");
0809 }
0810 
0811 BOOST_AUTO_TEST_CASE(DidNotConverge) {
0812   ACTS_INFO("*** Test: DidNotConverge -- Start");
0813 
0814   std::default_random_engine rng(42);
0815 
0816   ACTS_DEBUG("Create the detector");
0817   const std::size_t nSurfaces = 5;
0818   Detector detector;
0819   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0820 
0821   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0822   const auto parametersMeasurements = makeParameters();
0823   const auto startParametersFit = makeParameters(
0824       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0825 
0826   ACTS_DEBUG("Create the measurements");
0827   // simulation propagator
0828   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0829   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0830   const auto measurements =
0831       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0832                          resMapAllPixel, rng);
0833   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0834   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0835 
0836   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0837 
0838   ACTS_DEBUG("Set up the fitter");
0839   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0840 
0841   using RecoStepper = EigenStepper<>;
0842   const auto recoPropagator =
0843       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0844 
0845   using RecoPropagator = decltype(recoPropagator);
0846   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
0847   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
0848 
0849   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0850   extensions.calibrator
0851       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0852   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0853   extensions.surfaceAccessor
0854       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0855 
0856   // The relChi2changeCutOff = 0 prevents to stop the fitter after convergence,
0857   // therefore all updates will be done (even if the result does not change).
0858   // Since we didn't break due to convergence, we reach nUpdatesMax and
0859   // therefore fail the fit.
0860   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0861                                      PropagatorPlainOptions(geoCtx, magCtx),
0862                                      rSurface, false, false,
0863                                      FreeToBoundCorrection(false), 6, 0);
0864 
0865   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0866 
0867   ACTS_DEBUG("Fit the track");
0868   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0869   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0870   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0871                               startParametersFit, gx2fOptions, tracks);
0872 
0873   BOOST_REQUIRE(!res.ok());
0874   BOOST_CHECK_EQUAL(res.error(), GlobalChiSquareFitterError::DidNotConverge);
0875 
0876   ACTS_INFO("*** Test: DidNotConverge -- Finish");
0877 }
0878 
0879 BOOST_AUTO_TEST_CASE(NotEnoughMeasurements) {
0880   ACTS_INFO("*** Test: NotEnoughMeasurements -- Start");
0881 
0882   std::default_random_engine rng(42);
0883 
0884   ACTS_DEBUG("Create the detector");
0885   const std::size_t nSurfaces = 2;
0886   Detector detector;
0887   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0888 
0889   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0890   const auto parametersMeasurements = makeParameters();
0891   const auto startParametersFit = makeParameters(
0892       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0893 
0894   ACTS_DEBUG("Create the measurements");
0895   // simulation propagator
0896   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0897   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0898   const auto measurements =
0899       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0900                          resMapAllPixel, rng);
0901   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0902   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
0903 
0904   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
0905 
0906   ACTS_DEBUG("Set up the fitter");
0907   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0908 
0909   using RecoStepper = EigenStepper<>;
0910   const auto recoPropagator =
0911       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0912 
0913   using RecoPropagator = decltype(recoPropagator);
0914   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
0915   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
0916 
0917   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
0918   extensions.calibrator
0919       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
0920   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
0921   extensions.surfaceAccessor
0922       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
0923 
0924   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
0925                                      PropagatorPlainOptions(geoCtx, magCtx),
0926                                      rSurface, false, false,
0927                                      FreeToBoundCorrection(false), 6, 0);
0928 
0929   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
0930 
0931   ACTS_DEBUG("Fit the track");
0932   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
0933   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
0934   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
0935                               startParametersFit, gx2fOptions, tracks);
0936 
0937   BOOST_REQUIRE(!res.ok());
0938   BOOST_CHECK_EQUAL(res.error(),
0939                     GlobalChiSquareFitterError::NotEnoughMeasurements);
0940 
0941   ACTS_INFO("*** Test: NotEnoughMeasurements -- Finish");
0942 }
0943 
0944 BOOST_AUTO_TEST_CASE(FindHoles) {
0945   ACTS_INFO("*** Test: FindHoles -- Start");
0946 
0947   std::default_random_engine rng(42);
0948 
0949   ACTS_DEBUG("Create the detector");
0950   //  const std::size_t nSurfaces = 7;
0951   const std::size_t nSurfaces = 8;
0952   Detector detector;
0953   detector.geometry = makeToyDetector(geoCtx, nSurfaces);
0954 
0955   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
0956   const auto parametersMeasurements = makeParameters();
0957   const auto startParametersFit = makeParameters(
0958       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
0959 
0960   ACTS_DEBUG("Create the measurements");
0961   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
0962   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
0963   const auto measurements =
0964       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
0965                          resMapAllPixel, rng);
0966   auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
0967   ACTS_VERBOSE("sourceLinks.size() [before] = " << sourceLinks.size());
0968 
0969   // We remove the first measurement in the list. This does not create a hole.
0970   sourceLinks.erase(std::next(sourceLinks.begin(), 0));
0971   ACTS_VERBOSE(
0972       "sourceLinks.size() [after first erase] = " << sourceLinks.size());
0973 
0974   // We remove the last measurement in the list. This does not create a hole.
0975   sourceLinks.pop_back();
0976   ACTS_VERBOSE("sourceLinks.size() [after pop] = " << sourceLinks.size());
0977 
0978   // We remove the second to last measurement in the list. This effectively
0979   // creates a hole on that surface.
0980   const std::size_t indexHole = sourceLinks.size() - 2;
0981   ACTS_VERBOSE("Remove measurement " << indexHole);
0982   sourceLinks.erase(std::next(sourceLinks.begin(), indexHole));
0983   ACTS_VERBOSE("sourceLinks.size() [after second-to-last erase]= "
0984                << sourceLinks.size());
0985 
0986   // We removed 3 measurements
0987   //  const std::size_t nMeasurements = nSurfaces - 2;
0988   const std::size_t nMeasurements = nSurfaces - 3;
0989   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nMeasurements);
0990 
0991   ACTS_DEBUG("Set up the fitter");
0992   const Surface* rSurface = &parametersMeasurements.referenceSurface();
0993 
0994   using RecoStepper = EigenStepper<>;
0995   const auto recoPropagator =
0996       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
0997 
0998   using RecoPropagator = decltype(recoPropagator);
0999   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
1000   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
1001 
1002   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
1003   extensions.calibrator
1004       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
1005   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
1006   extensions.surfaceAccessor
1007       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
1008 
1009   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
1010                                      PropagatorPlainOptions(geoCtx, magCtx),
1011                                      rSurface, false, false,
1012                                      FreeToBoundCorrection(false), 20, 1e-5);
1013 
1014   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
1015 
1016   ACTS_DEBUG("Fit the track");
1017   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
1018   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
1019   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
1020                               startParametersFit, gx2fOptions, tracks);
1021 
1022   BOOST_REQUIRE(res.ok());
1023 
1024   const auto& track = *res;
1025 
1026   // -1, because the index starts at 0
1027   // -2, because the first and the last surface are not part of the track
1028   BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1 - 2);
1029   BOOST_CHECK(track.hasReferenceSurface());
1030 
1031   // Track quantities
1032   CHECK_CLOSE_ABS(track.chi2(), 6.5, 2.);
1033   BOOST_CHECK_EQUAL(track.nDoF(), 10u);
1034   BOOST_CHECK_EQUAL(track.nHoles(), 1u);
1035   BOOST_CHECK_EQUAL(track.nMeasurements(), nMeasurements);
1036   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
1037   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
1038 
1039   // Parameters
1040   // We need quite coarse checks here, since on different builds
1041   // the created measurements differ in the randomness
1042   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 7e0);
1043   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 6e0);
1044   BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-5, 1e3);
1045   BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], std::numbers::pi / 2,
1046                     1e-3);
1047   BOOST_CHECK_EQUAL(track.parameters()[eBoundQOverP], 1);
1048   BOOST_CHECK_CLOSE(track.parameters()[eBoundTime],
1049                     startParametersFit.parameters()[eBoundTime], 1e-6);
1050   BOOST_CHECK_CLOSE(track.covariance().determinant(), 4.7e-28, 2e0);
1051 
1052   ACTS_INFO("*** Test: FindHoles -- Finish");
1053 }
1054 
1055 BOOST_AUTO_TEST_CASE(Material) {
1056   ACTS_INFO("*** Test: Material -- Start");
1057 
1058   std::default_random_engine rng(42);
1059 
1060   ACTS_DEBUG("Create the detector");
1061   const std::size_t nSurfaces = 7;
1062   const std::set<std::size_t> surfaceIndexWithMaterial = {4};
1063   Detector detector;
1064   detector.geometry =
1065       makeToyDetector(geoCtx, nSurfaces, surfaceIndexWithMaterial);
1066 
1067   ACTS_DEBUG("Set the start parameters for measurement creation and fit");
1068   const auto parametersMeasurements = makeParameters();
1069   const auto startParametersFit = makeParameters(
1070       7_mm, 11_mm, 15_mm, 42_ns, 10_degree, 80_degree, 1_GeV, 1_e);
1071 
1072   ACTS_DEBUG("Create the measurements");
1073   using SimPropagator = Propagator<StraightLineStepper, Navigator>;
1074   const SimPropagator simPropagator = makeStraightPropagator(detector.geometry);
1075   auto measurements =
1076       createMeasurements(simPropagator, geoCtx, magCtx, parametersMeasurements,
1077                          resMapAllPixel, rng);
1078 
1079   const Vector2 scatterOffset = {100_mm, 100_mm};
1080   const std::size_t indexMaterialSurface = 3;
1081   for (std::size_t iMeas = indexMaterialSurface; iMeas < nSurfaces; iMeas++) {
1082     // This only works, because our detector is evenly spaced
1083     const std::size_t offsetFactor = iMeas - indexMaterialSurface;
1084 
1085     auto& sl = measurements.sourceLinks[iMeas];
1086     sl.parameters[0] += scatterOffset[0] * offsetFactor;
1087     sl.parameters[1] += scatterOffset[1] * offsetFactor;
1088   }
1089 
1090   const auto sourceLinks = prepareSourceLinks(measurements.sourceLinks);
1091   ACTS_VERBOSE("sourceLinks.size() = " << sourceLinks.size());
1092 
1093   BOOST_REQUIRE_EQUAL(sourceLinks.size(), nSurfaces);
1094 
1095   ACTS_DEBUG("Set up the fitter");
1096   const Surface* rSurface = &parametersMeasurements.referenceSurface();
1097 
1098   using RecoStepper = EigenStepper<>;
1099   const auto recoPropagator =
1100       makeConstantFieldPropagator<RecoStepper>(detector.geometry, 0_T);
1101 
1102   using RecoPropagator = decltype(recoPropagator);
1103   using Gx2Fitter = Gx2Fitter<RecoPropagator, VectorMultiTrajectory>;
1104   const Gx2Fitter fitter(recoPropagator, gx2fLogger->clone());
1105 
1106   Gx2FitterExtensions<VectorMultiTrajectory> extensions;
1107   extensions.calibrator
1108       .connect<&testSourceLinkCalibratorStrict<VectorMultiTrajectory>>();
1109   TestSourceLink::SurfaceAccessor surfaceAccessor{*detector.geometry};
1110   extensions.surfaceAccessor
1111       .connect<&TestSourceLink::SurfaceAccessor::operator()>(&surfaceAccessor);
1112 
1113   const Gx2FitterOptions gx2fOptions(geoCtx, magCtx, calCtx, extensions,
1114                                      PropagatorPlainOptions(geoCtx, magCtx),
1115                                      rSurface, true, false,
1116                                      FreeToBoundCorrection(false), 5, 0);
1117 
1118   TrackContainer tracks{VectorTrackContainer{}, VectorMultiTrajectory{}};
1119 
1120   ACTS_DEBUG("Fit the track");
1121   ACTS_VERBOSE("startParameter unsmeared:\n" << parametersMeasurements);
1122   ACTS_VERBOSE("startParameter fit:\n" << startParametersFit);
1123   const auto res = fitter.fit(sourceLinks.begin(), sourceLinks.end(),
1124                               startParametersFit, gx2fOptions, tracks);
1125 
1126   // Helper to visualise the detector
1127   {
1128     std::cout << "\n*** Create .obj of Detector ***\n" << std::endl;
1129     // Only need for obj
1130     ObjVisualization3D obj;
1131 
1132     bool triangulate = true;
1133     ViewConfig viewSensitive = {.color = {0, 180, 240}};
1134     viewSensitive.triangulate = triangulate;
1135     ViewConfig viewPassive = {.color = {240, 280, 0}};
1136     viewPassive.triangulate = triangulate;
1137     ViewConfig viewVolume = {.color = {220, 220, 0}};
1138     viewVolume.triangulate = triangulate;
1139     ViewConfig viewContainer = {.color = {220, 220, 0}};
1140     viewContainer.triangulate = triangulate;
1141     ViewConfig viewGrid = {.color = {220, 0, 0}};
1142     viewGrid.quarterSegments = 8;
1143     viewGrid.offset = 3.;
1144     viewGrid.triangulate = triangulate;
1145 
1146     std::string tag = "gx2f_toydet";
1147 
1148     const TrackingVolume& tgVolume =
1149         *(detector.geometry->highestTrackingVolume());
1150 
1151     GeometryView3D::drawTrackingVolume(obj, tgVolume, geoCtx, viewContainer,
1152                                        viewVolume, viewPassive, viewSensitive,
1153                                        viewGrid, true, tag);
1154   }
1155   // Helper to visualise the measurements
1156   {
1157     std::cout << "\n*** Create .obj of measurements ***\n" << std::endl;
1158     ObjVisualization3D obj;
1159 
1160     double localErrorScale = 10000000.;
1161     ViewConfig mcolor{.color = {255, 145, 48}};
1162     mcolor.offset = 2;
1163     //  mcolor.visible = true;
1164 
1165     drawMeasurements(obj, measurements, detector.geometry, geoCtx,
1166                      localErrorScale, mcolor);
1167 
1168     obj.write("meas");
1169   }
1170 
1171   BOOST_REQUIRE(res.ok());
1172 
1173   const auto& track = *res;
1174 
1175   BOOST_CHECK_EQUAL(track.tipIndex(), nSurfaces - 1);
1176   BOOST_CHECK(track.hasReferenceSurface());
1177 
1178   // TODO Add material handling to the gx2f, to pass the 6 commented tests
1179   // Track quantities
1180   //  CHECK_CLOSE_ABS(track.chi2(), 8., 2.);
1181   BOOST_CHECK_EQUAL(track.nDoF(), nSurfaces * 2);
1182   BOOST_CHECK_EQUAL(track.nHoles(), 0u);
1183   BOOST_CHECK_EQUAL(track.nMeasurements(), nSurfaces);
1184   BOOST_CHECK_EQUAL(track.nSharedHits(), 0u);
1185   BOOST_CHECK_EQUAL(track.nOutliers(), 0u);
1186 
1187   // Parameters
1188   // We need quite coarse checks here, since on different builds
1189   // the created measurements differ in the randomness
1190   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc0], -11., 26e0);
1191   BOOST_CHECK_CLOSE(track.parameters()[eBoundLoc1], -15., 15e0);
1192   BOOST_CHECK_CLOSE(track.parameters()[eBoundPhi], 1e-5, 1.1e3);
1193   BOOST_CHECK_CLOSE(track.parameters()[eBoundTheta], std::numbers::pi / 2,
1194                     2e-2);
1195   BOOST_CHECK_EQUAL(track.parameters()[eBoundQOverP], 1);
1196   BOOST_CHECK_CLOSE(track.parameters()[eBoundTime],
1197                     startParametersFit.parameters()[eBoundTime], 1e-6);
1198   BOOST_CHECK_CLOSE(track.covariance().determinant(), 3.5e-27, 1e1);
1199 
1200   // Convergence
1201   BOOST_CHECK_EQUAL(
1202       (track.template component<
1203           std::uint32_t, hashString(Gx2fConstants::gx2fnUpdateColumn)>()),
1204       4);
1205 
1206   ACTS_INFO("*** Test: Material -- Finish");
1207 }
1208 BOOST_AUTO_TEST_SUITE_END()
1209 
1210 }  // namespace ActsTests