File indexing completed on 2025-01-18 09:13:03
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/Tests/CommonHelpers/FloatComparisons.hpp"
0026 #include "Acts/Utilities/Result.hpp"
0027 #include "Acts/Vertexing/FullBilloirVertexFitter.hpp"
0028 #include "Acts/Vertexing/HelicalTrackLinearizer.hpp"
0029 #include "Acts/Vertexing/ImpactPointEstimator.hpp"
0030 #include "Acts/Vertexing/VertexingOptions.hpp"
0031 #include "Acts/Vertexing/ZScanVertexFinder.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::UnitLiterals;
0049
0050 namespace Acts::Test {
0051
0052 using Covariance = BoundSquareMatrix;
0053 using Propagator = Acts::Propagator<EigenStepper<>>;
0054
0055
0056 GeometryContext geoContext = GeometryContext();
0057 MagneticFieldContext magFieldContext = MagneticFieldContext();
0058
0059
0060 std::uniform_real_distribution<double> vXYDist(-0.1_mm, 0.1_mm);
0061
0062 std::uniform_real_distribution<double> vZDist(-20_mm, 20_mm);
0063
0064 std::uniform_real_distribution<double> d0Dist(-0.01_mm, 0.01_mm);
0065
0066 std::uniform_real_distribution<double> z0Dist(-0.2_mm, 0.2_mm);
0067
0068 std::uniform_real_distribution<double> pTDist(0.4_GeV, 10_GeV);
0069
0070 std::uniform_real_distribution<double> phiDist(-std::numbers::pi,
0071 std::numbers::pi);
0072
0073 std::uniform_real_distribution<double> thetaDist(1., std::numbers::pi - 1.);
0074
0075 std::uniform_real_distribution<double> qDist(-1, 1);
0076
0077 std::uniform_real_distribution<double> resIPDist(0., 100_um);
0078
0079 std::uniform_real_distribution<double> resAngDist(0., 0.1);
0080
0081 std::uniform_real_distribution<double> resQoPDist(-0.01, 0.01);
0082
0083
0084
0085
0086 BOOST_AUTO_TEST_CASE(zscan_finder_test) {
0087 unsigned int nTests = 50;
0088
0089 for (unsigned int iTest = 0; iTest < nTests; ++iTest) {
0090
0091 unsigned int nTracks = 30;
0092
0093
0094 int mySeed = 31415;
0095 std::mt19937 gen(mySeed);
0096
0097
0098 auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0099
0100
0101 EigenStepper<> stepper(bField);
0102
0103
0104 auto propagator = std::make_shared<Propagator>(stepper);
0105
0106
0107 std::shared_ptr<PerigeeSurface> perigeeSurface =
0108 Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0109
0110
0111 double x = vXYDist(gen);
0112 double y = vXYDist(gen);
0113 double z = vZDist(gen);
0114
0115
0116 double d0_v = std::hypot(x, y);
0117 double z0_v = z;
0118
0119
0120 std::vector<BoundTrackParameters> tracks;
0121
0122
0123
0124 for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
0125
0126 double q = qDist(gen) < 0 ? -1. : 1.;
0127
0128
0129 BoundVector paramVec = BoundVector::Zero();
0130 paramVec[eBoundLoc0] = d0_v + d0Dist(gen);
0131 paramVec[eBoundLoc1] = z0_v + z0Dist(gen);
0132 paramVec[eBoundPhi] = phiDist(gen);
0133 paramVec[eBoundTheta] = thetaDist(gen);
0134 paramVec[eBoundQOverP] = q / pTDist(gen);
0135
0136
0137 double resD0 = resIPDist(gen);
0138 double resZ0 = resIPDist(gen);
0139 double resPh = resAngDist(gen);
0140 double resTh = resAngDist(gen);
0141 double resQp = resQoPDist(gen);
0142
0143
0144 Covariance covMat;
0145
0146 covMat << resD0 * resD0, 0., 0., 0., 0., 0., 0., resZ0 * resZ0, 0., 0.,
0147 0., 0., 0., 0., resPh * resPh, 0., 0., 0., 0., 0., 0., resTh * resTh,
0148 0., 0., 0., 0., 0., 0., resQp * resQp, 0., 0., 0., 0., 0., 0., 1.;
0149
0150 tracks.emplace_back(perigeeSurface, paramVec, std::move(covMat),
0151 ParticleHypothesis::pion());
0152 }
0153
0154 std::vector<InputTrack> inputTracks;
0155 for (const auto& trk : tracks) {
0156 inputTracks.emplace_back(&trk);
0157 }
0158
0159 using VertexFinder = ZScanVertexFinder;
0160
0161 ImpactPointEstimator::Config ipEstimatorCfg(bField, propagator);
0162 ImpactPointEstimator ipEstimator(ipEstimatorCfg);
0163
0164 VertexFinder::Config cfg(ipEstimator);
0165 cfg.extractParameters.connect<&InputTrack::extractParameters>();
0166
0167 VertexFinder finder(cfg);
0168
0169 VertexingOptions vertexingOptions(geoContext, magFieldContext);
0170
0171 auto state = finder.makeState(magFieldContext);
0172 auto res = finder.find(inputTracks, vertexingOptions, state);
0173
0174 BOOST_CHECK(res.ok());
0175
0176 if (!res.ok()) {
0177 std::cout << res.error().message() << std::endl;
0178 }
0179
0180 if (res.ok()) {
0181 BOOST_CHECK(!(*res).empty());
0182 Vector3 result = (*res).back().position();
0183 CHECK_CLOSE_ABS(result[eZ], z, 1_mm);
0184 }
0185 }
0186 }
0187
0188
0189 struct InputTrackStub {
0190 InputTrackStub(const BoundTrackParameters& params) : m_parameters(params) {}
0191
0192 const BoundTrackParameters& parameters() const { return m_parameters; }
0193
0194
0195
0196 private:
0197 BoundTrackParameters m_parameters;
0198 };
0199
0200
0201
0202
0203 BOOST_AUTO_TEST_CASE(zscan_finder_usertrack_test) {
0204 unsigned int nTests = 50;
0205
0206 for (unsigned int iTest = 0; iTest < nTests; ++iTest) {
0207
0208 unsigned int nTracks = 30;
0209
0210
0211 int mySeed = 31415;
0212 std::mt19937 gen(mySeed);
0213
0214
0215 auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0216
0217
0218 EigenStepper<> stepper(bField);
0219
0220
0221 auto propagator = std::make_shared<Propagator>(stepper);
0222
0223
0224 std::shared_ptr<PerigeeSurface> perigeeSurface =
0225 Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0226
0227
0228 double x = vXYDist(gen);
0229 double y = vXYDist(gen);
0230 double z = vZDist(gen);
0231
0232
0233 double d0_v = std::hypot(x, y);
0234 double z0_v = z;
0235
0236
0237 std::vector<InputTrackStub> tracks;
0238
0239
0240
0241 for (unsigned int iTrack = 0; iTrack < nTracks; iTrack++) {
0242
0243 double q = qDist(gen) < 0 ? -1. : 1.;
0244
0245
0246 BoundVector paramVec;
0247 double z0track = z0_v + z0Dist(gen);
0248 paramVec << d0_v + d0Dist(gen), z0track, phiDist(gen), thetaDist(gen),
0249 q / pTDist(gen), 0.;
0250
0251
0252 double resD0 = resIPDist(gen);
0253 double resZ0 = resIPDist(gen);
0254 double resPh = resAngDist(gen);
0255 double resTh = resAngDist(gen);
0256 double resQp = resQoPDist(gen);
0257
0258
0259 Covariance covMat;
0260
0261 covMat << resD0 * resD0, 0., 0., 0., 0., 0., 0., resZ0 * resZ0, 0., 0.,
0262 0., 0., 0., 0., resPh * resPh, 0., 0., 0., 0., 0., 0., resTh * resTh,
0263 0., 0., 0., 0., 0., 0., resQp * resQp, 0., 0., 0., 0., 0., 0., 1.;
0264
0265 tracks.emplace_back(BoundTrackParameters(perigeeSurface, paramVec,
0266 std::move(covMat),
0267 ParticleHypothesis::pion()));
0268 }
0269
0270 std::vector<InputTrack> inputTracks;
0271 for (const auto& trk : tracks) {
0272 inputTracks.emplace_back(&trk);
0273 }
0274
0275 using VertexFinder = ZScanVertexFinder;
0276
0277 ImpactPointEstimator::Config ipEstimatorCfg(bField, propagator);
0278 ImpactPointEstimator ipEstimator(ipEstimatorCfg);
0279
0280 VertexFinder::Config cfg(ipEstimator);
0281
0282
0283
0284 auto extractParameters = [](const InputTrack& params) {
0285 return params.as<InputTrackStub>()->parameters();
0286 };
0287
0288 cfg.extractParameters.connect(extractParameters);
0289 VertexFinder finder(cfg);
0290 auto state = finder.makeState(magFieldContext);
0291
0292 VertexingOptions vertexingOptions(geoContext, magFieldContext);
0293
0294 auto res = finder.find(inputTracks, vertexingOptions, state);
0295
0296 BOOST_CHECK(res.ok());
0297
0298 if (!res.ok()) {
0299 std::cout << res.error().message() << std::endl;
0300 }
0301
0302 if (res.ok()) {
0303 BOOST_CHECK(!(*res).empty());
0304 Vector3 result = (*res).back().position();
0305 CHECK_CLOSE_ABS(result[eZ], z, 1_mm);
0306 }
0307 }
0308 }
0309
0310 }