Back to home page

EIC code displayed by LXR
 
 

    

File indexing completed on 2026-04-08 07:46:52 

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/Axis.hpp"
0012 #include "Acts/Utilities/IAxis.hpp"
0013 
0014 #include <array>
0015 #include <set>
0016 #include <tuple>
0017 #include <utility>
0018 
0019 #include <boost/container/small_vector.hpp>
0020 
0021 namespace Acts::detail {
0022 
0023 template <typename T>
0024 constexpr T ipow(T num, unsigned int pow) {
0025   return (pow >= sizeof(unsigned int) * 8) ? 0
0026          : pow == 0                        ? 1
0027                                            : num * ipow(num, pow - 1);
0028 }
0029 
0030 // This object can be iterated to produce the (ordered) set of global indices
0031 // associated with a neighborhood around a certain point on a grid.
0032 //
0033 // The goal is to emulate the effect of enumerating the global indices into
0034 // an std::set (or into an std::vector that gets subsequently sorted), without
0035 // paying the price of dynamic memory allocation in hot magnetic field
0036 // interpolation code.
0037 //
0038 template <std::size_t DIM>
0039 class GlobalNeighborHoodIndices {
0040  public:
0041   // You can get the local neighbor indices from
0042   // grid_helper_impl<DIM>::neighborHoodIndices and the number of bins in
0043   // each direction from grid_helper_impl<DIM>::getNBins.
0044   GlobalNeighborHoodIndices(
0045       std::array<NeighborHoodIndices, DIM>& neighborIndices,
0046       const std::array<std::size_t, DIM>& nBinsArray)
0047       : m_localIndices(neighborIndices) {
0048     if (DIM == 1) {
0049       return;
0050     }
0051     std::size_t globalStride = 1;
0052     for (long i = DIM - 2; i >= 0; --i) {
0053       globalStride *= (nBinsArray[i + 1] + 2);
0054       m_globalStrides[i] = globalStride;
0055     }
0056   }
0057 
0058   class iterator {
0059    public:
0060     using iterator_category = std::forward_iterator_tag;
0061     using value_type = std::size_t;
0062     using difference_type = std::ptrdiff_t;
0063     using pointer = std::size_t*;
0064     using reference = std::size_t&;
0065 
0066     iterator() = default;
0067 
0068     iterator(const GlobalNeighborHoodIndices& parent,
0069              std::array<NeighborHoodIndices::iterator, DIM>&& localIndicesIter)
0070         : m_localIndicesIter(std::move(localIndicesIter)), m_parent(&parent) {}
0071 
0072     std::size_t operator*() const {
0073       std::size_t globalIndex = *m_localIndicesIter[DIM - 1];
0074       if (DIM == 1) {
0075         return globalIndex;
0076       }
0077       for (std::size_t i = 0; i < DIM - 1; ++i) {
0078         globalIndex += m_parent->m_globalStrides[i] * (*m_localIndicesIter[i]);
0079       }
0080       return globalIndex;
0081     }
0082 
0083     iterator& operator++() {
0084       const auto& localIndices = m_parent->m_localIndices;
0085 
0086       // Go to the next global index via a lexicographic increment:
0087       // - Start by incrementing the last local index
0088       // - If it reaches the end, reset it and increment the previous one...
0089       for (long i = DIM - 1; i > 0; --i) {
0090         ++m_localIndicesIter[i];
0091         if (m_localIndicesIter[i] != localIndices[i].end()) {
0092           return *this;
0093         }
0094         m_localIndicesIter[i] = localIndices[i].begin();
0095       }
0096 
0097       // The first index should stay at the end value when it reaches it, so
0098       // that we know when we've reached the end of iteration.
0099       ++m_localIndicesIter[0];
0100       return *this;
0101     }
0102 
0103     iterator operator++(int) {
0104       iterator tmp = *this;
0105       ++(*this);
0106       return tmp;
0107     }
0108 
0109     bool isEqual(const iterator& b) const {
0110       if (b.m_parent == nullptr) {
0111         return m_localIndicesIter[0] == m_parent->m_localIndices[0].end();
0112       } else {
0113         return m_localIndicesIter == b.m_localIndicesIter;
0114       }
0115     }
0116 
0117     friend bool operator==(const iterator& a, const iterator& b) {
0118       return a.isEqual(b);
0119     }
0120 
0121    private:
0122     std::array<NeighborHoodIndices::iterator, DIM> m_localIndicesIter;
0123     const GlobalNeighborHoodIndices* m_parent = nullptr;
0124   };
0125 
0126   iterator begin() const {
0127     std::array<NeighborHoodIndices::iterator, DIM> localIndicesIter{};
0128     for (std::size_t i = 0; i < DIM; ++i) {
0129       localIndicesIter[i] = m_localIndices[i].begin();
0130     }
0131     return iterator(*this, std::move(localIndicesIter));
0132   }
0133 
0134   iterator end() const { return iterator(); }
0135 
0136   // Number of indices that will be produced if this sequence is iterated
0137   std::size_t size() const {
0138     std::size_t result = m_localIndices[0].size();
0139     for (std::size_t i = 1; i < DIM; ++i) {
0140       result *= m_localIndices[i].size();
0141     }
0142     return result;
0143   }
0144 
0145   // Collect the sequence of indices into a vector
0146   auto collect() const {
0147     boost::container::small_vector<std::size_t, ipow(3, DIM)> result;
0148     result.reserve(this->size());
0149     for (std::size_t idx : *this) {
0150       result.push_back(idx);
0151     }
0152     return result;
0153   }
0154 
0155   std::vector<std::size_t> collectVector() const {
0156     auto result = collect();
0157     return {result.begin(), result.end()};
0158   }
0159 
0160  private:
0161   std::array<NeighborHoodIndices, DIM> m_localIndices{};
0162   std::array<std::size_t, DIM - 1> m_globalStrides{};
0163 };
0164 
0165 /// @cond
0166 /// @brief helper struct to calculate number of bins inside a grid
0167 ///
0168 /// @tparam N number of axes to consider
0169 template <std::size_t N>
0170 struct grid_helper_impl;
0171 
0172 template <std::size_t N>
0173 struct grid_helper_impl {
0174   template <class... Axes>
0175   static void getBinCenter(
0176       std::array<double, sizeof...(Axes)>& center,
0177       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0178       const std::tuple<Axes...>& axes) {
0179     center.at(N) = std::get<N>(axes).getBinCenter(localIndices.at(N));
0180     grid_helper_impl<N - 1>::getBinCenter(center, localIndices, axes);
0181   }
0182 
0183   template <class... Axes>
0184   static void getGlobalBin(
0185       const std::array<std::size_t, sizeof...(Axes)>& localBins,
0186       const std::tuple<Axes...>& axes, std::size_t& bin, std::size_t& area) {
0187     const auto& thisAxis = std::get<N>(axes);
0188     bin += area * localBins.at(N);
0189     // make sure to account for under-/overflow bins
0190     area *= (thisAxis.getNBins() + 2);
0191     grid_helper_impl<N - 1>::getGlobalBin(localBins, axes, bin, area);
0192   }
0193 
0194   template <class Point, class... Axes>
0195   static void getLocalBinIndices(
0196       const Point& point, const std::tuple<Axes...>& axes,
0197       std::array<std::size_t, sizeof...(Axes)>& indices) {
0198     const auto& thisAxis = std::get<N>(axes);
0199     indices.at(N) = static_cast<std::size_t>(thisAxis.getBin(point[N]));
0200     grid_helper_impl<N - 1>::getLocalBinIndices(point, axes, indices);
0201   }
0202 
0203   template <class... Axes>
0204   static void getLocalBinIndices(
0205       std::size_t& bin, const std::tuple<Axes...>& axes, std::size_t& area,
0206       std::array<std::size_t, sizeof...(Axes)>& indices) {
0207     const auto& thisAxis = std::get<N>(axes);
0208     // make sure to account for under-/overflow bins
0209     std::size_t new_area = area * (thisAxis.getNBins() + 2);
0210     grid_helper_impl<N - 1>::getLocalBinIndices(bin, axes, new_area, indices);
0211     indices.at(N) = bin / area;
0212     bin %= area;
0213   }
0214 
0215   template <class... Axes>
0216   static void getLowerLeftBinEdge(
0217       std::array<double, sizeof...(Axes)>& llEdge,
0218       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0219       const std::tuple<Axes...>& axes) {
0220     llEdge.at(N) = std::get<N>(axes).getBinLowerBound(localIndices.at(N));
0221     grid_helper_impl<N - 1>::getLowerLeftBinEdge(llEdge, localIndices, axes);
0222   }
0223 
0224   template <class... Axes>
0225   static void getLowerLeftBinIndices(
0226       std::array<std::size_t, sizeof...(Axes)>& localIndices,
0227       const std::tuple<Axes...>& axes) {
0228     localIndices.at(N) = std::get<N>(axes).wrapBin(localIndices.at(N) - 1);
0229     grid_helper_impl<N - 1>::getLowerLeftBinIndices(localIndices, axes);
0230   }
0231 
0232   template <class... Axes>
0233   static void getNBins(const std::tuple<Axes...>& axes,
0234                        std::array<std::size_t, sizeof...(Axes)>& nBinsArray) {
0235     // by convention getNBins does not include under-/overflow bins
0236     nBinsArray[N] = std::get<N>(axes).getNBins();
0237     grid_helper_impl<N - 1>::getNBins(axes, nBinsArray);
0238   }
0239 
0240   template <class... Axes>
0241   static void getAxes(const std::tuple<Axes...>& axes,
0242                       std::array<const IAxis*, sizeof...(Axes)>& axesArr) {
0243     axesArr[N] = static_cast<const IAxis*>(&std::get<N>(axes));
0244     grid_helper_impl<N - 1>::getAxes(axes, axesArr);
0245   }
0246 
0247   template <class... Axes>
0248   static void getUpperRightBinEdge(
0249       std::array<double, sizeof...(Axes)>& urEdge,
0250       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0251       const std::tuple<Axes...>& axes) {
0252     urEdge.at(N) = std::get<N>(axes).getBinUpperBound(localIndices.at(N));
0253     grid_helper_impl<N - 1>::getUpperRightBinEdge(urEdge, localIndices, axes);
0254   }
0255 
0256   template <class... Axes>
0257   static void getUpperRightBinIndices(
0258       std::array<std::size_t, sizeof...(Axes)>& localIndices,
0259       const std::tuple<Axes...>& axes) {
0260     localIndices.at(N) = std::get<N>(axes).wrapBin(localIndices.at(N) + 1);
0261     grid_helper_impl<N - 1>::getUpperRightBinIndices(localIndices, axes);
0262   }
0263 
0264   template <class... Axes>
0265   static void getMin(const std::tuple<Axes...>& axes,
0266                      std::array<double, sizeof...(Axes)>& minArray) {
0267     minArray[N] = std::get<N>(axes).getMin();
0268     grid_helper_impl<N - 1>::getMin(axes, minArray);
0269   }
0270 
0271   template <class... Axes>
0272   static void getMax(const std::tuple<Axes...>& axes,
0273                      std::array<double, sizeof...(Axes)>& maxArray) {
0274     maxArray[N] = std::get<N>(axes).getMax();
0275     grid_helper_impl<N - 1>::getMax(axes, maxArray);
0276   }
0277 
0278   template <class... Axes>
0279   static void getWidth(const std::tuple<Axes...>& axes,
0280                        std::array<double, sizeof...(Axes)>& widthArray) {
0281     widthArray[N] = std::get<N>(axes).getBinWidth();
0282     grid_helper_impl<N - 1>::getWidth(axes, widthArray);
0283   }
0284 
0285   template <class Point, class... Axes>
0286   static bool isInside(const Point& position, const std::tuple<Axes...>& axes) {
0287     bool insideThisAxis = std::get<N>(axes).isInside(position[N]);
0288     return insideThisAxis && grid_helper_impl<N - 1>::isInside(position, axes);
0289   }
0290 
0291   template <class... Axes>
0292   static void neighborHoodIndices(
0293       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0294       std::pair<int, int> sizes, const std::tuple<Axes...>& axes,
0295       std::array<NeighborHoodIndices, sizeof...(Axes)>& neighborIndices) {
0296     // ask n-th axis
0297     std::size_t locIdx = localIndices.at(N);
0298     NeighborHoodIndices locNeighbors =
0299         std::get<N>(axes).neighborHoodIndices(locIdx, sizes);
0300     neighborIndices.at(N) = locNeighbors;
0301 
0302     grid_helper_impl<N - 1>::neighborHoodIndices(localIndices, sizes, axes,
0303                                                  neighborIndices);
0304   }
0305 
0306   template <class... Axes>
0307   static void neighborHoodIndices(
0308       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0309       std::array<std::pair<int, int>, sizeof...(Axes)> sizes,
0310       const std::tuple<Axes...>& axes,
0311       std::array<NeighborHoodIndices, sizeof...(Axes)>& neighborIndices) {
0312     // ask n-th axis
0313     std::size_t locIdx = localIndices.at(N);
0314     NeighborHoodIndices locNeighbors =
0315         std::get<N>(axes).neighborHoodIndices(locIdx, sizes.at(N));
0316     neighborIndices.at(N) = locNeighbors;
0317 
0318     grid_helper_impl<N - 1>::neighborHoodIndices(localIndices, sizes, axes,
0319                                                  neighborIndices);
0320   }
0321 
0322   template <class... Axes>
0323   static void exteriorBinIndices(std::array<std::size_t, sizeof...(Axes)>& idx,
0324                                  std::array<bool, sizeof...(Axes)> isExterior,
0325                                  std::set<std::size_t>& combinations,
0326                                  const std::tuple<Axes...>& axes) {
0327     // iterate over this axis' bins, remembering which bins are exterior
0328     for (std::size_t i = 0; i < std::get<N>(axes).getNBins() + 2; ++i) {
0329       idx.at(N) = i;
0330       isExterior.at(N) = (i == 0) || (i == std::get<N>(axes).getNBins() + 1);
0331       // vary other axes recursively
0332       grid_helper_impl<N - 1>::exteriorBinIndices(idx, isExterior, combinations,
0333                                                   axes);
0334     }
0335   }
0336 };
0337 
0338 template <>
0339 struct grid_helper_impl<0u> {
0340   template <class... Axes>
0341   static void getBinCenter(
0342       std::array<double, sizeof...(Axes)>& center,
0343       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0344       const std::tuple<Axes...>& axes) {
0345     center.at(0u) = std::get<0u>(axes).getBinCenter(localIndices.at(0u));
0346   }
0347 
0348   template <class... Axes>
0349   static void getGlobalBin(
0350       const std::array<std::size_t, sizeof...(Axes)>& localBins,
0351       const std::tuple<Axes...>& /*axes*/, std::size_t& bin,
0352       std::size_t& area) {
0353     bin += area * localBins.at(0u);
0354   }
0355 
0356   template <class Point, class... Axes>
0357   static void getLocalBinIndices(
0358       const Point& point, const std::tuple<Axes...>& axes,
0359       std::array<std::size_t, sizeof...(Axes)>& indices) {
0360     const auto& thisAxis = std::get<0u>(axes);
0361     indices.at(0u) = thisAxis.getBin(point[0u]);
0362   }
0363 
0364   template <class... Axes>
0365   static void getLocalBinIndices(
0366       std::size_t& bin, const std::tuple<Axes...>& /*axes*/, std::size_t& area,
0367       std::array<std::size_t, sizeof...(Axes)>& indices) {
0368     // make sure to account for under-/overflow bins
0369     indices.at(0u) = bin / area;
0370     bin %= area;
0371   }
0372 
0373   template <class... Axes>
0374   static void getLowerLeftBinEdge(
0375       std::array<double, sizeof...(Axes)>& llEdge,
0376       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0377       const std::tuple<Axes...>& axes) {
0378     llEdge.at(0u) = std::get<0u>(axes).getBinLowerBound(localIndices.at(0u));
0379   }
0380 
0381   template <class... Axes>
0382   static void getLowerLeftBinIndices(
0383       std::array<std::size_t, sizeof...(Axes)>& localIndices,
0384       const std::tuple<Axes...>& axes) {
0385     localIndices.at(0u) = std::get<0u>(axes).wrapBin(localIndices.at(0u) - 1);
0386   }
0387 
0388   template <class... Axes>
0389   static void getNBins(const std::tuple<Axes...>& axes,
0390                        std::array<std::size_t, sizeof...(Axes)>& nBinsArray) {
0391     // by convention getNBins does not include under-/overflow bins
0392     nBinsArray[0u] = std::get<0u>(axes).getNBins();
0393   }
0394 
0395   template <class... Axes>
0396   static void getAxes(const std::tuple<Axes...>& axes,
0397                       std::array<const IAxis*, sizeof...(Axes)>& axesArr) {
0398     axesArr[0u] = static_cast<const IAxis*>(&std::get<0u>(axes));
0399   }
0400 
0401   template <class... Axes>
0402   static void getUpperRightBinEdge(
0403       std::array<double, sizeof...(Axes)>& urEdge,
0404       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0405       const std::tuple<Axes...>& axes) {
0406     urEdge.at(0u) = std::get<0u>(axes).getBinUpperBound(localIndices.at(0u));
0407   }
0408 
0409   template <class... Axes>
0410   static void getUpperRightBinIndices(
0411       std::array<std::size_t, sizeof...(Axes)>& localIndices,
0412       const std::tuple<Axes...>& axes) {
0413     localIndices.at(0u) = std::get<0u>(axes).wrapBin(localIndices.at(0u) + 1);
0414   }
0415 
0416   template <class... Axes>
0417   static void getMin(const std::tuple<Axes...>& axes,
0418                      std::array<double, sizeof...(Axes)>& minArray) {
0419     minArray[0u] = std::get<0u>(axes).getMin();
0420   }
0421 
0422   template <class... Axes>
0423   static void getMax(const std::tuple<Axes...>& axes,
0424                      std::array<double, sizeof...(Axes)>& maxArray) {
0425     maxArray[0u] = std::get<0u>(axes).getMax();
0426   }
0427 
0428   template <class... Axes>
0429   static void getWidth(const std::tuple<Axes...>& axes,
0430                        std::array<double, sizeof...(Axes)>& widthArray) {
0431     widthArray[0u] = std::get<0u>(axes).getBinWidth();
0432   }
0433 
0434   template <class Point, class... Axes>
0435   static bool isInside(const Point& position, const std::tuple<Axes...>& axes) {
0436     return std::get<0u>(axes).isInside(position[0u]);
0437   }
0438 
0439   template <class... Axes>
0440   static void neighborHoodIndices(
0441       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0442       std::pair<int, int> sizes, const std::tuple<Axes...>& axes,
0443       std::array<NeighborHoodIndices, sizeof...(Axes)>& neighborIndices) {
0444     // ask 0-th axis
0445     std::size_t locIdx = localIndices.at(0u);
0446     NeighborHoodIndices locNeighbors =
0447         std::get<0u>(axes).neighborHoodIndices(locIdx, sizes);
0448     neighborIndices.at(0u) = locNeighbors;
0449   }
0450 
0451   template <class... Axes>
0452   static void neighborHoodIndices(
0453       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0454       std::array<std::pair<int, int>, sizeof...(Axes)> sizes,
0455       const std::tuple<Axes...>& axes,
0456       std::array<NeighborHoodIndices, sizeof...(Axes)>& neighborIndices) {
0457     // ask 0-th axis
0458     std::size_t locIdx = localIndices.at(0u);
0459     NeighborHoodIndices locNeighbors =
0460         std::get<0u>(axes).neighborHoodIndices(locIdx, sizes.at(0u));
0461     neighborIndices.at(0u) = locNeighbors;
0462   }
0463 
0464   template <class... Axes>
0465   static void exteriorBinIndices(std::array<std::size_t, sizeof...(Axes)>& idx,
0466                                  std::array<bool, sizeof...(Axes)> isExterior,
0467                                  std::set<std::size_t>& combinations,
0468                                  const std::tuple<Axes...>& axes) {
0469     // For each exterior bin on this axis, we will do this
0470     auto recordExteriorBin = [&](std::size_t i) {
0471       idx.at(0u) = i;
0472       // at this point, combinations are complete: save the global bin
0473       std::size_t bin = 0, area = 1;
0474       grid_helper_impl<sizeof...(Axes) - 1>::getGlobalBin(idx, axes, bin, area);
0475       combinations.insert(bin);
0476     };
0477 
0478     // The first and last bins on this axis are exterior by definition
0479     for (std::size_t i :
0480          {static_cast<std::size_t>(0), std::get<0u>(axes).getNBins() + 1}) {
0481       recordExteriorBin(i);
0482     }
0483 
0484     // If no other axis is on an exterior index, stop here
0485     bool otherAxisExterior = false;
0486     for (std::size_t N = 1; N < sizeof...(Axes); ++N) {
0487       otherAxisExterior = otherAxisExterior | isExterior[N];
0488     }
0489     if (!otherAxisExterior) {
0490       return;
0491     }
0492 
0493     // Otherwise, we're on a grid border: iterate over all the other indices
0494     for (std::size_t i = 1; i <= std::get<0u>(axes).getNBins(); ++i) {
0495       recordExteriorBin(i);
0496     }
0497   }
0498 };
0499 /// @endcond
0500 
0501 /// @brief helper functions for grid-related operations
0502 struct grid_helper {
0503   /// @brief get the global indices for closest points on grid
0504   ///
0505   /// @tparam Axes parameter pack of axis types defining the grid
0506   /// @param  [in] bin  global bin index for bin of interest
0507   /// @param  [in] axes actual axis objects spanning the grid
0508   /// @return Sorted collection of global bin indices for bins whose
0509   ///         lower-left corners are the closest points on the grid to every
0510   ///         point in the given bin
0511   ///
0512   /// @note @c bin must be a valid bin index (excluding under-/overflow bins
0513   ///       along any axis).
0514   template <class... Axes>
0515   static GlobalNeighborHoodIndices<sizeof...(Axes)> closestPointsIndices(
0516       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0517       const std::tuple<Axes...>& axes) {
0518     // get neighboring bins, but only increment.
0519     return neighborHoodIndices(localIndices, std::make_pair(0, 1), axes);
0520   }
0521 
0522   /// @brief retrieve bin center from set of local bin indices
0523   ///
0524   /// @tparam Axes parameter pack of axis types defining the grid
0525   /// @param  [in] localIndices local bin indices along each axis
0526   /// @param  [in] axes         actual axis objects spanning the grid
0527   /// @return center position of bin
0528   ///
0529   /// @pre @c localIndices must only contain valid bin indices (i.e. excluding
0530   ///      under-/overflow bins).
0531   template <class... Axes>
0532   static std::array<double, sizeof...(Axes)> getBinCenter(
0533       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0534       const std::tuple<Axes...>& axes) {
0535     std::array<double, sizeof...(Axes)> center{};
0536     grid_helper_impl<sizeof...(Axes) - 1>::getBinCenter(center, localIndices,
0537                                                         axes);
0538 
0539     return center;
0540   }
0541 
0542   /// @brief determine global bin index from local indices along each axis
0543   ///
0544   /// @tparam Axes parameter pack of axis types defining the grid
0545   ///
0546   /// @param  [in] localBins local bin indices along each axis
0547   /// @param  [in] axes  actual axis objects spanning the grid
0548   /// @return global index for bin defined by the local bin indices
0549   ///
0550   /// @pre All local bin indices must be a valid index for the corresponding
0551   ///      axis (including the under-/overflow bin for this axis).
0552   template <class... Axes>
0553   static std::size_t getGlobalBin(
0554       const std::array<std::size_t, sizeof...(Axes)>& localBins,
0555       const std::tuple<Axes...>& axes) {
0556     std::size_t area = 1;
0557     std::size_t bin = 0;
0558 
0559     grid_helper_impl<sizeof...(Axes) - 1>::getGlobalBin(localBins, axes, bin,
0560                                                         area);
0561 
0562     return bin;
0563   }
0564 
0565   /// @brief determine local bin index for each axis from point
0566   ///
0567   /// @tparam Point any type with point semantics supporting component access
0568   ///               through @c operator[]
0569   /// @tparam Axes parameter pack of axis types defining the grid
0570   ///
0571   /// @param  [in] point point to look up in the grid
0572   /// @param  [in] axes  actual axis objects spanning the grid
0573   /// @return array with local bin indices along each axis (in same order as
0574   ///         given @c axes object)
0575   ///
0576   /// @pre The given @c Point type must represent a point in d (or higher)
0577   ///      dimensions where d is the number of axis objects in the tuple.
0578   /// @note This could be a under-/overflow bin along one or more axes.
0579   template <class Point, class... Axes>
0580   static std::array<std::size_t, sizeof...(Axes)> getLocalBinIndices(
0581       const Point& point, const std::tuple<Axes...>& axes) {
0582     std::array<std::size_t, sizeof...(Axes)> indices{};
0583 
0584     grid_helper_impl<sizeof...(Axes) - 1>::getLocalBinIndices(point, axes,
0585                                                               indices);
0586 
0587     return indices;
0588   }
0589 
0590   /// @brief determine local bin index for each axis from global bin index
0591   ///
0592   /// @tparam Axes parameter pack of axis types defining the grid
0593   ///
0594   /// @param  [in] bin  global bin index
0595   /// @param  [in] axes actual axis objects spanning the grid
0596   /// @return array with local bin indices along each axis (in same order as
0597   ///         given @c axes object)
0598   ///
0599   /// @note Local bin indices can contain under-/overflow bins along any axis.
0600   template <class... Axes>
0601   static std::array<std::size_t, sizeof...(Axes)> getLocalBinIndices(
0602       std::size_t bin, const std::tuple<Axes...>& axes) {
0603     std::size_t area = 1;
0604     std::array<std::size_t, sizeof...(Axes)> indices{};
0605 
0606     grid_helper_impl<sizeof...(Axes) - 1>::getLocalBinIndices(bin, axes, area,
0607                                                               indices);
0608 
0609     return indices;
0610   }
0611 
0612   /// @brief retrieve lower-left bin edge from set of local bin indices
0613   ///
0614   /// @tparam Axes parameter pack of axis types defining the grid
0615   /// @param  [in] localIndices local bin indices along each axis
0616   /// @param  [in] axes         actual axis objects spanning the grid
0617   /// @return generalized lower-left bin edge
0618   ///
0619   /// @pre @c localIndices must only contain valid bin indices (excluding
0620   ///      underflow bins).
0621   template <class... Axes>
0622   static std::array<double, sizeof...(Axes)> getLowerLeftBinEdge(
0623       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0624       const std::tuple<Axes...>& axes) {
0625     std::array<double, sizeof...(Axes)> llEdge{};
0626     grid_helper_impl<sizeof...(Axes) - 1>::getLowerLeftBinEdge(
0627         llEdge, localIndices, axes);
0628 
0629     return llEdge;
0630   }
0631 
0632   /// @brief get local bin indices for lower-left neighboring bin
0633   ///
0634   /// @tparam Axes parameter pack of axis types defining the grid
0635   /// @param  [in] localIndices local bin indices along each axis
0636   /// @param  [in] axes         actual axis objects spanning the grid
0637   /// @return array with local bin indices of lower-left neighbor bin
0638   ///
0639   /// @pre @c localIndices must only contain valid bin indices (excluding
0640   ///      underflow bins).
0641   ///
0642   /// This function returns the local bin indices for the generalized
0643   /// lower-left neighbor which simply means that all local bin indices are
0644   /// decremented by one.
0645   template <class... Axes>
0646   static std::array<std::size_t, sizeof...(Axes)> getLowerLeftBinIndices(
0647       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0648       const std::tuple<Axes...>& axes) {
0649     auto llIndices = localIndices;
0650     grid_helper_impl<sizeof...(Axes) - 1>::getLowerLeftBinIndices(llIndices,
0651                                                                   axes);
0652 
0653     return llIndices;
0654   }
0655 
0656   /// @brief calculate number of bins in a grid defined by a set of
0657   /// axes for each axis
0658   ///
0659   /// @tparam Axes parameter pack of axis types defining the grid
0660   /// @param  [in] axes actual axis objects spanning the grid
0661   /// @return array of number of bins for each axis of the grid
0662   ///
0663   /// @note This does not include under-/overflow bins along each axis.
0664   template <class... Axes>
0665   static std::array<std::size_t, sizeof...(Axes)> getNBins(
0666       const std::tuple<Axes...>& axes) {
0667     std::array<std::size_t, sizeof...(Axes)> nBinsArray{};
0668     grid_helper_impl<sizeof...(Axes) - 1>::getNBins(axes, nBinsArray);
0669     return nBinsArray;
0670   }
0671 
0672   /// @brief return an array with copies of the axes, converted
0673   /// to type AnyAxis
0674   ///
0675   /// @tparam Axes parameter pack of axis types defining the grid
0676   /// @param [in] axes actual axis objects spanning the grid
0677   /// @return array with copies of the axis
0678   template <class... Axes>
0679   static std::array<const IAxis*, sizeof...(Axes)> getAxes(
0680       const std::tuple<Axes...>& axes) {
0681     std::array<const IAxis*, sizeof...(Axes)> arr{};
0682     grid_helper_impl<sizeof...(Axes) - 1>::getAxes(axes, arr);
0683     return arr;
0684   }
0685 
0686   /// @brief retrieve upper-right bin edge from set of local bin indices
0687   ///
0688   /// @tparam Axes parameter pack of axis types defining the grid
0689   /// @param  [in] localIndices local bin indices along each axis
0690   /// @param  [in] axes         actual axis objects spanning the grid
0691   /// @return generalized upper-right bin edge
0692   ///
0693   /// @pre @c localIndices must only contain valid bin indices (excluding
0694   ///      overflow bins).
0695   template <class... Axes>
0696   static std::array<double, sizeof...(Axes)> getUpperRightBinEdge(
0697       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0698       const std::tuple<Axes...>& axes) {
0699     std::array<double, sizeof...(Axes)> urEdge{};
0700     grid_helper_impl<sizeof...(Axes) - 1>::getUpperRightBinEdge(
0701         urEdge, localIndices, axes);
0702 
0703     return urEdge;
0704   }
0705 
0706   /// @brief get local bin indices for upper-right neighboring bin
0707   ///
0708   /// @tparam Axes parameter pack of axis types defining the grid
0709   /// @param  [in] localIndices local bin indices along each axis
0710   /// @param  [in] axes         actual axis objects spanning the grid
0711   /// @return array with local bin indices of upper-right neighbor bin
0712   ///
0713   /// @pre @c localIndices must only contain valid bin indices (excluding
0714   ///      overflow bins).
0715   ///
0716   /// This function returns the local bin indices for the generalized
0717   /// upper-right neighbor which simply means that all local bin indices are
0718   /// incremented by one.
0719   template <class... Axes>
0720   static std::array<std::size_t, sizeof...(Axes)> getUpperRightBinIndices(
0721       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0722       const std::tuple<Axes...>& axes) {
0723     auto urIndices = localIndices;
0724     grid_helper_impl<sizeof...(Axes) - 1>::getUpperRightBinIndices(urIndices,
0725                                                                    axes);
0726 
0727     return urIndices;
0728   }
0729 
0730   /// @brief get the minimum value of all axes of one grid
0731   ///
0732   /// @tparam Axes parameter pack of axis types defining the grid
0733   /// @param  [in] axes actual axis objects spanning the grid
0734   /// @return array returning the minima of all given axes
0735   template <class... Axes>
0736   static std::array<double, sizeof...(Axes)> getMin(
0737       const std::tuple<Axes...>& axes) {
0738     std::array<double, sizeof...(Axes)> minArray{};
0739     grid_helper_impl<sizeof...(Axes) - 1>::getMin(axes, minArray);
0740     return minArray;
0741   }
0742 
0743   /// @brief get the maximum value of all axes of one grid
0744   ///
0745   /// @tparam Axes parameter pack of axis types defining the grid
0746   /// @param  [in] axes actual axis objects spanning the grid
0747   /// @return array returning the maxima of all given axes
0748   template <class... Axes>
0749   static std::array<double, sizeof...(Axes)> getMax(
0750       const std::tuple<Axes...>& axes) {
0751     std::array<double, sizeof...(Axes)> maxArray{};
0752     grid_helper_impl<sizeof...(Axes) - 1>::getMax(axes, maxArray);
0753     return maxArray;
0754   }
0755 
0756   /// @brief get the bin width of all axes of one grid
0757   ///
0758   /// @tparam Axes parameter pack of axis types defining the grid
0759   /// @param  [in] axes actual axis objects spanning the grid
0760   /// @return array returning the maxima of all given axes
0761   template <class... Axes>
0762   static std::array<double, sizeof...(Axes)> getWidth(
0763       const std::tuple<Axes...>& axes) {
0764     std::array<double, sizeof...(Axes)> widthArray{};
0765     grid_helper_impl<sizeof...(Axes) - 1>::getWidth(axes, widthArray);
0766     return widthArray;
0767   }
0768 
0769   /// @brief get global bin indices for bins in specified neighborhood
0770   ///
0771   /// @tparam Axes parameter pack of axis types defining the grid
0772   /// @param  [in] localIndices local bin indices along each axis
0773   /// @param  [in] size         size of neighborhood determining how many
0774   ///                           adjacent bins along each axis are considered
0775   /// @param  [in] axes         actual axis objects spanning the grid
0776   /// @return Sorted collection of global bin indices for all bins in
0777   ///         the neighborhood
0778   ///
0779   /// @note Over-/underflow bins are included in the neighborhood.
0780   /// @note The @c size parameter sets the range by how many units each local
0781   ///       bin index is allowed to be varied. All local bin indices are
0782   ///       varied independently, that is diagonal neighbors are included.
0783   ///       Ignoring the truncation of the neighborhood size reaching beyond
0784   ///       over-/underflow bins, the neighborhood is of size \f$2 \times
0785   ///       \text{size}+1\f$ along each dimension.
0786   /// @note The concrete bins which are returned depend on the WrappingTypes
0787   ///       of the contained axes
0788   ///
0789   template <class... Axes>
0790   static GlobalNeighborHoodIndices<sizeof...(Axes)> neighborHoodIndices(
0791       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0792       std::pair<std::size_t, std::size_t> sizes,
0793       const std::tuple<Axes...>& axes) {
0794     // length N array which contains local neighbors based on size par
0795     std::array<NeighborHoodIndices, sizeof...(Axes)> neighborIndices{};
0796     // get local bin indices for neighboring bins
0797     grid_helper_impl<sizeof...(Axes) - 1>::neighborHoodIndices(
0798         localIndices, sizes, axes, neighborIndices);
0799 
0800     // Query the number of bins
0801     std::array<std::size_t, sizeof...(Axes)> nBinsArray = getNBins(axes);
0802 
0803     // Produce iterator of global indices
0804     return GlobalNeighborHoodIndices(neighborIndices, nBinsArray);
0805   }
0806 
0807   template <class... Axes>
0808   static GlobalNeighborHoodIndices<sizeof...(Axes)> neighborHoodIndices(
0809       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0810       std::size_t size, const std::tuple<Axes...>& axes) {
0811     return neighborHoodIndices(
0812         localIndices, std::make_pair(-static_cast<int>(size), size), axes);
0813   }
0814 
0815   /// @brief get global bin indices for bins in specified neighborhood
0816   ///        for each axis
0817   ///
0818   /// @tparam Axes parameter pack of axis types defining the grid
0819   /// @param  [in] localIndices local bin indices along each axis
0820   /// @param  [in] size         size of neighborhood for each axis, which
0821   ///                           bins along each axis are considered
0822   /// @param  [in] axes         actual axis objects spanning the grid
0823   /// @return Sorted collection of global bin indices for all bins in
0824   ///         the neighborhood
0825   ///
0826   /// @note Over-/underflow bins are included in the neighborhood.
0827   /// @note The @c size parameter sets the range by how many units each local
0828   ///       bin index is allowed to be varied. All local bin indices are
0829   ///       varied independently, that is diagonal neighbors are included.
0830   ///       Ignoring the truncation of the neighborhood size reaching beyond
0831   ///       over-/underflow bins, the neighborhood is of size \f$2 \times
0832   ///       \text{size}+1\f$ along each dimension.
0833   /// @note The concrete bins which are returned depend on the WrappingTypes
0834   ///       of the contained axes
0835   ///
0836   template <class... Axes>
0837   static GlobalNeighborHoodIndices<sizeof...(Axes)> neighborHoodIndices(
0838       const std::array<std::size_t, sizeof...(Axes)>& localIndices,
0839       std::array<std::pair<int, int>, sizeof...(Axes)>& sizes,
0840       const std::tuple<Axes...>& axes) {
0841     // length N array which contains local neighbors based on size par
0842     std::array<NeighborHoodIndices, sizeof...(Axes)> neighborIndices{};
0843     // get local bin indices for neighboring bins
0844     grid_helper_impl<sizeof...(Axes) - 1>::neighborHoodIndices(
0845         localIndices, sizes, axes, neighborIndices);
0846 
0847     // Query the number of bins
0848     std::array<std::size_t, sizeof...(Axes)> nBinsArray = getNBins(axes);
0849 
0850     // Produce iterator of global indices
0851     return GlobalNeighborHoodIndices(neighborIndices, nBinsArray);
0852   }
0853 
0854   /// @brief get bin indices of all overflow and underflow bins
0855   ///
0856   /// @tparam Axes parameter pack of axis types defining the grid
0857   /// @param  [in] axes         actual axis objects spanning the grid
0858   /// @return set of global bin indices for all over- and underflow bins
0859   template <class... Axes>
0860   static std::set<std::size_t> exteriorBinIndices(
0861       const std::tuple<Axes...>& axes) {
0862     std::array<std::size_t, sizeof...(Axes)> idx{};
0863     std::array<bool, sizeof...(Axes)> isExterior{};
0864     std::set<std::size_t> combinations;
0865     grid_helper_impl<sizeof...(Axes) - 1>::exteriorBinIndices(
0866         idx, isExterior, combinations, axes);
0867 
0868     return combinations;
0869   }
0870 
0871   /// @brief check whether given point is inside axes limits
0872   ///
0873   /// @tparam Point any type with point semantics supporting component access
0874   ///               through @c operator[]
0875   /// @tparam Axes parameter pack of axis types defining the grid
0876   ///
0877   /// @param  [in] position point to look up in the grid
0878   /// @param  [in] axes     actual axis objects spanning the grid
0879   /// @return @c true if \f$\text{xmin_i} \le x_i < \text{xmax}_i \forall i=0,
0880   ///         \dots, d-1\f$, otherwise @c false
0881   ///
0882   /// @pre The given @c Point type must represent a point in d (or higher)
0883   ///      dimensions where d is the number of axis objects in the tuple.
0884   template <class Point, class... Axes>
0885   static bool isInside(const Point& position, const std::tuple<Axes...>& axes) {
0886     constexpr std::size_t MAX = sizeof...(Axes) - 1;
0887     return grid_helper_impl<MAX>::isInside(position, axes);
0888   }
0889 };
0890 
0891 }  // namespace Acts::detail
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!