Examples/Framework/ParticleDeltaRTests.cpp

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/Utilities/VectorHelpers.hpp"
0012 #include "ActsExamples/Validation/PatternRecognitionPerformanceCollector.hpp"
0013
0014 #include <array>
0015 #include <cmath>
0016 #include <limits>
0017 #include <memory>
0018 #include <numbers>
0019 #include <random>
0020 #include <vector>
0021
0022 using namespace Acts;
0023 using namespace ActsExamples;
0024
0025 namespace {
0026
0027 void checkAgainstBruteForce(const std::vector<std::array<double, 2>>& points,
0028                             double maxDeltaR = 0.3) {
0029   SimParticleContainer particles;
0030   ParticleTrackMatching matching;
0031   TrackParticleMatching trackMatching;
0032   auto trackBackend = std::make_shared<VectorTrackContainer>();
0033   auto stateBackend = std::make_shared<VectorMultiTrajectory>();
0034   ActsExamples::TrackContainer tracks(trackBackend, stateBackend);
0035   for (std::size_t i = 0; i < points.size(); ++i) {
0036     auto barcode = SimBarcode{}.withVertexPrimary(1).withParticle(i + 1);
0037     SimParticle particle(barcode, PdgParticle::eMuon);
0038     const auto [eta, phi] = points[i];
0039     particle.initialState().setDirection(
0040         std::isinf(eta)
0041             ? Vector3(0., 0., std::copysign(1., eta))
0042             : Vector3(std::cos(phi), std::sin(phi), std::sinh(eta)));
0043     particle.initialState().setAbsoluteMomentum(10.);
0044     particles.insert(particle);
0045     if (i % 2 == 0) {
0046       auto track = tracks.makeTrack();
0047       matching[barcode].track = track.index();
0048       trackMatching[track.index()] = {
0049           TrackMatchClassification::Matched, barcode, {}};
0050     }
0051   }
0052   ConstTrackContainer constTracks(
0053       std::make_shared<ConstVectorTrackContainer>(*trackBackend),
0054       std::make_shared<ConstVectorMultiTrajectory>(*stateBackend));
0055   PatternRecognitionPerformanceCollector::Config config;
0056   config.effPlotToolConfig.varBinning["DeltaR"] =
0057       EffPlotTool::BoostRegularAxis(100, 0., maxDeltaR, "deltaR");
0058   PatternRecognitionPerformanceCollector collector(
0059       config, getDefaultLogger("Test", Logging::WARNING));
0060   EffPlotTool reference(config.effPlotToolConfig, Logging::WARNING);
0061   auto gctx = GeometryContext::dangerouslyDefaultConstruct();
0062   collector.fill(gctx, constTracks, particles, trackMatching, matching, {});
0063   for (const auto& particle : particles) {
0064     double nearest = -1.;
0065     for (const auto& other : particles) {
0066       if (particle.particleId() == other.particleId()) {
0067         continue;
0068       }
0069       const double distance =
0070           VectorHelpers::deltaR(particle.direction(), other.direction());
0071       if (nearest == -1. || distance < nearest) {
0072         nearest = distance;
0073       }
0074     }
0075     reference.fill(gctx, particle.initialState(), nearest,
0076                    matching.contains(particle.particleId()));
0077   }
0078   for (const auto& [name, expected] : reference.efficiencies1D()) {
0079     BOOST_TEST_CONTEXT(name) {
0080       const auto& actual = collector.effPlotTool().efficiencies1D().at(name);
0081       for (int bin = -1; bin <= expected.totalHistogram().axis(0).size();
0082            ++bin) {
0083         BOOST_CHECK_EQUAL(actual.totalHistogram().at(bin).value(),
0084                           expected.totalHistogram().at(bin).value());
0085         BOOST_CHECK_EQUAL(actual.acceptedHistogram().at(bin).value(),
0086                           expected.acceptedHistogram().at(bin).value());
0087       }
0088     }
0089   }
0090 }
0091
0092 }  // namespace
0093
0094 BOOST_AUTO_TEST_CASE(EmptySingletonAndIsolatedParticles) {
0095   checkAgainstBruteForce({});
0096   checkAgainstBruteForce({{0., 0.}});
0097   checkAgainstBruteForce({{-5., -2.}, {0., 0.}, {5., 2.}});
0098 }
0099
0100 BOOST_AUTO_TEST_CASE(PeriodicPhiCoincidenceAndCutoff) {
0101   const double pi = std::numbers::pi;
0102   checkAgainstBruteForce({{0., pi - 0.001},
0103                           {0., -pi + 0.001},
0104                           {1., 0.},
0105                           {1., 0.},
0106                           {2., 0.},
0107                           {2., 0.3},
0108                           {3., 0.},
0109                           {3., std::nextafter(0.3, 0.)},
0110                           {4., 0.},
0111                           {4., std::nextafter(0.3, 1.)}});
0112 }
0113
0114 BOOST_AUTO_TEST_CASE(ConfiguredRangeAndRandomParticles) {
0115   std::mt19937 rng(42);
0116   std::uniform_real_distribution<double> eta(-6., 6.);
0117   std::uniform_real_distribution<double> phi(-std::numbers::pi,
0118                                              std::numbers::pi);
0119   std::vector<std::array<double, 2>> points;
0120   for (std::size_t i = 0; i < 512; ++i) {
0121     points.push_back({eta(rng), phi(rng)});
0122   }
0123   for (double maxDeltaR : {0.01, 0.3, 1., 10.}) {
0124     BOOST_TEST_CONTEXT(maxDeltaR) {
0125       checkAgainstBruteForce(points, maxDeltaR);
0126     }
0127   }
0128 }
0129
0130 BOOST_AUTO_TEST_CASE(BeamAxisParticles) {
0131   const double infinity = std::numeric_limits<double>::infinity();
0132   checkAgainstBruteForce({{infinity, 0.}, {0., 0.}, {0., 0.1}});
0133   checkAgainstBruteForce({{-infinity, 0.}, {0., 0.}});
0134 }