File indexing completed on 2025-10-18 08:22:59
0001
0002
0003
0004
0005
0006
0007
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
0057 GeometryContext geoContext = GeometryContext();
0058 MagneticFieldContext magFieldContext = MagneticFieldContext();
0059
0060
0061 std::uniform_real_distribution<double> vXYDist(-0.1_mm, 0.1_mm);
0062
0063 std::uniform_real_distribution<double> vZDist(-20_mm, 20_mm);
0064
0065 std::uniform_real_distribution<double> d0Dist(-0.01_mm, 0.01_mm);
0066
0067 std::uniform_real_distribution<double> z0Dist(-0.2_mm, 0.2_mm);
0068
0069 std::uniform_real_distribution<double> pTDist(0.4_GeV, 10_GeV);
0070
0071 std::uniform_real_distribution<double> phiDist(-std::numbers::pi,
0072 std::numbers::pi);
0073
0074 std::uniform_real_distribution<double> thetaDist(1., std::numbers::pi - 1.);
0075
0076 std::uniform_real_distribution<double> qDist(-1, 1);
0077
0078 std::uniform_real_distribution<double> resIPDist(0., 100_um);
0079
0080 std::uniform_real_distribution<double> resAngDist(0., 0.1);
0081
0082 std::uniform_real_distribution<double> resQoPDist(-0.01, 0.01);
0083
0084 BOOST_AUTO_TEST_SUITE(VertexingSuite)
0085
0086
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
0093 unsigned int nTracks = 30;
0094
0095
0096 int mySeed = 31415;
0097 std::mt19937 gen(mySeed);
0098
0099
0100 auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0101
0102
0103 EigenStepper<> stepper(bField);
0104
0105
0106 auto propagator = std::make_shared<Propagator>(stepper);
0107
0108
0109 std::shared_ptr<PerigeeSurface> perigeeSurface =
0110 Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0111
0112
0113 double x = vXYDist(gen);
0114 double y = vXYDist(gen);
0115 double z = vZDist(gen);
0116
0117
0118 double d0_v = std::hypot(x, y);
0119 double z0_v = z;
0120
0121
0122 std::vector<BoundTrackParameters> tracks;
0123
0124
0125
0126 for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
0127
0128 double q = qDist(gen) < 0 ? -1. : 1.;
0129
0130
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
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
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
0191 struct InputTrackStub {
0192 explicit InputTrackStub(const BoundTrackParameters& params)
0193 : m_parameters(params) {}
0194
0195 const BoundTrackParameters& parameters() const { return m_parameters; }
0196
0197
0198
0199 private:
0200 BoundTrackParameters m_parameters;
0201 };
0202
0203
0204
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
0211 unsigned int nTracks = 30;
0212
0213
0214 int mySeed = 31415;
0215 std::mt19937 gen(mySeed);
0216
0217
0218 auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0219
0220
0221 EigenStepper<> stepper(bField);
0222
0223
0224 auto propagator = std::make_shared<Propagator>(stepper);
0225
0226
0227 std::shared_ptr<PerigeeSurface> perigeeSurface =
0228 Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0229
0230
0231 double x = vXYDist(gen);
0232 double y = vXYDist(gen);
0233 double z = vZDist(gen);
0234
0235
0236 double d0_v = std::hypot(x, y);
0237 double z0_v = z;
0238
0239
0240 std::vector<InputTrackStub> tracks;
0241
0242
0243
0244 for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
0245
0246 double q = qDist(gen) < 0 ? -1. : 1.;
0247
0248
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
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
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
0286
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 }