Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-10-18 08:22:59

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