diff --git a/Core/include/Acts/EventData/SubspaceHelpers.hpp b/Core/include/Acts/EventData/SubspaceHelpers.hpp index 737c12211a5..35bde575830 100644 --- a/Core/include/Acts/EventData/SubspaceHelpers.hpp +++ b/Core/include/Acts/EventData/SubspaceHelpers.hpp @@ -28,7 +28,7 @@ inline static bool checkSubspaceIndices(const Container& container, if (subspaceSize > fullSize) { return false; } - if (container.size() != subspaceSize) { + if (static_cast(container.size()) != subspaceSize) { return false; } for (auto it = container.begin(); it != container.end();) { @@ -49,6 +49,7 @@ class SubspaceHelperBase { public: static constexpr std::size_t kFullSize = FullSize; + using FullVector = ActsVector; using FullSquareMatrix = ActsSquareMatrix; bool empty() const { return self().empty(); } @@ -59,6 +60,19 @@ class SubspaceHelperBase { auto begin() const { return self().begin(); } auto end() const { return self().end(); } + template + bool contains(indices_t i) const { + return std::find(begin(), end(), static_cast(i)) != + end(); + } + + template + std::size_t indexOf(indices_t i) const { + auto it = std::find(begin(), end(), static_cast(i)); + assert(it != end()); + return std::distance(begin(), it); + } + FullSquareMatrix fullProjector() const { FullSquareMatrix result = FullSquareMatrix::Zero(); for (auto [i, index] : enumerate(*this)) { @@ -81,6 +95,33 @@ class SubspaceHelperBase { return matrixToBitset(fullProjector()).to_ullong(); } + template + FullVector expandVector( + const Eigen::DenseBase& subspaceVector) const { + assert(static_cast(subspaceVector.size()) == size() && + "Invalid subspace vector size"); + FullVector result = FullVector::Zero(); + for (auto [i, index] : enumerate(*this)) { + result(index) = subspaceVector(i); + } + return result; + } + + template + FullSquareMatrix expandMatrix( + const Eigen::DenseBase& subspaceMatrix) const { + assert(static_cast(subspaceMatrix.rows()) == size() && + static_cast(subspaceMatrix.cols()) == size() && + "Invalid subspace matrix size"); + FullSquareMatrix result = FullSquareMatrix::Zero(); + for (auto [i, indexI] : enumerate(*this)) { + for (auto [j, indexJ] : enumerate(*this)) { + result(indexI, indexJ) = subspaceMatrix(i, j); + } + } + return result; + } + private: const Derived& self() const { return static_cast(*this); } }; @@ -112,6 +153,8 @@ class VariableSubspaceHelper auto begin() const { return m_indices.begin(); } auto end() const { return m_indices.end(); } + const Container& indices() const { return m_indices; } + private: Container m_indices; }; @@ -154,6 +197,8 @@ class FixedSubspaceHelper auto begin() const { return m_indices.begin(); } auto end() const { return m_indices.end(); } + const Container& indices() const { return m_indices; } + Projector projector() const { Projector result = Projector::Zero(); for (auto [i, index] : enumerate(*this)) { diff --git a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/MeasurementCreation.hpp b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/MeasurementCreation.hpp index 3d9476815c9..e168858cfd1 100644 --- a/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/MeasurementCreation.hpp +++ b/Examples/Algorithms/Digitization/include/ActsExamples/Digitization/MeasurementCreation.hpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2021 CERN for the benefit of the Acts project +// Copyright (C) 2021-2024 CERN for the benefit of the Acts project // // This Source Code Form is subject to the terms of the Mozilla Public // License, v. 2.0. If a copy of the MPL was not distributed with this @@ -40,9 +40,10 @@ struct DigitizedParameters { /// /// To be used also by the e I/O system /// -/// @return a variant measurement -Measurement createMeasurement(const DigitizedParameters& dParams, - const IndexSourceLink& isl) noexcept(false); +/// @return the measurement proxy +ActsExamples::VariableBoundMeasurementProxy createMeasurement( + MeasurementContainer& container, const DigitizedParameters& dParams, + const IndexSourceLink& isl) noexcept(false); /// Construct the constituents of a measurement. /// diff --git a/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp b/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp index f3b04ac0aff..3ac84f45bac 100644 --- a/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp +++ b/Examples/Algorithms/Digitization/src/DigitizationAlgorithm.cpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2021 CERN for the benefit of the Acts project +// Copyright (C) 2021-2024 CERN for the benefit of the Acts project // // This Source Code Form is subject to the terms of the Mozilla Public // License, v. 2.0. If a copy of the MPL was not distributed with this @@ -268,8 +268,7 @@ ActsExamples::ProcessCode ActsExamples::DigitizationAlgorithm::execute( // be added at the end. sourceLinks.insert(sourceLinks.end(), sourceLink); - measurements.emplace_back( - createMeasurement(dParameters, sourceLink)); + createMeasurement(measurements, dParameters, sourceLink); 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. diff --git a/Examples/Algorithms/Digitization/src/MeasurementCreation.cpp b/Examples/Algorithms/Digitization/src/MeasurementCreation.cpp index 7338b431f8c..88ba53e1009 100644 --- a/Examples/Algorithms/Digitization/src/MeasurementCreation.cpp +++ b/Examples/Algorithms/Digitization/src/MeasurementCreation.cpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2021 CERN for the benefit of the Acts project +// Copyright (C) 2021-2024 CERN for the benefit of the Acts project // // This Source Code Form is subject to the terms of the Mozilla Public // License, v. 2.0. If a copy of the MPL was not distributed with this @@ -8,36 +8,34 @@ #include "ActsExamples/Digitization/MeasurementCreation.hpp" +#include "Acts/EventData/MeasurementHelpers.hpp" #include "Acts/EventData/SourceLink.hpp" #include "ActsExamples/EventData/IndexSourceLink.hpp" +#include "ActsExamples/EventData/Measurement.hpp" #include #include #include -ActsExamples::Measurement ActsExamples::createMeasurement( - const DigitizedParameters& dParams, const IndexSourceLink& isl) { +ActsExamples::VariableBoundMeasurementProxy ActsExamples::createMeasurement( + MeasurementContainer& container, const DigitizedParameters& dParams, + const IndexSourceLink& isl) { Acts::SourceLink sl{isl}; - switch (dParams.indices.size()) { - case 1u: { - auto [indices, par, cov] = measurementConstituents<1>(dParams); - return ActsExamples::Measurement(std::move(sl), indices, par, cov); - } - case 2u: { - auto [indices, par, cov] = measurementConstituents<2>(dParams); - return ActsExamples::Measurement(std::move(sl), indices, par, cov); - }; - case 3u: { - auto [indices, par, cov] = measurementConstituents<3>(dParams); - return ActsExamples::Measurement(std::move(sl), indices, par, cov); - }; - case 4u: { - auto [indices, par, cov] = measurementConstituents<4>(dParams); - return ActsExamples::Measurement(std::move(sl), indices, par, cov); - }; - default: - std::string errorMsg = "Invalid/mismatching measurement dimension: " + - std::to_string(dParams.indices.size()); - throw std::runtime_error(errorMsg.c_str()); + + if (dParams.indices.size() > 4u) { + std::string errorMsg = "Invalid/mismatching measurement dimension: " + + std::to_string(dParams.indices.size()); + throw std::runtime_error(errorMsg.c_str()); } + + return Acts::visit_measurement( + dParams.indices.size(), [&](auto dim) -> VariableBoundMeasurementProxy { + auto [indices, par, cov] = measurementConstituents(dParams); + FixedBoundMeasurementProxy measurement = + container.makeMeasurement(); + measurement.setSubspaceIndices(indices); + measurement.parameters() = par; + measurement.covariance() = cov; + return measurement; + }); } diff --git a/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp b/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp index 12764a360d5..63a5d502211 100644 --- a/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp +++ b/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp @@ -525,7 +525,8 @@ void ActsExamples::HoughTransformSeeder::addMeasurements( // are transformed to the bound space where we do know their location. // if the local parameters are not measured, this results in a // zero location, which is a reasonable default fall-back. - const auto& measurement = measurements[sourceLink.index()]; + const ConstVariableBoundMeasurementProxy measurement = + measurements.getMeasurement(sourceLink.index()); assert(measurement.contains(Acts::eBoundLoc0) && "Measurement does not contain the required bound loc0"); diff --git a/Examples/Algorithms/TrackFinding/src/SpacePointMaker.cpp b/Examples/Algorithms/TrackFinding/src/SpacePointMaker.cpp index 49802b57087..0605e1b4bc4 100644 --- a/Examples/Algorithms/TrackFinding/src/SpacePointMaker.cpp +++ b/Examples/Algorithms/TrackFinding/src/SpacePointMaker.cpp @@ -126,7 +126,8 @@ ActsExamples::ProcessCode ActsExamples::SpacePointMaker::execute( spOpt.paramCovAccessor = [&measurements](Acts::SourceLink slink) { const auto islink = slink.get(); - const auto& meas = measurements[islink.index()]; + const ConstVariableBoundMeasurementProxy meas = + measurements.getMeasurement(islink.index()); return std::make_pair(meas.fullParameters(), meas.fullCovariance()); }; diff --git a/Examples/Framework/CMakeLists.txt b/Examples/Framework/CMakeLists.txt index bfe33f2e8e6..692db8f2f15 100644 --- a/Examples/Framework/CMakeLists.txt +++ b/Examples/Framework/CMakeLists.txt @@ -5,6 +5,7 @@ set(ActsExamplesFramework_SOURCES) add_library( ActsExamplesFramework SHARED + src/EventData/Measurement.cpp src/EventData/MeasurementCalibration.cpp src/EventData/ScalingCalibrator.cpp src/Framework/IAlgorithm.cpp diff --git a/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp b/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp index 87c4c546c08..71e9e1c12e5 100644 --- a/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp @@ -8,6 +8,7 @@ #pragma once +#include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/EventData/MeasurementHelpers.hpp" #include "Acts/EventData/SourceLink.hpp" @@ -28,253 +29,353 @@ namespace ActsExamples { -/// A measurement of a variable-size subspace of the full parameters. -/// -/// @tparam indices_t Parameter index type, determines the full parameter space -/// -/// The measurement intentionally does not store a pointer/reference to the -/// reference object in the geometry hierarchy, i.e. the surface or volume. The -/// reference object can already be identified via the geometry identifier -/// provided by the source link. Since a measurement **must** be anchored within -/// the geometry hierarchy, all measurement surfaces and volumes **must** -/// provide valid geometry identifiers. In all use-cases, e.g. Kalman filtering, -/// a pointer/reference to the reference object is available before the -/// measurement is accessed; e.g. the propagator provides the surface pointer -/// during navigation, which is then used to lookup possible measurements. -/// -/// The pointed-to geometry object would differ depending on the parameter type. -/// This means either, that there needs to be an additional variable type or -/// that a pointer to a base object is stored (requiring a `dynamic_cast` later -/// on). Both variants add additional complications. Since the geometry object -/// is not required anyway (as discussed above), not storing it removes all -/// these complications altogether. -template -class VariableSizeMeasurement { +template +class MeasurementProxyBase; +template +class FixedMeasurementProxy; +template +class VariableMeasurementProxy; + +template +using FixedBoundMeasurementProxy = + FixedMeasurementProxy; +template +using ConstFixedBoundMeasurementProxy = + FixedMeasurementProxy; +using VariableBoundMeasurementProxy = + VariableMeasurementProxy; +using ConstVariableBoundMeasurementProxy = + VariableMeasurementProxy; + +class MeasurementContainer { public: - static constexpr std::size_t kFullSize = - Acts::detail::kParametersSize; - - using Scalar = Acts::ActsScalar; - - using SubspaceIndex = std::uint8_t; - using SubspaceIndices = - boost::container::static_vector; - - /// Vector type containing for measured parameter values. - template - using ParametersVector = Eigen::Matrix; - template - using ParametersVectorMap = Eigen::Map>; - template - using ConstParametersVectorMap = Eigen::Map>; - using EffectiveParametersVector = Eigen::Matrix; - using EffectiveParametersVectorMap = Eigen::Map; - using ConstEffectiveParametersVectorMap = - Eigen::Map; - - /// Matrix type for the measurement covariance. - template - using CovarianceMatrix = Eigen::Matrix; - template - using CovarianceMatrixMap = Eigen::Map>; - template - using ConstCovarianceMatrixMap = Eigen::Map>; - using EffectiveCovarianceMatrix = - Eigen::Matrix; - using EffectiveCovarianceMatrixMap = Eigen::Map; - using ConstEffectiveCovarianceMatrixMap = - Eigen::Map; - - using FullParametersVector = Acts::ActsVector; - using FullCovarianceMatrix = Acts::ActsSquareMatrix; - - using ProjectionMatrix = Eigen::Matrix; - using ExpansionMatrix = Eigen::Matrix; - - /// Construct from source link, subset indices, and measured data. - /// - /// @tparam parameters_t Input parameters vector type - /// @tparam covariance_t Input covariance matrix type - /// @param source The link that connects to the underlying detector readout - /// @param subspaceIndices Which parameters are measured - /// @param params Measured parameters values - /// @param cov Measured parameters covariance - /// - /// @note The indices must be ordered and must describe/match the content - /// of parameters and covariance. - template - VariableSizeMeasurement( - Acts::SourceLink source, - const std::array& subspaceIndices, - const Eigen::MatrixBase& params, - const Eigen::MatrixBase& cov) - : m_source(std::move(source)) { - static_assert(kSize == parameters_t::RowsAtCompileTime, - "Parameter size mismatch"); - static_assert(kSize == covariance_t::RowsAtCompileTime, - "Covariance rows mismatch"); - static_assert(kSize == covariance_t::ColsAtCompileTime, - "Covariance cols mismatch"); - - m_subspaceIndices.resize(subspaceIndices.size()); - std::transform(subspaceIndices.begin(), subspaceIndices.end(), - m_subspaceIndices.begin(), [](auto index) { - return static_cast(index); - }); - - parameters() = params; - covariance() = cov; + template + using FixedProxy = FixedMeasurementProxy; + template + using ConstFixedProxy = FixedMeasurementProxy; + using VariableProxy = VariableMeasurementProxy; + using ConstVariableProxy = VariableMeasurementProxy; + + MeasurementContainer(); + + std::size_t size() const; + + void reserve(std::size_t size); + + std::size_t addMeasurement(std::uint8_t size); + + VariableProxy getMeasurement(std::size_t index); + ConstVariableProxy getMeasurement(std::size_t index) const; + + template + FixedProxy getMeasurement(std::size_t index) { + return FixedProxy{*this, index}; + } + template + ConstFixedProxy getMeasurement(std::size_t index) const { + return ConstFixedProxy{*this, index}; } - /// A measurement can only be constructed with valid parameters. - VariableSizeMeasurement() = delete; - VariableSizeMeasurement(const VariableSizeMeasurement&) = default; - VariableSizeMeasurement(VariableSizeMeasurement&&) = default; - ~VariableSizeMeasurement() = default; - VariableSizeMeasurement& operator=(const VariableSizeMeasurement&) = default; - VariableSizeMeasurement& operator=(VariableSizeMeasurement&&) = default; - /// Source link that connects to the underlying detector readout. - const Acts::SourceLink& sourceLink() const { return m_source; } + VariableProxy makeMeasurement(std::uint8_t size); + template + FixedProxy makeMeasurement() { + return getMeasurement(addMeasurement(Size)); + } + + template + class IteratorImpl { + public: + using value_type = + std::conditional_t; + using reference = value_type; + using pointer = value_type*; + using difference_type = std::ptrdiff_t; + using iterator_category = std::forward_iterator_tag; + + using Container = std::conditional_t; + + IteratorImpl(Container& container, std::size_t index) + : m_container(container), m_index(index) {} + + reference operator*() const { + if constexpr (Const) { + return m_container.getMeasurement(m_index); + } else { + return m_container.getMeasurement(m_index); + } + } + + pointer operator->() const { return &operator*(); } + + IteratorImpl& operator++() { + ++m_index; + return *this; + } + + IteratorImpl operator++(int) { + auto copy = *this; + ++*this; + return copy; + } + + bool operator==(const IteratorImpl& other) const { + return m_index == other.m_index; + } + + bool operator!=(const IteratorImpl& other) const { + return !(*this == other); + } + + private: + Container& m_container; + std::size_t m_index; + }; + + using Iterator = IteratorImpl; + using ConstIterator = IteratorImpl; + + Iterator begin(); + Iterator end(); + ConstIterator begin() const; + ConstIterator end() const; + + public: + struct MeasurementEntry { + std::size_t subspaceIndexOffset{}; + std::size_t parameterOffset{}; + std::size_t covarianceOffset{}; + std::uint8_t size{}; + }; + + std::vector m_entries; + + std::vector m_sourceLinks; + std::vector m_subspaceIndices; + std::vector m_parameters; + std::vector m_covariances; +}; + +template +class MeasurementProxyBase { + public: + using Index = std::uint8_t; + using Scalar = double; + + using FullVector = Acts::ActsVector; + using FullSquareMatrix = Acts::ActsSquareMatrix; - constexpr std::size_t size() const { return m_subspaceIndices.size(); } + using Container = std::conditional_t; - /// Check if a specific parameter is part of this measurement. + MeasurementProxyBase(Container& container_, std::size_t index_) + : m_container(&container_), m_index(index_) {} + template + MeasurementProxyBase( + const MeasurementProxyBase& other) + requires(ReadOnly == OtherReadOnly || ReadOnly) + : m_container(&other.container()), m_index(other.index()) {} + + Container& container() const { return *m_container; } + std::size_t index() const { return m_index; } + + std::size_t size() const { return container().m_entries.at(m_index).size; } + + template bool contains(indices_t i) const { - return std::find(m_subspaceIndices.begin(), m_subspaceIndices.end(), i) != - m_subspaceIndices.end(); + return self().subspaceHelper().contains(i); } + template std::size_t indexOf(indices_t i) const { - auto it = std::find(m_subspaceIndices.begin(), m_subspaceIndices.end(), i); - assert(it != m_subspaceIndices.end()); - return std::distance(m_subspaceIndices.begin(), it); + return self().subspaceHelper().indexOf(i); } - /// The measurement indices - const SubspaceIndices& subspaceIndices() const { return m_subspaceIndices; } + Acts::SourceLink& sourceLink() { + return container().m_sourceLinks.at(m_index); + } - template - Acts::SubspaceIndices subspaceIndices() const { - assert(dim == size()); - Acts::SubspaceIndices result; - std::copy(m_subspaceIndices.begin(), m_subspaceIndices.end(), - result.begin()); - return result; + const Acts::SourceLink& sourceLink() const { + return container().m_sourceLinks.at(m_index); } - Acts::BoundSubspaceIndices boundSubsetIndices() const - requires(std::is_same_v) + template + void setSubspaceIndices(const IndexContainer& indices) + requires(!ReadOnly) { - Acts::BoundSubspaceIndices result = Acts::kBoundSubspaceIndicesInvalid; - std::copy(m_subspaceIndices.begin(), m_subspaceIndices.end(), - result.begin()); - return result; + assert(checkSubspaceIndices(indices, FullSize, size()) && + "Invalid indices"); + std::transform(indices.begin(), indices.end(), + self().subspaceIndexVector().begin(), + [](auto index) { return static_cast(index); }); } - template - ConstParametersVectorMap parameters() const { - assert(dim == size()); - return ConstParametersVectorMap{m_params.data()}; + FullVector fullParameters() const { + return self().subspaceHelper().expandVector(self().parameters()); } - template - ParametersVectorMap parameters() { - assert(dim == size()); - return ParametersVectorMap{m_params.data()}; - } - ConstEffectiveParametersVectorMap parameters() const { - return ConstEffectiveParametersVectorMap{m_params.data(), - static_cast(size())}; + + FullSquareMatrix fullCovariance() const { + return self().subspaceHelper().expandMatrix(self().covariance()); } - EffectiveParametersVectorMap parameters() { - return EffectiveParametersVectorMap{m_params.data(), - static_cast(size())}; + + template + void copyFrom(const OtherDerived& other) + requires(!ReadOnly) + { + assert(size() == other.size() && "Size mismatch"); + sourceLink() = other.sourceLink(); + self().subspaceIndexVector() = other.subspaceIndexVector(); + self().parameters() = other.parameters(); + self().covariance() = other.covariance(); } - template - ConstCovarianceMatrixMap covariance() const { - assert(dim == size()); - return ConstCovarianceMatrixMap{m_cov.data()}; + protected: + Derived& self() { return static_cast(*this); } + const Derived& self() const { return static_cast(*this); } + + Container* m_container; + std::size_t m_index; +}; + +template +class FixedMeasurementProxy + : public MeasurementProxyBase< + FixedMeasurementProxy, FullSize, ReadOnly> { + public: + using Base = + MeasurementProxyBase, + FullSize, ReadOnly>; + using Index = Base::Index; + using Scalar = Base::Scalar; + using Container = Base::Container; + + using SubspaceHelper = Acts::FixedSubspaceHelper; + + using SubspaceVector = Eigen::Matrix; + using SubspaceVectorMap = + std::conditional_t, + Eigen::Map>; + + using ParametersVector = Eigen::Matrix; + using ParametersVectorMap = + std::conditional_t, + Eigen::Map>; + + using CovarianceMatrix = Eigen::Matrix; + using CovarianceMatrixMap = + std::conditional_t, + Eigen::Map>; + + FixedMeasurementProxy(Container& container_, std::size_t index_) + : Base(container_, index_) { + assert(container().m_entries.at(index()).size == Size && "Size mismatch"); } - template - CovarianceMatrixMap covariance() { - assert(dim == size()); - return CovarianceMatrixMap{m_cov.data()}; + template + FixedMeasurementProxy( + const MeasurementProxyBase& other) + requires(ReadOnly == OtherReadOnly || ReadOnly) + : Base(other) { + assert(container().m_entries.at(index()).size == Size && "Size mismatch"); } - ConstEffectiveCovarianceMatrixMap covariance() const { - return ConstEffectiveCovarianceMatrixMap{m_cov.data(), - static_cast(size()), - static_cast(size())}; + + using Base::container; + using Base::index; + + static constexpr std::size_t size() { return Size; } + + SubspaceHelper subspaceHelper() const { + return SubspaceHelper{subspaceIndexVector()}; } - EffectiveCovarianceMatrixMap covariance() { - return EffectiveCovarianceMatrixMap{m_cov.data(), - static_cast(size()), - static_cast(size())}; + + Acts::SubspaceIndices subspaceIndices() const { + return subspaceHelper().indices(); } - FullParametersVector fullParameters() const { - FullParametersVector result = FullParametersVector::Zero(); - for (std::size_t i = 0; i < size(); ++i) { - result[m_subspaceIndices[i]] = parameters()[i]; - } - return result; + SubspaceVectorMap subspaceIndexVector() const { + return SubspaceVectorMap{ + container().m_subspaceIndices.data() + + container().m_entries.at(index()).subspaceIndexOffset}; } - FullCovarianceMatrix fullCovariance() const { - FullCovarianceMatrix result = FullCovarianceMatrix::Zero(); - for (std::size_t i = 0; i < size(); ++i) { - for (std::size_t j = 0; j < size(); ++j) { - result(m_subspaceIndices[i], m_subspaceIndices[j]) = covariance()(i, j); - } - } - return result; + ParametersVectorMap parameters() const { + return ParametersVectorMap{ + container().m_parameters.data() + + container().m_entries.at(index()).parameterOffset}; } - private: - Acts::SourceLink m_source; - SubspaceIndices m_subspaceIndices; - std::array m_params{}; - std::array m_cov{}; + CovarianceMatrixMap covariance() const { + return CovarianceMatrixMap{ + container().m_covariances.data() + + container().m_entries.at(index()).covarianceOffset}; + } }; -/// Construct a variable-size measurement for the given indices. -/// -/// @tparam parameters_t Input parameters vector type -/// @tparam covariance_t Input covariance matrix type -/// @tparam indices_t Parameter index type, determines the full parameter space -/// @tparam tail_indices_t Helper types required to support variadic arguments; -/// all types must be convertibale to `indices_t`. -/// @param source The link that connects to the underlying detector readout -/// @param params Measured parameters values -/// @param cov Measured parameters covariance -/// @param index0 Required parameter index, measurement must be at least 1d -/// @param tailIndices Additional parameter indices for larger measurements -/// @return Variable-size measurement w/ the correct type and given inputs -/// -/// @note The indices must be ordered and must be consistent with the content of -/// parameters and covariance. -template -VariableSizeMeasurement makeVariableSizeMeasurement( - Acts::SourceLink source, const Eigen::MatrixBase& params, - const Eigen::MatrixBase& cov, indices_t index0, - tail_indices_t... tailIndices) { - using IndexContainer = std::array; - return {std::move(source), IndexContainer{index0, tailIndices...}, params, - cov}; -} - -/// Type that can hold all possible bound measurements. -using BoundVariableMeasurement = VariableSizeMeasurement; - -/// Variable measurement type that can contain all possible combinations. -using Measurement = BoundVariableMeasurement; - -/// Container of measurements. -/// -/// In contrast to the source links, the measurements themself must not be -/// orderable. The source links stored in the measurements are treated -/// as opaque here and no ordering is enforced on the stored measurements. -using MeasurementContainer = std::vector; +template +class VariableMeasurementProxy + : public MeasurementProxyBase, + FullSize, ReadOnly> { + public: + using Base = + MeasurementProxyBase, + FullSize, ReadOnly>; + using Index = Base::Index; + using Scalar = Base::Scalar; + using Container = Base::Container; + + using SubspaceHelper = Acts::VariableSubspaceHelper; + + using SubspaceVector = Eigen::Matrix; + using SubspaceVectorMap = + std::conditional_t, + Eigen::Map>; + + using ParametersVector = Eigen::Matrix; + using ParametersVectorMap = + std::conditional_t, + Eigen::Map>; + + using CovarianceMatrix = + Eigen::Matrix; + using CovarianceMatrixMap = + std::conditional_t, + Eigen::Map>; + + VariableMeasurementProxy(Container& container_, std::size_t index_) + : Base(container_, index_) {} + template + VariableMeasurementProxy( + const MeasurementProxyBase& other) + requires(ReadOnly == OtherReadOnly || ReadOnly) + : Base(other) {} + + using Base::container; + using Base::index; + + SubspaceHelper subspaceHelper() const { + return SubspaceHelper{subspaceIndexVector()}; + } + + SubspaceVectorMap subspaceIndexVector() const { + return SubspaceVectorMap{ + container().m_subspaceIndices.data() + + container().m_entries.at(index()).subspaceIndexOffset, + static_cast(this->size())}; + } + + ParametersVectorMap parameters() const { + return ParametersVectorMap{ + container().m_parameters.data() + + container().m_entries.at(index()).parameterOffset, + static_cast(this->size())}; + } + + CovarianceMatrixMap covariance() const { + const auto size = this->size(); + return CovarianceMatrixMap{ + container().m_covariances.data() + + container().m_entries.at(index()).covarianceOffset, + static_cast(size), static_cast(size)}; + } +}; } // namespace ActsExamples diff --git a/Examples/Framework/src/EventData/MeasurementCalibration.cpp b/Examples/Framework/src/EventData/MeasurementCalibration.cpp index 560593dbd1f..84b1c5a2d45 100644 --- a/Examples/Framework/src/EventData/MeasurementCalibration.cpp +++ b/Examples/Framework/src/EventData/MeasurementCalibration.cpp @@ -6,11 +6,12 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at http://mozilla.org/MPL/2.0/. +#include "ActsExamples/EventData/MeasurementCalibration.hpp" + #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/EventData/SourceLink.hpp" #include "ActsExamples/EventData/IndexSourceLink.hpp" #include "ActsExamples/EventData/Measurement.hpp" -#include #include #include @@ -31,18 +32,19 @@ void ActsExamples::PassThroughCalibrator::calibrate( assert((idxSourceLink.index() < measurements.size()) && "Source link index is outside the container bounds"); - const auto& measurement = measurements[idxSourceLink.index()]; + const ConstVariableBoundMeasurementProxy measurement = + measurements.getMeasurement(idxSourceLink.index()); Acts::visit_measurement(measurement.size(), [&](auto N) -> void { constexpr std::size_t kMeasurementSize = decltype(N)::value; + const ConstFixedBoundMeasurementProxy fixedMeasurement = + measurement; trackState.allocateCalibrated(kMeasurementSize); - trackState.calibrated() = - measurement.parameters(); + trackState.calibrated() = fixedMeasurement.parameters(); trackState.calibratedCovariance() = - measurement.covariance(); - trackState.setSubspaceIndices( - measurement.subspaceIndices()); + fixedMeasurement.covariance(); + trackState.setSubspaceIndices(fixedMeasurement.subspaceIndices()); }); } diff --git a/Examples/Framework/src/EventData/ScalingCalibrator.cpp b/Examples/Framework/src/EventData/ScalingCalibrator.cpp index 6681af67ee4..5f7cc896883 100644 --- a/Examples/Framework/src/EventData/ScalingCalibrator.cpp +++ b/Examples/Framework/src/EventData/ScalingCalibrator.cpp @@ -151,7 +151,8 @@ void ActsExamples::ScalingCalibrator::calibrate( const Cluster& cl = clusters->at(idxSourceLink.index()); ConstantTuple ct = m_calib_maps.at(mgid).at(cl.sizeLoc0, cl.sizeLoc1); - const auto& measurement = measurements[idxSourceLink.index()]; + const ConstVariableBoundMeasurementProxy measurement = + measurements.getMeasurement(idxSourceLink.index()); assert(measurement.contains(Acts::eBoundLoc0) && "Measurement does not contain the required bound loc0"); @@ -161,21 +162,24 @@ void ActsExamples::ScalingCalibrator::calibrate( auto boundLoc0 = measurement.indexOf(Acts::eBoundLoc0); auto boundLoc1 = measurement.indexOf(Acts::eBoundLoc1); - Measurement measurementCopy = measurement; - measurementCopy.parameters()[boundLoc0] += ct.x_offset; - measurementCopy.parameters()[boundLoc1] += ct.y_offset; - measurementCopy.covariance()(boundLoc0, boundLoc0) *= ct.x_scale; - measurementCopy.covariance()(boundLoc1, boundLoc1) *= ct.y_scale; - Acts::visit_measurement(measurement.size(), [&](auto N) -> void { constexpr std::size_t kMeasurementSize = decltype(N)::value; + const ConstFixedBoundMeasurementProxy fixedMeasurement = + measurement; + + Acts::ActsVector calibratedParameters = + fixedMeasurement.parameters(); + Acts::ActsSquareMatrix calibratedCovariance = + fixedMeasurement.covariance(); + + calibratedParameters[boundLoc0] += ct.x_offset; + calibratedParameters[boundLoc1] += ct.y_offset; + calibratedCovariance(boundLoc0, boundLoc0) *= ct.x_scale; + calibratedCovariance(boundLoc1, boundLoc1) *= ct.y_scale; trackState.allocateCalibrated(kMeasurementSize); - trackState.calibrated() = - measurementCopy.parameters(); - trackState.calibratedCovariance() = - measurementCopy.covariance(); - trackState.setSubspaceIndices( - measurementCopy.subspaceIndices()); + trackState.calibrated() = calibratedParameters; + trackState.calibratedCovariance() = calibratedCovariance; + trackState.setSubspaceIndices(fixedMeasurement.subspaceIndices()); }); } diff --git a/Examples/Io/Csv/src/CsvMeasurementReader.cpp b/Examples/Io/Csv/src/CsvMeasurementReader.cpp index e1c31d52b63..ee036ea86a9 100644 --- a/Examples/Io/Csv/src/CsvMeasurementReader.cpp +++ b/Examples/Io/Csv/src/CsvMeasurementReader.cpp @@ -193,11 +193,14 @@ ActsExamples::ProcessCode ActsExamples::CsvMeasurementReader::read( readMeasurementsByGeometryId(m_cfg.inputDir, ctx.eventNumber); // Prepare containers for the hit data using the framework event data types - GeometryIdMultimap orderedMeasurements; + MeasurementContainer tmpMeasurements; + GeometryIdMultimap orderedMeasurements; IndexMultimap measurementSimHitsMap; IndexSourceLinkContainer sourceLinks; // need list here for stable addresses std::list sourceLinkStorage; + + tmpMeasurements.reserve(measurementData.size()); orderedMeasurements.reserve(measurementData.size()); // Safe long as we have single particle to sim hit association measurementSimHitsMap.reserve(measurementData.size()); @@ -251,14 +254,15 @@ ActsExamples::ProcessCode ActsExamples::CsvMeasurementReader::read( // the measurement will be stored is known before adding it. const Index index = orderedMeasurements.size(); IndexSourceLink& sourceLink = sourceLinkStorage.emplace_back(geoId, index); - auto measurement = createMeasurement(dParameters, sourceLink); + auto measurement = + createMeasurement(tmpMeasurements, dParameters, sourceLink); // Due to the previous sorting of the raw hit data by geometry id, new // measurements should always end up at the end of the container. previous // elements were not touched; cluster indices remain stable and can // be used to identify the m. - auto inserted = orderedMeasurements.emplace_hint( - orderedMeasurements.end(), geoId, std::move(measurement)); + auto inserted = orderedMeasurements.emplace_hint(orderedMeasurements.end(), + geoId, measurement); if (std::next(inserted) != orderedMeasurements.end()) { ACTS_FATAL("Something went horribly wrong with the hit sorting"); return ProcessCode::ABORT; @@ -269,7 +273,8 @@ ActsExamples::ProcessCode ActsExamples::CsvMeasurementReader::read( MeasurementContainer measurements; for (auto& [_, meas] : orderedMeasurements) { - measurements.emplace_back(std::move(meas)); + auto measurement = measurements.makeMeasurement(meas.size()); + measurement.copyFrom(meas); } // Generate measurement-particles-map diff --git a/Examples/Io/Csv/src/CsvMeasurementWriter.cpp b/Examples/Io/Csv/src/CsvMeasurementWriter.cpp index 0f65fca0d20..f31c812ba46 100644 --- a/Examples/Io/Csv/src/CsvMeasurementWriter.cpp +++ b/Examples/Io/Csv/src/CsvMeasurementWriter.cpp @@ -14,6 +14,7 @@ #include "ActsExamples/EventData/Cluster.hpp" #include "ActsExamples/EventData/Index.hpp" #include "ActsExamples/EventData/IndexSourceLink.hpp" +#include "ActsExamples/EventData/Measurement.hpp" #include "ActsExamples/Framework/AlgorithmContext.hpp" #include "ActsExamples/Io/Csv/CsvInputOutput.hpp" #include "ActsExamples/Utilities/Paths.hpp" @@ -92,7 +93,8 @@ ActsExamples::ProcessCode ActsExamples::CsvMeasurementWriter::writeT( << " measurements in this event."); for (Index measIdx = 0u; measIdx < measurements.size(); ++measIdx) { - const auto& measurement = measurements[measIdx]; + const ConstVariableBoundMeasurementProxy measurement = + measurements.getMeasurement(measIdx); auto simHitIndices = makeRange(measurementSimHitsMap.equal_range(measIdx)); for (auto [_, simHitIdx] : simHitIndices) { diff --git a/Examples/Io/Root/src/RootAthenaDumpReader.cpp b/Examples/Io/Root/src/RootAthenaDumpReader.cpp index 8300237ef9e..001a8d90bdb 100644 --- a/Examples/Io/Root/src/RootAthenaDumpReader.cpp +++ b/Examples/Io/Root/src/RootAthenaDumpReader.cpp @@ -372,7 +372,7 @@ ActsExamples::ProcessCode ActsExamples::RootAthenaDumpReader::read( IndexSourceLink sl(Acts::GeometryIdentifier{CLmoduleID[im]}, im); - measurements.push_back(createMeasurement(digiPars, sl)); + createMeasurement(measurements, digiPars, sl); // Create measurement particles map and particles container for (const auto& [subevt, bc] : Acts::zip(CLparticleLink_eventIndex->at(im), diff --git a/Examples/Io/Root/src/RootMeasurementWriter.cpp b/Examples/Io/Root/src/RootMeasurementWriter.cpp index 6dc10e3b4b1..88b66e192f0 100644 --- a/Examples/Io/Root/src/RootMeasurementWriter.cpp +++ b/Examples/Io/Root/src/RootMeasurementWriter.cpp @@ -13,6 +13,7 @@ #include "ActsExamples/EventData/AverageSimHits.hpp" #include "ActsExamples/EventData/Index.hpp" #include "ActsExamples/EventData/IndexSourceLink.hpp" +#include "ActsExamples/EventData/Measurement.hpp" #include "ActsExamples/Framework/AlgorithmContext.hpp" #include "ActsExamples/Utilities/Range.hpp" @@ -152,9 +153,9 @@ struct RootMeasurementWriter::DigitizationTree { /// Convenience function to fill bound parameters /// /// @param m The measurement - void fillBoundMeasurement(const Measurement& m) { + void fillBoundMeasurement(const ConstVariableBoundMeasurementProxy& m) { for (unsigned int i = 0; i < m.size(); ++i) { - auto ib = m.subspaceIndices()[i]; + auto ib = m.subspaceIndexVector()[i]; recBound[ib] = m.parameters()[i]; varBound[ib] = m.covariance()(i, i); @@ -266,7 +267,8 @@ ProcessCode RootMeasurementWriter::writeT( std::lock_guard lock(m_writeMutex); for (Index hitIdx = 0u; hitIdx < measurements.size(); ++hitIdx) { - const auto& meas = measurements[hitIdx]; + const ConstVariableBoundMeasurementProxy meas = + measurements.getMeasurement(hitIdx); Acts::GeometryIdentifier geoId = meas.sourceLink().template get().geometryId(); diff --git a/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp b/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp index bab0b16e9e9..f76d09b791a 100644 --- a/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp +++ b/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp @@ -50,8 +50,15 @@ std::default_random_engine rng(123); BOOST_AUTO_TEST_SUITE(EventDataMeasurement) BOOST_DATA_TEST_CASE(VariableBoundOne, bd::make(boundIndices), index) { + MeasurementContainer container; + auto [params, cov] = generateParametersCovariance(rng); - auto meas = makeVariableSizeMeasurement(source, params, cov, index); + + FixedBoundMeasurementProxy<1> meas = container.makeMeasurement<1>(); + meas.sourceLink() = source; + meas.setSubspaceIndices(std::array{index}); + meas.parameters() = params; + meas.covariance() = cov; BOOST_CHECK_EQUAL(meas.size(), 1); for (auto i : boundIndices) { @@ -68,10 +75,17 @@ BOOST_DATA_TEST_CASE(VariableBoundOne, bd::make(boundIndices), index) { } BOOST_AUTO_TEST_CASE(VariableBoundAll) { + MeasurementContainer container; + auto [params, cov] = generateBoundParametersCovariance(rng); - auto meas = makeVariableSizeMeasurement(source, params, cov, eBoundLoc0, - eBoundLoc1, eBoundPhi, eBoundTheta, - eBoundQOverP, eBoundTime); + + FixedBoundMeasurementProxy meas = + container.makeMeasurement(); + meas.sourceLink() = source; + meas.setSubspaceIndices(std::array{eBoundLoc0, eBoundLoc1, eBoundTime, + eBoundPhi, eBoundTheta, eBoundQOverP}); + meas.parameters() = params; + meas.covariance() = cov; BOOST_CHECK_EQUAL(meas.size(), eBoundSize); for (auto i : boundIndices) { @@ -83,9 +97,17 @@ BOOST_AUTO_TEST_CASE(VariableBoundAll) { } BOOST_AUTO_TEST_CASE(VariableBoundReassign) { + MeasurementContainer container; + // generate w/ a single parameter auto [par1, cov1] = generateParametersCovariance(rng); - auto meas = makeVariableSizeMeasurement(source, par1, cov1, eBoundTheta); + + VariableBoundMeasurementProxy meas = container.makeMeasurement(1); + meas.sourceLink() = source; + meas.setSubspaceIndices(std::array{eBoundTheta}); + meas.parameters() = par1; + meas.covariance() = cov1; + BOOST_CHECK_EQUAL(meas.size(), 1); BOOST_CHECK(!meas.contains(eBoundLoc0)); BOOST_CHECK(!meas.contains(eBoundLoc1)); @@ -96,9 +118,14 @@ BOOST_AUTO_TEST_CASE(VariableBoundReassign) { // reassign w/ all parameters auto [parN, covN] = generateBoundParametersCovariance(rng); - meas = makeVariableSizeMeasurement(source, parN, covN, eBoundLoc0, eBoundLoc1, - eBoundPhi, eBoundTheta, eBoundQOverP, - eBoundTime); + + meas = container.makeMeasurement(eBoundSize); + meas.sourceLink() = source; + meas.setSubspaceIndices(std::array{eBoundLoc0, eBoundLoc1, eBoundTime, + eBoundPhi, eBoundTheta, eBoundQOverP}); + meas.parameters() = parN; + meas.covariance() = covN; + BOOST_CHECK_EQUAL(meas.size(), eBoundSize); BOOST_CHECK(meas.contains(eBoundLoc0)); BOOST_CHECK(meas.contains(eBoundLoc1)); diff --git a/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp b/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp index 4b4b6d78bbe..6725a8f4668 100644 --- a/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp +++ b/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp @@ -50,11 +50,11 @@ BOOST_AUTO_TEST_CASE(CsvMeasurementRoundTrip) { Acts::Vector2 p = Acts::Vector2::Random(); Acts::SquareMatrix2 c = Acts::SquareMatrix2::Random(); - BoundVariableMeasurement m(Acts::SourceLink{sl}, - std::array{Acts::eBoundLoc0, Acts::eBoundLoc1}, - p, c); - - measOriginal.push_back(m); + FixedBoundMeasurementProxy<2> m = measOriginal.makeMeasurement<2>(); + m.sourceLink() = Acts::SourceLink(sl); + m.setSubspaceIndices(std::array{Acts::eBoundLoc0, Acts::eBoundLoc1}); + m.parameters() = p; + m.covariance() = c; ActsExamples::Cluster cl;