Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-07-12 07:52:57

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 "Acts/EventData/SubspaceHelpers.hpp"
0010 #include "Acts/EventData/TrackContainer.hpp"
0011 #include "Acts/EventData/TrackStatePropMask.hpp"
0012 #include "Acts/EventData/TrackStateType.hpp"
0013 #include "Acts/EventData/VectorMultiTrajectory.hpp"
0014 #include "Acts/EventData/VectorTrackContainer.hpp"
0015 #include "Acts/EventData/detail/GenerateParameters.hpp"
0016 #include "Acts/Geometry/GeometryIdentifier.hpp"
0017 #include "Acts/Surfaces/PerigeeSurface.hpp"
0018 #include "Acts/Surfaces/PlaneSurface.hpp"
0019 #include "Acts/Surfaces/RectangleBounds.hpp"
0020 #include "Acts/Utilities/TrackHelpers.hpp"
0021 
0022 #include <iostream>
0023 #include <numeric>
0024 #include <random>
0025 #include <type_traits>
0026 
0027 using namespace Acts;
0028 
0029 class BenchmarkSourceLink final {
0030  public:
0031   using Index = std::uint32_t;
0032 
0033   /// Construct from geometry identifier and index.
0034   constexpr BenchmarkSourceLink(Acts::GeometryIdentifier gid, Index idx)
0035       : m_geometryId(gid), m_index(idx) {}
0036 
0037   BenchmarkSourceLink() = default;
0038   BenchmarkSourceLink(const BenchmarkSourceLink&) = default;
0039   BenchmarkSourceLink(BenchmarkSourceLink&&) = default;
0040   BenchmarkSourceLink& operator=(const BenchmarkSourceLink&) = default;
0041   BenchmarkSourceLink& operator=(BenchmarkSourceLink&&) = default;
0042 
0043   /// Access the index.
0044   constexpr Index index() const { return m_index; }
0045 
0046   Acts::GeometryIdentifier geometryId() const { return m_geometryId; }
0047 
0048  private:
0049   Acts::GeometryIdentifier m_geometryId;
0050   Index m_index = 0;
0051 
0052   friend bool operator==(const BenchmarkSourceLink& lhs,
0053                          const BenchmarkSourceLink& rhs) {
0054     return (lhs.geometryId() == rhs.geometryId()) &&
0055            (lhs.m_index == rhs.m_index);
0056   }
0057 };
0058 
0059 int main(int /*argc*/, char** /*argv[]*/) {
0060   std::size_t runs = 1000000;
0061   std::size_t nTracks = 10000;
0062 
0063   VectorMultiTrajectory mtj;
0064   VectorTrackContainer vtc;
0065   TrackContainer tc{vtc, mtj};
0066 
0067   VectorMultiTrajectory mtjOut;
0068   VectorTrackContainer vtcOut;
0069   TrackContainer output{vtcOut, mtjOut};
0070 
0071   auto gid = GeometryIdentifier().withVolume(5).withLayer(3).withSensitive(1);
0072 
0073   static_assert(sizeof(BenchmarkSourceLink) <= ACTS_SOURCELINK_SBO_SIZE);
0074 
0075   static_assert(std::is_trivially_move_constructible_v<BenchmarkSourceLink>);
0076 
0077   std::uniform_int_distribution<> nStatesDist(1, 20);
0078   std::uniform_int_distribution<> measDimDist(1, 3);
0079   std::uniform_real_distribution<> typeDist(0, 1);
0080   std::uniform_real_distribution<> copyDist(0, 1);
0081   std::mt19937 rng{42};
0082 
0083   std::vector<std::shared_ptr<Surface>> surfaces;
0084   std::vector<std::pair<BoundVector, BoundMatrix>> parametersVector;
0085   for (std::size_t s = 0; s < 50; ++s) {
0086     surfaces.push_back(Surface::makeShared<PlaneSurface>(
0087         Transform3::Identity(), std::make_shared<RectangleBounds>(50, 50)));
0088 
0089     parametersVector.push_back(
0090         detail::Test::generateBoundParametersCovariance(rng, {}));
0091   }
0092 
0093   std::size_t nSurface = 0;
0094   auto surface = [&]() {
0095     nSurface++;
0096     return surfaces.at(nSurface % surfaces.size());
0097   };
0098 
0099   std::size_t nParams = 0;
0100   auto parameters = [&]() -> const std::pair<BoundVector, BoundMatrix>& {
0101     nParams++;
0102     return parametersVector.at(nParams % parametersVector.size());
0103   };
0104 
0105   auto perigee = Surface::makeShared<PerigeeSurface>(Vector3::Zero());
0106 
0107   std::cout << "Creating " << nTracks << " tracks x " << runs << " runs"
0108             << std::endl;
0109   for (std::size_t r = 0; r < runs; ++r) {
0110     tc.clear();
0111     output.clear();
0112 
0113     for (std::size_t i = 0; i < nTracks; ++i) {
0114       auto track = tc.makeTrack();
0115 
0116       std::size_t nStates = nStatesDist(rng);
0117 
0118       for (std::size_t j = 0; j < nStates; ++j) {
0119         auto trackState = track.appendTrackState(TrackStatePropMask::All);
0120         trackState.setReferenceSurface(surface());
0121 
0122         trackState.jacobian().setZero();
0123         trackState.jacobian().row(j % eBoundSize).setOnes();
0124 
0125         double crit = typeDist(rng);
0126 
0127         if (crit < 0.1) {
0128           // hole
0129           trackState.typeFlags().set(TrackStateFlag::HoleFlag);
0130         } else if (crit < 0.2) {
0131           // material
0132           trackState.typeFlags().set(TrackStateFlag::MaterialFlag);
0133         } else {
0134           BenchmarkSourceLink bsl{gid, 123};
0135           std::size_t measdim = measDimDist(rng);
0136 
0137           const auto& [predicted, covariance] = parameters();
0138           trackState.predicted() = predicted;
0139           trackState.predictedCovariance() = covariance;
0140 
0141           visit_measurement(
0142               measdim,
0143               [&]<std::size_t N>(std::integral_constant<std::size_t, N> /*d*/) {
0144                 trackState.allocateCalibrated(ActsVector<N>::Ones(),
0145                                               ActsSquareMatrix<N>::Identity());
0146 
0147                 std::array<std::uint8_t, eBoundSize> indices{0};
0148                 std::iota(indices.begin(), indices.end(), 0);
0149                 trackState.setProjectorSubspaceIndices(indices);
0150               });
0151 
0152           trackState.typeFlags().set(TrackStateFlag::MeasurementFlag);
0153           if (crit < 0.4) {
0154             // outlier
0155             trackState.typeFlags().set(TrackStateFlag::OutlierFlag);
0156           }
0157         }
0158       }
0159 
0160       track.setReferenceSurface(perigee);
0161 
0162       const auto& [ref, cov] = parameters();
0163       track.parameters() = ref;
0164       track.covariance() = cov;
0165 
0166       track.linkForward();
0167 
0168       calculateTrackQuantities(track);
0169     }
0170 
0171     for (const auto& track : tc) {
0172       if (copyDist(rng) > 0.1) {
0173         // copy only 10% of tracks
0174         continue;
0175       }
0176 
0177       auto target = output.makeTrack();
0178       target.copyFrom(track);
0179     }
0180   }
0181 
0182   return 0;
0183 }