From e0d235c170a15c72893d1f7e53b2efe711c67676 Mon Sep 17 00:00:00 2001 From: Andreas Stefl Date: Wed, 21 Aug 2024 09:36:54 +0200 Subject: [PATCH] Revert "refactor: Rework projector (#3453)" This reverts commit a28acc306ccd11f064c81a68e8eb7f3cc7526c2c. --- .../Acts/EventData/SubspaceHelpers.hpp | 253 -------------- .../Acts/EventData/TrackStateProxy.hpp | 85 +---- .../Acts/EventData/TrackStateProxy.ipp | 7 +- Core/include/Acts/EventData/Types.hpp | 11 +- .../Acts/EventData/VectorMultiTrajectory.hpp | 5 +- .../Acts/EventData/detail/TestSourceLink.hpp | 13 +- Core/include/Acts/Propagator/Navigator.hpp | 2 +- .../Acts/TrackFinding/MeasurementSelector.hpp | 5 +- .../Acts/TrackFinding/MeasurementSelector.ipp | 10 +- .../Acts/TrackFitting/GainMatrixUpdater.hpp | 6 +- .../detail/GainMatrixUpdaterImpl.hpp | 11 +- .../Acts/Utilities/detail/Subspace.hpp | 327 ++++++++++++++++++ Core/src/EventData/VectorMultiTrajectory.cpp | 6 +- Core/src/TrackFinding/MeasurementSelector.cpp | 21 +- .../TrackFinding/src/HoughTransformSeeder.cpp | 8 +- .../TrackFitting/src/RefittingCalibrator.cpp | 2 +- .../Framework/ML/src/NeuralCalibrator.cpp | 25 +- .../ActsExamples/EventData/Measurement.hpp | 90 ++--- .../src/EventData/MeasurementCalibration.cpp | 5 +- .../src/EventData/ScalingCalibrator.cpp | 17 +- .../Io/Root/src/RootMeasurementWriter.cpp | 6 +- Examples/Python/tests/root_file_hashes.txt | 32 +- Plugins/Podio/edm.yml | 2 +- .../Podio/PodioTrackStateContainer.hpp | 7 +- .../Core/TrackFitting/MbfSmootherTests.cpp | 11 +- Tests/UnitTests/Core/Utilities/CMakeLists.txt | 1 + .../Core/Utilities/SubspaceTests.cpp | 182 ++++++++++ .../Examples/EventData/MeasurementTests.cpp | 8 +- .../Io/Csv/MeasurementReaderWriterTests.cpp | 2 +- 29 files changed, 647 insertions(+), 513 deletions(-) delete mode 100644 Core/include/Acts/EventData/SubspaceHelpers.hpp create mode 100644 Core/include/Acts/Utilities/detail/Subspace.hpp create mode 100644 Tests/UnitTests/Core/Utilities/SubspaceTests.cpp diff --git a/Core/include/Acts/EventData/SubspaceHelpers.hpp b/Core/include/Acts/EventData/SubspaceHelpers.hpp deleted file mode 100644 index 737c12211a5..00000000000 --- a/Core/include/Acts/EventData/SubspaceHelpers.hpp +++ /dev/null @@ -1,253 +0,0 @@ -// This file is part of the Acts project. -// -// Copyright (C) 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 -// file, You can obtain one at http://mozilla.org/MPL/2.0/. - -#pragma once - -#include "Acts/Definitions/Algebra.hpp" -#include "Acts/Definitions/TrackParametrization.hpp" -#include "Acts/EventData/Types.hpp" -#include "Acts/Utilities/AlgebraHelpers.hpp" -#include "Acts/Utilities/Enumerate.hpp" - -#include -#include - -#include - -namespace Acts { - -template -inline static bool checkSubspaceIndices(const Container& container, - std::size_t fullSize, - std::size_t subspaceSize) { - if (subspaceSize > fullSize) { - return false; - } - if (container.size() != subspaceSize) { - return false; - } - for (auto it = container.begin(); it != container.end();) { - auto index = *it; - if (index >= fullSize) { - return false; - } - ++it; - if (std::find(it, container.end(), index) != container.end()) { - return false; - } - } - return true; -} - -template -class SubspaceHelperBase { - public: - static constexpr std::size_t kFullSize = FullSize; - - using FullSquareMatrix = ActsSquareMatrix; - - bool empty() const { return self().empty(); } - std::size_t size() const { return self().size(); } - - auto operator[](std::size_t i) const { return self()[i]; } - - auto begin() const { return self().begin(); } - auto end() const { return self().end(); } - - FullSquareMatrix fullProjector() const { - FullSquareMatrix result = FullSquareMatrix::Zero(); - for (auto [i, index] : enumerate(*this)) { - result(i, index) = 1; - } - return result; - } - - FullSquareMatrix fullExpander() const { - FullSquareMatrix result = FullSquareMatrix::Zero(); - for (auto [i, index] : enumerate(*this)) { - result(index, i) = 1; - } - return result; - } - - ProjectorBitset projectorBitset() const - requires(kFullSize <= 8) - { - return matrixToBitset(fullProjector()).to_ullong(); - } - - private: - const Derived& self() const { return static_cast(*this); } -}; - -template -class VariableSubspaceHelper - : public SubspaceHelperBase, - FullSize> { - public: - static constexpr std::size_t kFullSize = FullSize; - - using IndexType = index_t; - using Container = boost::container::static_vector; - - template - explicit VariableSubspaceHelper(const OtherContainer& indices) { - assert(checkSubspaceIndices(indices, kFullSize, indices.size()) && - "Invalid indices"); - m_indices.resize(indices.size()); - std::transform(indices.begin(), indices.end(), m_indices.begin(), - [](auto index) { return static_cast(index); }); - } - - bool empty() const { return m_indices.empty(); } - std::size_t size() const { return m_indices.size(); } - - IndexType operator[](std::size_t i) const { return m_indices[i]; } - - auto begin() const { return m_indices.begin(); } - auto end() const { return m_indices.end(); } - - private: - Container m_indices; -}; - -template -class FixedSubspaceHelper - : public SubspaceHelperBase< - FixedSubspaceHelper, FullSize> { - public: - static constexpr std::size_t kFullSize = FullSize; - static constexpr std::size_t kSubspaceSize = SubspaceSize; - static_assert(kSubspaceSize <= kFullSize, "Invalid subspace size"); - - using Projector = ActsMatrix; - using Expander = ActsMatrix; - using Vector = ActsVector; - using SquareMatrix = ActsSquareMatrix; - template - using ApplyLeftResult = ActsMatrix; - template - using ApplyRightResult = ActsMatrix; - - using IndexType = index_t; - using Container = std::array; - - template - explicit FixedSubspaceHelper(const OtherContainer& indices) { - assert(checkSubspaceIndices(indices, kFullSize, kSubspaceSize) && - "Invalid indices"); - std::transform(indices.begin(), indices.end(), m_indices.begin(), - [](auto index) { return static_cast(index); }); - } - - bool empty() const { return m_indices.empty(); } - std::size_t size() const { return m_indices.size(); } - - IndexType operator[](std::uint32_t i) const { return m_indices[i]; } - - auto begin() const { return m_indices.begin(); } - auto end() const { return m_indices.end(); } - - Projector projector() const { - Projector result = Projector::Zero(); - for (auto [i, index] : enumerate(*this)) { - result(i, index) = 1; - } - return result; - } - - Expander expander() const { - Expander result = Expander::Zero(); - for (auto [i, index] : enumerate(*this)) { - result(index, i) = 1; - } - return result; - } - - template - ApplyLeftResult applyLeftOf( - const Eigen::DenseBase& matrix) const { - assert(matrix.rows() == kFullSize && "Invalid matrix size"); - ApplyLeftResult result = ApplyLeftResult::Zero(); - for (auto [i, indexI] : enumerate(*this)) { - for (std::size_t j = 0; j < K; ++j) { - result(i, j) = matrix(indexI, j); - } - } - return result; - } - - template - ApplyRightResult applyRightOf( - const Eigen::DenseBase& matrix) const { - assert(matrix.cols() == kSubspaceSize && "Invalid matrix size"); - ApplyRightResult result = ApplyRightResult::Zero(); - for (std::size_t i = 0; i < N; ++i) { - for (auto [j, indexJ] : enumerate(*this)) { - result(i, j) = matrix(i, indexJ); - } - } - return result; - } - - template - Vector projectVector(const Eigen::DenseBase& fullVector) const { - assert(fullVector.size() == kFullSize && "Invalid full vector size"); - Vector result = Vector::Zero(); - for (auto [i, index] : enumerate(*this)) { - result(i) = fullVector(index); - } - return result; - } - - template - SquareMatrix projectMatrix( - const Eigen::DenseBase& fullMatrix) const { - assert(fullMatrix.rows() == kFullSize && fullMatrix.cols() == kFullSize && - "Invalid full matrix size"); - SquareMatrix result = SquareMatrix::Zero(); - for (auto [i, indexI] : enumerate(*this)) { - for (auto [j, indexJ] : enumerate(*this)) { - result(i, j) = fullMatrix(indexI, indexJ); - } - } - return result; - } - - private: - Container m_indices{}; -}; - -template -using FixedBoundSubspaceHelper = - FixedSubspaceHelper; -using VariableBoundSubspaceHelper = - VariableSubspaceHelper; - -template -SubspaceIndices projectorToSubspaceIndices( - const Eigen::DenseBase& projector) { - auto rows = static_cast(projector.rows()); - auto cols = static_cast(projector.cols()); - assert(cols == kFullSize && rows <= kFullSize && "Invalid projector size"); - SubspaceIndices result; - result.fill(kFullSize); - for (std::size_t i = 0; i < rows; ++i) { - for (std::size_t j = 0; j < cols; ++j) { - assert((projector(i, j) == 0 || projector(i, j) == 1) && - "Invalid projector value"); - if (projector(i, j) == 1) { - result[i] = j; - } - } - } - return result; -} - -} // namespace Acts diff --git a/Core/include/Acts/EventData/TrackStateProxy.hpp b/Core/include/Acts/EventData/TrackStateProxy.hpp index 3358a94cf62..1e68994fad8 100644 --- a/Core/include/Acts/EventData/TrackStateProxy.hpp +++ b/Core/include/Acts/EventData/TrackStateProxy.hpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2023-2024 CERN for the benefit of the Acts project +// Copyright (C) 2023 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 @@ -11,7 +11,6 @@ #include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/EventData/SourceLink.hpp" -#include "Acts/EventData/SubspaceHelpers.hpp" #include "Acts/EventData/TrackStatePropMask.hpp" #include "Acts/EventData/TrackStateProxyConcept.hpp" #include "Acts/EventData/TrackStateType.hpp" @@ -22,7 +21,6 @@ #include "Acts/Utilities/Helpers.hpp" #include -#include #include @@ -648,8 +646,7 @@ class TrackStateProxy { /// @param projector The projector in the form of a dense matrix /// @note @p projector is assumed to only have 0s or 1s as components. template - [[deprecated("use setProjector(span) instead")]] void setProjector( - const Eigen::MatrixBase& projector) + void setProjector(const Eigen::MatrixBase& projector) requires(!ReadOnly) { constexpr int rows = Eigen::MatrixBase::RowsAtCompileTime; @@ -672,8 +669,9 @@ class TrackStateProxy { fullProjector.template topLeftCorner() = projector; // convert to bitset before storing - ProjectorBitset projectorBitset = matrixToBitset(fullProjector).to_ulong(); - setProjectorBitset(projectorBitset); + auto projectorBitset = matrixToBitset(fullProjector); + component() = + projectorBitset.to_ullong(); } /// Directly get the projector bitset, a compressed form of a projection @@ -682,9 +680,9 @@ class TrackStateProxy { /// to another. Use the `projector` or `effectiveProjector` method if /// you want to access the matrix. /// @return The projector bitset - [[deprecated("use projector() instead")]] ProjectorBitset projectorBitset() - const { - return variableBoundSubspaceHelper().projectorBitset(); + ProjectorBitset projectorBitset() const { + assert(has()); + return component(); } /// Set the projector bitset, a compressed form of a projection matrix @@ -693,70 +691,11 @@ class TrackStateProxy { /// @note This is mainly to copy explicitly a projector from one state /// to another. If you have a projection matrix, set it with /// `setProjector`. - [[deprecated("use setProjector(span) instead")]] void setProjectorBitset( - ProjectorBitset proj) - requires(!ReadOnly) - { - BoundMatrix projMatrix = bitsetToMatrix(proj); - BoundSubspaceIndices boundSubspace = - projectorToSubspaceIndices(projMatrix); - setBoundSubspaceIndices(boundSubspace); - } - - BoundSubspaceIndices boundSubspaceIndices() const { - assert(has()); - return component(); - } - - template - SubspaceIndices subspaceIndices() const { - BoundSubspaceIndices boundSubspace = BoundSubspaceIndices(); - SubspaceIndices subspace; - std::copy(boundSubspace.begin(), boundSubspace.begin() + measdim, - subspace.begin()); - return subspace; - } - - void setBoundSubspaceIndices(BoundSubspaceIndices boundSubspace) + void setProjectorBitset(ProjectorBitset proj) requires(!ReadOnly) { assert(has()); - component() = boundSubspace; - } - - template - void setSubspaceIndices(SubspaceIndices proj) - requires(!ReadOnly) - { - assert(has()); - BoundSubspaceIndices& boundSubspace = - component(); - std::copy(proj.begin(), proj.end(), boundSubspace.begin()); - } - - template - void setSubspaceIndices(std::array subspaceIndices) - requires(!ReadOnly) - { - assert(has()); - BoundSubspaceIndices& boundSubspace = - component(); - std::transform(subspaceIndices.begin(), subspaceIndices.end(), - boundSubspace.begin(), - [](index_t i) { return static_cast(i); }); - } - - VariableBoundSubspaceHelper variableBoundSubspaceHelper() const { - BoundSubspaceIndices boundSubspace = boundSubspaceIndices(); - std::span validSubspaceIndices( - boundSubspace.begin(), boundSubspace.begin() + calibratedSize()); - return VariableBoundSubspaceHelper(validSubspaceIndices); - } - - template - FixedBoundSubspaceHelper fixedBoundSubspaceHelper() const { - SubspaceIndices subspace = subspaceIndices(); - return FixedBoundSubspaceHelper(subspace); + component() = proj; } /// Uncalibrated measurement in the form of a source link. Const version @@ -1021,7 +960,7 @@ class TrackStateProxy { other.template calibratedCovariance(); }); - setBoundSubspaceIndices(other.boundSubspaceIndices()); + setProjectorBitset(other.projectorBitset()); } } else { if (ACTS_CHECK_BIT(mask, PM::Predicted) && @@ -1068,7 +1007,7 @@ class TrackStateProxy { other.template calibratedCovariance(); }); - setBoundSubspaceIndices(other.boundSubspaceIndices()); + setProjectorBitset(other.projectorBitset()); } } diff --git a/Core/include/Acts/EventData/TrackStateProxy.ipp b/Core/include/Acts/EventData/TrackStateProxy.ipp index c12c34ac744..cfce0d9648c 100644 --- a/Core/include/Acts/EventData/TrackStateProxy.ipp +++ b/Core/include/Acts/EventData/TrackStateProxy.ipp @@ -1,13 +1,12 @@ // This file is part of the Acts project. // -// Copyright (C) 2023-2024 CERN for the benefit of the Acts project +// Copyright (C) 2023 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 // file, You can obtain one at http://mozilla.org/MPL/2.0/. namespace Acts { - template inline TrackStateProxy::TrackStateProxy( detail_lt::ConstIf, ReadOnly>& trajectory, @@ -63,7 +62,9 @@ inline auto TrackStateProxy::covariance() const template inline auto TrackStateProxy::projector() const -> Projector { - return variableBoundSubspaceHelper().fullProjector(); + assert(has()); + return bitsetToMatrix( + component()); } template diff --git a/Core/include/Acts/EventData/Types.hpp b/Core/include/Acts/EventData/Types.hpp index 3aabf64c1e1..c212e0bf5d0 100644 --- a/Core/include/Acts/EventData/Types.hpp +++ b/Core/include/Acts/EventData/Types.hpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2023-2024 CERN for the benefit of the Acts project +// Copyright (C) 2023 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,23 +8,14 @@ #pragma once -#include "Acts/Definitions/TrackParametrization.hpp" - #include #include namespace Acts { - using TrackIndexType = std::uint32_t; static constexpr TrackIndexType kTrackIndexInvalid = std::numeric_limits::max(); using ProjectorBitset = std::uint64_t; -template -using SubspaceIndices = std::array; -using BoundSubspaceIndices = SubspaceIndices; -static constexpr BoundSubspaceIndices kBoundSubspaceIndicesInvalid = { - eBoundSize, eBoundSize, eBoundSize, eBoundSize, eBoundSize, eBoundSize}; - } // namespace Acts diff --git a/Core/include/Acts/EventData/VectorMultiTrajectory.hpp b/Core/include/Acts/EventData/VectorMultiTrajectory.hpp index 17a7da758e4..8db631b9d33 100644 --- a/Core/include/Acts/EventData/VectorMultiTrajectory.hpp +++ b/Core/include/Acts/EventData/VectorMultiTrajectory.hpp @@ -13,7 +13,6 @@ #include "Acts/EventData/MultiTrajectoryBackendConcept.hpp" #include "Acts/EventData/SourceLink.hpp" #include "Acts/EventData/TrackStatePropMask.hpp" -#include "Acts/EventData/Types.hpp" #include "Acts/EventData/detail/DynamicColumn.hpp" #include "Acts/EventData/detail/DynamicKeyIterator.hpp" #include "Acts/Utilities/HashedString.hpp" @@ -133,7 +132,7 @@ class VectorMultiTrajectoryBase { h("meas", isMeas, weight(meas_size)); h("measCov", isMeas, weight(meas_cov_size)); h("sourceLinks", isMeas, weight(sizeof(const SourceLink))); - h("projectors", isMeas, weight(sizeof(BoundSubspaceIndices))); + h("projectors", isMeas, weight(sizeof(ProjectorBitset))); } if (ts.hasJacobian() && @@ -327,7 +326,7 @@ class VectorMultiTrajectoryBase { std::vector::Covariance> m_jac; std::vector> m_sourceLinks; - std::vector m_projectors; + std::vector m_projectors; // owning vector of shared pointers to surfaces // diff --git a/Core/include/Acts/EventData/detail/TestSourceLink.hpp b/Core/include/Acts/EventData/detail/TestSourceLink.hpp index aadd8e4dc91..d0bedcce0ef 100644 --- a/Core/include/Acts/EventData/detail/TestSourceLink.hpp +++ b/Core/include/Acts/EventData/detail/TestSourceLink.hpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2020-2024 CERN for the benefit of the Acts project +// Copyright (C) 2020 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 @@ -17,6 +17,7 @@ #include "Acts/Geometry/GeometryIdentifier.hpp" #include "Acts/Geometry/TrackingGeometry.hpp" #include "Acts/Utilities/CalibrationContext.hpp" +#include "Acts/Utilities/detail/Subspace.hpp" #include #include @@ -138,20 +139,22 @@ void testSourceLinkCalibratorReturn( trackState.allocateCalibrated(2); trackState.template calibrated<2>() = sl.parameters; trackState.template calibratedCovariance<2>() = sl.covariance; - trackState.template setSubspaceIndices( - std::array{sl.indices[0], sl.indices[1]}); + trackState.setProjector(FixedSizeSubspace( + std::array{sl.indices[0], sl.indices[1]}) + .projector()); } else if (sl.indices[0] != Acts::eBoundSize) { trackState.allocateCalibrated(1); trackState.template calibrated<1>() = sl.parameters.head<1>(); trackState.template calibratedCovariance<1>() = sl.covariance.topLeftCorner<1, 1>(); - trackState.template setSubspaceIndices(std::array{sl.indices[0]}); + trackState.setProjector(FixedSizeSubspace( + std::array{sl.indices[0]}) + .projector()); } else { throw std::runtime_error( "Tried to extract measurement from invalid TestSourceLink"); } } - /// Extract the measurement from a TestSourceLink. /// /// @param gctx Unused diff --git a/Core/include/Acts/Propagator/Navigator.hpp b/Core/include/Acts/Propagator/Navigator.hpp index 3b2def43002..4e33c6b1bc0 100644 --- a/Core/include/Acts/Propagator/Navigator.hpp +++ b/Core/include/Acts/Propagator/Navigator.hpp @@ -988,7 +988,7 @@ class Navigator { << targetIntersection.position().x() << ", " << targetIntersection.position().y() << ", " << targetIntersection.position().z() << ")"); - // get the target volume from the intersection + /// get the target volume from the intersection auto tPosition = targetIntersection.position(); state.navigation.targetVolume = m_cfg.trackingGeometry->lowestTrackingVolume(state.geoContext, diff --git a/Core/include/Acts/TrackFinding/MeasurementSelector.hpp b/Core/include/Acts/TrackFinding/MeasurementSelector.hpp index 32d6bcf12c0..a05f67d7cc0 100644 --- a/Core/include/Acts/TrackFinding/MeasurementSelector.hpp +++ b/Core/include/Acts/TrackFinding/MeasurementSelector.hpp @@ -12,7 +12,6 @@ #include "Acts/EventData/MeasurementHelpers.hpp" #include "Acts/EventData/MultiTrajectory.hpp" #include "Acts/EventData/TrackParameters.hpp" -#include "Acts/EventData/Types.hpp" #include "Acts/Geometry/GeometryHierarchyMap.hpp" #include "Acts/Geometry/GeometryIdentifier.hpp" #include "Acts/TrackFinding/CombinatorialKalmanFilterError.hpp" @@ -106,7 +105,9 @@ class MeasurementSelector { false>::Parameters predicted, TrackStateTraits::Covariance predictedCovariance, - BoundSubspaceIndices projector, unsigned int calibratedSize) const; + TrackStateTraits::Projector projector, + unsigned int calibratedSize) const; Config m_config; }; diff --git a/Core/include/Acts/TrackFinding/MeasurementSelector.ipp b/Core/include/Acts/TrackFinding/MeasurementSelector.ipp index 49fefe3581b..49b03e4730e 100644 --- a/Core/include/Acts/TrackFinding/MeasurementSelector.ipp +++ b/Core/include/Acts/TrackFinding/MeasurementSelector.ipp @@ -67,11 +67,11 @@ MeasurementSelector::select( // with compile time size. That way the Eigen math operations are // still done with compile time size and no dynamic memory allocation // is needed. - double chi2 = calculateChi2( - trackState.effectiveCalibrated().data(), - trackState.effectiveCalibratedCovariance().data(), - trackState.predicted(), trackState.predictedCovariance(), - trackState.boundSubspaceIndices(), trackState.calibratedSize()); + double chi2 = + calculateChi2(trackState.effectiveCalibrated().data(), + trackState.effectiveCalibratedCovariance().data(), + trackState.predicted(), trackState.predictedCovariance(), + trackState.projector(), trackState.calibratedSize()); trackState.chi2() = chi2; if (chi2 < minChi2) { diff --git a/Core/include/Acts/TrackFitting/GainMatrixUpdater.hpp b/Core/include/Acts/TrackFitting/GainMatrixUpdater.hpp index 538c875709b..b5a49a6295d 100644 --- a/Core/include/Acts/TrackFitting/GainMatrixUpdater.hpp +++ b/Core/include/Acts/TrackFitting/GainMatrixUpdater.hpp @@ -10,7 +10,6 @@ #include "Acts/EventData/MeasurementHelpers.hpp" #include "Acts/EventData/MultiTrajectory.hpp" -#include "Acts/EventData/Types.hpp" #include "Acts/Geometry/GeometryContext.hpp" #include "Acts/TrackFitting/KalmanFitterError.hpp" #include "Acts/Utilities/Logger.hpp" @@ -29,7 +28,8 @@ class GainMatrixUpdater { // This is used to build a covariance matrix view in the .cpp file const double* calibrated; const double* calibratedCovariance; - BoundSubspaceIndices projector; + TrackStateTraits::Projector projector; TrackStateTraits::Parameters predicted; @@ -80,7 +80,7 @@ class GainMatrixUpdater { // shape later trackState.effectiveCalibrated().data(), trackState.effectiveCalibratedCovariance().data(), - trackState.boundSubspaceIndices(), + trackState.projector(), trackState.predicted(), trackState.predictedCovariance(), trackState.filtered(), diff --git a/Core/include/Acts/TrackFitting/detail/GainMatrixUpdaterImpl.hpp b/Core/include/Acts/TrackFitting/detail/GainMatrixUpdaterImpl.hpp index cc0f7b45596..80b7729b835 100644 --- a/Core/include/Acts/TrackFitting/detail/GainMatrixUpdaterImpl.hpp +++ b/Core/include/Acts/TrackFitting/detail/GainMatrixUpdaterImpl.hpp @@ -34,14 +34,9 @@ std::tuple GainMatrixUpdater::visitMeasurementImpl( ACTS_VERBOSE("Calibrated measurement: " << calibrated.transpose()); ACTS_VERBOSE("Calibrated measurement covariance:\n" << calibratedCovariance); - std::span validSubspaceIndices( - trackState.projector.begin(), - trackState.projector.begin() + kMeasurementSize); - FixedBoundSubspaceHelper subspaceHelper( - validSubspaceIndices); - - // TODO use subspace helper for projection instead - const auto H = subspaceHelper.projector(); + const auto H = trackState.projector + .template topLeftCorner() + .eval(); ACTS_VERBOSE("Measurement projector H:\n" << H); diff --git a/Core/include/Acts/Utilities/detail/Subspace.hpp b/Core/include/Acts/Utilities/detail/Subspace.hpp new file mode 100644 index 00000000000..fe355b44edb --- /dev/null +++ b/Core/include/Acts/Utilities/detail/Subspace.hpp @@ -0,0 +1,327 @@ +// This file is part of the Acts project. +// +// Copyright (C) 2020 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 +// file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#pragma once + +#include "Acts/Definitions/Algebra.hpp" + +#include +#include +#include +#include + +namespace Acts::detail { + +/// @defgroup subspace Linear subspace definitions +/// +/// All types in this group define a linear subspace of a larger vector space +/// with the following properties: the subspace is defined by a set of axis +/// indices in the full space, i.e. there is no need for a non-trivial +/// projection matrix, all subspace axis indices are unique and there are no +/// duplicated ones, and the set of axis indices must be ordered, i.e. the axis +/// order in the subspace is the same as in the full space. The last requirement +/// is not strictly required by the implementation, but is still added to +/// simplify reasoning. +/// +/// Only the size of the subspace are defined by the types at compile time. +/// Which axes comprise the subspace is defined at runtime. This allows to use +/// fixed-size computations (selected at compile time) but avoids the +/// combinatorial explosion of also having to handle all possible combination of +/// axes at compile time. This was tried previously and resulted in sizable +/// resource consumption at compile time without runtime benefits. +/// +/// All types are intentionally using `std::size_t` as their template values, +/// instead of the more specific index enums, to reduce the number of templates. +/// This is fully compatible as the index enums are required to be convertible +/// to `std::size_t`. +/// +/// All types intentionally only define the subspace but not how vectors +/// and matrices are stored to avoid unnecessary coupling between components, +/// i.e here between the pure definition and the storage. +/// +/// All types provide `.projectVector(...)` and `.exandVector(...)` methods to +/// convert to/from the subspace. They also provide `.projector()` and +/// `.expander()` methods to create projection and expansion matrices if they +/// are required explicitly. For the specific subspace requirements listed +/// above, the projection and expansion matrices are transpose to each other. In +/// the general case, this does not have to be the case and users are encouraged +/// to use `.projector()` and `.expander()` instead of e.g. +/// `.projector().transpose()`. +/// +/// @{ + +/// Fixed-size subspace representation. +/// +/// @tparam kFullSize Size of the full vector space +/// @tparam kSize Size of the subspace +template +class FixedSizeSubspace { + static_assert(kFullSize <= static_cast( + std::numeric_limits::max()), + "Full vector space size is larger than the supported range"); + static_assert(1u <= kSize, "Subspace size must be at least 1"); + static_assert(kSize <= kFullSize, + "Subspace can only be as large as the full space"); + + template + using SubspaceVectorFor = Eigen::Matrix; + template + using FullspaceVectorFor = + Eigen::Matrix; + template + using ProjectionMatrix = Eigen::Matrix; + template + using ExpansionMatrix = Eigen::Matrix; + + // the functionality could also be implemented using a std::bitset where each + // bit corresponds to an axis in the fullspace and set bits indicate which + // bits make up the subspace. this would be a more compact representation but + // complicates the implementation since we can not easily iterate over the + // indices of the subspace. storing the subspace indices directly requires a + // bit more memory but is easier to work with. for our typical use cases with + // n<=8, this still takes only 64bit of memory. + std::array m_axes; + + public: + /// Construct from a container of axis indices. + /// + /// @tparam index_t Input index type, must be convertible to std::uint8_t + /// @param indices Unique, ordered indices + template + constexpr explicit FixedSizeSubspace( + const std::array& indices) { + for (std::size_t i = 0u; i < kSize; ++i) { + assert((indices[i] < kFullSize) && + "Axis indices must be within the full space"); + if (0u < i) { + assert((indices[i - 1u] < indices[i]) && + "Axis indices must be unique and ordered"); + } + } + for (std::size_t i = 0; i < kSize; ++i) { + m_axes[i] = static_cast(indices[i]); + } + } + + /// Size of the subspace. + static constexpr std::size_t size() { return kSize; } + /// Size of the full vector space. + static constexpr std::size_t fullSize() { return kFullSize; } + + /// Access axis index by position. + /// + /// @param i Position in the subspace + /// @return Axis index in the full space + constexpr std::size_t operator[](std::size_t i) const { return m_axes[i]; } + + std::size_t indexOf(std::size_t axis) const { + for (std::size_t i = 0; i < kSize; ++i) { + if (m_axes[i] == axis) { + return i; + } + } + return kSize; + } + + /// Axis indices that comprise the subspace. + /// + /// The specific container and index type should be considered an + /// implementation detail. Users should treat the return type as a generic + /// container whose elements are convertible to `std::size_t`. + constexpr const std::array& indices() const { + return m_axes; + } + + /// Check if the given axis index in the full space is part of the subspace. + constexpr bool contains(std::size_t index) const { + bool isContained = false; + // always iterate over all elements to avoid branching and hope the compiler + // can optimise this for us. + for (std::uint8_t a : m_axes) { + isContained = isContained || (a == index); + } + return isContained; + } + + /// Project a full vector into the subspace. + /// + /// @tparam fullspace_vector_t Vector type in the full space + /// @param full Vector in the full space + /// @return Subspace vector w/ just the configured axis components + /// + /// @note Always returns a column vector regardless of the input + template + SubspaceVectorFor projectVector( + const Eigen::MatrixBase& full) const { + EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(fullspace_vector_t, kFullSize); + + SubspaceVectorFor sub; + for (std::size_t i = 0u; i < kSize; ++i) { + sub[i] = full[m_axes[i]]; + } + return sub; + } + + /// Expand a subspace vector into the full space. + /// + /// @tparam subspace_vector_t Subspace vector type + /// @param sub Subspace vector + /// @return Vector in the full space w/ zeros for non-subspace components + /// + /// @note Always returns a column vector regardless of the input + template + FullspaceVectorFor expandVector( + const Eigen::MatrixBase& sub) const { + EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(subspace_vector_t, kSize); + + FullspaceVectorFor full; + full.setZero(); + for (std::size_t i = 0u; i < kSize; ++i) { + full[m_axes[i]] = sub[i]; + } + return full; + } + + /// Projection matrix that maps from the full space into the subspace. + /// + /// @tparam scalar_t Scalar type for the projection matrix + template + ProjectionMatrix projector() const { + ProjectionMatrix proj; + proj.setZero(); + for (std::size_t i = 0u; i < kSize; ++i) { + proj(i, m_axes[i]) = 1; + } + return proj; + } + + /// Expansion matrix that maps from the subspace into the full space. + /// + /// @tparam scalar_t Scalar type of the generated expansion matrix + template + ExpansionMatrix expander() const { + ExpansionMatrix expn; + expn.setZero(); + for (std::size_t i = 0u; i < kSize; ++i) { + expn(m_axes[i], i) = 1; + } + return expn; + } +}; + +/// Variable-size subspace representation. +/// +/// @tparam kFullSize Size of the full vector space +/// @tparam kSize Size of the subspace +template +class VariableSizeSubspace { + static_assert(kFullSize <= static_cast(UINT8_MAX), + "Full vector space size is larger than the supported range"); + + template + using FullProjectionMatrix = Eigen::Matrix; + template + using FullExpansionMatrix = Eigen::Matrix; + + std::size_t m_size{}; + + // the functionality could also be implemented using a std::bitset where each + // bit corresponds to an axis in the fullspace and set bits indicate which + // bits make up the subspace. this would be a more compact representation but + // complicates the implementation since we can not easily iterate over the + // indices of the subspace. storing the subspace indices directly requires a + // bit more memory but is easier to work with. for our typical use cases with + // n<=8, this still takes only 64bit of memory. + std::array m_axes{}; + + public: + /// Construct from a container of axis indices. + /// + /// @tparam index_t Input index type, must be convertible to std::uint8_t + /// @param indices Unique, ordered indices + template + constexpr explicit VariableSizeSubspace( + const std::array& indices) { + m_size = kSize; + for (std::size_t i = 0u; i < kSize; ++i) { + assert((indices[i] < kFullSize) && + "Axis indices must be within the full space"); + if (0u < i) { + assert((indices[i - 1u] < indices[i]) && + "Axis indices must be unique and ordered"); + } + } + for (std::size_t i = 0; i < kSize; ++i) { + m_axes[i] = static_cast(indices[i]); + } + } + + /// Size of the subspace. + constexpr std::size_t size() const { return m_size; } + /// Size of the full vector space. + static constexpr std::size_t fullSize() { return kFullSize; } + + /// Access axis index by position. + /// + /// @param i Position in the subspace + /// @return Axis index in the full space + constexpr std::size_t operator[](std::size_t i) const { + assert(i < m_size); + return m_axes[i]; + } + + std::size_t indexOf(std::size_t axis) const { + for (std::size_t i = 0; i < m_size; ++i) { + if (m_axes[i] == axis) { + return i; + } + } + return m_size; + } + + /// Check if the given axis index in the full space is part of the subspace. + constexpr bool contains(std::size_t index) const { + bool isContained = false; + // always iterate over all elements to avoid branching and hope the compiler + // can optimise this for us. + for (std::size_t i = 0; i < kFullSize; ++i) { + isContained = isContained || ((i < m_size) && (m_axes[i] == index)); + } + return isContained; + } + + /// Projection matrix that maps from the full space into the subspace. + /// + /// @tparam scalar_t Scalar type for the projection matrix + template + FullProjectionMatrix fullProjector() const { + FullProjectionMatrix proj; + proj.setZero(); + for (std::size_t i = 0u; i < m_size; ++i) { + proj(i, m_axes[i]) = 1; + } + return proj; + } + + /// Expansion matrix that maps from the subspace into the full space. + /// + /// @tparam scalar_t Scalar type of the generated expansion matrix + template + FullExpansionMatrix fullExpander() const { + FullExpansionMatrix expn; + expn.setZero(); + for (std::size_t i = 0u; i < m_size; ++i) { + expn(m_axes[i], i) = 1; + } + return expn; + } +}; + +/// @} + +} // namespace Acts::detail diff --git a/Core/src/EventData/VectorMultiTrajectory.cpp b/Core/src/EventData/VectorMultiTrajectory.cpp index 23bc78550aa..e6c1fb56161 100644 --- a/Core/src/EventData/VectorMultiTrajectory.cpp +++ b/Core/src/EventData/VectorMultiTrajectory.cpp @@ -8,10 +8,8 @@ #include "Acts/EventData/VectorMultiTrajectory.hpp" -#include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/EventData/MultiTrajectory.hpp" #include "Acts/EventData/TrackStatePropMask.hpp" -#include "Acts/EventData/Types.hpp" #include "Acts/Utilities/Helpers.hpp" #include @@ -76,7 +74,7 @@ auto VectorMultiTrajectory::addTrackState_impl( m_sourceLinks.emplace_back(std::nullopt); p.icalibratedsourcelink = m_sourceLinks.size() - 1; - m_projectors.push_back(kBoundSubspaceIndicesInvalid); + m_projectors.emplace_back(); p.iprojector = m_projectors.size() - 1; } @@ -131,7 +129,7 @@ void VectorMultiTrajectory::addTrackStateComponents_impl( m_sourceLinks.emplace_back(std::nullopt); p.icalibratedsourcelink = m_sourceLinks.size() - 1; - m_projectors.push_back(kBoundSubspaceIndicesInvalid); + m_projectors.emplace_back(); p.iprojector = m_projectors.size() - 1; } diff --git a/Core/src/TrackFinding/MeasurementSelector.cpp b/Core/src/TrackFinding/MeasurementSelector.cpp index 278157e7682..946ca981949 100644 --- a/Core/src/TrackFinding/MeasurementSelector.cpp +++ b/Core/src/TrackFinding/MeasurementSelector.cpp @@ -9,10 +9,7 @@ #include "Acts/TrackFinding/MeasurementSelector.hpp" #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 @@ -33,7 +30,9 @@ double MeasurementSelector::calculateChi2( false>::Parameters predicted, TrackStateTraits::Covariance predictedCovariance, - BoundSubspaceIndices projector, unsigned int calibratedSize) const { + TrackStateTraits::Projector projector, + unsigned int calibratedSize) const { return visit_measurement( calibratedSize, [&fullCalibrated, &fullCalibratedCovariance, &predicted, @@ -48,19 +47,17 @@ double MeasurementSelector::calculateChi2( using ParametersVector = ActsVector; - std::span validSubspaceIndices( - projector.begin(), projector.begin() + kMeasurementSize); - FixedBoundSubspaceHelper subspaceHelper( - validSubspaceIndices); + // Take the projector (measurement mapping function) + const auto H = + projector.template topLeftCorner() + .eval(); // Get the residuals - ParametersVector res = - calibrated - subspaceHelper.projectVector(predicted); + ParametersVector res = calibrated - H * predicted; // Get the chi2 return (res.transpose() * - (calibratedCovariance + - subspaceHelper.projectMatrix(predictedCovariance)) + (calibratedCovariance + H * predictedCovariance * H.transpose()) .inverse() * res) .eval()(0, 0); diff --git a/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp b/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp index 12764a360d5..f458b33795a 100644 --- a/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp +++ b/Examples/Algorithms/TrackFinding/src/HoughTransformSeeder.cpp @@ -532,11 +532,11 @@ void ActsExamples::HoughTransformSeeder::addMeasurements( assert(measurement.contains(Acts::eBoundLoc1) && "Measurement does not contain the required bound loc1"); - auto boundLoc0 = measurement.indexOf(Acts::eBoundLoc0); - auto boundLoc1 = measurement.indexOf(Acts::eBoundLoc1); + auto boundLoc0 = measurement.subspace().indexOf(Acts::eBoundLoc0); + auto boundLoc1 = measurement.subspace().indexOf(Acts::eBoundLoc1); - Acts::Vector2 localPos{measurement.parameters()[boundLoc0], - measurement.parameters()[boundLoc1]}; + Acts::Vector2 localPos{measurement.effectiveParameters()[boundLoc0], + measurement.effectiveParameters()[boundLoc1]}; // transform local position to global coordinates Acts::Vector3 globalFakeMom(1, 1, 1); diff --git a/Examples/Algorithms/TrackFitting/src/RefittingCalibrator.cpp b/Examples/Algorithms/TrackFitting/src/RefittingCalibrator.cpp index 680d12fe044..f63a94df499 100644 --- a/Examples/Algorithms/TrackFitting/src/RefittingCalibrator.cpp +++ b/Examples/Algorithms/TrackFitting/src/RefittingCalibrator.cpp @@ -37,7 +37,7 @@ void RefittingCalibrator::calibrate(const Acts::GeometryContext& /*gctx*/, sl.state.template calibratedCovariance(); }); - trackState.setBoundSubspaceIndices(sl.state.boundSubspaceIndices()); + trackState.setProjectorBitset(sl.state.projectorBitset()); } } // namespace ActsExamples diff --git a/Examples/Framework/ML/src/NeuralCalibrator.cpp b/Examples/Framework/ML/src/NeuralCalibrator.cpp index 77af3920a6a..6dcbe1cc7cf 100644 --- a/Examples/Framework/ML/src/NeuralCalibrator.cpp +++ b/Examples/Framework/ML/src/NeuralCalibrator.cpp @@ -113,13 +113,14 @@ void ActsExamples::NeuralCalibrator::calibrate( assert(measurement.contains(Acts::eBoundLoc1) && "Measurement does not contain the required bound loc1"); - auto boundLoc0 = measurement.indexOf(Acts::eBoundLoc0); - auto boundLoc1 = measurement.indexOf(Acts::eBoundLoc1); + auto boundLoc0 = measurement.subspace().indexOf(Acts::eBoundLoc0); + auto boundLoc1 = measurement.subspace().indexOf(Acts::eBoundLoc1); - Acts::Vector2 localPosition{measurement.parameters()[boundLoc0], - measurement.parameters()[boundLoc1]}; - Acts::Vector2 localCovariance{measurement.covariance()(boundLoc0, boundLoc0), - measurement.covariance()(boundLoc1, boundLoc1)}; + Acts::Vector2 localPosition{measurement.effectiveParameters()[boundLoc0], + measurement.effectiveParameters()[boundLoc1]}; + Acts::Vector2 localCovariance{ + measurement.effectiveCovariance()(boundLoc0, boundLoc0), + measurement.effectiveCovariance()(boundLoc1, boundLoc1)}; Acts::Vector3 dir = Acts::makeDirectionFromPhiTheta( trackParameters[Acts::eBoundPhi], trackParameters[Acts::eBoundTheta]); @@ -172,10 +173,11 @@ void ActsExamples::NeuralCalibrator::calibrate( std::size_t iVar0 = 3 * m_nComponents + iMax * 2; Measurement measurementCopy = measurement; - measurementCopy.parameters()[boundLoc0] = output[iLoc0]; - measurementCopy.parameters()[boundLoc1] = output[iLoc0 + 1]; - measurementCopy.covariance()(boundLoc0, boundLoc0) = output[iVar0]; - measurementCopy.covariance()(boundLoc1, boundLoc1) = output[iVar0 + 1]; + measurementCopy.effectiveParameters()[boundLoc0] = output[iLoc0]; + measurementCopy.effectiveParameters()[boundLoc1] = output[iLoc0 + 1]; + measurementCopy.effectiveCovariance()(boundLoc0, boundLoc0) = output[iVar0]; + measurementCopy.effectiveCovariance()(boundLoc1, boundLoc1) = + output[iVar0 + 1]; Acts::visit_measurement(measurement.size(), [&](auto N) -> void { constexpr std::size_t kMeasurementSize = decltype(N)::value; @@ -185,7 +187,6 @@ void ActsExamples::NeuralCalibrator::calibrate( measurementCopy.parameters(); trackState.calibratedCovariance() = measurementCopy.covariance(); - trackState.setSubspaceIndices( - measurementCopy.subspaceIndices()); + trackState.setProjector(measurementCopy.subspace().fullProjector()); }); } diff --git a/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp b/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp index 87c4c546c08..ec423531cfa 100644 --- a/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp +++ b/Examples/Framework/include/ActsExamples/EventData/Measurement.hpp @@ -11,21 +11,17 @@ #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/EventData/MeasurementHelpers.hpp" #include "Acts/EventData/SourceLink.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/Utilities/detail/Subspace.hpp" #include #include #include -#include #include #include -#include - namespace ActsExamples { /// A measurement of a variable-size subspace of the full parameters. @@ -50,15 +46,13 @@ namespace ActsExamples { /// these complications altogether. template class VariableSizeMeasurement { - public: static constexpr std::size_t kFullSize = Acts::detail::kParametersSize; - using Scalar = Acts::ActsScalar; + using Subspace = Acts::detail::VariableSizeSubspace; - using SubspaceIndex = std::uint8_t; - using SubspaceIndices = - boost::container::static_vector; + public: + using Scalar = Acts::ActsScalar; /// Vector type containing for measured parameter values. template @@ -96,20 +90,18 @@ class VariableSizeMeasurement { /// @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 indices 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)) { + template + VariableSizeMeasurement(Acts::SourceLink source, + const std::array& indices, + const Eigen::MatrixBase& params, + const Eigen::MatrixBase& cov) + : m_source(std::move(source)), m_subspace(indices) { static_assert(kSize == parameters_t::RowsAtCompileTime, "Parameter size mismatch"); static_assert(kSize == covariance_t::RowsAtCompileTime, @@ -117,12 +109,6 @@ class VariableSizeMeasurement { 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; } @@ -134,79 +120,47 @@ class VariableSizeMeasurement { VariableSizeMeasurement& operator=(const VariableSizeMeasurement&) = default; VariableSizeMeasurement& operator=(VariableSizeMeasurement&&) = default; + constexpr std::size_t size() const { return m_subspace.size(); } + /// Source link that connects to the underlying detector readout. const Acts::SourceLink& sourceLink() const { return m_source; } - constexpr std::size_t size() const { return m_subspaceIndices.size(); } + const Subspace& subspace() const { return m_subspace; } /// Check if a specific parameter is part of this measurement. - bool contains(indices_t i) const { - return std::find(m_subspaceIndices.begin(), m_subspaceIndices.end(), i) != - m_subspaceIndices.end(); - } - - 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); - } - - /// The measurement indices - const SubspaceIndices& subspaceIndices() const { return m_subspaceIndices; } - - template - Acts::SubspaceIndices subspaceIndices() const { - assert(dim == size()); - Acts::SubspaceIndices result; - std::copy(m_subspaceIndices.begin(), m_subspaceIndices.end(), - result.begin()); - return result; - } - - Acts::BoundSubspaceIndices boundSubsetIndices() const - requires(std::is_same_v) - { - Acts::BoundSubspaceIndices result = Acts::kBoundSubspaceIndicesInvalid; - std::copy(m_subspaceIndices.begin(), m_subspaceIndices.end(), - result.begin()); - return result; - } + bool contains(indices_t i) const { return m_subspace.contains(i); } template ConstParametersVectorMap parameters() const { - assert(dim == size()); return ConstParametersVectorMap{m_params.data()}; } template ParametersVectorMap parameters() { - assert(dim == size()); return ParametersVectorMap{m_params.data()}; } - ConstEffectiveParametersVectorMap parameters() const { + ConstEffectiveParametersVectorMap effectiveParameters() const { return ConstEffectiveParametersVectorMap{m_params.data(), static_cast(size())}; } - EffectiveParametersVectorMap parameters() { + EffectiveParametersVectorMap effectiveParameters() { return EffectiveParametersVectorMap{m_params.data(), static_cast(size())}; } template ConstCovarianceMatrixMap covariance() const { - assert(dim == size()); return ConstCovarianceMatrixMap{m_cov.data()}; } template CovarianceMatrixMap covariance() { - assert(dim == size()); return CovarianceMatrixMap{m_cov.data()}; } - ConstEffectiveCovarianceMatrixMap covariance() const { + ConstEffectiveCovarianceMatrixMap effectiveCovariance() const { return ConstEffectiveCovarianceMatrixMap{m_cov.data(), static_cast(size()), static_cast(size())}; } - EffectiveCovarianceMatrixMap covariance() { + EffectiveCovarianceMatrixMap effectiveCovariance() { return EffectiveCovarianceMatrixMap{m_cov.data(), static_cast(size()), static_cast(size())}; @@ -215,7 +169,7 @@ class VariableSizeMeasurement { FullParametersVector fullParameters() const { FullParametersVector result = FullParametersVector::Zero(); for (std::size_t i = 0; i < size(); ++i) { - result[m_subspaceIndices[i]] = parameters()[i]; + result[m_subspace[i]] = effectiveParameters()[i]; } return result; } @@ -224,7 +178,7 @@ class VariableSizeMeasurement { 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); + result(m_subspace[i], m_subspace[j]) = effectiveCovariance()(i, j); } } return result; @@ -232,7 +186,7 @@ class VariableSizeMeasurement { private: Acts::SourceLink m_source; - SubspaceIndices m_subspaceIndices; + Subspace m_subspace; std::array m_params{}; std::array m_cov{}; }; diff --git a/Examples/Framework/src/EventData/MeasurementCalibration.cpp b/Examples/Framework/src/EventData/MeasurementCalibration.cpp index 560593dbd1f..444021c91ca 100644 --- a/Examples/Framework/src/EventData/MeasurementCalibration.cpp +++ b/Examples/Framework/src/EventData/MeasurementCalibration.cpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2023-2024 CERN for the benefit of the Acts project +// Copyright (C) 2023 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 @@ -41,8 +41,7 @@ void ActsExamples::PassThroughCalibrator::calibrate( measurement.parameters(); trackState.calibratedCovariance() = measurement.covariance(); - trackState.setSubspaceIndices( - measurement.subspaceIndices()); + trackState.setProjector(measurement.subspace().fullProjector()); }); } diff --git a/Examples/Framework/src/EventData/ScalingCalibrator.cpp b/Examples/Framework/src/EventData/ScalingCalibrator.cpp index 805d8cdef15..fd93c584c51 100644 --- a/Examples/Framework/src/EventData/ScalingCalibrator.cpp +++ b/Examples/Framework/src/EventData/ScalingCalibrator.cpp @@ -158,24 +158,23 @@ void ActsExamples::ScalingCalibrator::calibrate( assert(measurement.contains(Acts::eBoundLoc1) && "Measurement does not contain the required bound loc1"); - auto boundLoc0 = measurement.indexOf(Acts::eBoundLoc0); - auto boundLoc1 = measurement.indexOf(Acts::eBoundLoc1); + auto boundLoc0 = measurement.subspace().indexOf(Acts::eBoundLoc0); + auto boundLoc1 = measurement.subspace().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; + measurementCopy.effectiveParameters()[boundLoc0] += ct.x_offset; + measurementCopy.effectiveParameters()[boundLoc1] += ct.y_offset; + measurementCopy.effectiveCovariance()(boundLoc0, boundLoc0) *= ct.x_scale; + measurementCopy.effectiveCovariance()(boundLoc1, boundLoc1) *= ct.y_scale; Acts::visit_measurement(measurement.size(), [&](auto N) -> void { constexpr std::size_t kMeasurementSize = decltype(N)::value; trackState.allocateCalibrated(kMeasurementSize); trackState.calibrated() = - measurementCopy.parameters(); + measurement.parameters(); trackState.calibratedCovariance() = measurementCopy.covariance(); - trackState.setSubspaceIndices( - measurementCopy.subspaceIndices()); + trackState.setProjector(measurement.subspace().fullProjector()); }); } diff --git a/Examples/Io/Root/src/RootMeasurementWriter.cpp b/Examples/Io/Root/src/RootMeasurementWriter.cpp index 6dc10e3b4b1..0a476b031a1 100644 --- a/Examples/Io/Root/src/RootMeasurementWriter.cpp +++ b/Examples/Io/Root/src/RootMeasurementWriter.cpp @@ -154,10 +154,10 @@ struct RootMeasurementWriter::DigitizationTree { /// @param m The measurement void fillBoundMeasurement(const Measurement& m) { for (unsigned int i = 0; i < m.size(); ++i) { - auto ib = m.subspaceIndices()[i]; + auto ib = m.subspace()[i]; - recBound[ib] = m.parameters()[i]; - varBound[ib] = m.covariance()(i, i); + recBound[ib] = m.effectiveParameters()[i]; + varBound[ib] = m.effectiveCovariance()(i, i); residual[ib] = recBound[ib] - trueBound[ib]; pull[ib] = residual[ib] / std::sqrt(varBound[ib]); diff --git a/Examples/Python/tests/root_file_hashes.txt b/Examples/Python/tests/root_file_hashes.txt index f33c396938f..cf99adafefc 100644 --- a/Examples/Python/tests/root_file_hashes.txt +++ b/Examples/Python/tests/root_file_hashes.txt @@ -19,8 +19,8 @@ test_propagation__propagation_summary.root: 85d877c3a48a6235d3d5a398b47130d5d17d test_material_recording__geant4_material_tracks.root: c022b9362249b29f57a07926b20644e3ab4ab8ebcf03f773fbf46c446fc1a0a1 test_truth_tracking_gsf[generic]__trackstates_gsf.root: 1466671bee7b7cfffe90459d9aef316dcee104e6c32a2df4f2f11ea9f12ddc14 test_truth_tracking_gsf[generic]__tracksummary_gsf.root: b698e3d21eacc34fc8b0ce1d3fbe07405a4b8b549e07f0160573e64c3b401f04 -test_truth_tracking_gsf[odd]__trackstates_gsf.root: 847a60b97ef9c9328a73f27597085fea2434003e2c15d30f98ae7b3372eeb65f -test_truth_tracking_gsf[odd]__tracksummary_gsf.root: d8221a90ed24f9f136bf2001645ce50f3eacf18c73d7d515eec71962a555d14b +test_truth_tracking_gsf[odd]__trackstates_gsf.root: 194198c9c9657b02192e573ba9969230663f179104e32e180101ffcd677b4745 +test_truth_tracking_gsf[odd]__tracksummary_gsf.root: 80e60784d4d6eabe734f6fc1ffe04aefa0276cf8798430c590a06f40a07a484f test_particle_gun__particles.root: 5fe7dda2933ee6b9615b064d192322fe07831133cd998e5ed99a3b992b713a10 test_material_mapping__material-map_tracks.root: 938b1a855369e9304401cb10d2751df3fd7acf32a2573f2057eb1691cd94edf3 test_material_mapping__propagation-material.root: 5eacd0cb804d381171c8fb65d7b2d36e1d57db9f4cb2b58c0f24703479218cbf @@ -40,14 +40,14 @@ test_ckf_tracks_example[generic-truth_estimated]__tracksummary_ckf.root: 0df341f test_ckf_tracks_example[generic-truth_estimated]__performance_seeding.root: 1facb05c066221f6361b61f015cdf0918e94d9f3fce2269ec7b6a4dffeb2bc7e test_ckf_tracks_example[generic-truth_smeared]__trackstates_ckf.root: 25b985218418c91782696daab1b3fd3eaf7d15cd05681c7dfe2e4c3e44977e16 test_ckf_tracks_example[generic-truth_smeared]__tracksummary_ckf.root: be592112183bd6554ed5e00dcb62eb4a5a5f4a9bb7efa4b1e1ce699dec13157b -test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 8c31d51cef775b4f92f99644895beecfd854db2e7925d5f394deafdcc8f2e629 -test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: 3b27644bf25c20e5ed0d25f2273b68a470c63b07677a806cb176d4c91abd153e +test_ckf_tracks_example[odd-full_seeding]__trackstates_ckf.root: 935359233aa503a077faf2394ba988d0da97ea5db4f6d3294a95fb93c5625c68 +test_ckf_tracks_example[odd-full_seeding]__tracksummary_ckf.root: 59537c52d23a77dbed9230c7b92cf653f7592c72fb882e8a24a2b33f1d27019e test_ckf_tracks_example[odd-full_seeding]__performance_seeding_trees.root: 43c58577aafe07645e5660c4f43904efadf91d8cda45c5c04c248bbe0f59814f -test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: d980eebb7df9d7fd5a861228ca498f439bef4853a8dd9b8b590286781b69b8cb -test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: a5eecfb6907406286c56db4df19937539cb43e454307bc53bd1eda535526736c +test_ckf_tracks_example[odd-truth_estimated]__trackstates_ckf.root: 89fad58add1ead7c0d38a952391d3c4a720feda6e37c37ed53f994a99fe15f57 +test_ckf_tracks_example[odd-truth_estimated]__tracksummary_ckf.root: 3a092423e717b1fe76f78f6b95eb6df67648a196e8c4da32158061d4bc490578 test_ckf_tracks_example[odd-truth_estimated]__performance_seeding.root: 1a36b7017e59f1c08602ef3c2cb0483c51df248f112e3780c66594110719c575 -test_ckf_tracks_example[odd-truth_smeared]__trackstates_ckf.root: 92a1367d2c5d90ead26d177206d9f644e285e42de5dea2615c7bc8aeb321419b -test_ckf_tracks_example[odd-truth_smeared]__tracksummary_ckf.root: 3dd1d9b1075a9c2b2ac7f4143329e13b6cb63f6e1dfe450f8a913a0d8110aa72 +test_ckf_tracks_example[odd-truth_smeared]__trackstates_ckf.root: 268357818ee53cc13b26aa23af06d7f839d112e9dc2eb3ac4c5b6e6eda48e6af +test_ckf_tracks_example[odd-truth_smeared]__tracksummary_ckf.root: 98c9e97240f541c22dd41d05c2114d7e9455b117162d64baad2c3ecb5b7d52f7 test_vertex_fitting_reading[Truth-False-100]__performance_vertexing.root: 76ef6084d758dfdfc0151ddec2170e12d73394424e3dac4ffe46f0f339ec8293 test_vertex_fitting_reading[Iterative-False-100]__performance_vertexing.root: 60372210c830a04f95ceb78c6c68a9b0de217746ff59e8e73053750c837b57eb test_vertex_fitting_reading[Iterative-True-100]__performance_vertexing.root: e34f217d524a5051dbb04a811d3407df3ebe2cc4bb7f54f6bda0847dbd7b52c3 @@ -85,11 +85,11 @@ test_truth_tracking_kalman[generic-True-0.0]__trackstates_kf.root: a87b908ee14a3 test_truth_tracking_kalman[generic-True-0.0]__tracksummary_kf.root: bf46a89e429fa77a380d5ee48babb8af2044196eff872825e84c0d4401357114 test_truth_tracking_kalman[generic-True-1000.0]__trackstates_kf.root: 5b945782c1ff59c7c30d3d2a23992292735bc3da27da6424309da9c9e9660322 test_truth_tracking_kalman[generic-True-1000.0]__tracksummary_kf.root: 055a74d2747d360398cc846cc2791204491528896de78cca66b188e3ff530dc1 -test_truth_tracking_kalman[odd-False-0.0]__trackstates_kf.root: 3f442d5f46cde700ced08ce1720b6b843ca84d7f0cdc372d8571f06e005144c2 -test_truth_tracking_kalman[odd-False-0.0]__tracksummary_kf.root: 7f1bb68b39e52da7a77ea295819fa3776e947568455338738c3a6436c2d7599c -test_truth_tracking_kalman[odd-False-1000.0]__trackstates_kf.root: 12b07b688b7b4cbfc1aa51668c6b693cfbae565c9116f4ad084241805384fbd1 -test_truth_tracking_kalman[odd-False-1000.0]__tracksummary_kf.root: 87eaae3192ab29e2c2542c017071b6477c7237c5b8eaff107e84caed2a5e5b7a -test_truth_tracking_kalman[odd-True-0.0]__trackstates_kf.root: 3f442d5f46cde700ced08ce1720b6b843ca84d7f0cdc372d8571f06e005144c2 -test_truth_tracking_kalman[odd-True-0.0]__tracksummary_kf.root: 7f1bb68b39e52da7a77ea295819fa3776e947568455338738c3a6436c2d7599c -test_truth_tracking_kalman[odd-True-1000.0]__trackstates_kf.root: 12b07b688b7b4cbfc1aa51668c6b693cfbae565c9116f4ad084241805384fbd1 -test_truth_tracking_kalman[odd-True-1000.0]__tracksummary_kf.root: 87eaae3192ab29e2c2542c017071b6477c7237c5b8eaff107e84caed2a5e5b7a +test_truth_tracking_kalman[odd-False-0.0]__trackstates_kf.root: 6e4a34638a4d76607597db1e02525f6c6a386ca44727cae255a0aa0c2b588c7f +test_truth_tracking_kalman[odd-False-0.0]__tracksummary_kf.root: a45dbd0b6d221ec1b153f6e90d84601fe7a799134604c63be8e94bf8cd22af43 +test_truth_tracking_kalman[odd-False-1000.0]__trackstates_kf.root: 7e0b19d1f8c818269925b95817b9c92c68af3b841fdbd21560aa55688ddec786 +test_truth_tracking_kalman[odd-False-1000.0]__tracksummary_kf.root: 15a542955428dec0776c8317a33b7ef6c14858888ffc7406b00ac7044694b137 +test_truth_tracking_kalman[odd-True-0.0]__trackstates_kf.root: 6e4a34638a4d76607597db1e02525f6c6a386ca44727cae255a0aa0c2b588c7f +test_truth_tracking_kalman[odd-True-0.0]__tracksummary_kf.root: a45dbd0b6d221ec1b153f6e90d84601fe7a799134604c63be8e94bf8cd22af43 +test_truth_tracking_kalman[odd-True-1000.0]__trackstates_kf.root: 7e0b19d1f8c818269925b95817b9c92c68af3b841fdbd21560aa55688ddec786 +test_truth_tracking_kalman[odd-True-1000.0]__tracksummary_kf.root: 15a542955428dec0776c8317a33b7ef6c14858888ffc7406b00ac7044694b137 diff --git a/Plugins/Podio/edm.yml b/Plugins/Podio/edm.yml index ccbb65a7e94..0b5996c37c7 100644 --- a/Plugins/Podio/edm.yml +++ b/Plugins/Podio/edm.yml @@ -58,7 +58,7 @@ components: - uint32_t measdim; - uint32_t ijacobian - - std::array projector + - uint64_t projector - bool hasProjector - float chi2 diff --git a/Plugins/Podio/include/Acts/Plugins/Podio/PodioTrackStateContainer.hpp b/Plugins/Podio/include/Acts/Plugins/Podio/PodioTrackStateContainer.hpp index 0a45c73d3cd..60951942618 100644 --- a/Plugins/Podio/include/Acts/Plugins/Podio/PodioTrackStateContainer.hpp +++ b/Plugins/Podio/include/Acts/Plugins/Podio/PodioTrackStateContainer.hpp @@ -1,6 +1,6 @@ // This file is part of the Acts project. // -// Copyright (C) 2023-2024 CERN for the benefit of the Acts project +// Copyright (C) 2023 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 @@ -124,10 +124,7 @@ class PodioTrackStateContainerBase { case "smoothed"_hash: return &data.ismoothed; case "projector"_hash: - // workaround podio not allowing `std::uint8_t` as a type - return reinterpret_cast>( - &data.projector); + return &data.projector; case "measdim"_hash: return &data.measdim; case "chi2"_hash: diff --git a/Tests/UnitTests/Core/TrackFitting/MbfSmootherTests.cpp b/Tests/UnitTests/Core/TrackFitting/MbfSmootherTests.cpp index 0de58b88305..bdd6d3acab2 100644 --- a/Tests/UnitTests/Core/TrackFitting/MbfSmootherTests.cpp +++ b/Tests/UnitTests/Core/TrackFitting/MbfSmootherTests.cpp @@ -18,6 +18,7 @@ #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp" #include "Acts/TrackFitting/MbfSmoother.hpp" #include "Acts/Utilities/Result.hpp" +#include "Acts/Utilities/detail/Subspace.hpp" #include #include @@ -40,7 +41,9 @@ BOOST_AUTO_TEST_SUITE(TrackFittingMbfSmoother) BOOST_AUTO_TEST_CASE(Smooth) { VectorMultiTrajectory traj; - std::array projector{eBoundLoc0, eBoundLoc1}; + auto projector = detail::FixedSizeSubspace( + std::array{eBoundLoc0, eBoundLoc1}) + .projector(); // Make dummy track parameter CovarianceMatrix covTrk; @@ -55,10 +58,10 @@ BOOST_AUTO_TEST_CASE(Smooth) { ts.predicted() << 0.3, 0.5, 0.5 * M_PI, 0., 1 / 100., 0.; ts.predictedCovariance() = covTrk; + ts.setProjector(projector); ts.allocateCalibrated(2); ts.calibrated<2>() << 0.351, 0.473; ts.calibratedCovariance<2>() << 1e+8, 0., 0., 1e+8; - ts.setSubspaceIndices<2>(projector); ts.filtered() << 0.301, 0.503, 0.5 * M_PI, 0., 1 / 100., 0.; ts.filteredCovariance() = covTrk; @@ -72,10 +75,10 @@ BOOST_AUTO_TEST_CASE(Smooth) { ts.predicted() << 0.2, 0.5, 0.5 * M_PI, 0., 1 / 100., 0.; ts.predictedCovariance() = covTrk; + ts.setProjector(projector); ts.allocateCalibrated(2); ts.calibrated<2>() << 0.351, 0.473; ts.calibratedCovariance<2>() << 1e+8, 0., 0., 1e+8; - ts.setSubspaceIndices<2>(projector); ts.filtered() << 0.27, 0.53, 0.5 * M_PI, 0., 1 / 100., 0.; ts.filteredCovariance() = covTrk; @@ -89,10 +92,10 @@ BOOST_AUTO_TEST_CASE(Smooth) { ts.predicted() << 0.35, 0.49, 0.5 * M_PI, 0., 1 / 100., 0.; ts.predictedCovariance() = covTrk; + ts.setProjector(projector); ts.allocateCalibrated(2); ts.calibrated<2>() << 0.351, 0.473; ts.calibratedCovariance<2>() << 1e+8, 0., 0., 1e+8; - ts.setSubspaceIndices<2>(projector); ts.filtered() << 0.33, 0.43, 0.5 * M_PI, 0., 1 / 100., 0.; ts.filteredCovariance() = covTrk; diff --git a/Tests/UnitTests/Core/Utilities/CMakeLists.txt b/Tests/UnitTests/Core/Utilities/CMakeLists.txt index 06b4632aa97..30bfda33fc3 100644 --- a/Tests/UnitTests/Core/Utilities/CMakeLists.txt +++ b/Tests/UnitTests/Core/Utilities/CMakeLists.txt @@ -33,6 +33,7 @@ add_unittest(RangeXD RangeXDTests.cpp) add_unittest(Ray RayTest.cpp) add_unittest(RealQuadraticEquation RealQuadraticEquationTests.cpp) add_unittest(Result ResultTests.cpp) +add_unittest(Subspace SubspaceTests.cpp) add_unittest(TypeList TypeListTests.cpp) add_unittest(TypeTraits TypeTraitsTest.cpp) add_unittest(UnitVectors UnitVectorsTests.cpp) diff --git a/Tests/UnitTests/Core/Utilities/SubspaceTests.cpp b/Tests/UnitTests/Core/Utilities/SubspaceTests.cpp new file mode 100644 index 00000000000..2fbf0d2ffc4 --- /dev/null +++ b/Tests/UnitTests/Core/Utilities/SubspaceTests.cpp @@ -0,0 +1,182 @@ +// This file is part of the Acts project. +// +// Copyright (C) 2020 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 +// file, You can obtain one at http://mozilla.org/MPL/2.0/. + +#include + +#include "Acts/Definitions/Algebra.hpp" +#include "Acts/Definitions/TrackParametrization.hpp" +#include "Acts/Utilities/AlgebraHelpers.hpp" +#include "Acts/Utilities/detail/Subspace.hpp" + +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace { + +using namespace Acts; + +// meta-programming type list of scalar type + subspace type combinations +// clang-format off +using ScalarsAndFixedSizeSubspaces = std::tuple< + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple>, + std::tuple> +>; +// clang-format on + +/// Construct a random vector of the specified size. +template +Eigen::Matrix makeRandomVector() { + Eigen::Matrix vec; + vec.setRandom(); + return vec; +} + +/// Construct a vector w/ monotonically inceasing values starting at 0. +std::vector makeMonotonicIndices(std::size_t n) { + std::vector indices(n); + std::iota(indices.begin(), indices.end(), 0u); + return indices; +} + +/// Build a sorted array from the first kSize indices. +template +std::array selectFixedIndices( + const std::vector& fullIndices) { + std::array indices{}; + for (auto i = 0u; i < kSize; ++i) { + indices[i] = fullIndices[i]; + } + std::sort(indices.begin(), indices.end()); + return indices; +} + +} // namespace + +BOOST_AUTO_TEST_SUITE(UtilitiesSubspace) + +// all these cases should lead to compile-time errors +// BOOST_AUTO_TEST_CASE(ConstructFixedInvalid) { +// { +// using Subspace = detail::FixedSizeSubspace<6u, 7u>; +// Subspace subspace(0u, 1u, 2u, 3u, 4u, 5u, 6u); +// } +// { +// using Subspace = detail::FixedSizeSubspace<8u, 9u>; +// Subspace subspace(0u, 1u, 2u, 3u, 4u, 5u, 6u, 7u, 8u); +// } +// } + +BOOST_AUTO_TEST_CASE_TEMPLATE(FixedSizeSubspace, ScalarAndSubspace, + ScalarsAndFixedSizeSubspaces) { + // extract the test types + using Scalar = std::tuple_element_t<0, ScalarAndSubspace>; + using Subspace = std::tuple_element_t<1, ScalarAndSubspace>; + + auto x = makeRandomVector(); + auto fullIndices = makeMonotonicIndices(Subspace::fullSize()); + + // in principle, we would like to iterate over all possible ordered subsets + // from the full space indices with a size identical to the subspace. since i + // do not know how to do that in a simple manner, we are iterating over all + // permutations of the full indices and pick the first n elements. this should + // give a reasonable set of different subspace configurations. + do { + auto indices = selectFixedIndices(fullIndices); + Subspace subspace(indices); + + // verify projector/expander consistency + BOOST_CHECK_EQUAL(subspace.template projector().transpose(), + subspace.template expander()); + BOOST_CHECK_EQUAL(subspace.template expander().transpose(), + subspace.template projector()); + // project into the subspace + auto s0 = subspace.projectVector(x); + auto s1 = (subspace.template projector() * x).eval(); + for (auto i = 0u; i < subspace.size(); ++i) { + BOOST_TEST_INFO("Checking projected subspace component " << i); + BOOST_CHECK_EQUAL(s0[i], x[indices[i]]); + BOOST_CHECK_EQUAL(s1[i], x[indices[i]]); + } + // expand from the subspace back into the full space + auto y0 = subspace.expandVector(s1); + auto y1 = (subspace.template expander() * s0).eval(); + for (auto i = 0u; i < subspace.fullSize(); ++i) { + BOOST_TEST_INFO("Checking expanded fullspace component " << i); + BOOST_CHECK_EQUAL(y0[i], subspace.contains(i) ? x[i] : 0); + BOOST_CHECK_EQUAL(y1[i], subspace.contains(i) ? x[i] : 0); + } + } while (std::next_permutation(fullIndices.begin(), fullIndices.end())); +} + +BOOST_AUTO_TEST_CASE_TEMPLATE(VariableSizeSubspace, ScalarAndSubspace, + ScalarsAndFixedSizeSubspaces) { + // extract the test types + using Scalar = std::tuple_element_t<0, ScalarAndSubspace>; + using FixedSubspace = std::tuple_element_t<1, ScalarAndSubspace>; + using VariableSubspace = + detail::VariableSizeSubspace; + + auto fullIndices = makeMonotonicIndices(FixedSubspace::fullSize()); + + // in principle, we would like to iterate over all possible ordered subsets + // from the full space indices with a size identical to the subspace. since i + // do not know how to do that in a simple manner, we are iterating over all + // permutations of the full indices and pick the first n elements. this should + // give a reasonable set of different subspace configurations. + do { + auto indices = selectFixedIndices(fullIndices); + FixedSubspace fixedSubspace(indices); + VariableSubspace variableSubspace(indices); + + BOOST_CHECK_EQUAL(variableSubspace.size(), fixedSubspace.size()); + BOOST_CHECK_EQUAL(variableSubspace.fullSize(), fixedSubspace.fullSize()); + + auto fixedProjector = fixedSubspace.template projector(); + + Eigen::Matrix + fixedFullProjector; + fixedFullProjector.setZero(); + fixedFullProjector.template topLeftCorner() = + fixedProjector; + + auto variableFullProjector = + variableSubspace.template fullProjector(); + + BOOST_CHECK_EQUAL(variableFullProjector, fixedFullProjector); + } while (std::next_permutation(fullIndices.begin(), fullIndices.end())); +} + +BOOST_AUTO_TEST_SUITE_END() diff --git a/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp b/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp index bab0b16e9e9..cd9e8904a8b 100644 --- a/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp +++ b/Tests/UnitTests/Examples/EventData/MeasurementTests.cpp @@ -61,8 +61,8 @@ BOOST_DATA_TEST_CASE(VariableBoundOne, bd::make(boundIndices), index) { BOOST_CHECK(!meas.contains(i)); } } - BOOST_CHECK_EQUAL(meas.parameters(), params); - BOOST_CHECK_EQUAL(meas.covariance(), cov); + BOOST_CHECK_EQUAL(meas.effectiveParameters(), params); + BOOST_CHECK_EQUAL(meas.effectiveCovariance(), cov); BOOST_CHECK_EQUAL(meas.sourceLink().template get(), sourceOrig); } @@ -77,8 +77,8 @@ BOOST_AUTO_TEST_CASE(VariableBoundAll) { for (auto i : boundIndices) { BOOST_CHECK(meas.contains(i)); } - BOOST_CHECK_EQUAL(meas.parameters(), params); - BOOST_CHECK_EQUAL(meas.covariance(), cov); + BOOST_CHECK_EQUAL(meas.effectiveParameters(), params); + BOOST_CHECK_EQUAL(meas.effectiveCovariance(), cov); BOOST_CHECK_EQUAL(meas.sourceLink().get(), sourceOrig); } diff --git a/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp b/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp index 4b4b6d78bbe..c9e9fa4126e 100644 --- a/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp +++ b/Tests/UnitTests/Examples/Io/Csv/MeasurementReaderWriterTests.cpp @@ -135,7 +135,7 @@ BOOST_AUTO_TEST_CASE(CsvMeasurementRoundTrip) { BOOST_REQUIRE(measRead.size() == measOriginal.size()); for (const auto &[a, b] : Acts::zip(measRead, measOriginal)) { if (a.size() == b.size()) { - CHECK_CLOSE_REL(a.parameters(), b.parameters(), 1e-4); + CHECK_CLOSE_REL(a.effectiveParameters(), b.effectiveParameters(), 1e-4); } }