From 8009432fc99d63ed48d3193812da6101271e4eb4 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Fri, 22 Nov 2024 16:10:12 +0100 Subject: [PATCH 1/8] EigenStepperDefaultExtension --- .../Propagator/EigenStepperDefaultExtension.hpp | 14 +++++--------- 1 file changed, 5 insertions(+), 9 deletions(-) diff --git a/Core/include/Acts/Propagator/EigenStepperDefaultExtension.hpp b/Core/include/Acts/Propagator/EigenStepperDefaultExtension.hpp index 2810b11a2a2..cb7939534f5 100644 --- a/Core/include/Acts/Propagator/EigenStepperDefaultExtension.hpp +++ b/Core/include/Acts/Propagator/EigenStepperDefaultExtension.hpp @@ -16,14 +16,10 @@ namespace Acts { -/// @brief Default evaluater of the k_i's and elements of the transport matrix +/// @brief Default evaluator of the k_i's and elements of the transport matrix /// D of the RKN4 stepping. This is a pure implementation by textbook. struct EigenStepperDefaultExtension { - using Scalar = ActsScalar; - /// @brief Vector3 replacement for the custom scalar type - using ThisVector3 = Eigen::Matrix; - - /// @brief Evaluater of the k_i's of the RKN4. For the case of i = 0 this + /// @brief Evaluator of the k_i's of the RKN4. For the case of i = 0 this /// step sets up qop, too. /// /// @tparam i Index of the k_i, i = [0, 3] @@ -43,9 +39,9 @@ struct EigenStepperDefaultExtension { template bool k(const propagator_state_t& state, const stepper_t& stepper, - const navigator_t& /*navigator*/, ThisVector3& knew, - const Vector3& bField, std::array& kQoP, - const double h = 0., const ThisVector3& kprev = ThisVector3::Zero()) + const navigator_t& /*navigator*/, Vector3& knew, const Vector3& bField, + std::array& kQoP, const double h = 0., + const Vector3& kprev = Vector3::Zero()) requires(i >= 0 && i <= 3) { auto qop = stepper.qOverP(state.stepping); From f832715ed28e8620fc3b18bde39008a4b5d355d9 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Fri, 22 Nov 2024 16:16:36 +0100 Subject: [PATCH 2/8] EigenStepperDenseExtension --- .../Propagator/EigenStepperDenseExtension.hpp | 42 +++++++++---------- 1 file changed, 19 insertions(+), 23 deletions(-) diff --git a/Core/include/Acts/Propagator/EigenStepperDenseExtension.hpp b/Core/include/Acts/Propagator/EigenStepperDenseExtension.hpp index 70dbb1cb820..a3d88ea1997 100644 --- a/Core/include/Acts/Propagator/EigenStepperDenseExtension.hpp +++ b/Core/include/Acts/Propagator/EigenStepperDenseExtension.hpp @@ -21,44 +21,40 @@ namespace Acts { -/// @brief Evaluater of the k_i's and elements of the transport matrix +/// @brief Evaluator of the k_i's and elements of the transport matrix /// D of the RKN4 stepping. This implementation involves energy loss due to -/// ioninisation, bremsstrahlung, pair production and photonuclear interaction -/// in the propagation and the jacobian. These effects will only occur if the +/// ionisation, bremsstrahlung, pair production and photonuclear interaction +/// in the propagation and the Jacobian. These effects will only occur if the /// propagation is in a TrackingVolume with attached material. struct EigenStepperDenseExtension { - using Scalar = ActsScalar; - /// @brief Vector3 replacement for the custom scalar type - using ThisVector3 = Eigen::Matrix; - /// Fallback extension EigenStepperDefaultExtension defaultExtension; /// Momentum at a certain point - Scalar currentMomentum = 0.; + double currentMomentum = 0.; /// Particles momentum at k1 - Scalar initialMomentum = 0.; + double initialMomentum = 0.; /// Material that will be passed /// TODO : Might not be needed anymore Material material; /// Derivatives dLambda''dlambda at each sub-step point - std::array dLdl{}; + std::array dLdl{}; /// q/p at each sub-step - std::array qop{}; + std::array qop{}; /// Derivatives dPds at each sub-step - std::array dPds{}; + std::array dPds{}; /// Derivative d(dEds)d(q/p) evaluated at the initial point - Scalar dgdqopValue = 0.; + double dgdqopValue = 0.; /// Derivative dEds at the initial point - Scalar g = 0.; + double g = 0.; /// k_i equivalent for the time propagation - std::array tKi{}; + std::array tKi{}; /// Lambda''_i - std::array Lambdappi{}; + std::array Lambdappi{}; /// Energy at each sub-step - std::array energy{}; + std::array energy{}; - /// @brief Evaluater of the k_i's of the RKN4. For the case of i = 0 this + /// @brief Evaluator of the k_i's of the RKN4. For the case of i = 0 this /// step sets up member parameters, too. /// /// @tparam i Index of the k_i, i = [0, 3] @@ -79,9 +75,9 @@ struct EigenStepperDenseExtension { template bool k(const propagator_state_t& state, const stepper_t& stepper, - const navigator_t& navigator, ThisVector3& knew, const Vector3& bField, - std::array& kQoP, const double h = 0., - const ThisVector3& kprev = ThisVector3::Zero()) + const navigator_t& navigator, Vector3& knew, const Vector3& bField, + std::array& kQoP, const double h = 0., + const Vector3& kprev = Vector3::Zero()) requires(i >= 0 && i <= 3) { const auto* volumeMaterial = @@ -98,7 +94,7 @@ struct EigenStepperDenseExtension { // i = 0 is used for setup and evaluation of k if constexpr (i == 0) { // Set up for energy loss - ThisVector3 position = stepper.position(state.stepping); + Vector3 position = stepper.position(state.stepping); material = volumeMaterial->material(position.template cast()); initialMomentum = stepper.absoluteMomentum(state.stepping); currentMomentum = initialMomentum; @@ -217,7 +213,7 @@ struct EigenStepperDenseExtension { } private: - /// @brief Evaluates the transport matrix D for the jacobian + /// @brief Evaluates the transport matrix D for the Jacobian /// /// @tparam propagator_state_t Type of the state of the propagator /// @tparam stepper_t Type of the stepper From bee7017943e43f27e66cd13de3f445720222ebe5 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Fri, 22 Nov 2024 16:33:03 +0100 Subject: [PATCH 3/8] part 3 --- .../ParametricParticleGenerator.cpp | 2 +- .../ActsExamples/EventData/DriftCircle.hpp | 22 +++++----- .../ActsExamples/EventData/SimParticle.hpp | 12 ++--- .../ActsExamples/EventData/SimSpacePoint.hpp | 44 +++++++++---------- .../ActsExamples/EventData/SimVertex.hpp | 9 ++-- 5 files changed, 39 insertions(+), 50 deletions(-) diff --git a/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp b/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp index 34363410f42..6e244bcdf5f 100644 --- a/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp +++ b/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp @@ -93,7 +93,7 @@ ParametricParticleGenerator::operator()(RandomEngine& rng) { // create the primary vertex auto& primaryVertex = vertices.emplace_back( - SimVertexBarcode{0}, SimVertex::Vector4(0., 0., 0., 0.)); + SimVertexBarcode{0}, Acts::Vector4(0., 0., 0., 0.)); // counter will be reused as barcode particle number which must be non-zero. for (std::size_t ip = 1; ip <= m_cfg.numParticles; ++ip) { diff --git a/Examples/Framework/include/ActsExamples/EventData/DriftCircle.hpp b/Examples/Framework/include/ActsExamples/EventData/DriftCircle.hpp index 88cad84fffe..3d9ab8c241f 100644 --- a/Examples/Framework/include/ActsExamples/EventData/DriftCircle.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/DriftCircle.hpp @@ -23,8 +23,6 @@ namespace ActsExamples { /// representation of a drift circle measurement used for track finding class DriftCircle { - using Scalar = Acts::ActsScalar; - public: /// Construct the drift circle from the drift radius and tube location /// @@ -52,11 +50,11 @@ class DriftCircle { m_tubeLayer(tubeLayer), m_tube(tube) {} - constexpr Scalar x() const { return m_x; } - constexpr Scalar y() const { return m_y; } - constexpr Scalar z() const { return m_z; } - constexpr Scalar rDrift() const { return m_rho; } - constexpr Scalar rDriftError() const { return m_sigmaRho; } + constexpr double x() const { return m_x; } + constexpr double y() const { return m_y; } + constexpr double z() const { return m_z; } + constexpr double rDrift() const { return m_rho; } + constexpr double rDriftError() const { return m_sigmaRho; } constexpr int stationName() const { return m_stationName; } constexpr int stationEta() const { return m_stationEta; } constexpr int stationPhi() const { return m_stationPhi; } @@ -66,11 +64,11 @@ class DriftCircle { private: // Global position - Scalar m_x = 0.0f; - Scalar m_y = 0.0f; - Scalar m_z = 0.0f; - Scalar m_rho = 0.0f; - Scalar m_sigmaRho = 0.0f; + double m_x = 0.; + double m_y = 0.; + double m_z = 0.; + double m_rho = 0.; + double m_sigmaRho = 0.; int m_stationName = 0; int m_stationEta = 0; int m_stationPhi = 0; diff --git a/Examples/Framework/include/ActsExamples/EventData/SimParticle.hpp b/Examples/Framework/include/ActsExamples/EventData/SimParticle.hpp index cfcf8a5e25e..2bf9dcc8cbc 100644 --- a/Examples/Framework/include/ActsExamples/EventData/SimParticle.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/SimParticle.hpp @@ -23,10 +23,6 @@ using SimParticleState = ::ActsFatras::Particle; class SimParticle final { public: - using Scalar = Acts::ActsScalar; - using Vector3 = Acts::ActsVector<3>; - using Vector4 = Acts::ActsVector<4>; - /// Construct a default particle with invalid identity. SimParticle() = default; @@ -39,8 +35,8 @@ class SimParticle final { /// /// @warning It is the users responsibility that charge and mass match /// the PDG particle number. - SimParticle(SimBarcode particleId, Acts::PdgParticle pdg, Scalar charge, - Scalar mass) + SimParticle(SimBarcode particleId, Acts::PdgParticle pdg, double charge, + double mass) : m_initial(particleId, pdg, charge, mass), m_final(particleId, pdg, charge, mass) {} @@ -89,13 +85,13 @@ class SimParticle final { return *this; } /// Set the charge. - SimParticle& setCharge(Scalar charge) { + SimParticle& setCharge(double charge) { initial().setCharge(charge); final().setCharge(charge); return *this; } /// Set the mass. - SimParticle& setMass(Scalar mass) { + SimParticle& setMass(double mass) { initial().setMass(mass); final().setMass(mass); return *this; diff --git a/Examples/Framework/include/ActsExamples/EventData/SimSpacePoint.hpp b/Examples/Framework/include/ActsExamples/EventData/SimSpacePoint.hpp index 0359913199f..bf5780c5a76 100644 --- a/Examples/Framework/include/ActsExamples/EventData/SimSpacePoint.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/SimSpacePoint.hpp @@ -23,8 +23,6 @@ namespace ActsExamples { /// Space point representation of a measurement suitable for track seeding. class SimSpacePoint { - using Scalar = Acts::ActsScalar; - public: /// Construct the space point from global position and selected variances. /// @@ -44,10 +42,10 @@ class SimSpacePoint { /// @param validDoubleMeasurementDetails boolean to check if double measurements are valid template SimSpacePoint( - const Eigen::MatrixBase& pos, std::optional t, - Scalar varRho, Scalar varZ, std::optional varT, + const Eigen::MatrixBase& pos, std::optional t, + double varRho, double varZ, std::optional varT, boost::container::static_vector sourceLinks, - Scalar topHalfStripLength, Scalar bottomHalfStripLength, + double topHalfStripLength, double bottomHalfStripLength, const Acts::Vector3& topStripDirection, const Acts::Vector3& bottomStripDirection, const Acts::Vector3& stripCenterDistance, @@ -82,8 +80,8 @@ class SimSpacePoint { /// @param sourceLinks sourceLinks of the measurements template SimSpacePoint( - const Eigen::MatrixBase& pos, std::optional t, - Scalar varRho, Scalar varZ, std::optional varT, + const Eigen::MatrixBase& pos, std::optional t, + double varRho, double varZ, std::optional varT, boost::container::static_vector sourceLinks) : m_x(pos[Acts::ePos0]), m_y(pos[Acts::ePos1]), @@ -97,14 +95,14 @@ class SimSpacePoint { EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(position_t, 3); } - constexpr Scalar x() const { return m_x; } - constexpr Scalar y() const { return m_y; } - constexpr Scalar z() const { return m_z; } - constexpr std::optional t() const { return m_t; } - constexpr Scalar r() const { return m_rho; } - constexpr Scalar varianceR() const { return m_varianceRho; } - constexpr Scalar varianceZ() const { return m_varianceZ; } - constexpr std::optional varianceT() const { return m_varianceT; } + constexpr double x() const { return m_x; } + constexpr double y() const { return m_y; } + constexpr double z() const { return m_z; } + constexpr std::optional t() const { return m_t; } + constexpr double r() const { return m_rho; } + constexpr double varianceR() const { return m_varianceRho; } + constexpr double varianceZ() const { return m_varianceZ; } + constexpr std::optional varianceT() const { return m_varianceT; } const boost::container::static_vector& sourceLinks() const { @@ -127,15 +125,15 @@ class SimSpacePoint { private: // Global position - Scalar m_x; - Scalar m_y; - Scalar m_z; - std::optional m_t; - Scalar m_rho; + double m_x; + double m_y; + double m_z; + std::optional m_t; + double m_rho; // Variance in rho/z of the global coordinates - Scalar m_varianceRho; - Scalar m_varianceZ; - std::optional m_varianceT; + double m_varianceRho; + double m_varianceZ; + std::optional m_varianceT; // SourceLinks of the corresponding measurements. A Pixel (strip) SP has one // (two) sourceLink(s). boost::container::static_vector m_sourceLinks; diff --git a/Examples/Framework/include/ActsExamples/EventData/SimVertex.hpp b/Examples/Framework/include/ActsExamples/EventData/SimVertex.hpp index 0bdfbb4925c..db906790e8b 100644 --- a/Examples/Framework/include/ActsExamples/EventData/SimVertex.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/SimVertex.hpp @@ -77,13 +77,10 @@ class SimVertexBarcode { /// A simulated vertex e.g. from a physics process. struct SimVertex { - using Scalar = Acts::ActsScalar; - using Vector4 = Acts::ActsVector<4>; - /// The vertex ID SimVertexBarcode id; /// The vertex four-position - Vector4 position4 = Vector4::Zero(); + Acts::Vector4 position4 = Acts::Vector4::Zero(); /// The vertex process type ActsFatras::ProcessType process = ActsFatras::ProcessType::eUndefined; /// The incoming particles into the vertex @@ -99,7 +96,7 @@ struct SimVertex { /// Associated particles are left empty by default and must be filled by the /// user after construction. SimVertex( - SimVertexBarcode id_, const Vector4& position4_, + SimVertexBarcode id_, const Acts::Vector4& position4_, ActsFatras::ProcessType process_ = ActsFatras::ProcessType::eUndefined) : id(id_), position4(position4_), process(process_) {} // explicitly default rule-of-five. @@ -113,7 +110,7 @@ struct SimVertex { /// The vertex three-position. auto position() const { return position4.head<3>(); } /// The vertex time. - Scalar time() const { return position4[3]; } + double time() const { return position4[3]; } }; namespace detail { From 5cfde4a41a27fd6603af1332fdc0eba5d3fbe95e Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Fri, 22 Nov 2024 17:04:23 +0100 Subject: [PATCH 4/8] step 4 --- .../ParametricParticleGenerator.cpp | 4 +- Examples/Io/Csv/src/CsvDriftCircleReader.cpp | 7 +- Examples/Io/Csv/src/CsvMuonSimHitReader.cpp | 22 +++-- Examples/Io/Csv/src/CsvSimHitReader.cpp | 6 +- .../include/ActsFatras/Geant4/Geant4Decay.hpp | 18 ++-- Fatras/Geant4/src/Geant4Decay.cpp | 2 +- .../Digitization/DigitizationData.hpp | 5 +- .../Digitization/UncorrelatedHitSmearer.hpp | 2 - Fatras/include/ActsFatras/EventData/Hit.hpp | 30 +++---- .../include/ActsFatras/EventData/Particle.hpp | 86 +++++++++---------- .../ActsFatras/Kernel/InteractionList.hpp | 10 +-- .../include/ActsFatras/Kernel/Simulation.hpp | 4 +- .../ActsFatras/Kernel/SimulationResult.hpp | 7 +- .../Kernel/detail/SimulationActor.hpp | 4 +- .../ActsFatras/Physics/Decay/NoDecay.hpp | 4 +- .../Physics/ElectroMagnetic/BetheHeitler.hpp | 11 +-- .../NuclearInteraction/NuclearInteraction.hpp | 29 +++---- .../Selectors/ParticleSelectors.hpp | 8 +- Fatras/src/Physics/BetheHeitler.cpp | 12 +-- .../NuclearInteraction/NuclearInteraction.cpp | 10 +-- Fatras/src/Physics/PhotonConversion.cpp | 5 +- Tests/UnitTests/Fatras/EventData/HitTests.cpp | 64 +++++++------- .../Fatras/EventData/ParticleTests.cpp | 36 ++++---- .../Fatras/Kernel/InteractionListTests.cpp | 23 +++-- .../Fatras/Kernel/SimulationActorTests.cpp | 18 ++-- .../Fatras/Physics/PhotonConversionTests.cpp | 28 +++--- 26 files changed, 211 insertions(+), 244 deletions(-) diff --git a/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp b/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp index 6e244bcdf5f..00411fc64c3 100644 --- a/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp +++ b/Examples/Algorithms/Generators/ActsExamples/Generators/ParametricParticleGenerator.cpp @@ -92,8 +92,8 @@ ParametricParticleGenerator::operator()(RandomEngine& rng) { SimParticleContainer::sequence_type particles; // create the primary vertex - auto& primaryVertex = vertices.emplace_back( - SimVertexBarcode{0}, Acts::Vector4(0., 0., 0., 0.)); + auto& primaryVertex = + vertices.emplace_back(SimVertexBarcode{0}, Acts::Vector4(0., 0., 0., 0.)); // counter will be reused as barcode particle number which must be non-zero. for (std::size_t ip = 1; ip <= m_cfg.numParticles; ++ip) { diff --git a/Examples/Io/Csv/src/CsvDriftCircleReader.cpp b/Examples/Io/Csv/src/CsvDriftCircleReader.cpp index ce179d9d5be..ad2500637ed 100644 --- a/Examples/Io/Csv/src/CsvDriftCircleReader.cpp +++ b/Examples/Io/Csv/src/CsvDriftCircleReader.cpp @@ -61,10 +61,9 @@ ActsExamples::ProcessCode ActsExamples::CsvDriftCircleReader::read( MuonDriftCircleData data; while (reader.read(data)) { - ActsFatras::Hit::Vector3 tube_pos{ - data.tubePositionx * Acts::UnitConstants::mm, - data.tubePositiony * Acts::UnitConstants::mm, - data.tubePositionz * Acts::UnitConstants::mm}; + Acts::Vector3 tube_pos{data.tubePositionx * Acts::UnitConstants::mm, + data.tubePositiony * Acts::UnitConstants::mm, + data.tubePositionz * Acts::UnitConstants::mm}; DriftCircles.push_back(DriftCircle(std::move(tube_pos), data.driftRadius, 0.0f, data.stationName, data.stationEta, diff --git a/Examples/Io/Csv/src/CsvMuonSimHitReader.cpp b/Examples/Io/Csv/src/CsvMuonSimHitReader.cpp index 58240dcecce..28836dc48de 100644 --- a/Examples/Io/Csv/src/CsvMuonSimHitReader.cpp +++ b/Examples/Io/Csv/src/CsvMuonSimHitReader.cpp @@ -61,18 +61,16 @@ ActsExamples::ProcessCode ActsExamples::CsvMuonSimHitReader::read( SimHitContainer::sequence_type unordered; while (reader.read(data)) { - ActsFatras::Hit::Vector4 pos{ - data.LocalPositionExtrx * Acts::UnitConstants::mm, - data.LocalPositionExtry * Acts::UnitConstants::mm, - data.LocalPositionExtrz * Acts::UnitConstants::mm, 0}; - ActsFatras::Hit::Vector4 mom{ - data.LocalDirectionx * Acts::UnitConstants::GeV, - data.LocalDirectiony * Acts::UnitConstants::GeV, - data.LocalDirectionz * Acts::UnitConstants::GeV, - std::sqrt(data.LocalDirectionx * data.LocalDirectionx + - data.LocalDirectiony * data.LocalDirectiony + - data.LocalDirectionz * data.LocalDirectionz) * - Acts::UnitConstants::GeV}; + Acts::Vector4 pos{data.LocalPositionExtrx * Acts::UnitConstants::mm, + data.LocalPositionExtry * Acts::UnitConstants::mm, + data.LocalPositionExtrz * Acts::UnitConstants::mm, 0}; + Acts::Vector4 mom{data.LocalDirectionx * Acts::UnitConstants::GeV, + data.LocalDirectiony * Acts::UnitConstants::GeV, + data.LocalDirectionz * Acts::UnitConstants::GeV, + std::sqrt(data.LocalDirectionx * data.LocalDirectionx + + data.LocalDirectiony * data.LocalDirectiony + + data.LocalDirectionz * data.LocalDirectionz) * + Acts::UnitConstants::GeV}; MuonMdtIdentifierFields f; f.multilayer = 0; f.tube = 0; diff --git a/Examples/Io/Csv/src/CsvSimHitReader.cpp b/Examples/Io/Csv/src/CsvSimHitReader.cpp index 949b6d4ce5c..ca77c7e0aa4 100644 --- a/Examples/Io/Csv/src/CsvSimHitReader.cpp +++ b/Examples/Io/Csv/src/CsvSimHitReader.cpp @@ -66,19 +66,19 @@ ActsExamples::ProcessCode ActsExamples::CsvSimHitReader::read( // TODO validate geo id consistency const auto particleId = ActsFatras::Barcode(data.particle_id); - ActsFatras::Hit::Vector4 pos4{ + Acts::Vector4 pos4{ data.tx * Acts::UnitConstants::mm, data.ty * Acts::UnitConstants::mm, data.tz * Acts::UnitConstants::mm, data.tt * Acts::UnitConstants::mm, }; - ActsFatras::Hit::Vector4 mom4{ + Acts::Vector4 mom4{ data.tpx * Acts::UnitConstants::GeV, data.tpy * Acts::UnitConstants::GeV, data.tpz * Acts::UnitConstants::GeV, data.te * Acts::UnitConstants::GeV, }; - ActsFatras::Hit::Vector4 delta4{ + Acts::Vector4 delta4{ data.deltapx * Acts::UnitConstants::GeV, data.deltapy * Acts::UnitConstants::GeV, data.deltapz * Acts::UnitConstants::GeV, diff --git a/Fatras/Geant4/include/ActsFatras/Geant4/Geant4Decay.hpp b/Fatras/Geant4/include/ActsFatras/Geant4/Geant4Decay.hpp index 177e3ac6221..758df351071 100644 --- a/Fatras/Geant4/include/ActsFatras/Geant4/Geant4Decay.hpp +++ b/Fatras/Geant4/include/ActsFatras/Geant4/Geant4Decay.hpp @@ -27,8 +27,6 @@ namespace ActsFatras { /// Handle particle decays using the Geant4 decay models. class Geant4Decay { public: - using Scalar = Particle::Scalar; - /// Constructor Geant4Decay(); @@ -40,7 +38,7 @@ class Geant4Decay { /// /// @return Proper time limit of the particle template - Scalar generateProperTimeLimit(generator_t& generator, + double generateProperTimeLimit(generator_t& generator, const Particle& particle) const; /// Decay the particle and create the decay products. @@ -67,13 +65,13 @@ class Geant4Decay { }; template -Particle::Scalar Geant4Decay::generateProperTimeLimit( - generator_t& generator, const Particle& particle) const { +double Geant4Decay::generateProperTimeLimit(generator_t& generator, + const Particle& particle) const { // Get the particle properties const Acts::PdgParticle pdgCode = particle.pdg(); // Keep muons stable if (makeAbsolutePdgParticle(pdgCode) == Acts::PdgParticle::eMuon) { - return std::numeric_limits::infinity(); + return std::numeric_limits::infinity(); } // Get the Geant4 particle @@ -81,14 +79,14 @@ Particle::Scalar Geant4Decay::generateProperTimeLimit( // Fast exit if the particle is stable if (!pDef || pDef->GetPDGStable()) { - return std::numeric_limits::infinity(); + return std::numeric_limits::infinity(); } // Get average lifetime - constexpr Scalar convertTime = Acts::UnitConstants::mm / CLHEP::s; - const Scalar tau = pDef->GetPDGLifeTime() * convertTime; + constexpr double convertTime = Acts::UnitConstants::mm / CLHEP::s; + const double tau = pDef->GetPDGLifeTime() * convertTime; // Sample & return the lifetime - std::uniform_real_distribution uniformDistribution{0., 1.}; + std::uniform_real_distribution uniformDistribution{0., 1.}; return -tau * std::log(uniformDistribution(generator)); } diff --git a/Fatras/Geant4/src/Geant4Decay.cpp b/Fatras/Geant4/src/Geant4Decay.cpp index 41b5a3656b1..1afef4ca1d0 100644 --- a/Fatras/Geant4/src/Geant4Decay.cpp +++ b/Fatras/Geant4/src/Geant4Decay.cpp @@ -50,7 +50,7 @@ std::vector ActsFatras::Geant4Decay::decayParticle( } // Boost the decay products using the parents four-momentum - const Particle::Vector4 mom4 = parent.fourMomentum(); + const Acts::Vector4 mom4 = parent.fourMomentum(); products->Boost(mom4[Acts::eMom0] / mom4[Acts::eEnergy], mom4[Acts::eMom1] / mom4[Acts::eEnergy], mom4[Acts::eMom2] / mom4[Acts::eEnergy]); diff --git a/Fatras/include/ActsFatras/Digitization/DigitizationData.hpp b/Fatras/include/ActsFatras/Digitization/DigitizationData.hpp index 314cbcc5eb1..f2e808e00cb 100644 --- a/Fatras/include/ActsFatras/Digitization/DigitizationData.hpp +++ b/Fatras/include/ActsFatras/Digitization/DigitizationData.hpp @@ -17,7 +17,7 @@ namespace ActsFatras { /// A single cell definition: index, cell central value -using Cell = std::pair; +using Cell = std::pair; /// A channel definition: Cell identification, readout word, links /// @@ -35,7 +35,7 @@ struct Channel { /// Channel constructor /// - /// @param cellId_ The Cell idenficiation and position + /// @param cellId_ The Cell identification and position /// @param value_ The Cell value /// @param links_ The (optional) links to e.g. truth indices Channel(std::array cellId_, signal_t value_, @@ -51,7 +51,6 @@ struct Channel { /// @tparam kSize Number of cluster coordinates template struct Cluster { - using Scalar = Acts::ActsScalar; using ParametersVector = Acts::ActsVector; using CovarianceMatrix = Acts::ActsSquareMatrix; diff --git a/Fatras/include/ActsFatras/Digitization/UncorrelatedHitSmearer.hpp b/Fatras/include/ActsFatras/Digitization/UncorrelatedHitSmearer.hpp index 4420c8a6200..caf5efa8c61 100644 --- a/Fatras/include/ActsFatras/Digitization/UncorrelatedHitSmearer.hpp +++ b/Fatras/include/ActsFatras/Digitization/UncorrelatedHitSmearer.hpp @@ -41,7 +41,6 @@ using SingleParameterSmearFunction = /// vector and associated covariance matrix. template struct BoundParametersSmearer { - using Scalar = Acts::ActsScalar; using ParametersVector = Acts::ActsVector; using CovarianceMatrix = Acts::ActsSquareMatrix; using Result = Acts::Result>; @@ -116,7 +115,6 @@ struct BoundParametersSmearer { /// individually is not recommended template struct FreeParametersSmearer { - using Scalar = Acts::ActsScalar; using ParametersVector = Acts::ActsVector; using CovarianceMatrix = Acts::ActsSquareMatrix; using Result = Acts::Result>; diff --git a/Fatras/include/ActsFatras/EventData/Hit.hpp b/Fatras/include/ActsFatras/EventData/Hit.hpp index d5eb25662f8..b3838299ca6 100644 --- a/Fatras/include/ActsFatras/EventData/Hit.hpp +++ b/Fatras/include/ActsFatras/EventData/Hit.hpp @@ -26,10 +26,6 @@ namespace ActsFatras { /// thus stored as two separate four-vectors. class Hit { public: - using Scalar = Acts::ActsScalar; - using Vector3 = Acts::ActsVector<3>; - using Vector4 = Acts::ActsVector<4>; - /// Construct default hit with (mostly) invalid information. Hit() = default; /// Construct from four-position and four-momenta. @@ -45,8 +41,8 @@ class Hit { /// users responsibility to ensure that the position correspond to a /// position on the given surface. Hit(Acts::GeometryIdentifier geometryId, Barcode particleId, - const Vector4& pos4, const Vector4& before4, const Vector4& after4, - std::int32_t index_ = -1) + const Acts::Vector4& pos4, const Acts::Vector4& before4, + const Acts::Vector4& after4, std::int32_t index_ = -1) : m_geometryId(geometryId), m_particleId(particleId), m_index(index_), @@ -68,26 +64,26 @@ class Hit { constexpr std::int32_t index() const { return m_index; } /// Space-time position four-vector. - const Vector4& fourPosition() const { return m_pos4; } + const Acts::Vector4& fourPosition() const { return m_pos4; } /// Three-position, i.e. spatial coordinates without the time. auto position() const { return m_pos4.segment<3>(Acts::ePos0); } /// Time coordinate. - Scalar time() const { return m_pos4[Acts::eTime]; } + double time() const { return m_pos4[Acts::eTime]; } /// Particle four-momentum before the hit. - const Vector4& momentum4Before() const { return m_before4; } + const Acts::Vector4& momentum4Before() const { return m_before4; } /// Particle four-momentum after the hit. - const Vector4& momentum4After() const { return m_after4; } + const Acts::Vector4& momentum4After() const { return m_after4; } /// Normalized particle direction vector before the hit. - Vector3 directionBefore() const { + Acts::Vector3 directionBefore() const { return m_before4.segment<3>(Acts::eMom0).normalized(); } /// Normalized particle direction vector the hit. - Vector3 directionAfter() const { + Acts::Vector3 directionAfter() const { return m_after4.segment<3>(Acts::eMom0).normalized(); } /// Average normalized particle direction vector through the surface. - Vector3 direction() const { + Acts::Vector3 direction() const { auto dir0 = m_before4.segment<3>(Acts::eMom0).normalized(); auto dir1 = m_after4.segment<3>(Acts::eMom0).normalized(); return ((dir0 + dir1) / 2.).segment<3>(Acts::eMom0).normalized(); @@ -96,7 +92,7 @@ class Hit { /// /// @retval positive if the particle lost energy when it passed the surface /// @retval negative if magic was involved - Scalar depositedEnergy() const { + double depositedEnergy() const { return m_before4[Acts::eEnergy] - m_after4[Acts::eEnergy]; } @@ -108,11 +104,11 @@ class Hit { /// Index of the hit along the particle trajectory. std::int32_t m_index = -1; /// Global space-time position four-vector. - Vector4 m_pos4 = Vector4::Zero(); + Acts::Vector4 m_pos4 = Acts::Vector4::Zero(); /// Global particle energy-momentum four-vector before the hit. - Vector4 m_before4 = Vector4::Zero(); + Acts::Vector4 m_before4 = Acts::Vector4::Zero(); /// Global particle energy-momentum four-vector after the hit. - Vector4 m_after4 = Vector4::Zero(); + Acts::Vector4 m_after4 = Acts::Vector4::Zero(); }; } // namespace ActsFatras diff --git a/Fatras/include/ActsFatras/EventData/Particle.hpp b/Fatras/include/ActsFatras/EventData/Particle.hpp index 5b6eb2f7a61..fbf9244fcea 100644 --- a/Fatras/include/ActsFatras/EventData/Particle.hpp +++ b/Fatras/include/ActsFatras/EventData/Particle.hpp @@ -32,10 +32,6 @@ namespace ActsFatras { /// Also stores some simulation-specific properties. class Particle { public: - using Scalar = Acts::ActsScalar; - using Vector3 = Acts::ActsVector<3>; - using Vector4 = Acts::ActsVector<4>; - /// Construct a default particle with invalid identity. Particle() = default; /// Construct a particle at rest with explicit mass and charge. @@ -47,8 +43,8 @@ class Particle { /// /// @warning It is the users responsibility that charge and mass match /// the PDG particle number. - Particle(Barcode particleId, Acts::PdgParticle pdg, Scalar charge, - Scalar mass) + Particle(Barcode particleId, Acts::PdgParticle pdg, double charge, + double mass) : m_particleId(particleId), m_pdg(pdg), m_charge(charge), m_mass(mass) {} /// Construct a particle at rest from a PDG particle number. /// @@ -84,12 +80,12 @@ class Particle { return *this; } /// Set the charge. - Particle setCharge(Scalar charge) { + Particle setCharge(double charge) { m_charge = charge; return *this; } /// Set the mass. - Particle setMass(Scalar mass) { + Particle setMass(double mass) { m_mass = mass; return *this; } @@ -99,18 +95,18 @@ class Particle { return *this; } /// Set the space-time position four-vector. - Particle &setPosition4(const Vector4 &pos4) { + Particle &setPosition4(const Acts::Vector4 &pos4) { m_position4 = pos4; return *this; } /// Set the space-time position four-vector from three-position and time. - Particle &setPosition4(const Vector3 &position, Scalar time) { + Particle &setPosition4(const Acts::Vector3 &position, double time) { m_position4.segment<3>(Acts::ePos0) = position; m_position4[Acts::eTime] = time; return *this; } /// Set the space-time position four-vector from scalar components. - Particle &setPosition4(Scalar x, Scalar y, Scalar z, Scalar time) { + Particle &setPosition4(double x, double y, double z, double time) { m_position4[Acts::ePos0] = x; m_position4[Acts::ePos1] = y; m_position4[Acts::ePos2] = z; @@ -118,13 +114,13 @@ class Particle { return *this; } /// Set the direction three-vector - Particle &setDirection(const Vector3 &direction) { + Particle &setDirection(const Acts::Vector3 &direction) { m_direction = direction; m_direction.normalize(); return *this; } /// Set the direction three-vector from scalar components. - Particle &setDirection(Scalar dx, Scalar dy, Scalar dz) { + Particle &setDirection(double dx, double dy, double dz) { m_direction[Acts::ePos0] = dx; m_direction[Acts::ePos1] = dy; m_direction[Acts::ePos2] = dz; @@ -132,7 +128,7 @@ class Particle { return *this; } /// Set the absolute momentum. - Particle &setAbsoluteMomentum(Scalar absMomentum) { + Particle &setAbsoluteMomentum(double absMomentum) { m_absMomentum = absMomentum; return *this; } @@ -142,10 +138,10 @@ class Particle { /// Energy loss corresponds to a negative change. If the updated energy /// would result in an unphysical value, the particle is put to rest, i.e. /// its absolute momentum is set to zero. - Particle &correctEnergy(Scalar delta) { + Particle &correctEnergy(double delta) { const auto newEnergy = std::hypot(m_mass, m_absMomentum) + delta; if (newEnergy <= m_mass) { - m_absMomentum = Scalar{0}; + m_absMomentum = 0.; } else { m_absMomentum = std::sqrt(newEnergy * newEnergy - m_mass * m_mass); } @@ -163,30 +159,30 @@ class Particle { return Acts::makeAbsolutePdgParticle(pdg()); } /// Particle charge. - Scalar charge() const { return m_charge; } + double charge() const { return m_charge; } /// Particle absolute charge. - Scalar absoluteCharge() const { return std::abs(m_charge); } + double absoluteCharge() const { return std::abs(m_charge); } /// Particle mass. - Scalar mass() const { return m_mass; } + double mass() const { return m_mass; } /// Particle hypothesis. Acts::ParticleHypothesis hypothesis() const { return Acts::ParticleHypothesis(absolutePdg(), mass(), absoluteCharge()); } /// Particl qOverP. - Scalar qOverP() const { + double qOverP() const { return hypothesis().qOverP(absoluteMomentum(), charge()); } /// Space-time position four-vector. - const Vector4 &fourPosition() const { return m_position4; } + const Acts::Vector4 &fourPosition() const { return m_position4; } /// Three-position, i.e. spatial coordinates without the time. auto position() const { return m_position4.segment<3>(Acts::ePos0); } /// Time coordinate. - Scalar time() const { return m_position4[Acts::eTime]; } + double time() const { return m_position4[Acts::eTime]; } /// Energy-momentum four-vector. - Vector4 fourMomentum() const { - Vector4 mom4; + Acts::Vector4 fourMomentum() const { + Acts::Vector4 mom4; // stored direction is always normalized mom4[Acts::eMom0] = m_absMomentum * m_direction[Acts::ePos0]; mom4[Acts::eMom1] = m_absMomentum * m_direction[Acts::ePos1]; @@ -195,24 +191,24 @@ class Particle { return mom4; } /// Unit three-direction, i.e. the normalized momentum three-vector. - const Vector3 &direction() const { return m_direction; } + const Acts::Vector3 &direction() const { return m_direction; } /// Polar angle. - Scalar theta() const { return Acts::VectorHelpers::theta(direction()); } + double theta() const { return Acts::VectorHelpers::theta(direction()); } /// Azimuthal angle. - Scalar phi() const { return Acts::VectorHelpers::phi(direction()); } + double phi() const { return Acts::VectorHelpers::phi(direction()); } /// Absolute momentum in the x-y plane. - Scalar transverseMomentum() const { + double transverseMomentum() const { return m_absMomentum * m_direction.segment<2>(Acts::eMom0).norm(); } /// Absolute momentum. - Scalar absoluteMomentum() const { return m_absMomentum; } + double absoluteMomentum() const { return m_absMomentum; } /// Absolute momentum. - Vector3 momentum() const { return absoluteMomentum() * direction(); } + Acts::Vector3 momentum() const { return absoluteMomentum() * direction(); } /// Total energy, i.e. norm of the four-momentum. - Scalar energy() const { return std::hypot(m_mass, m_absMomentum); } + double energy() const { return std::hypot(m_mass, m_absMomentum); } /// Check if the particle is alive, i.e. is not at rest. - bool isAlive() const { return Scalar{0} < m_absMomentum; } + bool isAlive() const { return 0. < m_absMomentum; } /// Check if this is a secondary particle. bool isSecondary() const { @@ -225,26 +221,26 @@ class Particle { /// Set the proper time in the particle rest frame. /// /// @param properTime passed proper time in the rest frame - Particle &setProperTime(Scalar properTime) { + Particle &setProperTime(double properTime) { m_properTime = properTime; return *this; } /// Proper time in the particle rest frame. - Scalar properTime() const { return m_properTime; } + double properTime() const { return m_properTime; } /// Set the accumulated material measured in radiation/interaction lengths. /// /// @param pathInX0 accumulated material measured in radiation lengths /// @param pathInL0 accumulated material measured in interaction lengths - Particle &setMaterialPassed(Scalar pathInX0, Scalar pathInL0) { + Particle &setMaterialPassed(double pathInX0, double pathInL0) { m_pathInX0 = pathInX0; m_pathInL0 = pathInL0; return *this; } /// Accumulated path within material measured in radiation lengths. - Scalar pathInX0() const { return m_pathInX0; } + double pathInX0() const { return m_pathInX0; } /// Accumulated path within material measured in interaction lengths. - Scalar pathInL0() const { return m_pathInL0; } + double pathInL0() const { return m_pathInL0; } /// Set the reference surface. /// @@ -314,17 +310,17 @@ class Particle { /// PDG particle number. Acts::PdgParticle m_pdg = Acts::PdgParticle::eInvalid; // Particle charge and mass. - Scalar m_charge = Scalar{0}; - Scalar m_mass = Scalar{0}; + double m_charge = 0.; + double m_mass = 0.; // kinematics, i.e. things that change over the particle lifetime. - Vector3 m_direction = Vector3::UnitZ(); - Scalar m_absMomentum = Scalar{0}; - Vector4 m_position4 = Vector4::Zero(); + Acts::Vector3 m_direction = Acts::Vector3::UnitZ(); + double m_absMomentum = 0.; + Acts::Vector4 m_position4 = Acts::Vector4::Zero(); /// proper time in the particle rest frame - Scalar m_properTime = Scalar{0}; + double m_properTime = 0.; // accumulated material - Scalar m_pathInX0 = Scalar{0}; - Scalar m_pathInL0 = Scalar{0}; + double m_pathInX0 = 0.; + double m_pathInL0 = 0.; /// number of hits std::uint32_t m_numberOfHits = 0; /// reference surface diff --git a/Fatras/include/ActsFatras/Kernel/InteractionList.hpp b/Fatras/include/ActsFatras/Kernel/InteractionList.hpp index eb6816c490a..b8faaf7f0ad 100644 --- a/Fatras/include/ActsFatras/Kernel/InteractionList.hpp +++ b/Fatras/include/ActsFatras/Kernel/InteractionList.hpp @@ -83,9 +83,7 @@ template concept PointLikeProcessConcept = requires( const process_t& p, std::uniform_int_distribution& rng, const Particle& prt) { - { - p.generatePathLimits(rng, prt) - } -> std::same_as>; + { p.generatePathLimits(rng, prt) } -> std::same_as>; }; template @@ -177,10 +175,8 @@ class InteractionList { public: /// Point-like interaction selection. struct Selection { - Particle::Scalar x0Limit = - std::numeric_limits::infinity(); - Particle::Scalar l0Limit = - std::numeric_limits::infinity(); + double x0Limit = std::numeric_limits::infinity(); + double l0Limit = std::numeric_limits::infinity(); std::size_t x0Process = std::numeric_limits::max(); std::size_t l0Process = std::numeric_limits::max(); }; diff --git a/Fatras/include/ActsFatras/Kernel/Simulation.hpp b/Fatras/include/ActsFatras/Kernel/Simulation.hpp index b6975c34172..07c03bb5d40 100644 --- a/Fatras/include/ActsFatras/Kernel/Simulation.hpp +++ b/Fatras/include/ActsFatras/Kernel/Simulation.hpp @@ -224,7 +224,7 @@ struct Simulation { // only need to switch between charged/neutral. SingleParticleSimulationResult result = SingleParticleSimulationResult::success({}); - if (initialParticle.charge() != Particle::Scalar{0}) { + if (initialParticle.charge() != 0.) { result = charged.simulate(geoCtx, magCtx, generator, initialParticle); } else { result = neutral.simulate(geoCtx, magCtx, generator, initialParticle); @@ -267,7 +267,7 @@ struct Simulation { private: /// Select if the particle should be simulated at all. bool selectParticle(const Particle &particle) const { - if (particle.charge() != Particle::Scalar{0}) { + if (particle.charge() != 0.) { return selectCharged(particle); } else { return selectNeutral(particle); diff --git a/Fatras/include/ActsFatras/Kernel/SimulationResult.hpp b/Fatras/include/ActsFatras/Kernel/SimulationResult.hpp index 2750e01b616..fa6cf211673 100644 --- a/Fatras/include/ActsFatras/Kernel/SimulationResult.hpp +++ b/Fatras/include/ActsFatras/Kernel/SimulationResult.hpp @@ -38,11 +38,10 @@ struct SimulationResult { // Whether the particle is still alive and the simulation should continue bool isAlive = true; // Proper time limit before decay. - Particle::Scalar properTimeLimit = - std::numeric_limits::quiet_NaN(); + double properTimeLimit = std::numeric_limits::quiet_NaN(); // Accumulated radiation/interaction length limit before next interaction. - Particle::Scalar x0Limit = std::numeric_limits::quiet_NaN(); - Particle::Scalar l0Limit = std::numeric_limits::quiet_NaN(); + double x0Limit = std::numeric_limits::quiet_NaN(); + double l0Limit = std::numeric_limits::quiet_NaN(); // Process selection for the next interaction. std::size_t x0Process = std::numeric_limits::max(); std::size_t l0Process = std::numeric_limits::max(); diff --git a/Fatras/include/ActsFatras/Kernel/detail/SimulationActor.hpp b/Fatras/include/ActsFatras/Kernel/detail/SimulationActor.hpp index acbaaaa713d..bee570485b0 100644 --- a/Fatras/include/ActsFatras/Kernel/detail/SimulationActor.hpp +++ b/Fatras/include/ActsFatras/Kernel/detail/SimulationActor.hpp @@ -53,7 +53,7 @@ struct SimulationActor { Particle initialParticle; /// Relative tolerance of the particles proper time limit - Particle::Scalar properTimeRelativeTolerance = 1e-3; + double properTimeRelativeTolerance = 1e-3; /// Simulate the interaction with a single surface. /// @@ -180,7 +180,7 @@ struct SimulationActor { result.hits.emplace_back( surface.geometryId(), before.particleId(), // the interaction could potentially modify the particle position - Hit::Scalar{0.5} * (before.fourPosition() + after.fourPosition()), + 0.5 * (before.fourPosition() + after.fourPosition()), before.fourMomentum(), after.fourMomentum(), result.hits.size()); after.setNumberOfHits(result.hits.size()); diff --git a/Fatras/include/ActsFatras/Physics/Decay/NoDecay.hpp b/Fatras/include/ActsFatras/Physics/Decay/NoDecay.hpp index 17aefd7341b..11889f3af06 100644 --- a/Fatras/include/ActsFatras/Physics/Decay/NoDecay.hpp +++ b/Fatras/include/ActsFatras/Physics/Decay/NoDecay.hpp @@ -21,9 +21,9 @@ struct NoDecay { /// /// @returns Always returns infinity as limit. template - constexpr Particle::Scalar generateProperTimeLimit( + constexpr double generateProperTimeLimit( generator_t& /* rng */, const Particle& /* particle */) const { - return std::numeric_limits::infinity(); + return std::numeric_limits::infinity(); } /// Decay the particle without generating any descendant particles. template diff --git a/Fatras/include/ActsFatras/Physics/ElectroMagnetic/BetheHeitler.hpp b/Fatras/include/ActsFatras/Physics/ElectroMagnetic/BetheHeitler.hpp index d0529837b7f..6cfbce087c9 100644 --- a/Fatras/include/ActsFatras/Physics/ElectroMagnetic/BetheHeitler.hpp +++ b/Fatras/include/ActsFatras/Physics/ElectroMagnetic/BetheHeitler.hpp @@ -25,9 +25,6 @@ namespace ActsFatras { /// "A Gaussian-mixture approximation of the Bethe–Heitler model of electron /// energy loss by bremsstrahlung" R. Frühwirth struct BetheHeitler { - using Scalar = Particle::Scalar; - using Vector3 = Particle::Vector3; - /// A scaling factor to double scaleFactor = 1.; @@ -42,9 +39,9 @@ struct BetheHeitler { /// @param [in] rndTheta1 Random number for the polar angle /// @param [in] rndTheta2 Random number for the polar angle /// @param [in] rndTheta3 Random number for the polar angle - Particle bremPhoton(const Particle &particle, Scalar gammaE, Scalar rndPsi, - Scalar rndTheta1, Scalar rndTheta2, - Scalar rndTheta3) const; + Particle bremPhoton(const Particle &particle, double gammaE, double rndPsi, + double rndTheta1, double rndTheta2, + double rndTheta3) const; /// Simulate energy loss and update the particle parameters. /// @@ -67,7 +64,7 @@ struct BetheHeitler { const auto sampledEnergyLoss = std::abs(scaleFactor * particle.energy() * (z - 1.)); - std::uniform_real_distribution uDist(0., 1.); + std::uniform_real_distribution uDist(0., 1.); // Build the produced photon Particle photon = bremPhoton(particle, sampledEnergyLoss, uDist(generator), diff --git a/Fatras/include/ActsFatras/Physics/NuclearInteraction/NuclearInteraction.hpp b/Fatras/include/ActsFatras/Physics/NuclearInteraction/NuclearInteraction.hpp index 89072e78536..d372cf8ffc9 100644 --- a/Fatras/include/ActsFatras/Physics/NuclearInteraction/NuclearInteraction.hpp +++ b/Fatras/include/ActsFatras/Physics/NuclearInteraction/NuclearInteraction.hpp @@ -38,7 +38,6 @@ namespace ActsFatras { /// interaction. Either the initial particle survives (soft) or it gets /// destroyed (hard) by this process. struct NuclearInteraction { - using Scalar = Particle::Scalar; /// The storage of the parameterisation detail::MultiParticleNuclearInteractionParametrisation multiParticleParameterisation; @@ -55,12 +54,12 @@ struct NuclearInteraction { /// /// @return valid X0 limit and no limit on L0 template - std::pair generatePathLimits(generator_t& generator, + std::pair generatePathLimits(generator_t& generator, const Particle& particle) const { // Fast exit: No parameterisation provided if (multiParticleParameterisation.empty()) { - return std::make_pair(std::numeric_limits::infinity(), - std::numeric_limits::infinity()); + return std::make_pair(std::numeric_limits::infinity(), + std::numeric_limits::infinity()); } // Find the parametrisation that corresponds to the particle type for (const auto& particleParametrisation : multiParticleParameterisation) { @@ -79,14 +78,14 @@ struct NuclearInteraction { const auto& distribution = parametrisation.nuclearInteractionProbability; auto limits = - std::make_pair(std::numeric_limits::infinity(), + std::make_pair(std::numeric_limits::infinity(), sampleContinuousValues( uniformDistribution(generator), distribution)); return limits; } } - return std::make_pair(std::numeric_limits::infinity(), - std::numeric_limits::infinity()); + return std::make_pair(std::numeric_limits::infinity(), + std::numeric_limits::infinity()); } /// This method performs a nuclear interaction. @@ -292,10 +291,8 @@ struct NuclearInteraction { /// /// @return Azimuthal and polar angle of the second particle in the global /// coordinate system - std::pair - globalAngle(ActsFatras::Particle::Scalar phi1, - ActsFatras::Particle::Scalar theta1, float phi2, - float theta2) const; + std::pair globalAngle(double phi1, double theta1, float phi2, + float theta2) const; /// Converter from sampled numbers to a vector of particles /// @@ -337,7 +334,7 @@ struct NuclearInteraction { /// neighbouring bins should be performed instead of a bin lookup /// /// @return The sampled value - Scalar sampleContinuousValues( + double sampleContinuousValues( double rnd, const detail::NuclearInteractionParameters::CumulativeDistribution& distribution, @@ -416,8 +413,8 @@ Acts::ActsDynamicVector NuclearInteraction::sampleInvariantMasses( // Sample in the eigenspace for (unsigned int i = 0; i < size; i++) { float variance = parametrisation.eigenvaluesInvariantMass[i]; - std::normal_distribution dist{ - parametrisation.meanInvariantMass[i], std::sqrt(variance)}; + std::normal_distribution dist{parametrisation.meanInvariantMass[i], + std::sqrt(variance)}; parameters[i] = dist(generator); } // Transform to multivariate normal distribution @@ -446,8 +443,8 @@ Acts::ActsDynamicVector NuclearInteraction::sampleMomenta( // Sample in the eigenspace for (unsigned int i = 0; i < size; i++) { float variance = parametrisation.eigenvaluesMomentum[i]; - std::normal_distribution dist{ - parametrisation.meanMomentum[i], std::sqrt(variance)}; + std::normal_distribution dist{parametrisation.meanMomentum[i], + std::sqrt(variance)}; parameters[i] = dist(generator); } diff --git a/Fatras/include/ActsFatras/Selectors/ParticleSelectors.hpp b/Fatras/include/ActsFatras/Selectors/ParticleSelectors.hpp index 4deb110969e..d784048abbd 100644 --- a/Fatras/include/ActsFatras/Selectors/ParticleSelectors.hpp +++ b/Fatras/include/ActsFatras/Selectors/ParticleSelectors.hpp @@ -21,28 +21,28 @@ struct EveryParticle { /// Select neutral particles. struct NeutralSelector { bool operator()(const Particle &particle) const { - return (particle.charge() == Particle::Scalar{0}); + return (particle.charge() == 0.); } }; /// Select all charged particles. struct ChargedSelector { bool operator()(const Particle &particle) const { - return (particle.charge() != Particle::Scalar{0}); + return (particle.charge() != 0.); } }; /// Select positively charged particles. struct PositiveSelector { bool operator()(const Particle &particle) const { - return (Particle::Scalar{0} < particle.charge()); + return (0. < particle.charge()); } }; /// Select negatively charged particles. struct NegativeSelector { bool operator()(const Particle &particle) const { - return (particle.charge() < Particle::Scalar{0}); + return (particle.charge() < 0.); } }; diff --git a/Fatras/src/Physics/BetheHeitler.cpp b/Fatras/src/Physics/BetheHeitler.cpp index af8cbd2711f..6191d832c38 100644 --- a/Fatras/src/Physics/BetheHeitler.cpp +++ b/Fatras/src/Physics/BetheHeitler.cpp @@ -20,8 +20,8 @@ #include ActsFatras::Particle ActsFatras::BetheHeitler::bremPhoton( - const Particle &particle, Scalar gammaE, Scalar rndPsi, Scalar rndTheta1, - Scalar rndTheta2, Scalar rndTheta3) const { + const Particle &particle, double gammaE, double rndPsi, double rndTheta1, + double rndTheta2, double rndTheta3) const { // ------------------------------------------------------ // simple approach // (a) simulate theta uniform within the opening angle of the relativistic @@ -31,10 +31,10 @@ ActsFatras::Particle ActsFatras::BetheHeitler::bremPhoton( // later // the azimutal angle - Scalar psi = 2. * std::numbers::pi * rndPsi; + double psi = 2. * std::numbers::pi * rndPsi; // the start of the equation - Scalar theta = 0.; + double theta = 0.; if (uniformHertzDipoleAngle) { // the simplest simulation theta = particle.mass() / particle.energy() * rndTheta1; @@ -42,8 +42,8 @@ ActsFatras::Particle ActsFatras::BetheHeitler::bremPhoton( // -----> theta = particle.mass() / particle.energy(); // follow - constexpr Scalar a = 0.625; // 5/8 - Scalar u = -log(rndTheta2 * rndTheta3) / a; + constexpr double a = 0.625; // 5/8 + double u = -log(rndTheta2 * rndTheta3) / a; theta *= (rndTheta1 < 0.25) ? u : u / 3.; // 9./(9.+27) = 0.25 } diff --git a/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp b/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp index 1515fc3c626..497d530663b 100644 --- a/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp +++ b/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp @@ -67,7 +67,7 @@ unsigned int NuclearInteraction::sampleDiscreteValues( return static_cast(distribution.first[iBin]); } -Particle::Scalar NuclearInteraction::sampleContinuousValues( +double NuclearInteraction::sampleContinuousValues( double rnd, const detail::NuclearInteractionParameters::CumulativeDistribution& distribution, @@ -111,10 +111,10 @@ unsigned int NuclearInteraction::finalStateMultiplicity( return sampleDiscreteValues(rnd, distribution); } -std::pair -NuclearInteraction::globalAngle(ActsFatras::Particle::Scalar phi1, - ActsFatras::Particle::Scalar theta1, float phi2, - float theta2) const { +std::pair NuclearInteraction::globalAngle(double phi1, + double theta1, + float phi2, + float theta2) const { // Rotation around the global y-axis Acts::SquareMatrix3 rotY = Acts::SquareMatrix3::Zero(); rotY(0, 0) = std::cos(theta1); diff --git a/Fatras/src/Physics/PhotonConversion.cpp b/Fatras/src/Physics/PhotonConversion.cpp index 73bb83ad413..bf71af4cae2 100644 --- a/Fatras/src/Physics/PhotonConversion.cpp +++ b/Fatras/src/Physics/PhotonConversion.cpp @@ -10,6 +10,5 @@ #include "Acts/Definitions/ParticleData.hpp" -const ActsFatras::PhotonConversion::Scalar - ActsFatras::PhotonConversion::kElectronMass = - Acts::findMass(Acts::PdgParticle::eElectron).value(); +const double ActsFatras::PhotonConversion::kElectronMass = + Acts::findMass(Acts::PdgParticle::eElectron).value(); diff --git a/Tests/UnitTests/Fatras/EventData/HitTests.cpp b/Tests/UnitTests/Fatras/EventData/HitTests.cpp index 765db8378a8..4fc89f8b686 100644 --- a/Tests/UnitTests/Fatras/EventData/HitTests.cpp +++ b/Tests/UnitTests/Fatras/EventData/HitTests.cpp @@ -18,7 +18,7 @@ using namespace ActsFatras; namespace { -constexpr auto eps = std::numeric_limits::epsilon(); +constexpr auto eps = std::numeric_limits::epsilon(); const auto pid = Barcode().setVertexPrimary(12).setParticle(23); const auto gid = Acts::GeometryIdentifier().setVolume(1).setLayer(2).setSensitive(3); @@ -28,93 +28,95 @@ BOOST_AUTO_TEST_SUITE(FatrasHit) BOOST_AUTO_TEST_CASE(WithoutInteraction) { // some hit position - auto p4 = Hit::Vector4(1, 2, 3, 4); + auto p4 = Acts::Vector4(1, 2, 3, 4); // before/after four-momenta are the same - auto m4 = Hit::Vector4(1, 1, 1, 4); + auto m4 = Acts::Vector4(1, 1, 1, 4); auto h = Hit(gid, pid, p4, m4, m4, 12u); BOOST_CHECK_EQUAL(h.geometryId(), gid); BOOST_CHECK_EQUAL(h.particleId(), pid); BOOST_CHECK_EQUAL(h.index(), 12u); CHECK_CLOSE_REL(h.fourPosition(), p4, eps); - CHECK_CLOSE_REL(h.position(), Hit::Vector3(1, 2, 3), eps); + CHECK_CLOSE_REL(h.position(), Acts::Vector3(1, 2, 3), eps); CHECK_CLOSE_REL(h.time(), 4, eps); CHECK_CLOSE_REL(h.momentum4Before(), m4, eps); CHECK_CLOSE_REL(h.momentum4After(), m4, eps); - CHECK_CLOSE_REL(h.directionBefore(), Hit::Vector3(1, 1, 1).normalized(), eps); - CHECK_CLOSE_REL(h.directionAfter(), Hit::Vector3(1, 1, 1).normalized(), eps); - CHECK_CLOSE_REL(h.direction(), Hit::Vector3(1, 1, 1).normalized(), eps); + CHECK_CLOSE_REL(h.directionBefore(), Acts::Vector3(1, 1, 1).normalized(), + eps); + CHECK_CLOSE_REL(h.directionAfter(), Acts::Vector3(1, 1, 1).normalized(), eps); + CHECK_CLOSE_REL(h.direction(), Acts::Vector3(1, 1, 1).normalized(), eps); CHECK_SMALL(h.depositedEnergy(), eps); } BOOST_AUTO_TEST_CASE(WithEnergyLoss) { // some hit position - auto p4 = Hit::Vector4(1, 2, 3, 4); + auto p4 = Acts::Vector4(1, 2, 3, 4); // before/after four-momenta differ by energy loss, use zero mass to simplify - auto m40 = Hit::Vector4(2, 0, 0, 2); - auto m41 = Hit::Vector4(1.5, 0, 0, 1.5); + auto m40 = Acts::Vector4(2, 0, 0, 2); + auto m41 = Acts::Vector4(1.5, 0, 0, 1.5); auto h = Hit(gid, pid, p4, m40, m41, 13u); BOOST_CHECK_EQUAL(h.geometryId(), gid); BOOST_CHECK_EQUAL(h.particleId(), pid); BOOST_CHECK_EQUAL(h.index(), 13u); CHECK_CLOSE_REL(h.fourPosition(), p4, eps); - CHECK_CLOSE_REL(h.position(), Hit::Vector3(1, 2, 3), eps); + CHECK_CLOSE_REL(h.position(), Acts::Vector3(1, 2, 3), eps); CHECK_CLOSE_REL(h.time(), 4, eps); CHECK_CLOSE_OR_SMALL(h.momentum4Before(), m40, eps, eps); CHECK_CLOSE_OR_SMALL(h.momentum4After(), m41, eps, eps); - CHECK_CLOSE_OR_SMALL(h.directionBefore(), Hit::Vector3(1, 0, 0), eps, eps); - CHECK_CLOSE_OR_SMALL(h.directionAfter(), Hit::Vector3(1, 0, 0), eps, eps); - CHECK_CLOSE_OR_SMALL(h.direction(), Hit::Vector3(1, 0, 0), eps, eps); + CHECK_CLOSE_OR_SMALL(h.directionBefore(), Acts::Vector3(1, 0, 0), eps, eps); + CHECK_CLOSE_OR_SMALL(h.directionAfter(), Acts::Vector3(1, 0, 0), eps, eps); + CHECK_CLOSE_OR_SMALL(h.direction(), Acts::Vector3(1, 0, 0), eps, eps); CHECK_CLOSE_REL(h.depositedEnergy(), 0.5, eps); } BOOST_AUTO_TEST_CASE(WithScattering) { // some hit position - auto p4 = Hit::Vector4(1, 2, 3, 4); + auto p4 = Acts::Vector4(1, 2, 3, 4); // before/after four-momenta differ only by direction - auto m40 = Hit::Vector4(2, 0, 2, 5); - auto m41 = Hit::Vector4(0, -2, 2, 5); + auto m40 = Acts::Vector4(2, 0, 2, 5); + auto m41 = Acts::Vector4(0, -2, 2, 5); auto h = Hit(gid, pid, p4, m40, m41, 42u); BOOST_CHECK_EQUAL(h.geometryId(), gid); BOOST_CHECK_EQUAL(h.particleId(), pid); BOOST_CHECK_EQUAL(h.index(), 42u); CHECK_CLOSE_REL(h.fourPosition(), p4, eps); - CHECK_CLOSE_REL(h.position(), Hit::Vector3(1, 2, 3), eps); + CHECK_CLOSE_REL(h.position(), Acts::Vector3(1, 2, 3), eps); CHECK_CLOSE_REL(h.time(), 4, eps); CHECK_CLOSE_OR_SMALL(h.momentum4Before(), m40, eps, eps); CHECK_CLOSE_OR_SMALL(h.momentum4After(), m41, eps, eps); - CHECK_CLOSE_OR_SMALL(h.directionBefore(), Hit::Vector3(1, 0, 1).normalized(), + CHECK_CLOSE_OR_SMALL(h.directionBefore(), Acts::Vector3(1, 0, 1).normalized(), eps, eps); - CHECK_CLOSE_OR_SMALL(h.directionAfter(), Hit::Vector3(0, -1, 1).normalized(), + CHECK_CLOSE_OR_SMALL(h.directionAfter(), Acts::Vector3(0, -1, 1).normalized(), eps, eps); - CHECK_CLOSE_REL(h.direction(), Hit::Vector3(1, -1, 2).normalized(), eps); + CHECK_CLOSE_REL(h.direction(), Acts::Vector3(1, -1, 2).normalized(), eps); CHECK_SMALL(h.depositedEnergy(), eps); } BOOST_AUTO_TEST_CASE(WithEverything) { // some hit position - auto p4 = Hit::Vector4(1, 2, 3, 4); + auto p4 = Acts::Vector4(1, 2, 3, 4); // before/after four-momenta differ by direction and norm - auto m40 = Hit::Vector4(3, 2, 2, 5); - auto m41 = Hit::Vector4(2, 1, 2, 4); + auto m40 = Acts::Vector4(3, 2, 2, 5); + auto m41 = Acts::Vector4(2, 1, 2, 4); auto h = Hit(gid, pid, p4, m40, m41, 1u); BOOST_CHECK_EQUAL(h.geometryId(), gid); BOOST_CHECK_EQUAL(h.particleId(), pid); BOOST_CHECK_EQUAL(h.index(), 1u); CHECK_CLOSE_REL(h.fourPosition(), p4, eps); - CHECK_CLOSE_REL(h.position(), Hit::Vector3(1, 2, 3), eps); + CHECK_CLOSE_REL(h.position(), Acts::Vector3(1, 2, 3), eps); CHECK_CLOSE_REL(h.time(), 4, eps); CHECK_CLOSE_OR_SMALL(h.momentum4Before(), m40, eps, eps); CHECK_CLOSE_OR_SMALL(h.momentum4After(), m41, eps, eps); - CHECK_CLOSE_REL(h.directionBefore(), Hit::Vector3(3, 2, 2).normalized(), eps); - CHECK_CLOSE_REL(h.directionAfter(), Hit::Vector3(2, 1, 2).normalized(), eps); - CHECK_CLOSE_REL( - h.direction(), - Hit::Vector3(0.7023994590205035, 0.41229136135810396, 0.5802161953247991), - eps); + CHECK_CLOSE_REL(h.directionBefore(), Acts::Vector3(3, 2, 2).normalized(), + eps); + CHECK_CLOSE_REL(h.directionAfter(), Acts::Vector3(2, 1, 2).normalized(), eps); + CHECK_CLOSE_REL(h.direction(), + Acts::Vector3(0.7023994590205035, 0.41229136135810396, + 0.5802161953247991), + eps); CHECK_CLOSE_REL(h.depositedEnergy(), 1, eps); } diff --git a/Tests/UnitTests/Fatras/EventData/ParticleTests.cpp b/Tests/UnitTests/Fatras/EventData/ParticleTests.cpp index a6858dacc74..e3fce039dab 100644 --- a/Tests/UnitTests/Fatras/EventData/ParticleTests.cpp +++ b/Tests/UnitTests/Fatras/EventData/ParticleTests.cpp @@ -23,7 +23,7 @@ using ActsFatras::Particle; using namespace Acts::UnitLiterals; namespace { -constexpr auto eps = std::numeric_limits::epsilon(); +constexpr auto eps = std::numeric_limits::epsilon(); } BOOST_AUTO_TEST_SUITE(FatrasParticle) @@ -35,17 +35,17 @@ BOOST_AUTO_TEST_CASE(Construct) { BOOST_CHECK_EQUAL(particle.particleId(), pid); BOOST_CHECK_EQUAL(particle.pdg(), PdgParticle::eProton); // particle is at rest at the origin - BOOST_CHECK_EQUAL(particle.fourPosition(), Particle::Vector4::Zero()); - BOOST_CHECK_EQUAL(particle.position(), Particle::Vector3::Zero()); - BOOST_CHECK_EQUAL(particle.time(), Particle::Scalar{0}); + BOOST_CHECK_EQUAL(particle.fourPosition(), Acts::Vector4::Zero()); + BOOST_CHECK_EQUAL(particle.position(), Acts::Vector3::Zero()); + BOOST_CHECK_EQUAL(particle.time(), 0.); BOOST_CHECK_EQUAL(particle.fourPosition().x(), particle.position().x()); BOOST_CHECK_EQUAL(particle.fourPosition().y(), particle.position().y()); BOOST_CHECK_EQUAL(particle.fourPosition().z(), particle.position().z()); BOOST_CHECK_EQUAL(particle.fourPosition().w(), particle.time()); // particle direction is undefined, but must be normalized CHECK_CLOSE_REL(particle.direction().norm(), 1, eps); - BOOST_CHECK_EQUAL(particle.transverseMomentum(), Particle::Scalar{0}); - BOOST_CHECK_EQUAL(particle.absoluteMomentum(), Particle::Scalar{0}); + BOOST_CHECK_EQUAL(particle.transverseMomentum(), 0.); + BOOST_CHECK_EQUAL(particle.absoluteMomentum(), 0.); // particle is created at rest and thus not alive BOOST_CHECK(!particle.isAlive()); } @@ -53,7 +53,7 @@ BOOST_AUTO_TEST_CASE(Construct) { BOOST_AUTO_TEST_CASE(CorrectEnergy) { const auto pid = Barcode().setVertexPrimary(1).setParticle(42); auto particle = Particle(pid, PdgParticle::eProton, 1_e, 1_GeV) - .setDirection(Particle::Vector3::UnitX()) + .setDirection(Acts::Vector3::UnitX()) .setAbsoluteMomentum(2_GeV); BOOST_CHECK_EQUAL(particle.mass(), 1_GeV); @@ -68,39 +68,37 @@ BOOST_AUTO_TEST_CASE(CorrectEnergy) { // particle direction must be normalized CHECK_CLOSE_REL(particle.direction().norm(), 1, eps); - // loose some energy + // lose some energy particle.correctEnergy(-100_MeV); BOOST_CHECK_LT(particle.transverseMomentum(), 2_GeV); BOOST_CHECK_LT(particle.absoluteMomentum(), 2_GeV); - BOOST_CHECK_EQUAL(particle.energy(), - Particle::Scalar{std::hypot(1_GeV, 2_GeV) - 100_MeV}); + BOOST_CHECK_EQUAL(particle.energy(), std::hypot(1_GeV, 2_GeV) - 100_MeV); CHECK_CLOSE_REL(particle.direction().norm(), 1, eps); // particle is still alive BOOST_CHECK(particle.isAlive()); - // loose some more energy + // lose some more energy particle.correctEnergy(-200_MeV); BOOST_CHECK_LT(particle.transverseMomentum(), 2_GeV); BOOST_CHECK_LT(particle.absoluteMomentum(), 2_GeV); - BOOST_CHECK_EQUAL(particle.energy(), - Particle::Scalar{std::hypot(1_GeV, 2_GeV) - 300_MeV}); + BOOST_CHECK_EQUAL(particle.energy(), std::hypot(1_GeV, 2_GeV) - 300_MeV); CHECK_CLOSE_REL(particle.direction().norm(), 1, eps); // particle is still alive BOOST_CHECK(particle.isAlive()); - // loose a lot of energy + // lose a lot of energy particle.correctEnergy(-3_GeV); - BOOST_CHECK_EQUAL(particle.transverseMomentum(), Particle::Scalar{0}); - BOOST_CHECK_EQUAL(particle.absoluteMomentum(), Particle::Scalar{0}); + BOOST_CHECK_EQUAL(particle.transverseMomentum(), 0.); + BOOST_CHECK_EQUAL(particle.absoluteMomentum(), 0.); BOOST_CHECK_EQUAL(particle.energy(), particle.mass()); CHECK_CLOSE_REL(particle.direction().norm(), 1, eps); // particle is not alive anymore BOOST_CHECK(!particle.isAlive()); - // lossing even more energy does nothing + // losing even more energy does nothing particle.correctEnergy(-10_GeV); - BOOST_CHECK_EQUAL(particle.transverseMomentum(), Particle::Scalar{0}); - BOOST_CHECK_EQUAL(particle.absoluteMomentum(), Particle::Scalar{0}); + BOOST_CHECK_EQUAL(particle.transverseMomentum(), 0.); + BOOST_CHECK_EQUAL(particle.absoluteMomentum(), 0.); BOOST_CHECK_EQUAL(particle.energy(), particle.mass()); CHECK_CLOSE_REL(particle.direction().norm(), 1, eps); // particle is still not alive diff --git a/Tests/UnitTests/Fatras/Kernel/InteractionListTests.cpp b/Tests/UnitTests/Fatras/Kernel/InteractionListTests.cpp index 167025ad739..0364fe810d6 100644 --- a/Tests/UnitTests/Fatras/Kernel/InteractionListTests.cpp +++ b/Tests/UnitTests/Fatras/Kernel/InteractionListTests.cpp @@ -25,7 +25,6 @@ using namespace Acts::UnitLiterals; using namespace ActsFatras; using Acts::MaterialSlab; -using Scalar = Particle::Scalar; namespace { @@ -63,9 +62,9 @@ static_assert(!detail::PointLikeProcessConcept, /// Each run call creates one descendant particle. struct X0PointLikeProcess { template - std::pair generatePathLimits( + std::pair generatePathLimits( generator_t & /*generator*/, const Particle & /*particle*/) const { - return {0.5, std::numeric_limits::infinity()}; + return {0.5, std::numeric_limits::infinity()}; } template @@ -87,9 +86,9 @@ static_assert(detail::PointLikeProcessConcept, /// Each run call creates two descendant particles. struct L0PointLikeProcess { template - std::pair generatePathLimits( + std::pair generatePathLimits( generator_t & /*generator*/, const Particle & /*particle*/) const { - return {std::numeric_limits::infinity(), 1.5}; + return {std::numeric_limits::infinity(), 1.5}; } template @@ -129,8 +128,8 @@ BOOST_AUTO_TEST_CASE(Empty) { // w/o processes there should be no selection auto sel = l.armPointLike(f.rng, f.incoming); - BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); - BOOST_CHECK_EQUAL(sel.l0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.l0Limit, std::numeric_limits::infinity()); BOOST_CHECK_EQUAL(sel.x0Process, std::numeric_limits::max()); BOOST_CHECK_EQUAL(sel.l0Process, std::numeric_limits::max()); @@ -179,7 +178,7 @@ BOOST_AUTO_TEST_CASE(PointLikeX0) { // w/o processes the list should never abort auto sel = l.armPointLike(f.rng, f.incoming); BOOST_CHECK_EQUAL(sel.x0Limit, 0.5); - BOOST_CHECK_EQUAL(sel.l0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.l0Limit, std::numeric_limits::infinity()); BOOST_CHECK_EQUAL(sel.x0Process, 0u); BOOST_CHECK_EQUAL(sel.l0Process, std::numeric_limits::max()); @@ -198,7 +197,7 @@ BOOST_AUTO_TEST_CASE(PointLikeL0) { // w/o processes the list should never abort auto sel = l.armPointLike(f.rng, f.incoming); - BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); BOOST_CHECK_EQUAL(sel.l0Limit, 1.5); BOOST_CHECK_EQUAL(sel.x0Process, std::numeric_limits::max()); BOOST_CHECK_EQUAL(sel.l0Process, 0u); @@ -253,7 +252,7 @@ BOOST_AUTO_TEST_CASE(Disable) { l.disable(); { auto sel = l.armPointLike(f.rng, f.incoming); - BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); BOOST_CHECK_EQUAL(sel.l0Limit, 1.5); BOOST_CHECK_EQUAL(sel.x0Process, std::numeric_limits::max()); BOOST_CHECK_EQUAL(sel.l0Process, 3u); @@ -271,8 +270,8 @@ BOOST_AUTO_TEST_CASE(Disable) { l.disable(); { auto sel = l.armPointLike(f.rng, f.incoming); - BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); - BOOST_CHECK_EQUAL(sel.l0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.x0Limit, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(sel.l0Limit, std::numeric_limits::infinity()); BOOST_CHECK_EQUAL(sel.x0Process, std::numeric_limits::max()); BOOST_CHECK_EQUAL(sel.l0Process, std::numeric_limits::max()); diff --git a/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp b/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp index 3306c5702d6..e5486f6a288 100644 --- a/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp +++ b/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp @@ -44,15 +44,15 @@ using namespace ActsFatras; namespace { -constexpr auto tol = 4 * std::numeric_limits::epsilon(); -constexpr auto inf = std::numeric_limits::infinity(); +constexpr auto tol = 4 * std::numeric_limits::epsilon(); +constexpr auto inf = std::numeric_limits::infinity(); struct MockDecay { - Particle::Scalar properTimeLimit = inf; + double properTimeLimit = inf; template - constexpr Particle::Scalar generateProperTimeLimit( - generator_t & /*generator*/, const Particle &particle) const { + constexpr double generateProperTimeLimit(generator_t & /*generator*/, + const Particle &particle) const { return particle.properTime() + properTimeLimit; } template @@ -181,10 +181,10 @@ struct Fixture { Barcode pid = Barcode().setVertexPrimary(12u).setParticle(3u); ProcessType proc = ProcessType::eUndefined; Acts::PdgParticle pdg = Acts::PdgParticle::eProton; - Particle::Scalar q = 1_e; - Particle::Scalar m = 1_GeV; - Particle::Scalar p = 1_GeV; - Particle::Scalar e; + double q = 1_e; + double m = 1_GeV; + double p = 1_GeV; + double e; Generator generator; std::shared_ptr surface; Actor actor; diff --git a/Tests/UnitTests/Fatras/Physics/PhotonConversionTests.cpp b/Tests/UnitTests/Fatras/Physics/PhotonConversionTests.cpp index f5340ab4cc7..3b7d6fe1b1c 100644 --- a/Tests/UnitTests/Fatras/Physics/PhotonConversionTests.cpp +++ b/Tests/UnitTests/Fatras/Physics/PhotonConversionTests.cpp @@ -32,7 +32,6 @@ BOOST_AUTO_TEST_SUITE(FatrasPhotonConversion) BOOST_DATA_TEST_CASE(NoPhoton, Dataset::parametersPhotonConversion, phi, theta, seed) { - using Scalar = ActsFatras::PhotonConversion::Scalar; using namespace Acts::UnitLiterals; Generator gen(seed); @@ -45,10 +44,10 @@ BOOST_DATA_TEST_CASE(NoPhoton, Dataset::parametersPhotonConversion, phi, theta, ActsFatras::PhotonConversion pc; // No limits should be set - std::pair limits; + std::pair limits; limits = pc.generatePathLimits(gen, particle); - BOOST_CHECK_EQUAL(limits.first, std::numeric_limits::infinity()); - BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(limits.first, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); // No particles should be generated std::vector generated; @@ -65,7 +64,6 @@ BOOST_DATA_TEST_CASE(NoPhoton, Dataset::parametersPhotonConversion, phi, theta, BOOST_DATA_TEST_CASE(DeadPhoton, Dataset::parametersPhotonConversion, phi, theta, seed) { - using Scalar = ActsFatras::PhotonConversion::Scalar; using namespace Acts::UnitLiterals; Generator gen(seed); @@ -78,9 +76,9 @@ BOOST_DATA_TEST_CASE(DeadPhoton, Dataset::parametersPhotonConversion, phi, ActsFatras::PhotonConversion pc; // No limits should be set - momentum too low - std::pair limits = pc.generatePathLimits(gen, particle); - BOOST_CHECK_EQUAL(limits.first, std::numeric_limits::infinity()); - BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); + std::pair limits = pc.generatePathLimits(gen, particle); + BOOST_CHECK_EQUAL(limits.first, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); // No particles should be generated - momentum too low std::vector generated; @@ -97,7 +95,6 @@ BOOST_DATA_TEST_CASE(DeadPhoton, Dataset::parametersPhotonConversion, phi, BOOST_DATA_TEST_CASE(LowMomentumPhoton, Dataset::parametersPhotonConversion, phi, theta, seed) { - using Scalar = ActsFatras::PhotonConversion::Scalar; using namespace Acts::UnitLiterals; Generator gen(seed); @@ -110,9 +107,9 @@ BOOST_DATA_TEST_CASE(LowMomentumPhoton, Dataset::parametersPhotonConversion, ActsFatras::PhotonConversion pc; // No limits should be set - momentum too low - std::pair limits = pc.generatePathLimits(gen, particle); - BOOST_CHECK_EQUAL(limits.first, std::numeric_limits::infinity()); - BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); + std::pair limits = pc.generatePathLimits(gen, particle); + BOOST_CHECK_EQUAL(limits.first, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); // No particles should be generated - momentum too low std::vector generated; @@ -129,7 +126,6 @@ BOOST_DATA_TEST_CASE(LowMomentumPhoton, Dataset::parametersPhotonConversion, BOOST_DATA_TEST_CASE(HighMomentumPhoton, Dataset::parametersPhotonConversion, phi, theta, seed) { - using Scalar = ActsFatras::PhotonConversion::Scalar; using namespace Acts::UnitLiterals; Generator gen(seed); @@ -142,9 +138,9 @@ BOOST_DATA_TEST_CASE(HighMomentumPhoton, Dataset::parametersPhotonConversion, ActsFatras::PhotonConversion pc; // No limits should be set - momentum too low - std::pair limits = pc.generatePathLimits(gen, particle); - BOOST_CHECK_NE(limits.first, std::numeric_limits::infinity()); - BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); + std::pair limits = pc.generatePathLimits(gen, particle); + BOOST_CHECK_NE(limits.first, std::numeric_limits::infinity()); + BOOST_CHECK_EQUAL(limits.second, std::numeric_limits::infinity()); // No particles should be generated - momentum too low std::vector generated; From c0f427277df96003c72ce04c51c6ba98cd626ec5 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Fri, 22 Nov 2024 17:09:44 +0100 Subject: [PATCH 5/8] step 5 --- .../ElectroMagnetic/PhotonConversion.hpp | 136 +++++++++--------- Fatras/src/Physics/BetheHeitler.cpp | 4 +- .../NuclearInteraction/NuclearInteraction.cpp | 4 +- 3 files changed, 70 insertions(+), 74 deletions(-) diff --git a/Fatras/include/ActsFatras/Physics/ElectroMagnetic/PhotonConversion.hpp b/Fatras/include/ActsFatras/Physics/ElectroMagnetic/PhotonConversion.hpp index d54b17e7407..af45481258e 100644 --- a/Fatras/include/ActsFatras/Physics/ElectroMagnetic/PhotonConversion.hpp +++ b/Fatras/include/ActsFatras/Physics/ElectroMagnetic/PhotonConversion.hpp @@ -34,12 +34,10 @@ namespace ActsFatras { /// interaction itself. class PhotonConversion { public: - using Scalar = ActsFatras::Particle::Scalar; - /// Scaling factor of children energy - Scalar childEnergyScaleFactor = 2.; + double childEnergyScaleFactor = 2.; /// Scaling factor for photon conversion probability - Scalar conversionProbScaleFactor = 0.98; + double conversionProbScaleFactor = 0.98; /// Method for evaluating the distance after which the photon /// conversion will occur. @@ -50,7 +48,7 @@ class PhotonConversion { /// /// @return valid X0 limit and no limit on L0 template - std::pair generatePathLimits(generator_t& generator, + std::pair generatePathLimits(generator_t& generator, const Particle& particle) const; /// This method evaluates the final state due to the photon conversion. @@ -74,8 +72,8 @@ class PhotonConversion { /// /// @return Array containing the produced leptons std::array generateChildren( - const Particle& photon, Scalar childEnergy, - const Particle::Vector3& childDirection) const; + const Particle& photon, double childEnergy, + const Acts::Vector3& childDirection) const; /// Generate the energy fraction of the first child particle. /// @@ -85,8 +83,8 @@ class PhotonConversion { /// /// @return The energy of the child particle template - Scalar generateFirstChildEnergyFraction(generator_t& generator, - Scalar gammaMom) const; + double generateFirstChildEnergyFraction(generator_t& generator, + double gammaMom) const; /// Generate the direction of the child particles. /// @@ -96,27 +94,27 @@ class PhotonConversion { /// /// @return The direction vector of the child particle template - Particle::Vector3 generateChildDirection(generator_t& generator, - const Particle& particle) const; + Acts::Vector3 generateChildDirection(generator_t& generator, + const Particle& particle) const; /// Helper methods for momentum evaluation /// @note These methods are taken from the Geant4 class /// G4PairProductionRelModel - Scalar screenFunction1(Scalar delta) const; - Scalar screenFunction2(Scalar delta) const; + double screenFunction1(double delta) const; + double screenFunction2(double delta) const; /// Electron mass. This is an static constant and not a member variable so the /// struct has no internal state. Otherwise, the interaction list breaks. - static const Scalar kElectronMass; + static const double kElectronMass; }; -inline Particle::Scalar PhotonConversion::screenFunction1(Scalar delta) const { +inline double PhotonConversion::screenFunction1(double delta) const { // Compute the value of the screening function 3*PHI1(delta) - PHI2(delta) return (delta > 1.4) ? 42.038 - 8.29 * std::log(delta + 0.958) : 42.184 - delta * (7.444 - 1.623 * delta); } -inline Particle::Scalar PhotonConversion::screenFunction2(Scalar delta) const { +inline double PhotonConversion::screenFunction2(double delta) const { // Compute the value of the screening function 1.5*PHI1(delta) // +0.5*PHI2(delta) return (delta > 1.4) ? 42.038 - 8.29 * std::log(delta + 0.958) @@ -124,16 +122,15 @@ inline Particle::Scalar PhotonConversion::screenFunction2(Scalar delta) const { } template -std::pair -PhotonConversion::generatePathLimits(generator_t& generator, - const Particle& particle) const { +std::pair PhotonConversion::generatePathLimits( + generator_t& generator, const Particle& particle) const { /// This method is based upon the Athena class PhotonConversionTool // Fast exit if not a photon or the energy is too low if (particle.pdg() != Acts::PdgParticle::eGamma || particle.absoluteMomentum() < (2 * kElectronMass)) { - return std::make_pair(std::numeric_limits::infinity(), - std::numeric_limits::infinity()); + return std::make_pair(std::numeric_limits::infinity(), + std::numeric_limits::infinity()); } // Use for the moment only Al data - Yung Tsai - Rev.Mod.Particle Physics Vol. @@ -149,73 +146,73 @@ PhotonConversion::generatePathLimits(generator_t& generator, // 1 p0 -7.01612e-03 8.43478e-01 1.62766e-04 1.11914e-05 // 2 p1 7.69040e-02 1.00059e+00 8.90718e-05 -8.41167e-07 // 3 p2 -6.07682e-01 5.13256e+00 6.07228e-04 -9.44448e-07 - constexpr Scalar p0 = -7.01612e-03; - constexpr Scalar p1 = 7.69040e-02; - constexpr Scalar p2 = -6.07682e-01; + constexpr double p0 = -7.01612e-03; + constexpr double p1 = 7.69040e-02; + constexpr double p2 = -6.07682e-01; // Calculate xi - const Scalar xi = p0 + p1 * std::pow(particle.absoluteMomentum(), p2); + const double xi = p0 + p1 * std::pow(particle.absoluteMomentum(), p2); - std::uniform_real_distribution uniformDistribution{0., 1.}; + std::uniform_real_distribution uniformDistribution{0., 1.}; // This is a transformation of eq. 3.75 return std::make_pair(-9. / 7. * std::log(conversionProbScaleFactor * (1 - uniformDistribution(generator))) / (1. - xi), - std::numeric_limits::infinity()); + std::numeric_limits::infinity()); } template -Particle::Scalar PhotonConversion::generateFirstChildEnergyFraction( - generator_t& generator, Scalar gammaMom) const { +double PhotonConversion::generateFirstChildEnergyFraction( + generator_t& generator, double gammaMom) const { /// This method is based upon the Geant4 class G4PairProductionRelModel /// @note This method is from the Geant4 class G4Element // // Compute Coulomb correction factor (Phys Rev. D50 3-1 (1994) page 1254) - constexpr Scalar k1 = 0.0083; - constexpr Scalar k2 = 0.20206; - constexpr Scalar k3 = 0.0020; // This term is missing in Athena - constexpr Scalar k4 = 0.0369; - constexpr Scalar alphaEM = 1. / 137.; - constexpr Scalar m_Z = 13.; // Aluminium - constexpr Scalar az2 = (alphaEM * m_Z) * (alphaEM * m_Z); - constexpr Scalar az4 = az2 * az2; - constexpr Scalar coulombFactor = + constexpr double k1 = 0.0083; + constexpr double k2 = 0.20206; + constexpr double k3 = 0.0020; // This term is missing in Athena + constexpr double k4 = 0.0369; + constexpr double alphaEM = 1. / 137.; + constexpr double m_Z = 13.; // Aluminium + constexpr double az2 = (alphaEM * m_Z) * (alphaEM * m_Z); + constexpr double az4 = az2 * az2; + constexpr double coulombFactor = (k1 * az4 + k2 + 1. / (1. + az2)) * az2 - (k3 * az4 + k4) * az4; - const Scalar logZ13 = std::log(m_Z) * 1. / 3.; - const Scalar FZ = 8. * (logZ13 + coulombFactor); - const Scalar deltaMax = exp((42.038 - FZ) * 0.1206) - 0.958; + const double logZ13 = std::log(m_Z) * 1. / 3.; + const double FZ = 8. * (logZ13 + coulombFactor); + const double deltaMax = std::exp((42.038 - FZ) * 0.1206) - 0.958; - const Scalar deltaPreFactor = 136. / std::pow(m_Z, 1. / 3.); - const Scalar eps0 = kElectronMass / gammaMom; - const Scalar deltaFactor = deltaPreFactor * eps0; - const Scalar deltaMin = 4. * deltaFactor; + const double deltaPreFactor = 136. / std::pow(m_Z, 1. / 3.); + const double eps0 = kElectronMass / gammaMom; + const double deltaFactor = deltaPreFactor * eps0; + const double deltaMin = 4. * deltaFactor; // Compute the limits of eps - const Scalar epsMin = + const double epsMin = std::max(eps0, 0.5 - 0.5 * std::sqrt(1. - deltaMin / deltaMax)); - const Scalar epsRange = 0.5 - epsMin; + const double epsRange = 0.5 - epsMin; // Sample the energy rate (eps) of the created electron (or positron) - const Scalar F10 = screenFunction1(deltaMin) - FZ; - const Scalar F20 = screenFunction2(deltaMin) - FZ; - const Scalar NormF1 = F10 * epsRange * epsRange; - const Scalar NormF2 = 1.5 * F20; + const double F10 = screenFunction1(deltaMin) - FZ; + const double F20 = screenFunction2(deltaMin) - FZ; + const double NormF1 = F10 * epsRange * epsRange; + const double NormF2 = 1.5 * F20; // We will need 3 uniform random number for each trial of sampling - Scalar greject = 0.; - Scalar eps = 0.; - std::uniform_real_distribution rndmEngine; + double greject = 0.; + double eps = 0.; + std::uniform_real_distribution rndmEngine; do { if (NormF1 > rndmEngine(generator) * (NormF1 + NormF2)) { eps = 0.5 - epsRange * std::pow(rndmEngine(generator), 1. / 3.); - const Scalar delta = deltaFactor / (eps * (1. - eps)); + const double delta = deltaFactor / (eps * (1. - eps)); greject = (screenFunction1(delta) - FZ) / F10; } else { eps = epsMin + epsRange * rndmEngine(generator); - const Scalar delta = deltaFactor / (eps * (1. - eps)); + const double delta = deltaFactor / (eps * (1. - eps)); greject = (screenFunction2(delta) - FZ) / F20; } } while (greject < rndmEngine(generator)); @@ -224,16 +221,16 @@ Particle::Scalar PhotonConversion::generateFirstChildEnergyFraction( } template -Particle::Vector3 PhotonConversion::generateChildDirection( +Acts::Vector3 PhotonConversion::generateChildDirection( generator_t& generator, const Particle& particle) const { /// This method is based upon the Athena class PhotonConversionTool // Following the Geant4 approximation from L. Urban // the azimutal angle - Scalar theta = kElectronMass / particle.energy(); + double theta = kElectronMass / particle.energy(); - std::uniform_real_distribution uniformDistribution{0., 1.}; - const Scalar u = -std::log(uniformDistribution(generator) * + std::uniform_real_distribution uniformDistribution{0., 1.}; + const double u = -std::log(uniformDistribution(generator) * uniformDistribution(generator)) * 1.6; @@ -257,20 +254,20 @@ Particle::Vector3 PhotonConversion::generateChildDirection( } inline std::array PhotonConversion::generateChildren( - const Particle& photon, Scalar childEnergy, - const Particle::Vector3& childDirection) const { + const Particle& photon, double childEnergy, + const Acts::Vector3& childDirection) const { using namespace Acts::UnitLiterals; // Calculate the child momentum - const Scalar massChild = kElectronMass; - const Scalar momentum1 = + const double massChild = kElectronMass; + const double momentum1 = sqrt(childEnergy * childEnergy - massChild * massChild); // Use energy-momentum conservation for the other child - const Particle::Vector3 vtmp = + const Acts::Vector3 vtmp = photon.fourMomentum().template segment<3>(Acts::eMom0) - momentum1 * childDirection; - const Scalar momentum2 = vtmp.norm(); + const double momentum2 = vtmp.norm(); // The daughter particles are created with the explicit electron mass used in // the calculations for consistency. Using the full Particle constructor with @@ -304,17 +301,16 @@ bool PhotonConversion::run(generator_t& generator, Particle& particle, } // Fast exit if momentum is too low - const Scalar p = particle.absoluteMomentum(); + const double p = particle.absoluteMomentum(); if (p < (2 * kElectronMass)) { return false; } // Get one child energy - const Scalar childEnergy = p * generateFirstChildEnergyFraction(generator, p); + const double childEnergy = p * generateFirstChildEnergyFraction(generator, p); // Now get the deflection - const Particle::Vector3 childDir = - generateChildDirection(generator, particle); + const Acts::Vector3 childDir = generateChildDirection(generator, particle); // Produce the final state const std::array finalState = diff --git a/Fatras/src/Physics/BetheHeitler.cpp b/Fatras/src/Physics/BetheHeitler.cpp index 6191d832c38..bb635ecb934 100644 --- a/Fatras/src/Physics/BetheHeitler.cpp +++ b/Fatras/src/Physics/BetheHeitler.cpp @@ -47,8 +47,8 @@ ActsFatras::Particle ActsFatras::BetheHeitler::bremPhoton( theta *= (rndTheta1 < 0.25) ? u : u / 3.; // 9./(9.+27) = 0.25 } - Vector3 particleDirection = particle.direction(); - Vector3 photonDirection = particleDirection; + Acts::Vector3 particleDirection = particle.direction(); + Acts::Vector3 photonDirection = particleDirection; // construct the combined rotation to the scattered direction Acts::RotationMatrix3 rotation( diff --git a/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp b/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp index 497d530663b..84ab93921fe 100644 --- a/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp +++ b/Fatras/src/Physics/NuclearInteraction/NuclearInteraction.cpp @@ -74,7 +74,7 @@ double NuclearInteraction::sampleContinuousValues( bool interpolate) const { // Fast exit if (distribution.second.empty()) { - return std::numeric_limits::infinity(); + return std::numeric_limits::infinity(); } // Find the bin @@ -82,7 +82,7 @@ double NuclearInteraction::sampleContinuousValues( std::numeric_limits::max() * rnd); // Fast exit for non-normalised CDFs like interaction probability if (int_rnd > distribution.second.back()) { - return std::numeric_limits::infinity(); + return std::numeric_limits::infinity(); } const auto it = std::upper_bound(distribution.second.begin(), distribution.second.end(), int_rnd); From ba4e47a6b662fb4d7679ce95719ac2fcf4cce543 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Fri, 22 Nov 2024 17:37:46 +0100 Subject: [PATCH 6/8] fixes --- .../Generators/Pythia8ProcessGenerator.cpp | 3 +-- Fatras/Geant4/src/Geant4Decay.cpp | 2 +- .../UnitTests/Core/Seeding/HoughTransformTest.cpp | 15 +++++++-------- 3 files changed, 9 insertions(+), 11 deletions(-) diff --git a/Examples/Algorithms/GeneratorsPythia8/ActsExamples/Generators/Pythia8ProcessGenerator.cpp b/Examples/Algorithms/GeneratorsPythia8/ActsExamples/Generators/Pythia8ProcessGenerator.cpp index cd80d64bf08..1b774a7bdf0 100644 --- a/Examples/Algorithms/GeneratorsPythia8/ActsExamples/Generators/Pythia8ProcessGenerator.cpp +++ b/Examples/Algorithms/GeneratorsPythia8/ActsExamples/Generators/Pythia8ProcessGenerator.cpp @@ -121,8 +121,7 @@ Pythia8Generator::operator()(RandomEngine& rng) { } // create the primary vertex - vertices.emplace_back(SimVertexBarcode{0}, - SimVertex::Vector4(0., 0., 0., 0.)); + vertices.emplace_back(SimVertexBarcode{0}, Acts::Vector4(0., 0., 0., 0.)); // convert generated final state particles into internal format for (int ip = 0; ip < m_pythia8->event.size(); ++ip) { diff --git a/Fatras/Geant4/src/Geant4Decay.cpp b/Fatras/Geant4/src/Geant4Decay.cpp index 1afef4ca1d0..bf23566aa05 100644 --- a/Fatras/Geant4/src/Geant4Decay.cpp +++ b/Fatras/Geant4/src/Geant4Decay.cpp @@ -64,7 +64,7 @@ std::vector ActsFatras::Geant4Decay::decayParticle( // Convert the decay product from Geant4 to Acts const G4ThreeVector& mom = prod->GetMomentum(); - constexpr Scalar convertEnergy = Acts::UnitConstants::GeV / CLHEP::GeV; + constexpr double convertEnergy = Acts::UnitConstants::GeV / CLHEP::GeV; Acts::Vector3 amgMom(mom.x(), mom.y(), mom.z()); amgMom *= convertEnergy; const std::int32_t pdg = prod->GetPDGcode(); diff --git a/Tests/UnitTests/Core/Seeding/HoughTransformTest.cpp b/Tests/UnitTests/Core/Seeding/HoughTransformTest.cpp index 30e3164f62c..b469c6d7c1c 100644 --- a/Tests/UnitTests/Core/Seeding/HoughTransformTest.cpp +++ b/Tests/UnitTests/Core/Seeding/HoughTransformTest.cpp @@ -19,17 +19,16 @@ namespace Acts::Test { -using Scalar = Acts::ActsScalar; auto logger = Acts::getDefaultLogger("UnitTests", Acts::Logging::VERBOSE); struct DriftCircle { - Scalar y{0.}; - Scalar z{0.}; - Scalar rDrift{0.}; - Scalar rDriftError{0.}; + double y{0.}; + double z{0.}; + double rDrift{0.}; + double rDriftError{0.}; - DriftCircle(const Scalar _y, const Scalar _z, const Scalar _r, - const Scalar _rUncert) + DriftCircle(const double _y, const double _z, const double _r, + const double _rUncert) : y{_y}, z{_z}, rDrift{_r}, rDriftError{_rUncert} {} }; @@ -92,7 +91,7 @@ BOOST_AUTO_TEST_CASE(hough_transform_seeder) { Acts::HoughTransformUtils::HoughPlane houghPlane(planeCfg); - // also insantiate the peak finder + // also instantiate the peak finder Acts::HoughTransformUtils::PeakFinders::IslandsAroundMax< Acts::GeometryIdentifier::Value> peakFinder(peakFinderCfg); From cf5c818fdf5fa81a1523b067102c8b88d5a3c811 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Sat, 23 Nov 2024 10:04:31 +0100 Subject: [PATCH 7/8] boundingbox test --- .../Core/Utilities/BoundingBoxTest.cpp | 129 +++++++++--------- 1 file changed, 62 insertions(+), 67 deletions(-) diff --git a/Tests/UnitTests/Core/Utilities/BoundingBoxTest.cpp b/Tests/UnitTests/Core/Utilities/BoundingBoxTest.cpp index 675d7d49320..b65d533acfa 100644 --- a/Tests/UnitTests/Core/Utilities/BoundingBoxTest.cpp +++ b/Tests/UnitTests/Core/Utilities/BoundingBoxTest.cpp @@ -37,13 +37,11 @@ namespace Acts::Test { struct Object {}; -using BoundingBoxScalar = ActsScalar; +using ObjectBBox = Acts::AxisAlignedBoundingBox; -using ObjectBBox = Acts::AxisAlignedBoundingBox; - -using Vector2F = Eigen::Matrix; -using Vector3F = Eigen::Matrix; -using AngleAxis3F = Eigen::AngleAxis; +using Vector2F = Eigen::Matrix; +using Vector3F = Eigen::Matrix; +using AngleAxis3F = Eigen::AngleAxis; std::filesystem::path tmp_path = []() { auto tmpPath = std::filesystem::temp_directory_path() / "acts_unit_tests"; @@ -59,11 +57,11 @@ std::ofstream tmp(const std::string& path) { BOOST_AUTO_TEST_CASE(box_construction) { BOOST_TEST_CONTEXT("2D") { Object o; - using Box = Acts::AxisAlignedBoundingBox; + using Box = Acts::AxisAlignedBoundingBox; Box bb(&o, {-1, -1}, {2, 2}); typename Box::transform_type rot; - rot = Eigen::Rotation2D(std::numbers::pi / 7.); + rot = Eigen::Rotation2D(std::numbers::pi / 7.); Box bb_rot = bb.transformed(rot); CHECK_CLOSE_ABS(bb_rot.min(), Vector2F(-1.76874, -1.33485), 1e-4); @@ -72,7 +70,7 @@ BOOST_AUTO_TEST_CASE(box_construction) { BOOST_TEST_CONTEXT("3D") { Object o; - using Box = Acts::AxisAlignedBoundingBox; + using Box = Acts::AxisAlignedBoundingBox; Box bb(&o, {-1, -1, -1}, {2, 2, 2}); typename Box::transform_type rot; @@ -141,14 +139,14 @@ BOOST_AUTO_TEST_CASE(intersect_points) { BOOST_AUTO_TEST_CASE(intersect_rays) { BOOST_TEST_CONTEXT("2D") { - using Box = AxisAlignedBoundingBox; + using Box = AxisAlignedBoundingBox; Object o; Box bb(&o, {-1, -1}, {1, 1}); // ray in positive x direction - Ray ray({-2, 0}, {1, 0}); + Ray ray({-2, 0}, {1, 0}); BOOST_CHECK(bb.intersect(ray)); ray = {{-2, 2}, {1, 0}}; @@ -273,10 +271,10 @@ BOOST_AUTO_TEST_CASE(intersect_rays) { // let's make sure it also works in 3d ObjectBBox bb3(&o, {-1, -1, -1}, {1, 1, 1}); - Ray ray3({0, 0, -2}, {0, 0, 1}); + Ray ray3({0, 0, -2}, {0, 0, 1}); BOOST_CHECK(bb3.intersect(ray3)); - PlyVisualization3D ply; + PlyVisualization3D ply; ray3.draw(ply); auto os = tmp("ray3d.ply"); @@ -290,7 +288,7 @@ BOOST_AUTO_TEST_CASE(intersect_rays) { // let's make sure it also works in 3d ObjectBBox bb3(&o, {-1, -1, -1}, {1, 1, 1}); - Ray ray3({0, 0, -2}, {0, 0, 1}); + Ray ray3({0, 0, -2}, {0, 0, 1}); BOOST_CHECK(bb3.intersect(ray3)); // facing away from box @@ -451,7 +449,7 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { return os; }; - using Frustum2 = Frustum; + using Frustum2 = Frustum; std::ofstream os; @@ -460,19 +458,18 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { const std::size_t nSteps = 10; // Visualise the parameters - const BoundingBoxScalar min = -20; - const BoundingBoxScalar max = 20; + const double min = -20; + const double max = 20; os = make_svg("frust2d_parameters.svg", svgWidth, svgHeight); - const BoundingBoxScalar step = (max - min) / nSteps; + const double step = (max - min) / nSteps; for (std::size_t i = 0; i <= nSteps; i++) { for (std::size_t j = 0; j <= nSteps; j++) { Vector2F dir(1, 0); Vector2F origin(min + step * i, min + step * j); origin.x() *= 1.10; - Eigen::Rotation2D rot(2 * std::numbers::pi / nSteps * - i); - BoundingBoxScalar angle = std::numbers::pi / 2. / nSteps * j; + Eigen::Rotation2D rot(2 * std::numbers::pi / nSteps * i); + double angle = std::numbers::pi / 2. / nSteps * j; Frustum2 fr(origin, rot * dir, angle); fr.svg(os, svgWidth, svgHeight, 2); } @@ -481,18 +478,18 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { os << ""; os.close(); - const BoundingBoxScalar unit = 20; + const double unit = 20.; - using Box = AxisAlignedBoundingBox; + using Box = AxisAlignedBoundingBox; Object o; - Box::Size size(Eigen::Matrix(2, 2)); + Box::Size size(Eigen::Matrix(2, 2)); - const BoundingBoxScalar minX = -20; - const BoundingBoxScalar minY = -20; - const BoundingBoxScalar maxX = 20; - const BoundingBoxScalar maxY = 20; - const BoundingBoxScalar stepX = (maxX - minX) / nSteps; - const BoundingBoxScalar stepY = (maxY - minY) / nSteps; + const double minX = -20.; + const double minY = -20.; + const double maxX = 20.; + const double maxY = 20.; + const double stepX = (maxX - minX) / nSteps; + const double stepY = (maxY - minY) / nSteps; std::set idxsAct; @@ -579,10 +576,10 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { boxes.reserve((nSteps + 1) * (nSteps + 1)); for (std::size_t i = 0; i <= nSteps; i++) { for (std::size_t j = 0; j <= nSteps; j++) { - boxes.emplace_back(&o, - Eigen::Matrix{ - minX + i * stepX, minY + j * stepY}, - size); + boxes.emplace_back( + &o, + Eigen::Matrix{minX + i * stepX, minY + j * stepY}, + size); std::stringstream st; st << boxes.size() - 1; @@ -605,9 +602,9 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { } } - PlyVisualization3D helper; + PlyVisualization3D helper; BOOST_TEST_CONTEXT("3D - 3 Sides") { - using Frustum3 = Frustum; + using Frustum3 = Frustum; std::ofstream os; std::size_t n = 10; const std::size_t nSteps = 5; @@ -619,7 +616,7 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { auto make = [&](double angle, const Vector3F& origin, std::ofstream& osTmp) { helper.clear(); - BoundingBoxScalar far = 1; + double far = 1.; Frustum3 fr(origin, {0, 0, 1}, angle); fr.draw(helper, far); fr = Frustum3(origin, {0, 0, -1}, angle); @@ -655,7 +652,7 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { for (std::size_t i = 0; i <= n; i++) { Vector3F origin(i * 4, 0, 0); AngleAxis3F rot(std::numbers::pi / n * i, Vector3F::UnitY()); - BoundingBoxScalar angle = (std::numbers::pi / 2.) / n * (1 + i); + double angle = (std::numbers::pi / 2.) / n * (1 + i); Vector3F dir(1, 0, 0); Frustum3 fr(origin, rot * dir, angle); fr.draw(helper, 2); @@ -850,16 +847,16 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { step = (max - min) / n; Object o; - using Box = AxisAlignedBoundingBox; - Box::Size size(Eigen::Matrix(2, 2, 2)); + using Box = AxisAlignedBoundingBox; + Box::Size size(Eigen::Matrix(2, 2, 2)); std::size_t idx = 0; for (std::size_t i = 0; i <= n; i++) { for (std::size_t j = 0; j <= n; j++) { for (std::size_t k = 0; k <= n; k++) { - Eigen::Matrix pos( - min + i * step, min + j * step, min + k * step); + Eigen::Matrix pos(min + i * step, min + j * step, + min + k * step); Box bb(&o, pos, size); Color color = {255, 0, 0}; @@ -883,7 +880,7 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { } BOOST_TEST_CONTEXT("3D - 4 Sides") { - using Frustum34 = Frustum; + using Frustum34 = Frustum; const std::size_t n = 10; double min = -10; double max = 10; @@ -902,7 +899,7 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { Vector3F origin(min + i * step, min + j * step, min + k * step); Vector3F dir(1, 0, 0); - BoundingBoxScalar piOverS = std::numbers::pi_v / s; + double piOverS = std::numbers::pi / s; AngleAxis3F rot; rot = AngleAxis3F(piOverS * i, Vector3F::UnitX()) * AngleAxis3F(piOverS * j, Vector3F::UnitY()) * @@ -1104,15 +1101,15 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { fr.draw(helper, 50); Object o; - using Box = AxisAlignedBoundingBox; - Box::Size size(Eigen::Matrix(2, 2, 2)); + using Box = AxisAlignedBoundingBox; + Box::Size size(Eigen::Matrix(2, 2, 2)); std::size_t idx = 0; for (std::size_t i = 0; i <= n; i++) { for (std::size_t j = 0; j <= n; j++) { for (std::size_t k = 0; k <= n; k++) { - Eigen::Matrix pos( - min + i * step, min + j * step, min + k * step); + Eigen::Matrix pos(min + i * step, min + j * step, + min + k * step); Box bb(&o, pos, size); Color color = {255, 0, 0}; @@ -1136,13 +1133,13 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { } BOOST_TEST_CONTEXT("3D - 5 Sides") { - using Frustum = Frustum; - using Box = AxisAlignedBoundingBox; - Box::Size size(Eigen::Matrix(2, 2, 2)); + using Frustum = Frustum; + using Box = AxisAlignedBoundingBox; + Box::Size size(Acts::Vector3(2, 2, 2)); Object o; - PlyVisualization3D ply; + PlyVisualization3D ply; Frustum fr({0, 0, 0}, {0, 0, 1}, std::numbers::pi / 8.); fr.draw(ply, 10); @@ -1158,17 +1155,16 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { } BOOST_TEST_CONTEXT("3D - 10 Sides") { - using Frustum = Frustum; - using Box = AxisAlignedBoundingBox; - using vec3 = Eigen::Matrix; - Box::Size size(vec3(2, 2, 2)); + using Frustum = Frustum; + using Box = AxisAlignedBoundingBox; + Box::Size size(Acts::Vector3(2, 2, 2)); Object o; - PlyVisualization3D ply; + PlyVisualization3D ply; - vec3 pos = {-12.4205, 29.3578, 44.6207}; - vec3 dir = {-0.656862, 0.48138, 0.58035}; + Acts::Vector3 pos = {-12.4205, 29.3578, 44.6207}; + Acts::Vector3 dir = {-0.656862, 0.48138, 0.58035}; Frustum fr(pos, dir, 0.972419); fr.draw(ply, 10); @@ -1183,20 +1179,19 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { } BOOST_TEST_CONTEXT("3D - 4 Sides - Big box") { - using Frustum = Frustum; - using Box = AxisAlignedBoundingBox; - using vec3 = Eigen::Matrix; + using Frustum = Frustum; + using Box = AxisAlignedBoundingBox; Object o; - PlyVisualization3D ply; + PlyVisualization3D ply; - vec3 pos = {0, 0, 0}; - vec3 dir = {0, 0, 1}; + Acts::Vector3 pos = {0, 0, 0}; + Acts::Vector3 dir = {0, 0, 1}; Frustum fr(pos, dir, 0.972419); fr.draw(ply, 10); - Box::Size size(vec3(100, 100, 2)); + Box::Size size(Acts::Vector3(100, 100, 2)); Box bb(&o, pos + dir * 7, size); bb.draw(ply); @@ -1210,7 +1205,7 @@ BOOST_AUTO_TEST_CASE(frustum_intersect) { BOOST_AUTO_TEST_CASE(ostream_operator) { Object o; - using Box = Acts::AxisAlignedBoundingBox; + using Box = Acts::AxisAlignedBoundingBox; Box bb(&o, {-1, -1}, {2, 2}); std::stringstream ss; From cb206c1c4ad576647366990d95ee9d7bd04e7c82 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Sat, 23 Nov 2024 23:08:08 +0100 Subject: [PATCH 8/8] remaining cases --- .../EventData/GenericBoundTrackParameters.hpp | 21 ++++++----- .../GenericCurvilinearTrackParameters.hpp | 7 ++-- .../EventData/GenericFreeTrackParameters.hpp | 27 +++++++------- .../MultiComponentTrackParameters.hpp | 21 ++++++----- .../Acts/EventData/TrackParametersConcept.hpp | 1 - .../Acts/EventData/TrackStateProxy.hpp | 22 +++++------- .../Acts/Propagator/ConstrainedStep.hpp | 36 +++++++++---------- .../Acts/Plugins/Hashing/HashingAnnoy.ipp | 22 ++++++------ .../Acts/Plugins/Hashing/HashingTraining.ipp | 18 +++++----- .../Fatras/Kernel/SimulationActorTests.cpp | 17 ++++----- 10 files changed, 85 insertions(+), 107 deletions(-) diff --git a/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp b/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp index 38f297764b1..22567d8adf5 100644 --- a/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp +++ b/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp @@ -36,7 +36,6 @@ namespace Acts { template class GenericBoundTrackParameters { public: - using Scalar = ActsScalar; using ParametersVector = BoundVector; using CovarianceMatrix = BoundSquareMatrix; using ParticleHypothesis = particle_hypothesis_t; @@ -92,10 +91,10 @@ class GenericBoundTrackParameters { /// successfully be converted to on-surface parameters. static Result create( std::shared_ptr surface, const GeometryContext& geoCtx, - const Vector4& pos4, const Vector3& dir, Scalar qOverP, + const Vector4& pos4, const Vector3& dir, double qOverP, std::optional cov, ParticleHypothesis particleHypothesis, - ActsScalar tolerance = s_onSurfaceTolerance) { + double tolerance = s_onSurfaceTolerance) { Result bound = transformFreeToBoundParameters(pos4.segment<3>(ePos0), pos4[eTime], dir, qOverP, *surface, geoCtx, tolerance); @@ -168,7 +167,7 @@ class GenericBoundTrackParameters { /// /// @tparam kIndex Track parameter index template - Scalar get() const { + double get() const { return m_params[kIndex]; } @@ -203,14 +202,14 @@ class GenericBoundTrackParameters { return m_surface->localToGlobal(geoCtx, localPosition(), direction()); } /// Time coordinate. - Scalar time() const { return m_params[eBoundTime]; } + double time() const { return m_params[eBoundTime]; } /// Phi direction. - Scalar phi() const { return m_params[eBoundPhi]; } + double phi() const { return m_params[eBoundPhi]; } /// Theta direction. - Scalar theta() const { return m_params[eBoundTheta]; } + double theta() const { return m_params[eBoundTheta]; } /// Charge over momentum. - Scalar qOverP() const { return m_params[eBoundQOverP]; } + double qOverP() const { return m_params[eBoundQOverP]; } /// Unit direction three-vector, i.e. the normalized momentum /// three-vector. @@ -219,18 +218,18 @@ class GenericBoundTrackParameters { m_params[eBoundTheta]); } /// Absolute momentum. - Scalar absoluteMomentum() const { + double absoluteMomentum() const { return m_particleHypothesis.extractMomentum(m_params[eBoundQOverP]); } /// Transverse momentum. - Scalar transverseMomentum() const { + double transverseMomentum() const { return std::sin(m_params[eBoundTheta]) * absoluteMomentum(); } /// Momentum three-vector. Vector3 momentum() const { return absoluteMomentum() * direction(); } /// Particle electric charge. - Scalar charge() const { + double charge() const { return m_particleHypothesis.extractCharge(get()); } diff --git a/Core/include/Acts/EventData/GenericCurvilinearTrackParameters.hpp b/Core/include/Acts/EventData/GenericCurvilinearTrackParameters.hpp index 6e4ed9bae7a..9764ad68cfc 100644 --- a/Core/include/Acts/EventData/GenericCurvilinearTrackParameters.hpp +++ b/Core/include/Acts/EventData/GenericCurvilinearTrackParameters.hpp @@ -30,7 +30,6 @@ class GenericCurvilinearTrackParameters using Base = GenericBoundTrackParameters; public: - using Scalar = ActsScalar; using ParametersVector = BoundVector; using CovarianceMatrix = BoundSquareMatrix; using ParticleHypothesis = particle_hypothesis_t; @@ -43,7 +42,7 @@ class GenericCurvilinearTrackParameters /// @param cov Curvilinear bound parameters covariance matrix /// @param particleHypothesis Particle hypothesis GenericCurvilinearTrackParameters(const Vector4& pos4, const Vector3& dir, - Scalar qOverP, + double qOverP, std::optional cov, ParticleHypothesis particleHypothesis) : Base(CurvilinearSurface(pos4.segment<3>(ePos0), dir).surface(), @@ -58,8 +57,8 @@ class GenericCurvilinearTrackParameters /// @param qOverP Charge over momentum /// @param cov Curvilinear bound parameters covariance matrix /// @param particleHypothesis Particle hypothesis - GenericCurvilinearTrackParameters(const Vector4& pos4, Scalar phi, - Scalar theta, Scalar qOverP, + GenericCurvilinearTrackParameters(const Vector4& pos4, double phi, + double theta, double qOverP, std::optional cov, ParticleHypothesis particleHypothesis) : Base(CurvilinearSurface(pos4.segment<3>(ePos0), diff --git a/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp b/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp index 7a2b4b46522..b0499dd3137 100644 --- a/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp +++ b/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp @@ -33,7 +33,6 @@ namespace Acts { template class GenericFreeTrackParameters { public: - using Scalar = ActsScalar; using ParametersVector = FreeVector; using CovarianceMatrix = FreeSquareMatrix; using ParticleHypothesis = particle_hypothesis_t; @@ -66,7 +65,7 @@ class GenericFreeTrackParameters { /// @param cov Free parameters covariance matrix /// @param particleHypothesis Particle hypothesis GenericFreeTrackParameters(const Vector4& pos4, const Vector3& dir, - Scalar qOverP, std::optional cov, + double qOverP, std::optional cov, ParticleHypothesis particleHypothesis) : m_params(FreeVector::Zero()), m_cov(std::move(cov)), @@ -91,8 +90,8 @@ class GenericFreeTrackParameters { /// @param qOverP Charge over momentum /// @param cov Free parameters covariance matrix /// @param particleHypothesis Particle hypothesis - GenericFreeTrackParameters(const Vector4& pos4, Scalar phi, Scalar theta, - Scalar qOverP, std::optional cov, + GenericFreeTrackParameters(const Vector4& pos4, double phi, double theta, + double qOverP, std::optional cov, ParticleHypothesis particleHypothesis) : m_params(FreeVector::Zero()), m_cov(std::move(cov)), @@ -146,7 +145,7 @@ class GenericFreeTrackParameters { /// /// @tparam kIndex Track parameter index template - Scalar get() const { + double get() const { return m_params[kIndex]; } @@ -162,33 +161,33 @@ class GenericFreeTrackParameters { /// Spatial position three-vector. Vector3 position() const { return m_params.segment<3>(eFreePos0); } /// Time coordinate. - Scalar time() const { return m_params[eFreeTime]; } + double time() const { return m_params[eFreeTime]; } /// Phi direction. - Scalar phi() const { return VectorHelpers::phi(direction()); } + double phi() const { return VectorHelpers::phi(direction()); } /// Theta direction. - Scalar theta() const { return VectorHelpers::theta(direction()); } + double theta() const { return VectorHelpers::theta(direction()); } /// Charge over momentum. - Scalar qOverP() const { return m_params[eFreeQOverP]; } + double qOverP() const { return m_params[eFreeQOverP]; } /// Unit direction three-vector, i.e. the normalized momentum three-vector. Vector3 direction() const { return m_params.segment<3>(eFreeDir0).normalized(); } /// Absolute momentum. - Scalar absoluteMomentum() const { + double absoluteMomentum() const { return m_particleHypothesis.extractMomentum(m_params[eFreeQOverP]); } /// Transverse momentum. - Scalar transverseMomentum() const { + double transverseMomentum() const { // direction vector w/ arbitrary normalization can be parametrized as // [f*sin(theta)*cos(phi), f*sin(theta)*sin(phi), f*cos(theta)] // w/ f,sin(theta) positive, the transverse magnitude is then // sqrt(f^2*sin^2(theta)) = f*sin(theta) - Scalar transverseMagnitude2 = + double transverseMagnitude2 = square(m_params[eFreeDir0]) + square(m_params[eFreeDir1]); // absolute magnitude is f by construction - Scalar magnitude2 = transverseMagnitude2 + square(m_params[eFreeDir2]); + double magnitude2 = transverseMagnitude2 + square(m_params[eFreeDir2]); // such that we can extract sin(theta) = f*sin(theta) / f return std::sqrt(transverseMagnitude2 / magnitude2) * absoluteMomentum(); } @@ -196,7 +195,7 @@ class GenericFreeTrackParameters { Vector3 momentum() const { return absoluteMomentum() * direction(); } /// Particle electric charge. - Scalar charge() const { + double charge() const { return m_particleHypothesis.extractCharge(get()); } diff --git a/Core/include/Acts/EventData/MultiComponentTrackParameters.hpp b/Core/include/Acts/EventData/MultiComponentTrackParameters.hpp index 70ad87f8d77..cf932fc0470 100644 --- a/Core/include/Acts/EventData/MultiComponentTrackParameters.hpp +++ b/Core/include/Acts/EventData/MultiComponentTrackParameters.hpp @@ -34,7 +34,6 @@ class MultiComponentBoundTrackParameters { public: using Parameters = BoundTrackParameters; using ParticleHypothesis = Parameters::ParticleHypothesis; - using Scalar = typename Parameters::Scalar; using ParametersVector = typename Parameters::ParametersVector; using CovarianceMatrix = typename Parameters::CovarianceMatrix; @@ -159,7 +158,7 @@ class MultiComponentBoundTrackParameters { /// /// @tparam kIndex Track parameter index template - Scalar get() const { + double get() const { return reduce([&](const Parameters& p) { return p.get(); }); } @@ -180,7 +179,7 @@ class MultiComponentBoundTrackParameters { } /// Time coordinate. - Scalar time() const { + double time() const { return reduce([](const Parameters& p) { return p.time(); }); } @@ -192,21 +191,21 @@ class MultiComponentBoundTrackParameters { } /// Phi direction. - Scalar phi() const { return VectorHelpers::phi(direction()); } + double phi() const { return VectorHelpers::phi(direction()); } /// Theta direction. - Scalar theta() const { return VectorHelpers::theta(direction()); } + double theta() const { return VectorHelpers::theta(direction()); } /// Charge over momentum. - Scalar qOverP() const { return get(); } + double qOverP() const { return get(); } /// Absolute momentum. - Scalar absoluteMomentum() const { + double absoluteMomentum() const { return reduce([](const Parameters& p) { return p.absoluteMomentum(); }); } /// Transverse momentum. - Scalar transverseMomentum() const { + double transverseMomentum() const { return reduce([](const Parameters& p) { return p.transverseMomentum(); }); } @@ -216,7 +215,7 @@ class MultiComponentBoundTrackParameters { } /// Particle electric charge. - Scalar charge() const { + double charge() const { return reduce([](const Parameters& p) { return p.charge(); }); } @@ -238,8 +237,8 @@ class MultiComponentCurvilinearTrackParameters using covariance_t = BoundSquareMatrix; public: - using ConstructionTuple = std::tuple; + using ConstructionTuple = + std::tuple; private: using Base = MultiComponentBoundTrackParameters; diff --git a/Core/include/Acts/EventData/TrackParametersConcept.hpp b/Core/include/Acts/EventData/TrackParametersConcept.hpp index 5e4982545a4..4369fe57141 100644 --- a/Core/include/Acts/EventData/TrackParametersConcept.hpp +++ b/Core/include/Acts/EventData/TrackParametersConcept.hpp @@ -21,7 +21,6 @@ class Surface; namespace Concepts { template concept BasicTrackParameters = requires { - typename Parameters::Scalar; typename Parameters::ParametersVector; typename Parameters::CovarianceMatrix; diff --git a/Core/include/Acts/EventData/TrackStateProxy.hpp b/Core/include/Acts/EventData/TrackStateProxy.hpp index 741bfa5b511..b7b991fa2ee 100644 --- a/Core/include/Acts/EventData/TrackStateProxy.hpp +++ b/Core/include/Acts/EventData/TrackStateProxy.hpp @@ -82,17 +82,15 @@ template struct FixedSizeTypes { constexpr static auto Flags = Eigen::ColMajor | Eigen::AutoAlign; - using Scalar = ActsScalar; - // single items - using Coefficients = Eigen::Matrix; - using Covariance = Eigen::Matrix; + using Coefficients = Eigen::Matrix; + using Covariance = Eigen::Matrix; using CoefficientsMap = Eigen::Map>; using CovarianceMap = Eigen::Map>; - using DynamicCoefficients = Eigen::Matrix; + using DynamicCoefficients = Eigen::Matrix; using DynamicCovariance = - Eigen::Matrix; + Eigen::Matrix; using DynamicCoefficientsMap = Eigen::Map>; using DynamicCovarianceMap = @@ -105,11 +103,9 @@ template struct DynamicSizeTypes { constexpr static auto Flags = Eigen::ColMajor | Eigen::AutoAlign; - using Scalar = ActsScalar; - - using Coefficients = Eigen::Matrix; + using Coefficients = Eigen::Matrix; using Covariance = - Eigen::Matrix; + Eigen::Matrix; using CoefficientsMap = Eigen::Map>; using CovarianceMap = Eigen::Map>; }; @@ -119,8 +115,6 @@ struct DynamicSizeTypes { // This is public template struct TrackStateTraits { - using Scalar = ActsScalar; - using Parameters = typename detail_lt::FixedSizeTypes::CoefficientsMap; using Covariance = @@ -135,8 +129,8 @@ struct TrackStateTraits { typename detail_lt::DynamicSizeTypes::CovarianceMap; constexpr static auto ProjectorFlags = Eigen::RowMajor | Eigen::AutoAlign; - using Projector = Eigen::Matrix; - using EffectiveProjector = Eigen::Matrix; + using EffectiveProjector = Eigen::Matrix; }; diff --git a/Core/include/Acts/Propagator/ConstrainedStep.hpp b/Core/include/Acts/Propagator/ConstrainedStep.hpp index e957b001281..412f308bac9 100644 --- a/Core/include/Acts/Propagator/ConstrainedStep.hpp +++ b/Core/include/Acts/Propagator/ConstrainedStep.hpp @@ -40,11 +40,9 @@ namespace Acts { /// The hierarchy is: /// - Overstepping resolution / backpropagation /// - Convergence -/// - Step into the void with `std::numeric_limits::max()` +/// - Step into the void with `std::numeric_limits::max()` class ConstrainedStep { public: - using Scalar = ActsScalar; - /// the types of constraints /// from actor - this would be a typical navigation step /// from aborter - this would be a target condition @@ -53,26 +51,26 @@ class ConstrainedStep { constexpr ConstrainedStep() = default; - /// constructor from Scalar + /// constructor /// @param value is the user given initial value - constexpr explicit ConstrainedStep(Scalar value) { setUser(value); } + constexpr explicit ConstrainedStep(double value) { setUser(value); } - /// set accuracy by one Scalar + /// set accuracy /// /// this will set only the accuracy, as this is the most /// exposed to the Propagator /// /// @param value is the new accuracy value - constexpr void setAccuracy(Scalar value) { + constexpr void setAccuracy(double value) { assert(value > 0 && "ConstrainedStep accuracy must be > 0."); // set the accuracy value m_accuracy = value; } - /// set user by one Scalar + /// set user /// /// @param value is the new user value - constexpr void setUser(Scalar value) { + constexpr void setUser(double value) { // TODO enable assert; see https://github.com/acts-project/acts/issues/2543 // assert(value != 0 && "ConstrainedStep user must be != 0."); // set the user value @@ -80,20 +78,20 @@ class ConstrainedStep { } /// returns the min step size - constexpr Scalar value() const { - Scalar min = *std::min_element(m_values.begin(), m_values.end()); + constexpr double value() const { + double min = *std::min_element(m_values.begin(), m_values.end()); // accuracy is always positive and therefore handled separately - Scalar result = std::min(std::abs(min), m_accuracy); + double result = std::min(std::abs(min), m_accuracy); return std::signbit(min) ? -result : result; } /// Access a specific value /// /// @param type is the requested parameter type - constexpr Scalar value(Type type) const { return m_values[type]; } + constexpr double value(Type type) const { return m_values[type]; } /// Access the accuracy value - constexpr Scalar accuracy() const { return m_accuracy; } + constexpr double accuracy() const { return m_accuracy; } /// release a certain constraint value /// @@ -111,7 +109,7 @@ class ConstrainedStep { /// @param value is the new value to be updated /// @param type is the constraint type /// @param releaseStep Allow step size to increase again - constexpr void update(Scalar value, Type type, bool releaseStep = false) { + constexpr void update(double value, Type type, bool releaseStep = false) { if (releaseStep) { release(type); } @@ -127,7 +125,7 @@ class ConstrainedStep { std::ostream& toStream(std::ostream& os) const { // Helper method to avoid unreadable screen output - auto streamValue = [&](Scalar val) { + auto streamValue = [&](double val) { os << std::setw(5); if (std::abs(val) == kNotSet) { os << (val > 0 ? "+∞" : "-∞"); @@ -156,12 +154,12 @@ class ConstrainedStep { } private: - static constexpr auto kNotSet = std::numeric_limits::max(); + static constexpr auto kNotSet = std::numeric_limits::max(); /// the step size tuple - std::array m_values = {kNotSet, kNotSet, kNotSet}; + std::array m_values = {kNotSet, kNotSet, kNotSet}; /// the accuracy value - this can vary up and down given a good step estimator - Scalar m_accuracy = kNotSet; + double m_accuracy = kNotSet; }; inline std::ostream& operator<<(std::ostream& os, const ConstrainedStep& step) { diff --git a/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingAnnoy.ipp b/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingAnnoy.ipp index d5e4f28f3a5..223fc3c4573 100644 --- a/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingAnnoy.ipp +++ b/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingAnnoy.ipp @@ -28,8 +28,6 @@ void HashingAnnoy:: const unsigned int phiBins, const double layerRMin, const double layerRMax, const double layerZMin, const double layerZMax) { - using Scalar = Acts::ActsScalar; - static thread_local std::vector> bucketsSetSPMap; bucketsSetSPMap.clear(); @@ -61,21 +59,21 @@ void HashingAnnoy:: }; // Functions to get the bin index - auto getBinIndexZ = [&zBins, &layerZMin, &layerZMax](Scalar z) { - Scalar binSize = (layerZMax - layerZMin) / zBins; + auto getBinIndexZ = [&zBins, &layerZMin, &layerZMax](double z) { + double binSize = (layerZMax - layerZMin) / zBins; auto binIndex = static_cast((z - layerZMin + 0.5 * binSize) / binSize); return binIndex; }; - auto getBinIndexPhi = [&phiBins](Scalar phi) { - Scalar binSize = 2 * std::numbers::pi / phiBins; + auto getBinIndexPhi = [&phiBins](double phi) { + double binSize = 2 * std::numbers::pi / phiBins; auto binIndex = static_cast((phi + std::numbers::pi) / binSize); return binIndex; }; // Function pointers to a unified bin index function for z and phi auto getBinIndex = [&zBins, &phiBins, &getBinIndexZ, &getBinIndexPhi]( - Scalar z, Scalar phi) -> int { + double z, double phi) -> int { if (zBins > 0) { return getBinIndexZ(z); } else if (phiBins > 0) { @@ -89,16 +87,16 @@ void HashingAnnoy:: for (unsigned int spacePointIndex = 0; spacePointIndex < spacePoints.size(); spacePointIndex++) { external_spacepoint_t spacePoint = spacePoints[spacePointIndex]; - Scalar x = spacePoint->x() / Acts::UnitConstants::mm; - Scalar y = spacePoint->y() / Acts::UnitConstants::mm; - Scalar z = spacePoint->z() / Acts::UnitConstants::mm; + double x = spacePoint->x() / Acts::UnitConstants::mm; + double y = spacePoint->y() / Acts::UnitConstants::mm; + double z = spacePoint->z() / Acts::UnitConstants::mm; // Helix transform - if (Scalar r2 = x * x + y * y; !layerSelection(r2, z)) { + if (double r2 = x * x + y * y; !layerSelection(r2, z)) { continue; } - Scalar phi = atan2(y, x); + double phi = atan2(y, x); int binIndex = getBinIndex(z, phi); if (binIndex < 0 || static_cast(binIndex) >= nBins) { diff --git a/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingTraining.ipp b/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingTraining.ipp index 3e24e828177..b3abefa2a11 100644 --- a/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingTraining.ipp +++ b/Plugins/Hashing/include/Acts/Plugins/Hashing/HashingTraining.ipp @@ -40,18 +40,16 @@ AnnoyModel HashingTrainingAlgorithm::execute( auto annoyModel = AnnoyModel(f); - using Scalar = Acts::ActsScalar; - annoyModel.set_seed(annoySeed); unsigned int spacePointIndex = 0; // Add spacePoints parameters to Annoy for (const auto& spacePoint : spacePoints) { - Scalar x = spacePoint->x() / Acts::UnitConstants::mm; - Scalar y = spacePoint->y() / Acts::UnitConstants::mm; + double x = spacePoint->x() / Acts::UnitConstants::mm; + double y = spacePoint->y() / Acts::UnitConstants::mm; // Helix transform - Scalar phi = std::atan2(y, x); + double phi = std::atan2(y, x); std::vector vec(f); // Avoid potential null pointer dereference @@ -59,11 +57,11 @@ AnnoyModel HashingTrainingAlgorithm::execute( vec[0] = phi; } if (f >= 2) { - Scalar z = spacePoint->z() / Acts::UnitConstants::mm; - Scalar r2 = x * x + y * y; - Scalar rho = std::sqrt(r2 + z * z); - Scalar theta = std::acos(z / rho); - Scalar eta = Acts::AngleHelpers::etaFromTheta(theta); + double z = spacePoint->z() / Acts::UnitConstants::mm; + double r2 = x * x + y * y; + double rho = std::sqrt(r2 + z * z); + double theta = std::acos(z / rho); + double eta = Acts::AngleHelpers::etaFromTheta(theta); vec[1] = eta; } diff --git a/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp b/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp index e5486f6a288..c31f33547e0 100644 --- a/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp +++ b/Tests/UnitTests/Fatras/Kernel/SimulationActorTests.cpp @@ -97,26 +97,21 @@ struct MockInteractionList { }; struct MockStepperState { - using Scalar = Acts::ActsScalar; - using Vector3 = Acts::ActsVector<3>; - - Vector3 pos = Vector3::Zero(); - Scalar time = 0; - Vector3 dir = Vector3::Zero(); - Scalar p = 0; + Acts::Vector3 pos = Acts::Vector3::Zero(); + double time = 0; + Acts::Vector3 dir = Acts::Vector3::Zero(); + double p = 0; }; struct MockStepper { using State = MockStepperState; - using Scalar = MockStepperState::Scalar; - using Vector3 = MockStepperState::Vector3; auto position(const State &state) const { return state.pos; } auto time(const State &state) const { return state.time; } auto direction(const State &state) const { return state.dir; } auto absoluteMomentum(const State &state) const { return state.p; } - void update(State &state, const Vector3 &pos, const Vector3 &dir, Scalar qop, - Scalar time) { + void update(State &state, const Acts::Vector3 &pos, const Acts::Vector3 &dir, + double qop, double time) { state.pos = pos; state.time = time; state.dir = dir;