Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-10-22 07:53:36

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/Direction.hpp"
0013 #include "Acts/Definitions/TrackParametrization.hpp"
0014 #include "Acts/Definitions/Units.hpp"
0015 #include "Acts/EventData/GenericBoundTrackParameters.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/HelicalTrackLinearizer.hpp"
0027 #include "Acts/Vertexing/ImpactPointEstimator.hpp"
0028 #include "Acts/Vertexing/KalmanVertexUpdater.hpp"
0029 #include "Acts/Vertexing/LinearizedTrack.hpp"
0030 #include "Acts/Vertexing/TrackAtVertex.hpp"
0031 #include "Acts/Vertexing/Vertex.hpp"
0032 
0033 #include <algorithm>
0034 #include <array>
0035 #include <cmath>
0036 #include <iostream>
0037 #include <memory>
0038 #include <numbers>
0039 #include <optional>
0040 #include <random>
0041 #include <tuple>
0042 #include <type_traits>
0043 #include <utility>
0044 #include <vector>
0045 
0046 using namespace Acts;
0047 using namespace Acts::UnitLiterals;
0048 
0049 namespace ActsTests {
0050 
0051 using Covariance = BoundSquareMatrix;
0052 using Propagator = Acts::Propagator<EigenStepper<>>;
0053 using Linearizer = HelicalTrackLinearizer;
0054 
0055 // Create a test context
0056 GeometryContext geoContext = GeometryContext();
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 // Number of vertices per test event distribution
0083 
0084 BOOST_AUTO_TEST_SUITE(VertexingSuite)
0085 ///
0086 /// @brief Unit test for KalmanVertexUpdater
0087 ///
0088 BOOST_AUTO_TEST_CASE(Kalman_Vertex_Updater) {
0089   bool debug = false;
0090 
0091   // Number of tests
0092   unsigned int nTests = 10;
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   // Linearizer for BoundTrackParameters type test
0108   Linearizer::Config ltConfig;
0109   ltConfig.bField = bField;
0110   ltConfig.propagator = propagator;
0111   Linearizer linearizer(ltConfig);
0112   auto fieldCache = bField->makeCache(magFieldContext);
0113 
0114   // Create perigee surface at origin
0115   std::shared_ptr<PerigeeSurface> perigeeSurface =
0116       Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0117 
0118   // Creates a random tracks around origin and a random vertex.
0119   // VertexUpdater adds track to vertex and updates the position
0120   // which should afterwards be closer to the origin/track
0121   for (unsigned int i = 0; i < nTests; ++i) {
0122     if (debug) {
0123       std::cout << "Test " << i + 1 << std::endl;
0124     }
0125     // Construct positive or negative charge randomly
0126     double q = qDist(gen) < 0 ? -1. : 1.;
0127 
0128     // Construct random track parameters around origin
0129     BoundTrackParameters::ParametersVector paramVec;
0130 
0131     paramVec << d0Dist(gen), z0Dist(gen), phiDist(gen), thetaDist(gen),
0132         q / pTDist(gen), 0.;
0133 
0134     if (debug) {
0135       std::cout << "Creating track parameters: " << paramVec << std::endl;
0136     }
0137 
0138     // Fill vector of track objects with simple covariance matrix
0139     Covariance covMat;
0140 
0141     // Resolutions
0142     double res_d0 = resIPDist(gen);
0143     double res_z0 = resIPDist(gen);
0144     double res_ph = resAngDist(gen);
0145     double res_th = resAngDist(gen);
0146     double res_qp = resQoPDist(gen);
0147 
0148     covMat << res_d0 * res_d0, 0., 0., 0., 0., 0., 0., res_z0 * res_z0, 0., 0.,
0149         0., 0., 0., 0., res_ph * res_ph, 0., 0., 0., 0., 0., 0.,
0150         res_th * res_th, 0., 0., 0., 0., 0., 0., res_qp * res_qp, 0., 0., 0.,
0151         0., 0., 0., 1.;
0152     BoundTrackParameters params(perigeeSurface, paramVec, std::move(covMat),
0153                                 ParticleHypothesis::pion());
0154 
0155     std::shared_ptr<PerigeeSurface> perigee =
0156         Surface::makeShared<PerigeeSurface>(Vector3::Zero());
0157 
0158     // Linearized state of the track
0159     LinearizedTrack linTrack =
0160         linearizer
0161             .linearizeTrack(params, 0, *perigee, geoContext, magFieldContext,
0162                             fieldCache)
0163             .value();
0164 
0165     // Create TrackAtVertex
0166     TrackAtVertex trkAtVtx(0., params, InputTrack{&params});
0167 
0168     // Set linearized state of trackAtVertex
0169     trkAtVtx.linearizedState = linTrack;
0170     trkAtVtx.isLinearized = true;
0171 
0172     // Create a vertex
0173     Vector3 vtxPos(vXYDist(gen), vXYDist(gen), vZDist(gen));
0174     Vertex vtx(vtxPos);
0175     vtx.setFullCovariance(SquareMatrix4::Identity() * 0.01);
0176 
0177     // Update trkAtVertex with assumption of originating from vtx
0178     KalmanVertexUpdater::updateVertexWithTrack(vtx, trkAtVtx, 3);
0179 
0180     if (debug) {
0181       std::cout << "Old vertex position: " << vtxPos << std::endl;
0182       std::cout << "New vertex position: " << vtx.position() << std::endl;
0183     }
0184 
0185     double oldDistance = vtxPos.norm();
0186     double newDistance = vtx.position().norm();
0187 
0188     if (debug) {
0189       std::cout << "Old distance: " << oldDistance << std::endl;
0190       std::cout << "New distance: " << newDistance << std::endl;
0191     }
0192 
0193     // After update, vertex should be closer to the track
0194     BOOST_CHECK_LT(newDistance, oldDistance);
0195 
0196     // Note: KalmanVertexUpdater updates the vertex w.r.t. the
0197     // newly given track, but does NOT add the track to the
0198     // TrackAtVertex list. Has to be done manually after calling
0199     // the update method.
0200     BOOST_CHECK(vtx.tracks().empty());
0201 
0202   }  // end for loop
0203 
0204 }  // end test case
0205 
0206 ///
0207 /// @brief Unit test for KalmanVertexTrackUpdater
0208 ///
0209 BOOST_AUTO_TEST_CASE(Kalman_Vertex_TrackUpdater) {
0210   bool debug = true;
0211 
0212   // Number of tests
0213   unsigned int nTests = 10;
0214 
0215   // Set up RNG
0216   int mySeed = 31415;
0217   std::mt19937 gen(mySeed);
0218 
0219   // Set up constant B-Field
0220   auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0221 
0222   // Set up Eigenstepper
0223   EigenStepper<> stepper(bField);
0224 
0225   // Set up propagator with void navigator
0226   auto propagator = std::make_shared<Propagator>(stepper);
0227 
0228   // Set up ImpactPointEstimator, used for comparisons later
0229   ImpactPointEstimator::Config ip3dEstConfig(bField, propagator);
0230   ImpactPointEstimator ip3dEst(ip3dEstConfig);
0231   ImpactPointEstimator::State state{bField->makeCache(magFieldContext)};
0232 
0233   // Set up HelicalTrackLinearizer, needed for linearizing the tracks
0234   // Linearizer for BoundTrackParameters type test
0235   Linearizer::Config ltConfig;
0236   ltConfig.bField = bField;
0237   ltConfig.propagator = propagator;
0238   Linearizer linearizer(ltConfig);
0239   auto fieldCache = bField->makeCache(magFieldContext);
0240 
0241   // Create perigee surface at origin
0242   std::shared_ptr<PerigeeSurface> perigeeSurface =
0243       Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0244 
0245   // Create random tracks around origin and a random vertex.
0246   // Update tracks with the assumption that they originate from
0247   // the vertex position and check if they are closer to the
0248   // vertex after the update process
0249   for (unsigned int i = 0; i < nTests; ++i) {
0250     // Construct positive or negative charge randomly
0251     double q = qDist(gen) < 0 ? -1. : 1.;
0252 
0253     // Construct random track parameters
0254     BoundTrackParameters::ParametersVector paramVec;
0255 
0256     paramVec << d0Dist(gen), z0Dist(gen), phiDist(gen), thetaDist(gen),
0257         q / pTDist(gen), 0.;
0258 
0259     if (debug) {
0260       std::cout << "Creating track parameters: " << paramVec << std::endl;
0261     }
0262 
0263     // Fill vector of track objects with simple covariance matrix
0264     Covariance covMat;
0265 
0266     // Resolutions
0267     double res_d0 = resIPDist(gen);
0268     double res_z0 = resIPDist(gen);
0269     double res_ph = resAngDist(gen);
0270     double res_th = resAngDist(gen);
0271     double res_qp = resQoPDist(gen);
0272 
0273     covMat << res_d0 * res_d0, 0., 0., 0., 0., 0., 0., res_z0 * res_z0, 0., 0.,
0274         0., 0., 0., 0., res_ph * res_ph, 0., 0., 0., 0., 0., 0.,
0275         res_th * res_th, 0., 0., 0., 0., 0., 0., res_qp * res_qp, 0., 0., 0.,
0276         0., 0., 0., 1.;
0277     BoundTrackParameters params(perigeeSurface, paramVec, std::move(covMat),
0278                                 ParticleHypothesis::pion());
0279 
0280     std::shared_ptr<PerigeeSurface> perigee =
0281         Surface::makeShared<PerigeeSurface>(Vector3::Zero());
0282 
0283     // Linearized state of the track
0284     LinearizedTrack linTrack =
0285         linearizer
0286             .linearizeTrack(params, 0, *perigee, geoContext, magFieldContext,
0287                             fieldCache)
0288             .value();
0289 
0290     // Create TrackAtVertex
0291     TrackAtVertex trkAtVtx(0., params, InputTrack{&params});
0292 
0293     // Set linearized state of trackAtVertex
0294     trkAtVtx.linearizedState = linTrack;
0295     trkAtVtx.isLinearized = true;
0296 
0297     // Copy parameters for later comparison of old and new version
0298     auto fittedParamsCopy = trkAtVtx.fittedParams;
0299 
0300     // Create a vertex
0301     Vector3 vtxPos(vXYDist(gen), vXYDist(gen), vZDist(gen));
0302     Vertex vtx(vtxPos);
0303 
0304     // Update trkAtVertex with assumption of originating from vtx
0305     KalmanVertexUpdater::updateTrackWithVertex(trkAtVtx, vtx, 3);
0306 
0307     // The old distance
0308     double oldDistance =
0309         ip3dEst.calculateDistance(geoContext, fittedParamsCopy, vtxPos, state)
0310             .value();
0311 
0312     // The new distance after update
0313     double newDistance =
0314         ip3dEst
0315             .calculateDistance(geoContext, trkAtVtx.fittedParams, vtxPos, state)
0316             .value();
0317     if (debug) {
0318       std::cout << "Old distance: " << oldDistance << std::endl;
0319       std::cout << "New distance: " << newDistance << std::endl;
0320     }
0321 
0322     // Parameters should have changed
0323     BOOST_CHECK_NE(fittedParamsCopy, trkAtVtx.fittedParams);
0324 
0325     // After update, track should be closer to the vertex
0326     BOOST_CHECK_LT(newDistance, oldDistance);
0327 
0328   }  // end for loop
0329 
0330 }  // end test case
0331 
0332 BOOST_AUTO_TEST_SUITE_END()
0333 
0334 }  // namespace ActsTests