From 4f0557f8c0a5bb6d94c54fb5bb995c0e870a508a Mon Sep 17 00:00:00 2001 From: Andreas Stefl Date: Sun, 10 Nov 2024 11:00:24 +0100 Subject: [PATCH 1/2] refactor: Relate measurements to truth in Examples --- .../Digitization/DigitizationAlgorithm.hpp | 28 +---- .../Digitization/DigitizationConfig.hpp | 16 +-- .../src/DigitizationAlgorithm.cpp | 78 ++++--------- .../ActsExamples/EventData/Cluster.hpp | 2 +- .../ActsExamples/EventData/Measurement.hpp | 109 ++++++++++++++++-- .../EventData/MeasurementConcept.hpp | 6 +- .../Framework/src/EventData/Measurement.cpp | 8 +- Examples/Python/src/Digitization.cpp | 4 - 8 files changed, 142 insertions(+), 109 deletions(-) diff --git a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp index fd8937bd300..2fa2392693b 100644 --- a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp +++ b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp @@ -8,16 +8,12 @@ #pragma once -#include "Acts/Definitions/Algebra.hpp" -#include "Acts/Definitions/TrackParametrization.hpp" -#include "Acts/Geometry/GeometryContext.hpp" #include "Acts/Geometry/GeometryHierarchyMap.hpp" #include "Acts/Utilities/Logger.hpp" #include "ActsExamples/Digitization/DigitizationConfig.hpp" #include "ActsExamples/Digitization/MeasurementCreation.hpp" #include "ActsExamples/Digitization/SmearingConfig.hpp" #include "ActsExamples/EventData/Cluster.hpp" -#include "ActsExamples/EventData/Index.hpp" #include "ActsExamples/EventData/Measurement.hpp" #include "ActsExamples/EventData/SimHit.hpp" #include "ActsExamples/Framework/DataHandle.hpp" @@ -29,17 +25,10 @@ #include "ActsFatras/Digitization/UncorrelatedHitSmearer.hpp" #include -#include #include -#include -#include #include #include -namespace ActsFatras { -class Barcode; -} // namespace ActsFatras - namespace ActsExamples { struct AlgorithmContext; @@ -98,17 +87,12 @@ class DigitizationAlgorithm final : public IAlgorithm { using CellsMap = std::map>; - ReadDataHandle m_simContainerReadHandle{this, - "SimHitContainer"}; - - WriteDataHandle m_measurementWriteHandle{ - this, "Measurements"}; - WriteDataHandle m_cellsWriteHandle{this, "Cells"}; - WriteDataHandle m_clusterWriteHandle{this, "Clusters"}; - WriteDataHandle> - m_measurementParticlesMapWriteHandle{this, "MeasurementParticlesMap"}; - WriteDataHandle> m_measurementSimHitsMapWriteHandle{ - this, "MeasurementSimHitsMap"}; + ReadDataHandle m_inputHits{this, "InputHits"}; + + WriteDataHandle m_outputMeasurements{ + this, "OutputMeasurements"}; + WriteDataHandle m_outputCells{this, "OutputCells"}; + WriteDataHandle m_outputClusters{this, "OutputClusters"}; /// Construct a fixed-size smearer from a configuration. /// diff --git a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp index 6a933f50075..53fbdc21d89 100644 --- a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp +++ b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp @@ -10,13 +10,9 @@ #include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/TrackParametrization.hpp" -#include "Acts/Definitions/Units.hpp" #include "Acts/Geometry/GeometryHierarchyMap.hpp" #include "Acts/Utilities/BinUtility.hpp" -#include "Acts/Utilities/BinningType.hpp" -#include "Acts/Utilities/Logger.hpp" #include "Acts/Utilities/Result.hpp" -#include "ActsExamples/Digitization/DigitizationConfig.hpp" #include "ActsExamples/Digitization/Smearers.hpp" #include "ActsExamples/Digitization/SmearingConfig.hpp" #include "ActsExamples/Framework/RandomNumbers.hpp" @@ -24,7 +20,6 @@ #include #include -#include #include #include #include @@ -141,10 +136,7 @@ class DigitizationConfig { std::string outputCells = "cells"; /// Output cluster collection. std::string outputClusters = "clusters"; - /// Output collection to map measured hits to contributing particles. - std::string outputMeasurementParticlesMap = "measurement_particles_map"; - /// Output collection to map measured hits to simulated hits. - std::string outputMeasurementSimHitsMap = "measurement_simhits_map"; + /// Map of surface by identifier to allow local - to global std::unordered_map surfaceByIdentifier; @@ -157,11 +149,11 @@ class DigitizationConfig { /// clusters, measurements, and sim-hit-maps are output. bool doClusterization = true; /// Do we merge hits or not - bool doMerge; + bool doMerge = false; /// How close do parameters have to be to consider merged - const double mergeNsigma; + double mergeNsigma = 3.0; /// Consider clusters that share a corner as merged (8-cell connectivity) - const bool mergeCommonCorner; + bool mergeCommonCorner = false; /// Energy deposit threshold for accepting a hit /// For a generic readout frontend we assume 1000 e/h pairs, in Si each /// e/h-pair requiers on average an energy of 3.65 eV (PDG review 2023, diff --git a/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp b/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp index 72f3905ab57..79f4a3b1a3f 100644 --- a/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp +++ b/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp @@ -11,27 +11,15 @@ #include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/Geometry/GeometryIdentifier.hpp" -#include "Acts/Geometry/TrackingGeometry.hpp" -#include "Acts/Surfaces/Surface.hpp" #include "Acts/Utilities/BinUtility.hpp" -#include "Acts/Utilities/Result.hpp" #include "ActsExamples/Digitization/ModuleClusters.hpp" -#include "ActsExamples/EventData/GeometryContainers.hpp" -#include "ActsExamples/EventData/Index.hpp" -#include "ActsExamples/EventData/SimHit.hpp" #include "ActsExamples/Framework/AlgorithmContext.hpp" -#include "ActsExamples/Utilities/GroupBy.hpp" -#include "ActsExamples/Utilities/Range.hpp" #include "ActsFatras/EventData/Barcode.hpp" -#include "ActsFatras/EventData/Hit.hpp" #include #include -#include -#include #include #include -#include #include #include #include @@ -60,21 +48,9 @@ ActsExamples::DigitizationAlgorithm::DigitizationAlgorithm( if (m_cfg.outputClusters.empty()) { throw std::invalid_argument("Missing cluster output collection"); } - if (m_cfg.outputMeasurementParticlesMap.empty()) { - throw std::invalid_argument( - "Missing hit-to-particles map output collection"); - } - if (m_cfg.outputMeasurementSimHitsMap.empty()) { - throw std::invalid_argument( - "Missing hit-to-simulated-hits map output collection"); - } - m_measurementWriteHandle.initialize(m_cfg.outputMeasurements); - m_clusterWriteHandle.initialize(m_cfg.outputClusters); - m_measurementParticlesMapWriteHandle.initialize( - m_cfg.outputMeasurementParticlesMap); - m_measurementSimHitsMapWriteHandle.initialize( - m_cfg.outputMeasurementSimHitsMap); + m_outputMeasurements.initialize(m_cfg.outputMeasurements); + m_outputClusters.initialize(m_cfg.outputClusters); } if (m_cfg.doOutputCells) { @@ -82,10 +58,10 @@ ActsExamples::DigitizationAlgorithm::DigitizationAlgorithm( throw std::invalid_argument("Missing cell output collection"); } - m_cellsWriteHandle.initialize(m_cfg.outputCells); + m_outputCells.initialize(m_cfg.outputCells); } - m_simContainerReadHandle.initialize(m_cfg.inputSimHits); + m_inputHits.initialize(m_cfg.inputSimHits); // Create the digitizers from the configuration std::vector> digitizerInput; @@ -142,18 +118,13 @@ ActsExamples::DigitizationAlgorithm::DigitizationAlgorithm( ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( const AlgorithmContext& ctx) const { // Retrieve input - const auto& simHits = m_simContainerReadHandle(ctx); + const auto& simHits = m_inputHits(ctx); ACTS_DEBUG("Loaded " << simHits.size() << " sim hits"); // Prepare output containers // need list here for stable addresses MeasurementContainer measurements; ClusterContainer clusters; - IndexMultimap measurementParticlesMap; - IndexMultimap measurementSimHitsMap; - measurements.reserve(simHits.size()); - measurementParticlesMap.reserve(simHits.size()); - measurementSimHitsMap.reserve(simHits.size()); // Setup random number generator auto rng = m_cfg.randomNumbers->spawnGenerator(ctx); @@ -271,7 +242,7 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( // Store the cell data into a map. if (m_cfg.doOutputCells) { std::vector cells; - for (const auto& [dParameters, simhits] : + for (const auto& [dParameters, simHitsIdxs] : digitizeParametersResult) { for (const auto& cell : dParameters.cluster.channels) { cells.push_back(cell); @@ -281,21 +252,17 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( } if (m_cfg.doClusterization) { - for (auto& [dParameters, simhits] : digitizeParametersResult) { - // The measurement container is unordered and the index under - // which the measurement will be stored is known before adding it. - Index measurementIdx = measurements.size(); - - createMeasurement(measurements, moduleGeoId, dParameters); + for (auto& [dParameters, simHitsIdxs] : digitizeParametersResult) { + auto measurement = + createMeasurement(measurements, moduleGeoId, dParameters); clusters.emplace_back(std::move(dParameters.cluster)); - // this digitization does hit merging so there can be more than - // one mapping entry for each digitized hit. - for (auto simHitIdx : simhits) { - measurementParticlesMap.emplace_hint( - measurementParticlesMap.end(), measurementIdx, - simHits.nth(simHitIdx)->particleId()); - measurementSimHitsMap.emplace_hint(measurementSimHitsMap.end(), - measurementIdx, simHitIdx); + + measurement.allocateHitIndices(simHitsIdxs.size()); + measurement.allocateParticleIds(simHitsIdxs.size()); + for (auto [i, simHitIdx] : Acts::enumerate(simHitsIdxs)) { + measurement.hitIndices()[i] = simHitIdx; + measurement.particleIds()[i] = + simHits.nth(simHitIdx)->particleId(); } } } @@ -309,16 +276,13 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( << " skipped in Digitization. Enable DEBUG mode to see more details."); } - if (m_cfg.doOutputCells) { - m_cellsWriteHandle(ctx, std::move(cellsMap)); + if (m_cfg.doClusterization) { + m_outputMeasurements(ctx, std::move(measurements)); + m_outputClusters(ctx, std::move(clusters)); } - if (m_cfg.doClusterization) { - m_measurementWriteHandle(ctx, std::move(measurements)); - m_clusterWriteHandle(ctx, std::move(clusters)); - m_measurementParticlesMapWriteHandle(ctx, - std::move(measurementParticlesMap)); - m_measurementSimHitsMapWriteHandle(ctx, std::move(measurementSimHitsMap)); + if (m_cfg.doOutputCells) { + m_outputCells(ctx, std::move(cellsMap)); } return ProcessCode::SUCCESS; diff --git a/Examples/Framework/include/ActsExamples/EventData/Cluster.hpp b/Examples/Framework/include/ActsExamples/EventData/Cluster.hpp index 3eb850b1214..bfd81c9f130 100644 --- a/Examples/Framework/include/ActsExamples/EventData/Cluster.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/Cluster.hpp @@ -8,10 +8,10 @@ #pragma once +#include "Acts/Definitions/Algebra.hpp" #include "ActsFatras/Digitization/Segmentizer.hpp" #include -#include #include namespace ActsExamples { diff --git a/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp b/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp index a0bf4559509..820507564c4 100644 --- a/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp @@ -10,26 +10,17 @@ #include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/TrackParametrization.hpp" -#include "Acts/EventData/MeasurementHelpers.hpp" #include "Acts/EventData/SubspaceHelpers.hpp" #include "Acts/EventData/Types.hpp" -#include "Acts/EventData/detail/CalculateResiduals.hpp" -#include "Acts/EventData/detail/ParameterTraits.hpp" -#include "Acts/EventData/detail/PrintParameters.hpp" #include "Acts/Geometry/GeometryIdentifier.hpp" #include "Acts/Utilities/Iterator.hpp" #include "ActsExamples/EventData/GeometryContainers.hpp" #include "ActsExamples/EventData/IndexSourceLink.hpp" #include "ActsExamples/EventData/MeasurementConcept.hpp" +#include "ActsExamples/EventData/SimParticle.hpp" -#include -#include -#include #include -#include -#include #include -#include #include #include @@ -168,6 +159,8 @@ class MeasurementContainer { std::size_t subspaceIndexOffset{}; std::size_t parameterOffset{}; std::size_t covarianceOffset{}; + std::pair hitIndexRange{}; + std::pair particleIdRange{}; std::uint8_t size{}; }; @@ -177,6 +170,8 @@ class MeasurementContainer { std::vector m_subspaceIndices; std::vector m_parameters; std::vector m_covariances; + std::vector m_hitIndices; + std::vector m_particleIds; OrderedIndices m_orderedIndices; }; @@ -268,6 +263,90 @@ class MeasurementProxyBase { return self().subspaceHelper().expandMatrix(self().covariance()); } + bool hasHitIndices() const { + const auto& entry = container().m_entries.at(m_index); + return entry.hitIndexRange.first != entry.hitIndexRange.second; + } + + std::span allocateHitIndices(std::size_t size) + requires(!ReadOnly) + { + if (hasHitIndices()) { + throw std::runtime_error("Hit indices already allocated"); + } + + auto& entry = container().m_entries.at(m_index); + entry.hitIndexRange.first = container().m_hitIndices.size(); + entry.hitIndexRange.second = size; + container().m_hitIndices.resize(entry.hitIndexRange.first + size); + + return hitIndices(); + } + + std::span hitIndices() + requires(!ReadOnly) + { + if (!hasHitIndices()) { + throw std::runtime_error("Hit indices not allocated"); + } + + const auto& entry = container().m_entries.at(m_index); + return {container().m_hitIndices.data() + entry.hitIndexRange.first, + entry.hitIndexRange.second}; + } + + std::span hitIndices() const { + if (!hasHitIndices()) { + throw std::runtime_error("Hit indices not allocated"); + } + + const auto& entry = container().m_entries.at(m_index); + return {container().m_hitIndices.data() + entry.hitIndexRange.first, + entry.hitIndexRange.second}; + } + + bool hasParticleIds() const { + const auto& entry = container().m_entries.at(m_index); + return entry.particleIdRange.first != entry.particleIdRange.second; + } + + std::span allocateParticleIds(std::size_t size) + requires(!ReadOnly) + { + if (hasParticleIds()) { + throw std::runtime_error("Particle IDs already allocated"); + } + + auto& entry = container().m_entries.at(m_index); + entry.particleIdRange.first = container().m_particleIds.size(); + entry.particleIdRange.second = size; + container().m_particleIds.resize(entry.particleIdRange.first + size); + + return particleIds(); + } + + std::span particleIds() + requires(!ReadOnly) + { + if (!hasParticleIds()) { + throw std::runtime_error("Particle IDs not allocated"); + } + + const auto& entry = container().m_entries.at(m_index); + return {container().m_particleIds.data() + entry.particleIdRange.first, + entry.particleIdRange.second}; + } + + std::span particleIds() const { + if (!hasParticleIds()) { + throw std::runtime_error("Particle IDs not allocated"); + } + + const auto& entry = container().m_entries.at(m_index); + return {container().m_particleIds.data() + entry.particleIdRange.first, + entry.particleIdRange.second}; + } + /// @brief Construct the measurement from a subspace vector, /// parameters, and covariance. /// @@ -292,6 +371,16 @@ class MeasurementProxyBase { { assert(size() == other.size() && "Size mismatch"); fill(other.subspaceIndexVector(), other.parameters(), other.covariance()); + if (other.hasHitIndices()) { + allocateHitIndices(other.hitIndices().size()); + std::copy(other.hitIndices().begin(), other.hitIndices().end(), + hitIndices().begin()); + } + if (other.hasParticleIds()) { + allocateParticleIds(other.particleIds().size()); + std::copy(other.particleIds().begin(), other.particleIds().end(), + particleIds().begin()); + } } /// @brief Copy the data from another measurement diff --git a/Examples/Framework/include/ActsExamples/EventData/MeasurementConcept.hpp b/Examples/Framework/include/ActsExamples/EventData/MeasurementConcept.hpp index 14627a78c53..e959b78bb90 100644 --- a/Examples/Framework/include/ActsExamples/EventData/MeasurementConcept.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/MeasurementConcept.hpp @@ -8,7 +8,6 @@ #pragma once -#include "Acts/EventData/SourceLink.hpp" #include "Acts/Geometry/GeometryIdentifier.hpp" #include @@ -22,5 +21,10 @@ concept MeasurementConcept = requires(const T& m) { { m.subspaceIndexVector() }; { m.parameters() }; { m.covariance() }; + { m.hasHitIndices() } -> std::same_as; + { m.hitIndices() }; + { m.hasParticleIds() } -> std::same_as; + { m.particleIds() }; }; + } // namespace ActsExamples diff --git a/Examples/Framework/src/EventData/Measurement.cpp b/Examples/Framework/src/EventData/Measurement.cpp index ab2dad07489..e01b417e43d 100644 --- a/Examples/Framework/src/EventData/Measurement.cpp +++ b/Examples/Framework/src/EventData/Measurement.cpp @@ -28,8 +28,12 @@ void MeasurementContainer::reserve(std::size_t size) { std::size_t MeasurementContainer::addMeasurement( std::uint8_t size, Acts::GeometryIdentifier geometryId) { - m_entries.push_back({m_subspaceIndices.size(), m_parameters.size(), - m_covariances.size(), size}); + m_entries.push_back({m_subspaceIndices.size(), + m_parameters.size(), + m_covariances.size(), + {}, + {}, + size}); m_geometryIds.emplace_back(geometryId); m_subspaceIndices.resize(m_subspaceIndices.size() + size); m_parameters.resize(m_parameters.size() + size); diff --git a/Examples/Python/src/Digitization.cpp b/Examples/Python/src/Digitization.cpp index 948bcc9593b..38dd8573a8a 100644 --- a/Examples/Python/src/Digitization.cpp +++ b/Examples/Python/src/Digitization.cpp @@ -6,7 +6,6 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at https://mozilla.org/MPL/2.0/. -#include "Acts/Definitions/Algebra.hpp" #include "Acts/Geometry/GeometryHierarchyMap.hpp" #include "Acts/Plugins/Python/Utilities.hpp" #include "Acts/Utilities/Logger.hpp" @@ -14,7 +13,6 @@ #include "ActsExamples/Digitization/DigitizationConfig.hpp" #include "ActsExamples/Digitization/DigitizationConfigurator.hpp" #include "ActsExamples/Digitization/DigitizationCoordinatesConverter.hpp" -#include "ActsExamples/Framework/AlgorithmContext.hpp" #include "ActsExamples/Io/Json/JsonDigitizationConfig.hpp" #include @@ -65,8 +63,6 @@ void addDigitization(Context& ctx) { ACTS_PYTHON_MEMBER(inputSimHits); ACTS_PYTHON_MEMBER(outputMeasurements); ACTS_PYTHON_MEMBER(outputClusters); - ACTS_PYTHON_MEMBER(outputMeasurementParticlesMap); - ACTS_PYTHON_MEMBER(outputMeasurementSimHitsMap); ACTS_PYTHON_MEMBER(surfaceByIdentifier); ACTS_PYTHON_MEMBER(randomNumbers); ACTS_PYTHON_MEMBER(doOutputCells); From f48e9b7911a59bed7118a67d257f3d5193556b89 Mon Sep 17 00:00:00 2001 From: Andreas Stefl Date: Sun, 10 Nov 2024 11:32:27 +0100 Subject: [PATCH 2/2] revert to original wiring --- .../Digitization/DigitizationAlgorithm.hpp | 6 +++++ .../Digitization/DigitizationConfig.hpp | 4 +++ .../src/DigitizationAlgorithm.cpp | 27 +++++++++++++++++++ 3 files changed, 37 insertions(+) diff --git a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp index 2fa2392693b..2144239a62c 100644 --- a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp +++ b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationAlgorithm.hpp @@ -16,6 +16,7 @@ #include "ActsExamples/EventData/Cluster.hpp" #include "ActsExamples/EventData/Measurement.hpp" #include "ActsExamples/EventData/SimHit.hpp" +#include "ActsExamples/EventData/SimParticle.hpp" #include "ActsExamples/Framework/DataHandle.hpp" #include "ActsExamples/Framework/IAlgorithm.hpp" #include "ActsExamples/Framework/ProcessCode.hpp" @@ -94,6 +95,11 @@ class DigitizationAlgorithm final : public IAlgorithm { WriteDataHandle m_outputCells{this, "OutputCells"}; WriteDataHandle m_outputClusters{this, "OutputClusters"}; + WriteDataHandle> m_outputMeasurementParticlesMap{ + this, "OutputMeasurementParticlesMap"}; + WriteDataHandle> m_outputMeasurementSimHitsMap{ + this, "OutputMeasurementSimHitsMap"}; + /// Construct a fixed-size smearer from a configuration. /// /// It's templated on the smearing dimension given by @tparam kSmearDIM diff --git a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp index 53fbdc21d89..4772aa40405 100644 --- a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp +++ b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/DigitizationConfig.hpp @@ -136,6 +136,10 @@ class DigitizationConfig { std::string outputCells = "cells"; /// Output cluster collection. std::string outputClusters = "clusters"; + /// Output collection to map measured hits to contributing particles. + std::string outputMeasurementParticlesMap = "measurement_particles_map"; + /// Output collection to map measured hits to simulated hits. + std::string outputMeasurementSimHitsMap = "measurement_simhits_map"; /// Map of surface by identifier to allow local - to global std::unordered_map diff --git a/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp b/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp index 79f4a3b1a3f..b0c44a88e4a 100644 --- a/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp +++ b/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp @@ -48,9 +48,20 @@ ActsExamples::DigitizationAlgorithm::DigitizationAlgorithm( if (m_cfg.outputClusters.empty()) { throw std::invalid_argument("Missing cluster output collection"); } + if (m_cfg.outputMeasurementParticlesMap.empty()) { + throw std::invalid_argument( + "Missing hit-to-particles map output collection"); + } + if (m_cfg.outputMeasurementSimHitsMap.empty()) { + throw std::invalid_argument( + "Missing hit-to-simulated-hits map output collection"); + } m_outputMeasurements.initialize(m_cfg.outputMeasurements); m_outputClusters.initialize(m_cfg.outputClusters); + m_outputMeasurementParticlesMap.initialize( + m_cfg.outputMeasurementParticlesMap); + m_outputMeasurementSimHitsMap.initialize(m_cfg.outputMeasurementSimHitsMap); } if (m_cfg.doOutputCells) { @@ -126,6 +137,12 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( MeasurementContainer measurements; ClusterContainer clusters; + IndexMultimap measurementParticlesMap; + IndexMultimap measurementSimHitsMap; + measurements.reserve(simHits.size()); + measurementParticlesMap.reserve(simHits.size()); + measurementSimHitsMap.reserve(simHits.size()); + // Setup random number generator auto rng = m_cfg.randomNumbers->spawnGenerator(ctx); @@ -263,6 +280,13 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( measurement.hitIndices()[i] = simHitIdx; measurement.particleIds()[i] = simHits.nth(simHitIdx)->particleId(); + + measurementParticlesMap.emplace_hint( + measurementParticlesMap.end(), measurement.index(), + simHits.nth(simHitIdx)->particleId()); + measurementSimHitsMap.emplace_hint(measurementSimHitsMap.end(), + measurement.index(), + simHitIdx); } } } @@ -279,6 +303,9 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( if (m_cfg.doClusterization) { m_outputMeasurements(ctx, std::move(measurements)); m_outputClusters(ctx, std::move(clusters)); + + m_outputMeasurementParticlesMap(ctx, std::move(measurementParticlesMap)); + m_outputMeasurementSimHitsMap(ctx, std::move(measurementSimHitsMap)); } if (m_cfg.doOutputCells) {