File indexing completed on 2025-01-18 09:13:02
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/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::UnitLiterals;
0047
0048 namespace Acts::Test {
0049
0050 using Covariance = BoundSquareMatrix;
0051 using Propagator = Acts::Propagator<EigenStepper<>>;
0052 using Linearizer = HelicalTrackLinearizer;
0053
0054
0055 GeometryContext geoContext = GeometryContext();
0056 MagneticFieldContext magFieldContext = MagneticFieldContext();
0057
0058
0059 std::uniform_real_distribution<double> vXYDist(-0.1_mm, 0.1_mm);
0060
0061 std::uniform_real_distribution<double> vZDist(-20_mm, 20_mm);
0062
0063 std::uniform_real_distribution<double> d0Dist(-0.01_mm, 0.01_mm);
0064
0065 std::uniform_real_distribution<double> z0Dist(-0.2_mm, 0.2_mm);
0066
0067 std::uniform_real_distribution<double> pTDist(0.4_GeV, 10_GeV);
0068
0069 std::uniform_real_distribution<double> phiDist(-std::numbers::pi,
0070 std::numbers::pi);
0071
0072 std::uniform_real_distribution<double> thetaDist(1., std::numbers::pi - 1.);
0073
0074 std::uniform_real_distribution<double> qDist(-1, 1);
0075
0076 std::uniform_real_distribution<double> resIPDist(0., 100_um);
0077
0078 std::uniform_real_distribution<double> resAngDist(0., 0.1);
0079
0080 std::uniform_real_distribution<double> resQoPDist(-0.01, 0.01);
0081
0082
0083
0084
0085
0086 BOOST_AUTO_TEST_CASE(Kalman_Vertex_Updater) {
0087 bool debug = false;
0088
0089
0090 unsigned int nTests = 10;
0091
0092
0093 int mySeed = 31415;
0094 std::mt19937 gen(mySeed);
0095
0096
0097 auto bField = std::make_shared<ConstantBField>(Vector3{0.0, 0.0, 1_T});
0098
0099
0100 EigenStepper<> stepper(bField);
0101
0102
0103 auto propagator = std::make_shared<Propagator>(stepper);
0104
0105
0106 Linearizer::Config ltConfig;
0107 ltConfig.bField = bField;
0108 ltConfig.propagator = propagator;
0109 Linearizer linearizer(ltConfig);
0110 auto fieldCache = bField->makeCache(magFieldContext);
0111
0112
0113 std::shared_ptr<PerigeeSurface> perigeeSurface =
0114 Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0115
0116
0117
0118
0119 for (unsigned int i = 0; i < nTests; ++i) {
0120 if (debug) {
0121 std::cout << "Test " << i + 1 << std::endl;
0122 }
0123
0124 double q = qDist(gen) < 0 ? -1. : 1.;
0125
0126
0127 BoundTrackParameters::ParametersVector paramVec;
0128
0129 paramVec << d0Dist(gen), z0Dist(gen), phiDist(gen), thetaDist(gen),
0130 q / pTDist(gen), 0.;
0131
0132 if (debug) {
0133 std::cout << "Creating track parameters: " << paramVec << std::endl;
0134 }
0135
0136
0137 Covariance covMat;
0138
0139
0140 double res_d0 = resIPDist(gen);
0141 double res_z0 = resIPDist(gen);
0142 double res_ph = resAngDist(gen);
0143 double res_th = resAngDist(gen);
0144 double res_qp = resQoPDist(gen);
0145
0146 covMat << res_d0 * res_d0, 0., 0., 0., 0., 0., 0., res_z0 * res_z0, 0., 0.,
0147 0., 0., 0., 0., res_ph * res_ph, 0., 0., 0., 0., 0., 0.,
0148 res_th * res_th, 0., 0., 0., 0., 0., 0., res_qp * res_qp, 0., 0., 0.,
0149 0., 0., 0., 1.;
0150 BoundTrackParameters params(perigeeSurface, paramVec, std::move(covMat),
0151 ParticleHypothesis::pion());
0152
0153 std::shared_ptr<PerigeeSurface> perigee =
0154 Surface::makeShared<PerigeeSurface>(Vector3::Zero());
0155
0156
0157 LinearizedTrack linTrack =
0158 linearizer
0159 .linearizeTrack(params, 0, *perigee, geoContext, magFieldContext,
0160 fieldCache)
0161 .value();
0162
0163
0164 TrackAtVertex trkAtVtx(0., params, InputTrack{¶ms});
0165
0166
0167 trkAtVtx.linearizedState = linTrack;
0168 trkAtVtx.isLinearized = true;
0169
0170
0171 Vector3 vtxPos(vXYDist(gen), vXYDist(gen), vZDist(gen));
0172 Vertex vtx(vtxPos);
0173 vtx.setFullCovariance(SquareMatrix4::Identity() * 0.01);
0174
0175
0176 KalmanVertexUpdater::updateVertexWithTrack(vtx, trkAtVtx, 3);
0177
0178 if (debug) {
0179 std::cout << "Old vertex position: " << vtxPos << std::endl;
0180 std::cout << "New vertex position: " << vtx.position() << std::endl;
0181 }
0182
0183 double oldDistance = vtxPos.norm();
0184 double newDistance = vtx.position().norm();
0185
0186 if (debug) {
0187 std::cout << "Old distance: " << oldDistance << std::endl;
0188 std::cout << "New distance: " << newDistance << std::endl;
0189 }
0190
0191
0192 BOOST_CHECK_LT(newDistance, oldDistance);
0193
0194
0195
0196
0197
0198 BOOST_CHECK(vtx.tracks().empty());
0199
0200 }
0201
0202 }
0203
0204
0205
0206
0207 BOOST_AUTO_TEST_CASE(Kalman_Vertex_TrackUpdater) {
0208 bool debug = true;
0209
0210
0211 unsigned int nTests = 10;
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 ImpactPointEstimator::Config ip3dEstConfig(bField, propagator);
0228 ImpactPointEstimator ip3dEst(ip3dEstConfig);
0229 ImpactPointEstimator::State state{bField->makeCache(magFieldContext)};
0230
0231
0232
0233 Linearizer::Config ltConfig;
0234 ltConfig.bField = bField;
0235 ltConfig.propagator = propagator;
0236 Linearizer linearizer(ltConfig);
0237 auto fieldCache = bField->makeCache(magFieldContext);
0238
0239
0240 std::shared_ptr<PerigeeSurface> perigeeSurface =
0241 Surface::makeShared<PerigeeSurface>(Vector3(0., 0., 0.));
0242
0243
0244
0245
0246
0247 for (unsigned int i = 0; i < nTests; ++i) {
0248
0249 double q = qDist(gen) < 0 ? -1. : 1.;
0250
0251
0252 BoundTrackParameters::ParametersVector paramVec;
0253
0254 paramVec << d0Dist(gen), z0Dist(gen), phiDist(gen), thetaDist(gen),
0255 q / pTDist(gen), 0.;
0256
0257 if (debug) {
0258 std::cout << "Creating track parameters: " << paramVec << std::endl;
0259 }
0260
0261
0262 Covariance covMat;
0263
0264
0265 double res_d0 = resIPDist(gen);
0266 double res_z0 = resIPDist(gen);
0267 double res_ph = resAngDist(gen);
0268 double res_th = resAngDist(gen);
0269 double res_qp = resQoPDist(gen);
0270
0271 covMat << res_d0 * res_d0, 0., 0., 0., 0., 0., 0., res_z0 * res_z0, 0., 0.,
0272 0., 0., 0., 0., res_ph * res_ph, 0., 0., 0., 0., 0., 0.,
0273 res_th * res_th, 0., 0., 0., 0., 0., 0., res_qp * res_qp, 0., 0., 0.,
0274 0., 0., 0., 1.;
0275 BoundTrackParameters params(perigeeSurface, paramVec, std::move(covMat),
0276 ParticleHypothesis::pion());
0277
0278 std::shared_ptr<PerigeeSurface> perigee =
0279 Surface::makeShared<PerigeeSurface>(Vector3::Zero());
0280
0281
0282 LinearizedTrack linTrack =
0283 linearizer
0284 .linearizeTrack(params, 0, *perigee, geoContext, magFieldContext,
0285 fieldCache)
0286 .value();
0287
0288
0289 TrackAtVertex trkAtVtx(0., params, InputTrack{¶ms});
0290
0291
0292 trkAtVtx.linearizedState = linTrack;
0293 trkAtVtx.isLinearized = true;
0294
0295
0296 auto fittedParamsCopy = trkAtVtx.fittedParams;
0297
0298
0299 Vector3 vtxPos(vXYDist(gen), vXYDist(gen), vZDist(gen));
0300 Vertex vtx(vtxPos);
0301
0302
0303 KalmanVertexUpdater::updateTrackWithVertex(trkAtVtx, vtx, 3);
0304
0305
0306 double oldDistance =
0307 ip3dEst.calculateDistance(geoContext, fittedParamsCopy, vtxPos, state)
0308 .value();
0309
0310
0311 double newDistance =
0312 ip3dEst
0313 .calculateDistance(geoContext, trkAtVtx.fittedParams, vtxPos, state)
0314 .value();
0315 if (debug) {
0316 std::cout << "Old distance: " << oldDistance << std::endl;
0317 std::cout << "New distance: " << newDistance << std::endl;
0318 }
0319
0320
0321 BOOST_CHECK_NE(fittedParamsCopy, trkAtVtx.fittedParams);
0322
0323
0324 BOOST_CHECK_LT(newDistance, oldDistance);
0325
0326 }
0327
0328 }
0329
0330 }