Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-05-17 07:36:34

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/Common.hpp"
0013 #include "Acts/Definitions/Direction.hpp"
0014 #include "Acts/Definitions/TrackParametrization.hpp"
0015 #include "Acts/Definitions/Units.hpp"
0016 #include "Acts/Geometry/GeometryContext.hpp"
0017 #include "Acts/Geometry/GeometryIdentifier.hpp"
0018 #include "Acts/MagneticField/ConstantBField.hpp"
0019 #include "Acts/MagneticField/MagneticFieldContext.hpp"
0020 #include "Acts/Propagator/EigenStepper.hpp"
0021 #include "Acts/Propagator/Propagator.hpp"
0022 #include "Acts/Surfaces/PerigeeSurface.hpp"
0023 #include "Acts/Surfaces/Surface.hpp"
0024 #include "Acts/Utilities/Result.hpp"
0025 #include "Acts/Vertexing/FullBilloirVertexFitter.hpp"
0026 #include "Acts/Vertexing/HelicalTrackLinearizer.hpp"
0027 #include "Acts/Vertexing/ImpactPointEstimator.hpp"
0028 #include "Acts/Vertexing/VertexingOptions.hpp"
0029 #include "Acts/Vertexing/ZScanVertexFinder.hpp"
0030 #include "ActsTests/CommonHelpers/FloatComparisons.hpp"
0031 
0032 #include <algorithm>
0033 #include <array>
0034 #include <cmath>
0035 #include <functional>
0036 #include <iostream>
0037 #include <memory>
0038 #include <numbers>
0039 #include <random>
0040 #include <string>
0041 #include <system_error>
0042 #include <tuple>
0043 #include <type_traits>
0044 #include <utility>
0045 #include <vector>
0046 
0047 using namespace Acts;
0048 using namespace Acts::UnitLiterals;
0049 
0050 namespace ActsTests {
0051 
0052 using Covariance = BoundMatrix;
0053 using Propagator = Acts::Propagator<EigenStepper<>>;
0054 
0055 // Create a test context
0056 GeometryContext geoContext = GeometryContext::dangerouslyDefaultConstruct();
0057 MagneticFieldContext magFieldContext = MagneticFieldContext();
0058 
0059 // Vertex x/y position distribution
0060 std::uniform_real_distribution<double> vXYDist(-0.1_mm, 0.1_mm);
0061 // Vertex z position distribution
0062 std::uniform_real_distribution<double> vZDist(-20_mm, 20_mm);
0063 // Track d0 distribution
0064 std::uniform_real_distribution<double> d0Dist(-0.01_mm, 0.01_mm);
0065 // Track z0 distribution
0066 std::uniform_real_distribution<double> z0Dist(-0.2_mm, 0.2_mm);
0067 // Track pT distribution
0068 std::uniform_real_distribution<double> pTDist(0.4_GeV, 10_GeV);
0069 // Track phi distribution
0070 std::uniform_real_distribution<double> phiDist(-std::numbers::pi,
0071                                                std::numbers::pi);
0072 // Track theta distribution
0073 std::uniform_real_distribution<double> thetaDist(1., std::numbers::pi - 1.);
0074 // Track charge helper distribution
0075 std::uniform_real_distribution<double> qDist(-1, 1);
0076 // Track IP resolution distribution
0077 std::uniform_real_distribution<double> resIPDist(0., 100_um);
0078 // Track angular distribution
0079 std::uniform_real_distribution<double> resAngDist(0., 0.1);
0080 // Track q/p resolution distribution
0081 std::uniform_real_distribution<double> resQoPDist(-0.01, 0.01);
0082 
0083 BOOST_AUTO_TEST_SUITE(VertexingSuite)
0084 ///
0085 /// @brief Unit test for ZScanVertexFinder
0086 ///
0087 BOOST_AUTO_TEST_CASE(zscan_finder_test) {
0088   unsigned int nTests = 50;
0089 
0090   for (unsigned int iTest = 0; iTest < nTests; ++iTest) {
0091     // Number of tracks
0092     unsigned int nTracks = 30;
0093 
0094     // Set up RNG
0095     int mySeed = 31415;
0096     std::mt19937 gen(mySeed);
0097 
0098     // Set up constant B-Field
0099     auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0100 
0101     // Set up Eigenstepper
0102     EigenStepper<> stepper(bField);
0103 
0104     // Set up propagator with void navigator
0105     auto propagator = std::make_shared<Propagator>(stepper);
0106 
0107     // Create perigee surface
0108     std::shared_ptr<PerigeeSurface> perigeeSurface =
0109         Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0110 
0111     // Create position of vertex and perigee surface
0112     double x = vXYDist(gen);
0113     double y = vXYDist(gen);
0114     double z = vZDist(gen);
0115 
0116     // Calculate d0 and z0 corresponding to vertex position
0117     double d0_v = std::hypot(x, y);
0118     double z0_v = z;
0119 
0120     // Start constructing nTracks tracks in the following
0121     std::vector<BoundTrackParameters> tracks;
0122 
0123     // Construct random track emerging from vicinity of vertex position
0124     // Vector to store track objects used for vertex fit
0125     for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
0126       // Construct positive or negative charge randomly
0127       double q = std::copysign(1., qDist(gen));
0128 
0129       // Construct random track parameters
0130       BoundVector paramVec = BoundVector::Zero();
0131       paramVec[eBoundLoc0] = d0_v + d0Dist(gen);
0132       paramVec[eBoundLoc1] = z0_v + z0Dist(gen);
0133       paramVec[eBoundPhi] = phiDist(gen);
0134       paramVec[eBoundTheta] = thetaDist(gen);
0135       paramVec[eBoundQOverP] = q / pTDist(gen);
0136 
0137       // Resolutions
0138       double resD0 = resIPDist(gen);
0139       double resZ0 = resIPDist(gen);
0140       double resPh = resAngDist(gen);
0141       double resTh = resAngDist(gen);
0142       double resQp = resQoPDist(gen);
0143 
0144       // Fill vector of track objects with simple covariance matrix
0145       Covariance covMat;
0146 
0147       covMat << resD0 * resD0, 0., 0., 0., 0., 0., 0., resZ0 * resZ0, 0., 0.,
0148           0., 0., 0., 0., resPh * resPh, 0., 0., 0., 0., 0., 0., resTh * resTh,
0149           0., 0., 0., 0., 0., 0., resQp * resQp, 0., 0., 0., 0., 0., 0., 1.;
0150 
0151       tracks.emplace_back(perigeeSurface, paramVec, std::move(covMat),
0152                           ParticleHypothesis::pion());
0153     }
0154 
0155     std::vector<InputTrack> inputTracks;
0156     for (const auto& trk : tracks) {
0157       inputTracks.emplace_back(&trk);
0158     }
0159 
0160     using VertexFinder = ZScanVertexFinder;
0161 
0162     ImpactPointEstimator::Config ipEstimatorCfg(bField, propagator);
0163     ImpactPointEstimator ipEstimator(ipEstimatorCfg);
0164 
0165     VertexFinder::Config cfg(ipEstimator);
0166     cfg.extractParameters.connect<&InputTrack::extractParameters>();
0167 
0168     VertexFinder finder(cfg);
0169 
0170     VertexingOptions vertexingOptions(geoContext, magFieldContext);
0171 
0172     auto state = finder.makeState(magFieldContext);
0173     auto res = finder.find(inputTracks, vertexingOptions, state);
0174 
0175     BOOST_CHECK(res.ok());
0176 
0177     if (!res.ok()) {
0178       std::cout << res.error().message() << std::endl;
0179     }
0180 
0181     if (res.ok()) {
0182       BOOST_CHECK(!(*res).empty());
0183       Vector3 result = (*res).back().position();
0184       CHECK_CLOSE_ABS(result[eZ], z, 1_mm);
0185     }
0186   }
0187 }
0188 
0189 // Dummy user-defined InputTrackStub type
0190 struct InputTrackStub {
0191   explicit InputTrackStub(const BoundTrackParameters& params)
0192       : m_parameters(params) {}
0193 
0194   const BoundTrackParameters& parameters() const { return m_parameters; }
0195 
0196   // store e.g. link to original objects here
0197 
0198  private:
0199   BoundTrackParameters m_parameters;
0200 };
0201 
0202 ///
0203 /// @brief Unit test for ZScanVertexFinder with user-defined input track type
0204 ///
0205 BOOST_AUTO_TEST_CASE(zscan_finder_usertrack_test) {
0206   unsigned int nTests = 50;
0207 
0208   for (unsigned int iTest = 0; iTest < nTests; ++iTest) {
0209     // Number of tracks
0210     unsigned int nTracks = 30;
0211 
0212     // Set up RNG
0213     int mySeed = 31415;
0214     std::mt19937 gen(mySeed);
0215 
0216     // Set up constant B-Field
0217     auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0218 
0219     // Set up Eigenstepper
0220     EigenStepper<> stepper(bField);
0221 
0222     // Set up propagator with void navigator
0223     auto propagator = std::make_shared<Propagator>(stepper);
0224 
0225     // Create perigee surface
0226     std::shared_ptr<PerigeeSurface> perigeeSurface =
0227         Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0228 
0229     // Create position of vertex and perigee surface
0230     double x = vXYDist(gen);
0231     double y = vXYDist(gen);
0232     double z = vZDist(gen);
0233 
0234     // Calculate d0 and z0 corresponding to vertex position
0235     double d0_v = std::hypot(x, y);
0236     double z0_v = z;
0237 
0238     // Start constructing nTracks tracks in the following
0239     std::vector<InputTrackStub> tracks;
0240 
0241     // Construct random track emerging from vicinity of vertex position
0242     // Vector to store track objects used for vertex fit
0243     for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
0244       // Construct positive or negative charge randomly
0245       double q = std::copysign(1., qDist(gen));
0246 
0247       // Construct random track parameters
0248       BoundVector paramVec;
0249       double z0track = z0_v + z0Dist(gen);
0250       paramVec << d0_v + d0Dist(gen), z0track, phiDist(gen), thetaDist(gen),
0251           q / pTDist(gen), 0.;
0252 
0253       // Resolutions
0254       double resD0 = resIPDist(gen);
0255       double resZ0 = resIPDist(gen);
0256       double resPh = resAngDist(gen);
0257       double resTh = resAngDist(gen);
0258       double resQp = resQoPDist(gen);
0259 
0260       // Fill vector of track objects with simple covariance matrix
0261       Covariance covMat;
0262 
0263       covMat << resD0 * resD0, 0., 0., 0., 0., 0., 0., resZ0 * resZ0, 0., 0.,
0264           0., 0., 0., 0., resPh * resPh, 0., 0., 0., 0., 0., 0., resTh * resTh,
0265           0., 0., 0., 0., 0., 0., resQp * resQp, 0., 0., 0., 0., 0., 0., 1.;
0266 
0267       tracks.emplace_back(BoundTrackParameters(perigeeSurface, paramVec,
0268                                                std::move(covMat),
0269                                                ParticleHypothesis::pion()));
0270     }
0271 
0272     std::vector<InputTrack> inputTracks;
0273     for (const auto& trk : tracks) {
0274       inputTracks.emplace_back(&trk);
0275     }
0276 
0277     using VertexFinder = ZScanVertexFinder;
0278 
0279     ImpactPointEstimator::Config ipEstimatorCfg(bField, propagator);
0280     ImpactPointEstimator ipEstimator(ipEstimatorCfg);
0281 
0282     VertexFinder::Config cfg(ipEstimator);
0283 
0284     // Create a custom std::function to extract BoundTrackParameters from
0285     // user-defined InputTrackStub
0286     auto extractParameters = [](const InputTrack& params) {
0287       return params.as<InputTrackStub>()->parameters();
0288     };
0289 
0290     cfg.extractParameters.connect(extractParameters);
0291     VertexFinder finder(cfg);
0292     auto state = finder.makeState(magFieldContext);
0293 
0294     VertexingOptions vertexingOptions(geoContext, magFieldContext);
0295 
0296     auto res = finder.find(inputTracks, vertexingOptions, state);
0297 
0298     BOOST_CHECK(res.ok());
0299 
0300     if (!res.ok()) {
0301       std::cout << res.error().message() << std::endl;
0302     }
0303 
0304     if (res.ok()) {
0305       BOOST_CHECK(!(*res).empty());
0306       Vector3 result = (*res).back().position();
0307       CHECK_CLOSE_ABS(result[eZ], z, 1_mm);
0308     }
0309   }
0310 }
0311 
0312 BOOST_AUTO_TEST_SUITE_END()
0313 
0314 }  // namespace ActsTests