Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-10-26 07:54:53

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/.
0008 
0009 #pragma once
0010 
0011 #include "Acts/Utilities/GridIterator.hpp"
0012 #include "Acts/Utilities/IAxis.hpp"
0013 #include "Acts/Utilities/Interpolation.hpp"
0014 #include "Acts/Utilities/TypeTag.hpp"
0015 #include "Acts/Utilities/detail/grid_helper.hpp"
0016 #include "Acts/Utilities/detail/interpolation_impl.hpp"
0017 
0018 #include <any>
0019 #include <array>
0020 #include <tuple>
0021 #include <type_traits>
0022 #include <typeinfo>
0023 #include <utility>
0024 #include <vector>
0025 
0026 #include <boost/container/small_vector.hpp>
0027 
0028 namespace Acts {
0029 
0030 namespace detail {
0031 
0032 template <typename>
0033 class AnyGridView;
0034 template <typename>
0035 class AnyGridConstView;
0036 
0037 }  // namespace detail
0038 
0039 /// Base class for all grid types
0040 class IGrid {
0041  public:
0042   virtual ~IGrid() = default;
0043 
0044   /// Get a dynamically sized vector of axis objects for inspection
0045   /// @return a vector of axis pointers
0046   virtual boost::container::small_vector<const IAxis*, 3> axes() const = 0;
0047 
0048   /// @brief Get the number of dimensions of the grid
0049   /// @return The number of dimensions of the grid
0050   virtual std::size_t dimensions() const = 0;
0051 
0052   /// @brief Get the type of the values stored in the grid
0053   /// @return The type of the values stored in the grid
0054   virtual std::type_info const& valueType() const = 0;
0055 
0056   /// Type-erased interface to access the contents of the grid
0057   ///
0058   /// @note This interface has non-negligible runtime overhead due to packing
0059   ///       and unpacking from/to @c std::any and the dynamically sized index and
0060   ///       point types. **USE WITH CARE!**
0061   ///
0062   /// @{
0063   using AnyIndexType = boost::container::small_vector<std::size_t, 3>;
0064   /// Type alias for dynamic point type (coordinates as vector of doubles)
0065   using AnyPointType = boost::container::small_vector<double, 3>;
0066 
0067   /// @brief Get the lower left edge of a bin for a given set of indices
0068   /// @param indices The indices to get the lower left edge of the bin for
0069   /// @return The lower left edge of the bin
0070   virtual AnyPointType lowerLeftBinEdgeAny(AnyIndexType indices) const = 0;
0071 
0072   /// @brief Get the upper right edge of a bin for a given set of indices
0073   /// @param indices The indices to get the upper right edge of the bin for
0074   /// @return The upper right edge of the bin
0075   virtual AnyPointType upperRightBinEdgeAny(AnyIndexType indices) const = 0;
0076 
0077   /// @brief Get the center of a bin for a given set of indices
0078   /// @param indices The indices to get the center of the bin for
0079   /// @return The center of the bin
0080   virtual AnyPointType binCenterAny(AnyIndexType indices) const = 0;
0081 
0082   /// @brief Get the number of local bins for a given set of indices
0083   /// @return The number of local bins
0084   virtual AnyIndexType numLocalBinsAny() const = 0;
0085 
0086   /// @}
0087 
0088   /// Helper to print out the grid
0089   /// @param os the output stream
0090   /// @param grid the grid to print
0091   /// @return the output stream
0092   friend std::ostream& operator<<(std::ostream& os, const IGrid& grid) {
0093     grid.toStream(os);
0094     return os;
0095   }
0096 
0097   friend bool operator==(const IGrid& lhs, const IGrid& rhs) {
0098     auto lhsAxes = lhs.axes();
0099     auto rhsAxes = rhs.axes();
0100     return lhsAxes.size() == rhsAxes.size() &&
0101            std::equal(lhsAxes.begin(), lhsAxes.end(), rhsAxes.begin(),
0102                       [](const IAxis* a, const IAxis* b) { return *a == *b; });
0103   }
0104 
0105  protected:
0106   /// @param os Output stream to write grid representation to
0107   virtual void toStream(std::ostream& os) const = 0;
0108 
0109   /// @brief Get the value of a bin for a given set of indices
0110   /// @param indices The indices to get the value of the bin for
0111   /// @return The value of the bin: the @c std::any contains a const pointer to
0112   ///         the value
0113   virtual std::any atLocalBinsAny(AnyIndexType indices) const = 0;
0114 
0115   /// @brief Get the value of a bin for a given set of indices
0116   /// @param indices The indices to get the value of the bin for
0117   /// @return The value of the bin: the @c std::any contains a pointer to the
0118   ///         value
0119   virtual std::any atLocalBinsAny(AnyIndexType indices) = 0;
0120 
0121   template <typename>
0122   friend class AnyGridView;
0123   template <typename>
0124   friend class AnyGridConstView;
0125 };
0126 
0127 /// @brief class for describing a regular multi-dimensional grid
0128 ///
0129 /// @tparam T    type of values stored inside the bins of the grid
0130 /// @tparam Axes parameter pack of axis types defining the grid
0131 ///
0132 /// Class describing a multi-dimensional, regular grid which can store objects
0133 /// in its multi-dimensional bins. Bins are hyper-boxes and can be accessed
0134 /// either by global bin index, local bin indices or position.
0135 ///
0136 /// @note @c T must be default-constructible.
0137 /// @note @c T must not be @c bool, because @c std::vector<bool> is special
0138 ///          and does not return references to its elements.
0139 template <typename T, class... Axes>
0140   requires(std::is_default_constructible_v<T> && !std::is_same_v<T, bool>)
0141 class Grid final : public IGrid {
0142  public:
0143   /// number of dimensions of the grid
0144   static constexpr std::size_t DIM = sizeof...(Axes);
0145 
0146   /// type of values stored
0147   using value_type = T;
0148   /// reference type to values stored
0149   using reference = value_type&;
0150   /// constant reference type to values stored
0151   using const_reference = const value_type&;
0152   /// type for points in d-dimensional grid space
0153   using point_t = std::array<double, DIM>;
0154   /// index type using local bin indices along each axis
0155   using index_t = std::array<std::size_t, DIM>;
0156   /// global iterator type
0157   using global_iterator_t = GridGlobalIterator<T, Axes...>;
0158   /// local iterator type
0159   using local_iterator_t = GridLocalIterator<T, Axes...>;
0160 
0161   /// @brief Constructor from const axis tuple, this will allow
0162   /// creating a grid with a different value type from a template
0163   /// grid object.
0164   ///
0165   /// @param axes
0166   explicit Grid(const std::tuple<Axes...>& axes) : m_axes(axes) {
0167     m_values.resize(size());
0168   }
0169 
0170   /// @brief Move constructor from axis tuple
0171   /// @param axes
0172   explicit Grid(std::tuple<Axes...>&& axes) : m_axes(std::move(axes)) {
0173     m_values.resize(size());
0174   }
0175 
0176   /// @brief constructor from parameters pack of axes
0177   /// @param axes
0178   explicit Grid(Axes&&... axes) : m_axes(std::forward_as_tuple(axes...)) {
0179     m_values.resize(size());
0180   }
0181 
0182   /// @brief constructor from parameters pack of axes
0183   /// @param axes
0184   explicit Grid(const Axes&... axes) : m_axes(std::tuple(axes...)) {
0185     m_values.resize(size());
0186   }
0187 
0188   /// @brief constructor from parameters pack of axes and type tag
0189   /// @param axes
0190   explicit Grid(TypeTag<T> /*tag*/, Axes&&... axes)
0191       : m_axes(std::forward_as_tuple(axes...)) {
0192     m_values.resize(size());
0193   }
0194 
0195   /// @brief constructor from parameters pack of axes and type tag
0196   /// @param axes
0197   explicit Grid(TypeTag<T> /*tag*/, const Axes&... axes)
0198       : m_axes(std::tuple(axes...)) {
0199     m_values.resize(size());
0200   }
0201 
0202   // Grid(TypeTag<T> /*tag*/, Axes&... axes) = delete;
0203 
0204   /// @brief access value stored in bin for a given point
0205   ///
0206   /// @tparam Point any type with point semantics supporting component access
0207   ///               through @c operator[]
0208   /// @param [in] point point used to look up the corresponding bin in the
0209   ///                   grid
0210   /// @return reference to value stored in bin containing the given point
0211   ///
0212   /// @pre The given @c Point type must represent a point in d (or higher)
0213   ///      dimensions where d is dimensionality of the grid.
0214   ///
0215   /// @note The look-up considers under-/overflow bins along each axis.
0216   ///       Therefore, the look-up will never fail.
0217   //
0218   template <class Point>
0219   reference atPosition(const Point& point) {
0220     return m_values.at(globalBinFromPosition(point));
0221   }
0222 
0223   /// @brief access value stored in bin for a given point
0224   ///
0225   /// @tparam Point any type with point semantics supporting component access
0226   ///               through @c operator[]
0227   /// @param [in] point point used to look up the corresponding bin in the
0228   ///                   grid
0229   /// @return const-reference to value stored in bin containing the given
0230   ///         point
0231   ///
0232   /// @pre The given @c Point type must represent a point in d (or higher)
0233   ///      dimensions where d is dimensionality of the grid.
0234   ///
0235   /// @note The look-up considers under-/overflow bins along each axis.
0236   ///       Therefore, the look-up will never fail.
0237   template <class Point>
0238   const_reference atPosition(const Point& point) const {
0239     return m_values.at(globalBinFromPosition(point));
0240   }
0241 
0242   /// @brief access value stored in bin with given global bin number
0243   ///
0244   /// @param  [in] bin global bin number
0245   /// @return reference to value stored in bin containing the given
0246   ///         point
0247   reference at(std::size_t bin) { return m_values.at(bin); }
0248 
0249   /// @brief access value stored in bin with given global bin number
0250   ///
0251   /// @param  [in] bin global bin number
0252   /// @return const-reference to value stored in bin containing the given
0253   ///         point
0254   const_reference at(std::size_t bin) const { return m_values.at(bin); }
0255 
0256   /// @brief access value stored in bin with given local bin numbers
0257   ///
0258   /// @param  [in] localBins local bin indices along each axis
0259   /// @return reference to value stored in bin containing the given
0260   ///         point
0261   ///
0262   /// @pre All local bin indices must be a valid index for the corresponding
0263   ///      axis (including the under-/overflow bin for this axis).
0264   reference atLocalBins(const index_t& localBins) {
0265     return m_values.at(globalBinFromLocalBins(localBins));
0266   }
0267 
0268   /// @copydoc Acts::IGrid::atLocalBinsAny
0269   std::any atLocalBinsAny(AnyIndexType indices) const override {
0270     const_reference cref = atLocalBins(toIndexType(indices));
0271     return &cref;
0272   }
0273 
0274   /// @brief access value stored in bin with given local bin numbers
0275   ///
0276   /// @param  [in] localBins local bin indices along each axis
0277   /// @return const-reference to value stored in bin containing the given
0278   ///         point
0279   ///
0280   /// @pre All local bin indices must be a valid index for the corresponding
0281   ///      axis (including the under-/overflow bin for this axis).
0282   const_reference atLocalBins(const index_t& localBins) const {
0283     return m_values.at(globalBinFromLocalBins(localBins));
0284   }
0285 
0286   /// @copydoc Acts::IGrid::atLocalBinsAny
0287   std::any atLocalBinsAny(AnyIndexType indices) override {
0288     reference ref = atLocalBins(toIndexType(indices));
0289     return &ref;
0290   }
0291 
0292   /// @brief get global bin indices for closest points on grid
0293   ///
0294   /// @tparam Point any type with point semantics supporting component access
0295   ///               through @c operator[]
0296   /// @param [in] position point of interest
0297   /// @return Iterable thatemits the indices of bins whose lower-left corners
0298   ///         are the closest points on the grid to the input.
0299   ///
0300   /// @pre The given @c Point type must represent a point in d (or higher)
0301   ///      dimensions where d is dimensionality of the grid. It must lie
0302   ///      within the grid range (i.e. not within a under-/overflow bin).
0303   template <class Point>
0304   detail::GlobalNeighborHoodIndices<DIM> closestPointsIndices(
0305       const Point& position) const {
0306     return rawClosestPointsIndices(localBinsFromPosition(position));
0307   }
0308 
0309   /// @brief dimensionality of grid
0310   ///
0311   /// @return number of axes spanning the grid
0312   std::size_t dimensions() const override { return DIM; }
0313 
0314   /// @copydoc Acts::IGrid::valueType
0315   const std::type_info& valueType() const override { return typeid(T); }
0316 
0317   /// @brief get center position of bin with given local bin numbers
0318   ///
0319   /// @param  [in] localBins local bin indices along each axis
0320   /// @return center position of bin
0321   ///
0322   /// @pre All local bin indices must be a valid index for the corresponding
0323   ///      axis (excluding the under-/overflow bins for each axis).
0324   point_t binCenter(const index_t& localBins) const {
0325     return detail::grid_helper::getBinCenter(localBins, m_axes);
0326   }
0327 
0328   AnyPointType binCenterAny(AnyIndexType indices) const override {
0329     return toAnyPointType(binCenter(toIndexType(indices)));
0330   }
0331 
0332   /// @brief determine global index for bin containing the given point
0333   ///
0334   /// @tparam Point any type with point semantics supporting component access
0335   ///               through @c operator[]
0336   ///
0337   /// @param  [in] point point to look up in the grid
0338   /// @return global index for bin containing the given point
0339   ///
0340   /// @pre The given @c Point type must represent a point in d (or higher)
0341   ///      dimensions where d is dimensionality of the grid.
0342   /// @note This could be a under-/overflow bin along one or more axes.
0343   template <class Point>
0344   std::size_t globalBinFromPosition(const Point& point) const {
0345     return globalBinFromLocalBins(localBinsFromPosition(point));
0346   }
0347 
0348   /// @brief determine global bin index from local bin indices along each axis
0349   ///
0350   /// @param  [in] localBins local bin indices along each axis
0351   /// @return global index for bin defined by the local bin indices
0352   ///
0353   /// @pre All local bin indices must be a valid index for the corresponding
0354   ///      axis (including the under-/overflow bin for this axis).
0355   std::size_t globalBinFromLocalBins(const index_t& localBins) const {
0356     return detail::grid_helper::getGlobalBin(localBins, m_axes);
0357   }
0358 
0359   /// @brief  determine global bin index of the bin with the lower left edge
0360   ///         closest to the given point for each axis
0361   ///
0362   /// @tparam Point any type with point semantics supporting component access
0363   ///               through @c operator[]
0364   ///
0365   /// @param  [in] point point to look up in the grid
0366   /// @return global index for bin containing the given point
0367   ///
0368   /// @pre The given @c Point type must represent a point in d (or higher)
0369   ///      dimensions where d is dimensionality of the grid.
0370   /// @note This could be a under-/overflow bin along one or more axes.
0371   template <class Point>
0372   std::size_t globalBinFromFromLowerLeftEdge(const Point& point) const {
0373     return globalBinFromLocalBins(localBinsFromLowerLeftEdge(point));
0374   }
0375 
0376   /// @brief  determine local bin index for each axis from the given point
0377   ///
0378   /// @tparam Point any type with point semantics supporting component access
0379   ///               through @c operator[]
0380   ///
0381   /// @param  [in] point point to look up in the grid
0382   /// @return array with local bin indices along each axis (in same order as
0383   ///         given @c axes object)
0384   ///
0385   /// @pre The given @c Point type must represent a point in d (or higher)
0386   ///      dimensions where d is dimensionality of the grid.
0387   /// @note This could be a under-/overflow bin along one or more axes.
0388   template <class Point>
0389   index_t localBinsFromPosition(const Point& point) const {
0390     return detail::grid_helper::getLocalBinIndices(point, m_axes);
0391   }
0392 
0393   /// @brief determine local bin index for each axis from global bin index
0394   ///
0395   /// @param  [in] bin global bin index
0396   /// @return array with local bin indices along each axis (in same order as
0397   ///         given @c axes object)
0398   ///
0399   /// @note Local bin indices can contain under-/overflow bins along the
0400   ///       corresponding axis.
0401   index_t localBinsFromGlobalBin(std::size_t bin) const {
0402     return detail::grid_helper::getLocalBinIndices(bin, m_axes);
0403   }
0404 
0405   /// @brief  determine local bin index of the bin with the lower left edge
0406   ///         closest to the given point for each axis
0407   ///
0408   /// @tparam Point any type with point semantics supporting component access
0409   ///               through @c operator[]
0410   ///
0411   /// @param  [in] point point to look up in the grid
0412   /// @return array with local bin indices along each axis (in same order as
0413   ///         given @c axes object)
0414   ///
0415   /// @pre The given @c Point type must represent a point in d (or higher)
0416   ///      dimensions where d is dimensionality of the grid.
0417   /// @note This could be a under-/overflow bin along one or more axes.
0418   template <class Point>
0419   index_t localBinsFromLowerLeftEdge(const Point& point) const {
0420     Point shiftedPoint;
0421     point_t width = detail::grid_helper::getWidth(m_axes);
0422     for (std::size_t i = 0; i < DIM; i++) {
0423       shiftedPoint[i] = point[i] + width[i] / 2;
0424     }
0425     return detail::grid_helper::getLocalBinIndices(shiftedPoint, m_axes);
0426   }
0427 
0428   /// @brief retrieve lower-left bin edge from set of local bin indices
0429   ///
0430   /// @param  [in] localBins local bin indices along each axis
0431   /// @return generalized lower-left bin edge position
0432   ///
0433   /// @pre @c localBins must only contain valid bin indices (excluding
0434   ///      underflow bins).
0435   point_t lowerLeftBinEdge(const index_t& localBins) const {
0436     return detail::grid_helper::getLowerLeftBinEdge(localBins, m_axes);
0437   }
0438 
0439   /// @copydoc Acts::IGrid::lowerLeftBinEdgeAny
0440   AnyPointType lowerLeftBinEdgeAny(AnyIndexType indices) const override {
0441     return toAnyPointType(lowerLeftBinEdge(toIndexType(indices)));
0442   }
0443 
0444   /// @brief retrieve upper-right bin edge from set of local bin indices
0445   ///
0446   /// @param  [in] localBins local bin indices along each axis
0447   /// @return generalized upper-right bin edge position
0448   ///
0449   /// @pre @c localBins must only contain valid bin indices (excluding
0450   ///      overflow bins).
0451   point_t upperRightBinEdge(const index_t& localBins) const {
0452     return detail::grid_helper::getUpperRightBinEdge(localBins, m_axes);
0453   }
0454 
0455   /// @copydoc Acts::IGrid::upperRightBinEdgeAny
0456   AnyPointType upperRightBinEdgeAny(AnyIndexType indices) const override {
0457     return toAnyPointType(upperRightBinEdge(toIndexType(indices)));
0458   }
0459 
0460   /// @brief get bin width along each specific axis
0461   ///
0462   /// @return array giving the bin width alonf all axes
0463   point_t binWidth() const { return detail::grid_helper::getWidth(m_axes); }
0464 
0465   /// @brief get number of bins along each specific axis
0466   ///
0467   /// @return array giving the number of bins along all axes
0468   ///
0469   /// @note Not including under- and overflow bins
0470   index_t numLocalBins() const { return detail::grid_helper::getNBins(m_axes); }
0471 
0472   /// @copydoc Acts::IGrid::numLocalBinsAny
0473   AnyIndexType numLocalBinsAny() const override {
0474     return toAnyIndexType(numLocalBins());
0475   }
0476 
0477   /// @brief get the minimum value of all axes of one grid
0478   ///
0479   /// @return array returning the minima of all given axes
0480   point_t minPosition() const { return detail::grid_helper::getMin(m_axes); }
0481 
0482   /// @brief get the maximum value of all axes of one grid
0483   ///
0484   /// @return array returning the maxima of all given axes
0485   point_t maxPosition() const { return detail::grid_helper::getMax(m_axes); }
0486 
0487   /// @brief set all overflow and underflow bins to a certain value
0488   ///
0489   /// @param [in] value value to be inserted in every overflow and underflow
0490   ///                   bin of the grid.
0491   ///
0492   void setExteriorBins(const value_type& value) {
0493     for (std::size_t index : detail::grid_helper::exteriorBinIndices(m_axes)) {
0494       at(index) = value;
0495     }
0496   }
0497 
0498   /// @brief interpolate grid values to given position
0499   ///
0500   /// @tparam Point type specifying geometric positions
0501   /// @tparam U     dummy template parameter identical to @c T
0502   ///
0503   /// @param [in] point location to which to interpolate grid values. The
0504   ///                   position must be within the grid dimensions and not
0505   ///                   lie in an under-/overflow bin along any axis.
0506   ///
0507   /// @return interpolated value at given position
0508   ///
0509   /// @pre The given @c Point type must represent a point in d (or higher)
0510   ///      dimensions where d is dimensionality of the grid.
0511   ///
0512   /// @note This function is available only if the following conditions are
0513   /// fulfilled:
0514   /// - Given @c U and @c V of value type @c T as well as two @c double
0515   /// @c a and @c b, then the following must be a valid expression <tt>a * U + b
0516   /// * V</tt> yielding an object which is (implicitly) convertible to @c T.
0517   /// - @c Point must represent a d-dimensional position and support
0518   /// coordinate access using @c operator[] which should return a @c
0519   /// double (or a value which is implicitly convertible). Coordinate
0520   /// indices must start at 0.
0521   /// @note Bin values are interpreted as being the field values at the
0522   /// lower-left corner of the corresponding hyper-box.
0523   template <class Point>
0524   T interpolate(const Point& point) const
0525     requires(Concepts::interpolatable<T, Point, std::array<double, DIM>,
0526                                       std::array<double, DIM>>)
0527   {
0528     // there are 2^DIM corner points used during the interpolation
0529     constexpr std::size_t nCorners = 1 << DIM;
0530 
0531     // construct vector of pairs of adjacent bin centers and values
0532     std::array<value_type, nCorners> neighbors{};
0533 
0534     // get local indices for current bin
0535     // value of bin is interpreted as being the field value at its lower left
0536     // corner
0537     const auto& llIndices = localBinsFromPosition(point);
0538 
0539     // get global indices for all surrounding corner points
0540     const auto& closestIndices = rawClosestPointsIndices(llIndices);
0541 
0542     // get values on grid points
0543     std::size_t i = 0;
0544     for (std::size_t index : closestIndices) {
0545       neighbors.at(i++) = at(index);
0546     }
0547 
0548     return Acts::interpolate(point, lowerLeftBinEdge(llIndices),
0549                              upperRightBinEdge(llIndices), neighbors);
0550   }
0551 
0552   /// @brief check whether given point is inside grid limits
0553   ///
0554   /// @param position Point to check for inclusion within grid boundaries
0555   /// @return @c true if \f$\text{xmin_i} \le x_i < \text{xmax}_i \forall i=0,
0556   ///         \dots, d-1\f$, otherwise @c false
0557   ///
0558   /// @pre The given @c Point type must represent a point in d (or higher)
0559   ///      dimensions where d is dimensionality of the grid.
0560   ///
0561   /// @post If @c true is returned, the global bin containing the given point
0562   ///       is a valid bin, i.e. it is neither a underflow nor an overflow bin
0563   ///       along any axis.
0564   template <class Point>
0565   bool isInside(const Point& position) const {
0566     return detail::grid_helper::isInside(position, m_axes);
0567   }
0568 
0569   /// @brief get global bin indices for neighborhood
0570   ///
0571   /// @param [in] localBins center bin defined by local bin indices along each
0572   ///                       axis
0573   /// @param [in] size      size of neighborhood determining how many adjacent
0574   ///                       bins along each axis are considered
0575   /// @return set of global bin indices for all bins in neighborhood
0576   ///
0577   /// @note Over-/underflow bins are included in the neighborhood.
0578   /// @note The @c size parameter sets the range by how many units each local
0579   ///       bin index is allowed to be varied. All local bin indices are
0580   ///       varied independently, that is diagonal neighbors are included.
0581   ///       Ignoring the truncation of the neighborhood size reaching beyond
0582   ///       over-/underflow bins, the neighborhood is of size \f$2 \times
0583   ///       \text{size}+1\f$ along each dimension.
0584   detail::GlobalNeighborHoodIndices<DIM> neighborHoodIndices(
0585       const index_t& localBins, std::size_t size = 1u) const {
0586     return detail::grid_helper::neighborHoodIndices(localBins, size, m_axes);
0587   }
0588 
0589   /// @brief get global bin   indices for neighborhood
0590   ///
0591   /// @param [in] localBins   center bin defined by local bin indices along
0592   ///                         each axis. If size is negative, center bin
0593   ///                         is not returned.
0594   /// @param [in] sizePerAxis size of neighborhood for each axis, how many
0595   ///                         adjacent bins along each axis are considered
0596   /// @return set of global bin indices for all bins in neighborhood
0597   ///
0598   /// @note Over-/underflow bins are included in the neighborhood.
0599   /// @note The @c size parameter sets the range by how many units each local
0600   ///       bin index is allowed to be varied. All local bin indices are
0601   ///       varied independently, that is diagonal neighbors are included.
0602   ///       Ignoring the truncation of the neighborhood size reaching beyond
0603   ///       over-/underflow bins, the neighborhood is of size \f$2 \times
0604   ///       \text{size}+1\f$ along each dimension.
0605   detail::GlobalNeighborHoodIndices<DIM> neighborHoodIndices(
0606       const index_t& localBins,
0607       std::array<std::pair<int, int>, DIM>& sizePerAxis) const {
0608     return detail::grid_helper::neighborHoodIndices(localBins, sizePerAxis,
0609                                                     m_axes);
0610   }
0611 
0612   /// @brief total number of bins
0613   ///
0614   /// @param fullCounter Whether to include under-and overflow bins in the count
0615   /// @return total number of bins in the grid
0616   ///
0617   /// @note This number contains under-and overflow bins along all axes.
0618   std::size_t size(bool fullCounter = true) const {
0619     index_t nBinsArray = numLocalBins();
0620     std::size_t current_size = 1;
0621     // add under-and overflow bins for each axis and multiply all bins
0622     if (fullCounter) {
0623       for (const auto& value : nBinsArray) {
0624         current_size *= value + 2;
0625       }
0626     }
0627     // ignore under-and overflow bins for each axis and multiply all bins
0628     else {
0629       for (const auto& value : nBinsArray) {
0630         current_size *= value;
0631       }
0632     }
0633     return current_size;
0634   }
0635 
0636   /// @brief Convenience function to convert the type of the grid
0637   /// to hold another object type.
0638   ///
0639   /// @tparam U the new grid value type
0640   ///
0641   /// @return a new grid with the same axes and a different value type
0642   template <typename U>
0643   Grid<U, Axes...> convertType() const {
0644     Grid<U, Axes...> cGrid(m_axes);
0645     return cGrid;
0646   }
0647 
0648   /// @brief Convenience function to convert the type of the grid
0649   /// to hold another object type.
0650   ///
0651   /// @tparam converter_t the converter type
0652   ///
0653   /// This is designed to be most flexible with a converter object
0654   /// as a visitor. If needed, such a visitor could also use
0655   /// caching or other techniques to speed up the conversion.
0656   ///
0657   /// @param cVisitor the converter object as visitor
0658   ///
0659   /// @return a new grid with the same axes and a different value type
0660   template <typename converter_t>
0661   Grid<typename converter_t::value_type, Axes...> convertGrid(
0662       converter_t& cVisitor) const {
0663     Grid<typename converter_t::value_type, Axes...> cGrid(m_axes);
0664     // Loop through the values and convert them
0665     for (std::size_t i = 0; i < size(); i++) {
0666       cGrid.at(i) = cVisitor(at(i));
0667     }
0668     return cGrid;
0669   }
0670 
0671   /// @brief get the axes as a tuple
0672   /// @return Reference to the tuple containing all grid axes
0673   const std::tuple<Axes...>& axesTuple() const { return m_axes; }
0674 
0675   /// @brief get the axes as an array of IAxis pointers
0676   /// @return Vector containing pointers to all grid axes
0677   boost::container::small_vector<const IAxis*, 3> axes() const override {
0678     boost::container::small_vector<const IAxis*, 3> result;
0679     auto axes = detail::grid_helper::getAxes(m_axes);
0680     std::copy(axes.begin(), axes.end(), std::back_inserter(result));
0681     return result;
0682   }
0683 
0684   /// begin iterator for global bins
0685   /// @return Iterator pointing to the first global bin
0686   global_iterator_t begin() const { return global_iterator_t(*this, 0); }
0687 
0688   /// end iterator for global bins
0689   /// @return Iterator pointing one past the last global bin
0690   global_iterator_t end() const { return global_iterator_t(*this, size()); }
0691 
0692   /// @brief begin iterator for local bins
0693   ///
0694   /// @param navigator is local navigator for the grid
0695   /// @return Iterator pointing to the first local bin
0696   local_iterator_t begin(
0697       const std::array<std::vector<std::size_t>, DIM>& navigator) const {
0698     std::array<std::size_t, DIM> localBin{};
0699     return local_iterator_t(*this, std::move(localBin), navigator);
0700   }
0701 
0702   /// @brief end iterator for local bins
0703   ///
0704   /// @param navigator is local navigator for the grid
0705   /// @return Iterator pointing one past the last local bin
0706   local_iterator_t end(
0707       const std::array<std::vector<std::size_t>, DIM>& navigator) const {
0708     std::array<std::size_t, DIM> endline{};
0709     for (std::size_t i(0ul); i < DIM; ++i) {
0710       endline[i] = navigator[i].size();
0711     }
0712     return local_iterator_t(*this, std::move(endline), navigator);
0713   }
0714 
0715  protected:
0716   void toStream(std::ostream& os) const override {
0717     printAxes(os, std::make_index_sequence<sizeof...(Axes)>());
0718   }
0719 
0720  private:
0721   /// set of axis defining the multi-dimensional grid
0722   std::tuple<Axes...> m_axes;
0723   /// linear value store for each bin
0724   std::vector<T> m_values;
0725 
0726   // Part of closestPointsIndices that goes after local bins resolution.
0727   // Used as an interpolation performance optimization, but not exposed as it
0728   // doesn't make that much sense from an API design standpoint.
0729   detail::GlobalNeighborHoodIndices<DIM> rawClosestPointsIndices(
0730       const index_t& localBins) const {
0731     return detail::grid_helper::closestPointsIndices(localBins, m_axes);
0732   }
0733 
0734   template <std::size_t... Is>
0735   void printAxes(std::ostream& os, std::index_sequence<Is...> /*s*/) const {
0736     auto printOne = [&os, this]<std::size_t index>(
0737                         std::integral_constant<std::size_t, index>) {
0738       if constexpr (index > 0) {
0739         os << ", ";
0740       }
0741       os << std::get<index>(m_axes);
0742     };
0743     (printOne(std::integral_constant<std::size_t, Is>()), ...);
0744   }
0745 
0746   static AnyIndexType toAnyIndexType(const index_t& indices) {
0747     AnyIndexType anyIndices;
0748     anyIndices.reserve(indices.size());
0749     std::ranges::copy(indices, std::back_inserter(anyIndices));
0750     return anyIndices;
0751   }
0752 
0753   static AnyPointType toAnyPointType(const point_t& point) {
0754     AnyPointType anyPoint;
0755     anyPoint.reserve(point.size());
0756     std::ranges::copy(point, std::back_inserter(anyPoint));
0757     return anyPoint;
0758   }
0759 
0760   static index_t toIndexType(const AnyIndexType& indices) {
0761     if (indices.size() != DIM) {
0762       throw std::invalid_argument("Invalid number of indices");
0763     }
0764     index_t concrete;
0765     std::ranges::copy(indices, concrete.begin());
0766     return concrete;
0767   }
0768 };
0769 
0770 /// Deduction guide for Grid with rvalue reference axes
0771 /// @param axes Variable number of axes (rvalue references)
0772 template <typename T, class... Axes>
0773 Grid(TypeTag<T> /*type*/, Axes&&... axes) -> Grid<T, Axes...>;
0774 
0775 /// Deduction guide for Grid with lvalue reference axes
0776 /// @param axes Variable number of axes (lvalue references)
0777 template <typename T, class... Axes>
0778 Grid(TypeTag<T> /*type*/, Axes&... axes) -> Grid<T, Axes...>;
0779 
0780 }  // namespace Acts