Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-06 08:11:36

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2019-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/data/test_case.hpp>
0010 #include <boost/test/tools/output_test_stream.hpp>
0011 #include <boost/test/unit_test.hpp>
0012 
0013 #include "Acts/Definitions/Algebra.hpp"
0014 #include "Acts/Definitions/Direction.hpp"
0015 #include "Acts/Definitions/TrackParametrization.hpp"
0016 #include "Acts/Definitions/Units.hpp"
0017 #include "Acts/EventData/GenericBoundTrackParameters.hpp"
0018 #include "Acts/EventData/TrackParameters.hpp"
0019 #include "Acts/Geometry/GeometryContext.hpp"
0020 #include "Acts/Geometry/GeometryIdentifier.hpp"
0021 #include "Acts/MagneticField/ConstantBField.hpp"
0022 #include "Acts/MagneticField/MagneticFieldContext.hpp"
0023 #include "Acts/Propagator/EigenStepper.hpp"
0024 #include "Acts/Propagator/Propagator.hpp"
0025 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0026 #include "Acts/Utilities/AnnealingUtility.hpp"
0027 #include "Acts/Utilities/Helpers.hpp"
0028 #include "Acts/Utilities/Result.hpp"
0029 #include "Acts/Vertexing/AdaptiveGridDensityVertexFinder.hpp"
0030 #include "Acts/Vertexing/AdaptiveMultiVertexFinder.hpp"
0031 #include "Acts/Vertexing/AdaptiveMultiVertexFitter.hpp"
0032 #include "Acts/Vertexing/GaussianTrackDensity.hpp"
0033 #include "Acts/Vertexing/GridDensityVertexFinder.hpp"
0034 #include "Acts/Vertexing/HelicalTrackLinearizer.hpp"
0035 #include "Acts/Vertexing/IVertexFinder.hpp"
0036 #include "Acts/Vertexing/ImpactPointEstimator.hpp"
0037 #include "Acts/Vertexing/TrackAtVertex.hpp"
0038 #include "Acts/Vertexing/TrackDensityVertexFinder.hpp"
0039 #include "Acts/Vertexing/Vertex.hpp"
0040 #include "Acts/Vertexing/VertexingOptions.hpp"
0041 
0042 #include <algorithm>
0043 #include <array>
0044 #include <chrono>
0045 #include <cmath>
0046 #include <functional>
0047 #include <iostream>
0048 #include <map>
0049 #include <memory>
0050 #include <string>
0051 #include <system_error>
0052 #include <tuple>
0053 #include <utility>
0054 #include <vector>
0055 
0056 #include "VertexingDataHelper.hpp"
0057 
0058 namespace Acts::Test {
0059 
0060 using namespace Acts::UnitLiterals;
0061 
0062 using Covariance = BoundSquareMatrix;
0063 using Propagator = Acts::Propagator<EigenStepper<>>;
0064 using Linearizer = HelicalTrackLinearizer;
0065 
0066 // Create a test context
0067 GeometryContext geoContext = GeometryContext();
0068 MagneticFieldContext magFieldContext = MagneticFieldContext();
0069 
0070 const std::string toolString = "AMVF";
0071 
0072 /// @brief AMVF test with Gaussian seed finder
0073 BOOST_AUTO_TEST_CASE(adaptive_multi_vertex_finder_test) {
0074   // Set debug mode
0075   bool debugMode = false;
0076   // Set up constant B-Field
0077   auto bField = std::make_shared<ConstantBField>(Vector3(0., 0., 2_T));
0078 
0079   // Set up EigenStepper
0080   // EigenStepper<> stepper(bField);
0081   EigenStepper<> stepper(bField);
0082 
0083   // Set up propagator with void navigator
0084   auto propagator = std::make_shared<Propagator>(stepper);
0085 
0086   // IP 3D Estimator
0087   ImpactPointEstimator::Config ipEstimatorCfg(bField, propagator);
0088   ImpactPointEstimator ipEstimator(ipEstimatorCfg);
0089 
0090   std::vector<double> temperatures{8.0, 4.0, 2.0, 1.4142136, 1.2247449, 1.0};
0091   AnnealingUtility::Config annealingConfig;
0092   annealingConfig.setOfTemperatures = temperatures;
0093   AnnealingUtility annealingUtility(annealingConfig);
0094 
0095   using Fitter = AdaptiveMultiVertexFitter;
0096 
0097   Fitter::Config fitterCfg(ipEstimator);
0098 
0099   fitterCfg.annealingTool = annealingUtility;
0100 
0101   // Linearizer for BoundTrackParameters type test
0102   Linearizer::Config ltConfig;
0103   ltConfig.bField = bField;
0104   ltConfig.propagator = propagator;
0105   Linearizer linearizer(ltConfig);
0106 
0107   // Test smoothing
0108   fitterCfg.doSmoothing = true;
0109   fitterCfg.extractParameters.connect<&InputTrack::extractParameters>();
0110   fitterCfg.trackLinearizer.connect<&Linearizer::linearizeTrack>(&linearizer);
0111 
0112   Fitter fitter(fitterCfg);
0113 
0114   GaussianTrackDensity::Config densityCfg;
0115   densityCfg.extractParameters.connect<&InputTrack::extractParameters>();
0116   auto seedFinder = std::make_shared<TrackDensityVertexFinder>(
0117       TrackDensityVertexFinder::Config{densityCfg});
0118 
0119   AdaptiveMultiVertexFinder::Config finderConfig(std::move(fitter), seedFinder,
0120                                                  ipEstimator, bField);
0121   finderConfig.extractParameters.connect<&InputTrack::extractParameters>();
0122 
0123   AdaptiveMultiVertexFinder finder(std::move(finderConfig));
0124   IVertexFinder::State state = finder.makeState(magFieldContext);
0125 
0126   auto csvData = readTracksAndVertexCSV(toolString);
0127   std::vector<BoundTrackParameters> tracks = std::get<TracksData>(csvData);
0128 
0129   if (debugMode) {
0130     std::cout << "Number of tracks in event: " << tracks.size() << std::endl;
0131     int maxCout = 10;
0132     int count = 0;
0133     for (const auto& trk : tracks) {
0134       std::cout << count << ". track: " << std::endl;
0135       std::cout << "params: " << trk << std::endl;
0136       count++;
0137       if (count == maxCout) {
0138         break;
0139       }
0140     }
0141   }
0142 
0143   std::vector<InputTrack> inputTracks;
0144   for (const auto& trk : tracks) {
0145     inputTracks.emplace_back(&trk);
0146   }
0147 
0148   // TODO: test without using beam spot constraint
0149   Vertex bsConstr = std::get<BeamSpotData>(csvData);
0150   VertexingOptions vertexingOptions(geoContext, magFieldContext, bsConstr);
0151 
0152   auto t1 = std::chrono::system_clock::now();
0153   auto findResult = finder.find(inputTracks, vertexingOptions, state);
0154   auto t2 = std::chrono::system_clock::now();
0155 
0156   auto timediff =
0157       std::chrono::duration_cast<std::chrono::milliseconds>(t2 - t1).count();
0158 
0159   if (!findResult.ok()) {
0160     std::cout << findResult.error().message() << std::endl;
0161   }
0162 
0163   BOOST_CHECK(findResult.ok());
0164 
0165   std::vector<Vertex> allVertices = *findResult;
0166 
0167   if (debugMode) {
0168     std::cout << "Time needed: " << timediff << " ms." << std::endl;
0169     std::cout << "Number of vertices reconstructed: " << allVertices.size()
0170               << std::endl;
0171 
0172     int count = 0;
0173     for (const auto& vtx : allVertices) {
0174       count++;
0175       std::cout << count << ". Vertex at position: " << vtx.position()[0]
0176                 << ", " << vtx.position()[1] << ", " << vtx.position()[2]
0177                 << std::endl;
0178       std::cout << count << ". Vertex with cov: " << vtx.covariance()
0179                 << std::endl;
0180       std::cout << "\t with n tracks: " << vtx.tracks().size() << std::endl;
0181     }
0182   }
0183 
0184   // Test expected outcomes from athena implementation
0185   // Number of reconstructed vertices
0186   auto verticesInfo = std::get<VerticesData>(csvData);
0187   const int expNRecoVertices = verticesInfo.size();
0188 
0189   BOOST_CHECK_EQUAL(allVertices.size(), expNRecoVertices);
0190 
0191   double relTol = 1e-2;
0192   double small = 1e-3;
0193   for (int i = 0; i < expNRecoVertices; i++) {
0194     auto recoVtx = allVertices[i];
0195     auto expVtx = verticesInfo[i];
0196     CHECK_CLOSE_OR_SMALL(recoVtx.position(), expVtx.position, relTol, small);
0197     CHECK_CLOSE_OR_SMALL(recoVtx.covariance(), expVtx.covariance, relTol,
0198                          small);
0199     BOOST_CHECK_EQUAL(recoVtx.tracks().size(), expVtx.nTracks);
0200     CHECK_CLOSE_OR_SMALL(recoVtx.tracks()[0].trackWeight, expVtx.trk1Weight,
0201                          relTol, small);
0202     CHECK_CLOSE_OR_SMALL(recoVtx.tracks()[0].vertexCompatibility,
0203                          expVtx.trk1Comp, relTol, small);
0204   }
0205 }
0206 
0207 // Dummy user-defined InputTrackStub type
0208 struct InputTrackStub {
0209   InputTrackStub(const BoundTrackParameters& params, int id)
0210       : m_parameters(params), m_id(id) {}
0211 
0212   const BoundTrackParameters& parameters() const { return m_parameters; }
0213   // store e.g. link to original objects here
0214 
0215   int id() const { return m_id; }
0216 
0217  private:
0218   BoundTrackParameters m_parameters;
0219 
0220   // Some test track ID
0221   int m_id;
0222 };
0223 
0224 /// @brief AMVF test with user-defined input track type
0225 BOOST_AUTO_TEST_CASE(adaptive_multi_vertex_finder_usertype_test) {
0226   // Set debug mode
0227   bool debugMode = false;
0228   // Set up constant B-Field
0229   auto bField = std::make_shared<ConstantBField>(Vector3(0., 0., 2_T));
0230 
0231   // Set up EigenStepper
0232   // EigenStepper<> stepper(bField);
0233   EigenStepper<> stepper(bField);
0234 
0235   // Set up propagator with void navigator
0236   auto propagator = std::make_shared<Propagator>(stepper);
0237 
0238   // Create a custom std::function to extract BoundTrackParameters from
0239   // user-defined InputTrackStub
0240   auto extractParameters = [](const InputTrack& track) {
0241     return track.as<InputTrackStub>()->parameters();
0242   };
0243 
0244   // IP 3D Estimator
0245   ImpactPointEstimator::Config ipEstimatorCfg(bField, propagator);
0246   ImpactPointEstimator ipEstimator(ipEstimatorCfg);
0247 
0248   std::vector<double> temperatures{8.0, 4.0, 2.0, 1.4142136, 1.2247449, 1.0};
0249   AnnealingUtility::Config annealingConfig;
0250   annealingConfig.setOfTemperatures = temperatures;
0251   AnnealingUtility annealingUtility(annealingConfig);
0252 
0253   using Fitter = AdaptiveMultiVertexFitter;
0254 
0255   Fitter::Config fitterCfg(ipEstimator);
0256 
0257   fitterCfg.annealingTool = annealingUtility;
0258 
0259   // Linearizer
0260   Linearizer::Config ltConfig;
0261   ltConfig.bField = bField;
0262   ltConfig.propagator = propagator;
0263   Linearizer linearizer(ltConfig);
0264 
0265   // Test smoothing
0266   fitterCfg.doSmoothing = true;
0267   fitterCfg.extractParameters.connect(extractParameters);
0268   fitterCfg.trackLinearizer.connect<&Linearizer::linearizeTrack>(&linearizer);
0269 
0270   Fitter fitter(fitterCfg);
0271 
0272   GaussianTrackDensity::Config densityCfg;
0273   densityCfg.extractParameters.connect(extractParameters);
0274   auto seedFinder = std::make_shared<TrackDensityVertexFinder>(
0275       TrackDensityVertexFinder::Config{densityCfg});
0276 
0277   AdaptiveMultiVertexFinder::Config finderConfig(
0278       std::move(fitter), std::move(seedFinder), ipEstimator, bField);
0279   finderConfig.extractParameters.connect(extractParameters);
0280 
0281   AdaptiveMultiVertexFinder finder(std::move(finderConfig));
0282   IVertexFinder::State state = finder.makeState(magFieldContext);
0283 
0284   auto csvData = readTracksAndVertexCSV(toolString);
0285   auto tracks = std::get<TracksData>(csvData);
0286 
0287   std::vector<InputTrackStub> userTracks;
0288   int idCount = 0;
0289   for (const auto& trk : tracks) {
0290     userTracks.push_back(InputTrackStub(trk, idCount));
0291     idCount++;
0292   }
0293 
0294   if (debugMode) {
0295     std::cout << "Number of tracks in event: " << tracks.size() << std::endl;
0296     int maxCout = 10;
0297     int count = 0;
0298     for (const auto& trk : tracks) {
0299       std::cout << count << ". track: " << std::endl;
0300       std::cout << "params: " << trk << std::endl;
0301       count++;
0302       if (count == maxCout) {
0303         break;
0304       }
0305     }
0306   }
0307 
0308   std::vector<InputTrack> userInputTracks;
0309   for (const auto& trk : userTracks) {
0310     userInputTracks.emplace_back(&trk);
0311   }
0312 
0313   Vertex constraintVtx;
0314   constraintVtx.setPosition(std::get<BeamSpotData>(csvData).position());
0315   constraintVtx.setCovariance(std::get<BeamSpotData>(csvData).covariance());
0316 
0317   VertexingOptions vertexingOptions(geoContext, magFieldContext, constraintVtx);
0318 
0319   auto findResult = finder.find(userInputTracks, vertexingOptions, state);
0320 
0321   if (!findResult.ok()) {
0322     std::cout << findResult.error().message() << std::endl;
0323   }
0324 
0325   BOOST_CHECK(findResult.ok());
0326 
0327   std::vector<Vertex> allVertices = *findResult;
0328 
0329   if (debugMode) {
0330     std::cout << "Number of vertices reconstructed: " << allVertices.size()
0331               << std::endl;
0332 
0333     int count = 0;
0334     for (const auto& vtx : allVertices) {
0335       count++;
0336       std::cout << count << ". Vertex at position: " << vtx.position()[0]
0337                 << ", " << vtx.position()[1] << ", " << vtx.position()[2]
0338                 << std::endl;
0339       std::cout << count << ". Vertex with cov: " << vtx.covariance()
0340                 << std::endl;
0341       std::cout << "\t with n tracks: " << vtx.tracks().size() << std::endl;
0342     }
0343     for (auto& trk : allVertices[0].tracks()) {
0344       std::cout << "Track ID at first vertex: "
0345                 << trk.originalParams.as<InputTrackStub>()->id() << std::endl;
0346     }
0347   }
0348 
0349   auto verticesInfo = std::get<VerticesData>(csvData);
0350   const int expNRecoVertices = verticesInfo.size();
0351 
0352   BOOST_CHECK_EQUAL(allVertices.size(), expNRecoVertices);
0353 
0354   double relTol = 1e-2;
0355   double small = 1e-3;
0356   for (int i = 0; i < expNRecoVertices; i++) {
0357     auto recoVtx = allVertices[i];
0358     auto expVtx = verticesInfo[i];
0359     CHECK_CLOSE_OR_SMALL(recoVtx.position(), expVtx.position, relTol, small);
0360     CHECK_CLOSE_OR_SMALL(recoVtx.covariance(), expVtx.covariance, relTol,
0361                          small);
0362     BOOST_CHECK_EQUAL(recoVtx.tracks().size(), expVtx.nTracks);
0363     CHECK_CLOSE_OR_SMALL(recoVtx.tracks()[0].trackWeight, expVtx.trk1Weight,
0364                          relTol, small);
0365     CHECK_CLOSE_OR_SMALL(recoVtx.tracks()[0].vertexCompatibility,
0366                          expVtx.trk1Comp, relTol, small);
0367   }
0368 }
0369 
0370 /// @brief AMVF test with grid seed finder
0371 BOOST_AUTO_TEST_CASE(adaptive_multi_vertex_finder_grid_seed_finder_test) {
0372   // Set debug mode
0373   bool debugMode = false;
0374   if (debugMode) {
0375     std::cout << "Starting AMVF test with grid seed finder..." << std::endl;
0376   }
0377   // Set up constant B-Field
0378   auto bField = std::make_shared<ConstantBField>(Vector3(0., 0., 2_T));
0379 
0380   // Set up EigenStepper
0381   // EigenStepper<> stepper(bField);
0382   EigenStepper<> stepper(bField);
0383 
0384   // Set up propagator with void navigator
0385   auto propagator = std::make_shared<Propagator>(stepper);
0386 
0387   // IP Estimator
0388   ImpactPointEstimator::Config ipEstCfg(bField, propagator);
0389   ImpactPointEstimator ipEst(ipEstCfg);
0390 
0391   std::vector<double> temperatures{8.0, 4.0, 2.0, 1.4142136, 1.2247449, 1.0};
0392   AnnealingUtility::Config annealingConfig;
0393   annealingConfig.setOfTemperatures = temperatures;
0394   AnnealingUtility annealingUtility(annealingConfig);
0395 
0396   using Fitter = AdaptiveMultiVertexFitter;
0397 
0398   Fitter::Config fitterCfg(ipEst);
0399 
0400   fitterCfg.annealingTool = annealingUtility;
0401 
0402   // Linearizer for BoundTrackParameters type test
0403   Linearizer::Config ltConfig;
0404   ltConfig.bField = bField;
0405   ltConfig.propagator = propagator;
0406   Linearizer linearizer(ltConfig);
0407 
0408   // Test smoothing
0409   fitterCfg.doSmoothing = true;
0410   fitterCfg.extractParameters.connect<&InputTrack::extractParameters>();
0411   fitterCfg.trackLinearizer.connect<&Linearizer::linearizeTrack>(&linearizer);
0412 
0413   Fitter fitter(fitterCfg);
0414 
0415   using SeedFinder = GridDensityVertexFinder;
0416   SeedFinder::Config seedFinderCfg{{{250, 4000, 55}}};
0417   seedFinderCfg.cacheGridStateForTrackRemoval = true;
0418   seedFinderCfg.extractParameters.connect<&InputTrack::extractParameters>();
0419 
0420   auto seedFinder = std::make_shared<SeedFinder>(seedFinderCfg);
0421 
0422   AdaptiveMultiVertexFinder::Config finderConfig(
0423       std::move(fitter), std::move(seedFinder), ipEst, bField);
0424   finderConfig.extractParameters.connect<&InputTrack::extractParameters>();
0425 
0426   AdaptiveMultiVertexFinder finder(std::move(finderConfig));
0427   IVertexFinder::State state = finder.makeState(magFieldContext);
0428 
0429   auto csvData = readTracksAndVertexCSV(toolString);
0430   auto tracks = std::get<TracksData>(csvData);
0431 
0432   if (debugMode) {
0433     std::cout << "Number of tracks in event: " << tracks.size() << std::endl;
0434     int maxCout = 10;
0435     int count = 0;
0436     for (const auto& trk : tracks) {
0437       std::cout << count << ". track: " << std::endl;
0438       std::cout << "params: " << trk << std::endl;
0439       count++;
0440       if (count == maxCout) {
0441         break;
0442       }
0443     }
0444   }
0445 
0446   std::vector<InputTrack> inputTracks;
0447   for (const auto& trk : tracks) {
0448     inputTracks.emplace_back(&trk);
0449   }
0450 
0451   // TODO: test using beam spot constraint
0452   Vertex bsConstr = std::get<BeamSpotData>(csvData);
0453   VertexingOptions vertexingOptions(geoContext, magFieldContext, bsConstr);
0454 
0455   auto t1 = std::chrono::system_clock::now();
0456   auto findResult = finder.find(inputTracks, vertexingOptions, state);
0457   auto t2 = std::chrono::system_clock::now();
0458 
0459   auto timediff =
0460       std::chrono::duration_cast<std::chrono::milliseconds>(t2 - t1).count();
0461 
0462   if (!findResult.ok()) {
0463     std::cout << findResult.error().message() << std::endl;
0464   }
0465 
0466   BOOST_CHECK(findResult.ok());
0467 
0468   std::vector<Vertex> allVertices = *findResult;
0469 
0470   if (debugMode) {
0471     std::cout << "Time needed: " << timediff << " ms." << std::endl;
0472     std::cout << "Number of vertices reconstructed: " << allVertices.size()
0473               << std::endl;
0474 
0475     int count = 0;
0476     for (const auto& vtx : allVertices) {
0477       count++;
0478       std::cout << count << ". Vertex at position: " << vtx.position()[0]
0479                 << ", " << vtx.position()[1] << ", " << vtx.position()[2]
0480                 << std::endl;
0481       std::cout << count << ". Vertex with cov: " << vtx.covariance()
0482                 << std::endl;
0483       std::cout << "\t with n tracks: " << vtx.tracks().size() << std::endl;
0484     }
0485   }
0486   // Test expected outcomes from athena implementation
0487   // Number of reconstructed vertices
0488   auto verticesInfo = std::get<VerticesData>(csvData);
0489   const int expNRecoVertices = verticesInfo.size();
0490 
0491   BOOST_CHECK_EQUAL(allVertices.size(), expNRecoVertices);
0492   std::vector<bool> vtxFound(expNRecoVertices, false);
0493 
0494   for (const auto& vtx : allVertices) {
0495     double vtxZ = vtx.position()[2];
0496     double diffZ = 1e5;
0497     int foundVtxIdx = -1;
0498     for (int i = 0; i < expNRecoVertices; i++) {
0499       if (!vtxFound[i]) {
0500         if (std::abs(vtxZ - verticesInfo[i].position[2]) < diffZ) {
0501           diffZ = std::abs(vtxZ - verticesInfo[i].position[2]);
0502           foundVtxIdx = i;
0503         }
0504       }
0505     }
0506     if (diffZ < 0.5_mm) {
0507       vtxFound[foundVtxIdx] = true;
0508       CHECK_CLOSE_ABS(vtx.tracks().size(), verticesInfo[foundVtxIdx].nTracks,
0509                       1);
0510     }
0511   }
0512   for (bool found : vtxFound) {
0513     BOOST_CHECK_EQUAL(found, true);
0514   }
0515 }
0516 
0517 /// @brief AMVF test with adaptive grid seed finder
0518 BOOST_AUTO_TEST_CASE(
0519     adaptive_multi_vertex_finder_adaptive_grid_seed_finder_test) {
0520   // Set debug mode
0521   bool debugMode = false;
0522   if (debugMode) {
0523     std::cout << "Starting AMVF test with adaptive grid seed finder..."
0524               << std::endl;
0525   }
0526   // Set up constant B-Field
0527   auto bField = std::make_shared<ConstantBField>(Vector3(0., 0., 2_T));
0528 
0529   // Set up EigenStepper
0530   // EigenStepper<> stepper(bField);
0531   EigenStepper<> stepper(bField);
0532 
0533   // Set up propagator with void navigator
0534   auto propagator = std::make_shared<Propagator>(stepper);
0535 
0536   // IP Estimator
0537   ImpactPointEstimator::Config ipEstCfg(bField, propagator);
0538   ImpactPointEstimator ipEst(ipEstCfg);
0539 
0540   std::vector<double> temperatures{8.0, 4.0, 2.0, 1.4142136, 1.2247449, 1.0};
0541   AnnealingUtility::Config annealingConfig;
0542   annealingConfig.setOfTemperatures = temperatures;
0543   AnnealingUtility annealingUtility(annealingConfig);
0544 
0545   using Fitter = AdaptiveMultiVertexFitter;
0546 
0547   Fitter::Config fitterCfg(ipEst);
0548 
0549   fitterCfg.annealingTool = annealingUtility;
0550 
0551   // Linearizer for BoundTrackParameters type test
0552   Linearizer::Config ltConfig;
0553   ltConfig.bField = bField;
0554   ltConfig.propagator = propagator;
0555   Linearizer linearizer(ltConfig);
0556 
0557   // Test smoothing
0558   fitterCfg.doSmoothing = true;
0559   fitterCfg.extractParameters.connect<&InputTrack::extractParameters>();
0560   fitterCfg.trackLinearizer.connect<&Linearizer::linearizeTrack>(&linearizer);
0561 
0562   Fitter fitter(fitterCfg);
0563 
0564   // Grid density used during vertex seed finding
0565   AdaptiveGridTrackDensity::Config gridDensityCfg;
0566   // force track to have exactly spatialTrkGridSize spatial bins for testing
0567   // purposes
0568   gridDensityCfg.spatialTrkGridSizeRange = {55, 55};
0569   gridDensityCfg.spatialBinExtent = 0.05;
0570   AdaptiveGridTrackDensity gridDensity(gridDensityCfg);
0571 
0572   using SeedFinder = AdaptiveGridDensityVertexFinder;
0573   SeedFinder::Config seedFinderCfg(gridDensity);
0574   seedFinderCfg.cacheGridStateForTrackRemoval = true;
0575   seedFinderCfg.extractParameters.connect<&InputTrack::extractParameters>();
0576 
0577   auto seedFinder = std::make_shared<SeedFinder>(seedFinderCfg);
0578 
0579   AdaptiveMultiVertexFinder::Config finderConfig(
0580       std::move(fitter), std::move(seedFinder), ipEst, bField);
0581   finderConfig.extractParameters.connect<&InputTrack::extractParameters>();
0582 
0583   AdaptiveMultiVertexFinder finder(std::move(finderConfig));
0584   IVertexFinder::State state = finder.makeState(magFieldContext);
0585 
0586   auto csvData = readTracksAndVertexCSV(toolString);
0587   auto tracks = std::get<TracksData>(csvData);
0588 
0589   if (debugMode) {
0590     std::cout << "Number of tracks in event: " << tracks.size() << std::endl;
0591     int maxCout = 10;
0592     int count = 0;
0593     for (const auto& trk : tracks) {
0594       std::cout << count << ". track: " << std::endl;
0595       std::cout << "params: " << trk << std::endl;
0596       count++;
0597       if (count == maxCout) {
0598         break;
0599       }
0600     }
0601   }
0602 
0603   std::vector<InputTrack> inputTracks;
0604   for (const auto& trk : tracks) {
0605     inputTracks.emplace_back(&trk);
0606   }
0607 
0608   Vertex bsConstr = std::get<BeamSpotData>(csvData);
0609   VertexingOptions vertexingOptions(geoContext, magFieldContext, bsConstr);
0610 
0611   auto t1 = std::chrono::system_clock::now();
0612   auto findResult = finder.find(inputTracks, vertexingOptions, state);
0613   auto t2 = std::chrono::system_clock::now();
0614 
0615   auto timediff =
0616       std::chrono::duration_cast<std::chrono::milliseconds>(t2 - t1).count();
0617 
0618   if (!findResult.ok()) {
0619     std::cout << findResult.error().message() << std::endl;
0620   }
0621 
0622   BOOST_CHECK(findResult.ok());
0623 
0624   std::vector<Vertex> allVertices = *findResult;
0625 
0626   if (debugMode) {
0627     std::cout << "Time needed: " << timediff << " ms." << std::endl;
0628     std::cout << "Number of vertices reconstructed: " << allVertices.size()
0629               << std::endl;
0630 
0631     int count = 0;
0632     for (const auto& vtx : allVertices) {
0633       count++;
0634       std::cout << count << ". Vertex at position: " << vtx.position()[0]
0635                 << ", " << vtx.position()[1] << ", " << vtx.position()[2]
0636                 << std::endl;
0637       std::cout << count << ". Vertex with cov: " << vtx.covariance()
0638                 << std::endl;
0639       std::cout << "\t with n tracks: " << vtx.tracks().size() << std::endl;
0640     }
0641   }
0642   // Test expected outcomes from athena implementation
0643   // Number of reconstructed vertices
0644   auto verticesInfo = std::get<VerticesData>(csvData);
0645   const int expNRecoVertices = verticesInfo.size();
0646 
0647   BOOST_CHECK_EQUAL(allVertices.size(), expNRecoVertices);
0648   std::vector<bool> vtxFound(expNRecoVertices, false);
0649 
0650   for (const auto& vtx : allVertices) {
0651     double vtxZ = vtx.position()[2];
0652     double diffZ = 1e5;
0653     int foundVtxIdx = -1;
0654     for (int i = 0; i < expNRecoVertices; i++) {
0655       if (!vtxFound[i]) {
0656         if (std::abs(vtxZ - verticesInfo[i].position[2]) < diffZ) {
0657           diffZ = std::abs(vtxZ - verticesInfo[i].position[2]);
0658           foundVtxIdx = i;
0659         }
0660       }
0661     }
0662     if (diffZ < 0.5_mm) {
0663       vtxFound[foundVtxIdx] = true;
0664       CHECK_CLOSE_ABS(vtx.tracks().size(), verticesInfo[foundVtxIdx].nTracks,
0665                       2);
0666     }
0667   }
0668   for (bool found : vtxFound) {
0669     BOOST_CHECK_EQUAL(found, true);
0670   }
0671 }
0672 
0673 }  // namespace Acts::Test