Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-09-17 08:02:02

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