ActsExamples/EventData/Measurement.hpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2020-2024 CERN for the benefit of the Acts project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include "Acts/Definitions/TrackParametrization.hpp"
0012 #include "Acts/EventData/MeasurementHelpers.hpp"
0013 #include "Acts/EventData/SourceLink.hpp"
0014 #include "Acts/EventData/SubspaceHelpers.hpp"
0015 #include "Acts/EventData/Types.hpp"
0016 #include "Acts/EventData/detail/CalculateResiduals.hpp"
0017 #include "Acts/EventData/detail/ParameterTraits.hpp"
0018 #include "Acts/EventData/detail/PrintParameters.hpp"
0019
0020 #include <array>
0021 #include <cstddef>
0022 #include <iosfwd>
0023 #include <type_traits>
0024 #include <variant>
0025 #include <vector>
0026
0027 #include <boost/container/static_vector.hpp>
0028
0029 namespace ActsExamples {
0030
0031 /// A measurement of a variable-size subspace of the full parameters.
0032 ///
0033 /// @tparam indices_t Parameter index type, determines the full parameter space
0034 ///
0035 /// The measurement intentionally does not store a pointer/reference to the
0036 /// reference object in the geometry hierarchy, i.e. the surface or volume. The
0037 /// reference object can already be identified via the geometry identifier
0038 /// provided by the source link. Since a measurement **must** be anchored within
0039 /// the geometry hierarchy, all measurement surfaces and volumes **must**
0040 /// provide valid geometry identifiers. In all use-cases, e.g. Kalman filtering,
0041 /// a pointer/reference to the reference object is available before the
0042 /// measurement is accessed; e.g. the propagator provides the surface pointer
0043 /// during navigation, which is then used to lookup possible measurements.
0044 ///
0045 /// The pointed-to geometry object would differ depending on the parameter type.
0046 /// This means either, that there needs to be an additional variable type or
0047 /// that a pointer to a base object is stored (requiring a `dynamic_cast` later
0048 /// on). Both variants add additional complications. Since the geometry object
0049 /// is not required anyway (as discussed above), not storing it removes all
0050 /// these complications altogether.
0051 template <typename indices_t>
0052 class VariableSizeMeasurement {
0053  public:
0054   static constexpr std::size_t kFullSize =
0055       Acts::detail::kParametersSize<indices_t>;
0056
0057   using Scalar = Acts::ActsScalar;
0058
0059   using SubspaceIndex = std::uint8_t;
0060   using SubspaceIndices =
0061       boost::container::static_vector<SubspaceIndex, kFullSize>;
0062
0063   /// Vector type containing for measured parameter values.
0064   template <std::size_t dim>
0065   using ParametersVector = Eigen::Matrix<Scalar, dim, 1>;
0066   template <std::size_t dim>
0067   using ParametersVectorMap = Eigen::Map<ParametersVector<dim>>;
0068   template <std::size_t dim>
0069   using ConstParametersVectorMap = Eigen::Map<const ParametersVector<dim>>;
0070   using EffectiveParametersVector = Eigen::Matrix<Scalar, Eigen::Dynamic, 1>;
0071   using EffectiveParametersVectorMap = Eigen::Map<EffectiveParametersVector>;
0072   using ConstEffectiveParametersVectorMap =
0073       Eigen::Map<const EffectiveParametersVector>;
0074
0075   /// Matrix type for the measurement covariance.
0076   template <std::size_t dim>
0077   using CovarianceMatrix = Eigen::Matrix<Scalar, dim, dim>;
0078   template <std::size_t dim>
0079   using CovarianceMatrixMap = Eigen::Map<CovarianceMatrix<dim>>;
0080   template <std::size_t dim>
0081   using ConstCovarianceMatrixMap = Eigen::Map<const CovarianceMatrix<dim>>;
0082   using EffectiveCovarianceMatrix =
0083       Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>;
0084   using EffectiveCovarianceMatrixMap = Eigen::Map<EffectiveCovarianceMatrix>;
0085   using ConstEffectiveCovarianceMatrixMap =
0086       Eigen::Map<const EffectiveCovarianceMatrix>;
0087
0088   using FullParametersVector = Acts::ActsVector<kFullSize>;
0089   using FullCovarianceMatrix = Acts::ActsSquareMatrix<kFullSize>;
0090
0091   using ProjectionMatrix = Eigen::Matrix<Scalar, Eigen::Dynamic, kFullSize>;
0092   using ExpansionMatrix = Eigen::Matrix<Scalar, kFullSize, Eigen::Dynamic>;
0093
0094   /// Construct from source link, subset indices, and measured data.
0095   ///
0096   /// @tparam parameters_t Input parameters vector type
0097   /// @tparam covariance_t Input covariance matrix type
0098   /// @param source The link that connects to the underlying detector readout
0099   /// @param subspaceIndices Which parameters are measured
0100   /// @param params Measured parameters values
0101   /// @param cov Measured parameters covariance
0102   ///
0103   /// @note The indices must be ordered and must describe/match the content
0104   ///   of parameters and covariance.
0105   template <typename other_indices_t, std::size_t kSize, typename parameters_t,
0106             typename covariance_t>
0107   VariableSizeMeasurement(
0108       Acts::SourceLink source,
0109       const std::array<other_indices_t, kSize>& subspaceIndices,
0110       const Eigen::MatrixBase<parameters_t>& params,
0111       const Eigen::MatrixBase<covariance_t>& cov)
0112       : m_source(std::move(source)) {
0113     static_assert(kSize == parameters_t::RowsAtCompileTime,
0114                   "Parameter size mismatch");
0115     static_assert(kSize == covariance_t::RowsAtCompileTime,
0116                   "Covariance rows mismatch");
0117     static_assert(kSize == covariance_t::ColsAtCompileTime,
0118                   "Covariance cols mismatch");
0119
0120     m_subspaceIndices.resize(subspaceIndices.size());
0121     std::transform(subspaceIndices.begin(), subspaceIndices.end(),
0122                    m_subspaceIndices.begin(), [](auto index) {
0123                      return static_cast<SubspaceIndex>(index);
0124                    });
0125
0126     parameters<kSize>() = params;
0127     covariance<kSize>() = cov;
0128   }
0129   /// A measurement can only be constructed with valid parameters.
0130   VariableSizeMeasurement() = delete;
0131   VariableSizeMeasurement(const VariableSizeMeasurement&) = default;
0132   VariableSizeMeasurement(VariableSizeMeasurement&&) = default;
0133   ~VariableSizeMeasurement() = default;
0134   VariableSizeMeasurement& operator=(const VariableSizeMeasurement&) = default;
0135   VariableSizeMeasurement& operator=(VariableSizeMeasurement&&) = default;
0136
0137   /// Source link that connects to the underlying detector readout.
0138   const Acts::SourceLink& sourceLink() const { return m_source; }
0139
0140   constexpr std::size_t size() const { return m_subspaceIndices.size(); }
0141
0142   /// Check if a specific parameter is part of this measurement.
0143   bool contains(indices_t i) const {
0144     return std::find(m_subspaceIndices.begin(), m_subspaceIndices.end(), i) !=
0145            m_subspaceIndices.end();
0146   }
0147
0148   std::size_t indexOf(indices_t i) const {
0149     auto it = std::find(m_subspaceIndices.begin(), m_subspaceIndices.end(), i);
0150     assert(it != m_subspaceIndices.end());
0151     return std::distance(m_subspaceIndices.begin(), it);
0152   }
0153
0154   /// The measurement indices
0155   const SubspaceIndices& subspaceIndices() const { return m_subspaceIndices; }
0156
0157   template <std::size_t dim>
0158   Acts::SubspaceIndices<dim> subspaceIndices() const {
0159     assert(dim == size());
0160     Acts::SubspaceIndices<dim> result;
0161     std::copy(m_subspaceIndices.begin(), m_subspaceIndices.end(),
0162               result.begin());
0163     return result;
0164   }
0165
0166   Acts::BoundSubspaceIndices boundSubsetIndices() const
0167     requires(std::is_same_v<indices_t, Acts::BoundIndices>)
0168   {
0169     Acts::BoundSubspaceIndices result = Acts::kBoundSubspaceIndicesInvalid;
0170     std::copy(m_subspaceIndices.begin(), m_subspaceIndices.end(),
0171               result.begin());
0172     return result;
0173   }
0174
0175   template <std::size_t dim>
0176   ConstParametersVectorMap<dim> parameters() const {
0177     assert(dim == size());
0178     return ConstParametersVectorMap<dim>{m_params.data()};
0179   }
0180   template <std::size_t dim>
0181   ParametersVectorMap<dim> parameters() {
0182     assert(dim == size());
0183     return ParametersVectorMap<dim>{m_params.data()};
0184   }
0185   ConstEffectiveParametersVectorMap parameters() const {
0186     return ConstEffectiveParametersVectorMap{m_params.data(),
0187                                              static_cast<Eigen::Index>(size())};
0188   }
0189   EffectiveParametersVectorMap parameters() {
0190     return EffectiveParametersVectorMap{m_params.data(),
0191                                         static_cast<Eigen::Index>(size())};
0192   }
0193
0194   template <std::size_t dim>
0195   ConstCovarianceMatrixMap<dim> covariance() const {
0196     assert(dim == size());
0197     return ConstCovarianceMatrixMap<dim>{m_cov.data()};
0198   }
0199   template <std::size_t dim>
0200   CovarianceMatrixMap<dim> covariance() {
0201     assert(dim == size());
0202     return CovarianceMatrixMap<dim>{m_cov.data()};
0203   }
0204   ConstEffectiveCovarianceMatrixMap covariance() const {
0205     return ConstEffectiveCovarianceMatrixMap{m_cov.data(),
0206                                              static_cast<Eigen::Index>(size()),
0207                                              static_cast<Eigen::Index>(size())};
0208   }
0209   EffectiveCovarianceMatrixMap covariance() {
0210     return EffectiveCovarianceMatrixMap{m_cov.data(),
0211                                         static_cast<Eigen::Index>(size()),
0212                                         static_cast<Eigen::Index>(size())};
0213   }
0214
0215   FullParametersVector fullParameters() const {
0216     FullParametersVector result = FullParametersVector::Zero();
0217     for (std::size_t i = 0; i < size(); ++i) {
0218       result[m_subspaceIndices[i]] = parameters()[i];
0219     }
0220     return result;
0221   }
0222
0223   FullCovarianceMatrix fullCovariance() const {
0224     FullCovarianceMatrix result = FullCovarianceMatrix::Zero();
0225     for (std::size_t i = 0; i < size(); ++i) {
0226       for (std::size_t j = 0; j < size(); ++j) {
0227         result(m_subspaceIndices[i], m_subspaceIndices[j]) = covariance()(i, j);
0228       }
0229     }
0230     return result;
0231   }
0232
0233  private:
0234   Acts::SourceLink m_source;
0235   SubspaceIndices m_subspaceIndices;
0236   std::array<Scalar, kFullSize> m_params{};
0237   std::array<Scalar, kFullSize * kFullSize> m_cov{};
0238 };
0239
0240 /// Construct a variable-size measurement for the given indices.
0241 ///
0242 /// @tparam parameters_t Input parameters vector type
0243 /// @tparam covariance_t Input covariance matrix type
0244 /// @tparam indices_t Parameter index type, determines the full parameter space
0245 /// @tparam tail_indices_t Helper types required to support variadic arguments;
0246 ///   all types must be convertibale to `indices_t`.
0247 /// @param source The link that connects to the underlying detector readout
0248 /// @param params Measured parameters values
0249 /// @param cov Measured parameters covariance
0250 /// @param index0 Required parameter index, measurement must be at least 1d
0251 /// @param tailIndices Additional parameter indices for larger measurements
0252 /// @return Variable-size measurement w/ the correct type and given inputs
0253 ///
0254 /// @note The indices must be ordered and must be consistent with the content of
0255 /// parameters and covariance.
0256 template <typename parameters_t, typename covariance_t, typename indices_t,
0257           typename... tail_indices_t>
0258 VariableSizeMeasurement<indices_t> makeVariableSizeMeasurement(
0259     Acts::SourceLink source, const Eigen::MatrixBase<parameters_t>& params,
0260     const Eigen::MatrixBase<covariance_t>& cov, indices_t index0,
0261     tail_indices_t... tailIndices) {
0262   using IndexContainer = std::array<indices_t, 1u + sizeof...(tail_indices_t)>;
0263   return {std::move(source), IndexContainer{index0, tailIndices...}, params,
0264           cov};
0265 }
0266
0267 /// Type that can hold all possible bound measurements.
0268 using BoundVariableMeasurement = VariableSizeMeasurement<Acts::BoundIndices>;
0269
0270 /// Variable measurement type that can contain all possible combinations.
0271 using Measurement = BoundVariableMeasurement;
0272
0273 /// Container of measurements.
0274 ///
0275 /// In contrast to the source links, the measurements themself must not be
0276 /// orderable. The source links stored in the measurements are treated
0277 /// as opaque here and no ordering is enforced on the stored measurements.
0278 using MeasurementContainer = std::vector<Measurement>;
0279
0280 }  // namespace ActsExamples