EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
grid_helper.hpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file grid_helper.hpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2017-2018 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #pragma once
10 
11 #include "Acts/Utilities/IAxis.hpp"
12 #include "Acts/Utilities/detail/Axis.hpp"
13 
14 #include <array>
15 #include <set>
16 #include <tuple>
17 #include <utility>
18 
19 namespace Acts {
20 
21 namespace detail {
22 
23 // This object can be iterated to produce the (ordered) set of global indices
24 // associated with a neighborhood around a certain point on a grid.
25 //
26 // The goal is to emulate the effect of enumerating the global indices into
27 // an std::set (or into an std::vector that gets subsequently sorted), without
28 // paying the price of dynamic memory allocation in hot magnetic field
29 // interpolation code.
30 //
31 template <size_t DIM>
32 class GlobalNeighborHoodIndices {
33  public:
34  // You can get the local neighbor indices from
35  // grid_helper_impl<DIM>::neighborHoodIndices and the number of bins in
36  // each direction from grid_helper_impl<DIM>::getNBins.
37  GlobalNeighborHoodIndices(
38  std::array<NeighborHoodIndices, DIM>& neighborIndices,
39  const std::array<size_t, DIM>& nBinsArray)
40  : m_localIndices(neighborIndices) {
41  if (DIM == 1)
42  return;
43  size_t globalStride = 1;
44  for (long i = DIM - 2; i >= 0; --i) {
45  globalStride *= (nBinsArray[i + 1] + 2);
46  m_globalStrides[i] = globalStride;
47  }
48  }
49 
50  class iterator {
51  public:
52  iterator() = default;
53 
54  iterator(const GlobalNeighborHoodIndices& parent,
55  std::array<NeighborHoodIndices::iterator, DIM>&& localIndicesIter)
56  : m_localIndicesIter(std::move(localIndicesIter)), m_parent(&parent) {}
57 
58  size_t operator*() const {
59  size_t globalIndex = *m_localIndicesIter[DIM - 1];
60  if (DIM == 1)
61  return globalIndex;
62  for (size_t i = 0; i < DIM - 1; ++i) {
63  globalIndex += m_parent->m_globalStrides[i] * (*m_localIndicesIter[i]);
64  }
65  return globalIndex;
66  }
67 
68  iterator& operator++() {
69  const auto& localIndices = m_parent->m_localIndices;
70 
71  // Go to the next global index via a lexicographic increment:
72  // - Start by incrementing the last local index
73  // - If it reaches the end, reset it and increment the previous one...
74  for (long i = DIM - 1; i > 0; --i) {
75  ++m_localIndicesIter[i];
76  if (m_localIndicesIter[i] != localIndices[i].end())
77  return *this;
78  m_localIndicesIter[i] = localIndices[i].begin();
79  }
80 
81  // The first index should stay at the end value when it reaches it, so
82  // that we know when we've reached the end of iteration.
83  ++m_localIndicesIter[0];
84  return *this;
85  }
86 
87  bool operator==(const iterator& it) {
88  // We know when we've reached the end, so we don't need an end-iterator.
89  // Sadly, in C++, there has to be one. Therefore, we special-case it
90  // heavily so that it's super-efficient to create and compare to.
91  if (it.m_parent == nullptr) {
92  return m_localIndicesIter[0] == m_parent->m_localIndices[0].end();
93  } else {
94  return m_localIndicesIter == it.m_localIndicesIter;
95  }
96  }
97 
98  bool operator!=(const iterator& it) { return !(*this == it); }
99 
100  private:
101  std::array<NeighborHoodIndices::iterator, DIM> m_localIndicesIter;
102  const GlobalNeighborHoodIndices* m_parent = nullptr;
103  };
104 
105  iterator begin() const {
106  std::array<NeighborHoodIndices::iterator, DIM> localIndicesIter;
107  for (size_t i = 0; i < DIM; ++i) {
108  localIndicesIter[i] = m_localIndices[i].begin();
109  }
110  return iterator(*this, std::move(localIndicesIter));
111  }
112 
113  iterator end() const { return iterator(); }
114 
115  // Number of indices that will be produced if this sequence is iterated
116  size_t size() const {
117  size_t result = m_localIndices[0].size();
118  for (size_t i = 1; i < DIM; ++i) {
119  result *= m_localIndices[i].size();
120  }
121  return result;
122  }
123 
124  // Collect the sequence of indices into an std::vector
125  std::vector<size_t> collect() const {
126  std::vector<size_t> result;
127  result.reserve(this->size());
128  for (size_t idx : *this) {
129  result.push_back(idx);
130  }
131  return result;
132  }
133 
134  private:
135  std::array<NeighborHoodIndices, DIM> m_localIndices;
136  std::array<size_t, DIM - 1> m_globalStrides;
137 };
138 
143 template <size_t N>
144 struct grid_helper_impl;
145 
146 template <size_t N>
147 struct grid_helper_impl {
148  template <class... Axes>
149  static void getBinCenter(
150  std::array<double, sizeof...(Axes)>& center,
151  const std::array<size_t, sizeof...(Axes)>& localIndices,
152  const std::tuple<Axes...>& axes) {
153  center.at(N) = std::get<N>(axes).getBinCenter(localIndices.at(N));
154  grid_helper_impl<N - 1>::getBinCenter(center, localIndices, axes);
155  }
156 
157  template <class... Axes>
158  static void getGlobalBin(const std::array<size_t, sizeof...(Axes)>& localBins,
159  const std::tuple<Axes...>& axes, size_t& bin,
160  size_t& area) {
161  const auto& thisAxis = std::get<N>(axes);
162  bin += area * localBins.at(N);
163  // make sure to account for under-/overflow bins
164  area *= (thisAxis.getNBins() + 2);
165  grid_helper_impl<N - 1>::getGlobalBin(localBins, axes, bin, area);
166  }
167 
168  template <class Point, class... Axes>
169  static void getLocalBinIndices(const Point& point,
170  const std::tuple<Axes...>& axes,
171  std::array<size_t, sizeof...(Axes)>& indices) {
172  const auto& thisAxis = std::get<N>(axes);
173  indices.at(N) = thisAxis.getBin(point[N]);
174  grid_helper_impl<N - 1>::getLocalBinIndices(point, axes, indices);
175  }
176 
177  template <class... Axes>
178  static void getLocalBinIndices(size_t& bin, const std::tuple<Axes...>& axes,
179  size_t& area,
180  std::array<size_t, sizeof...(Axes)>& indices) {
181  const auto& thisAxis = std::get<N>(axes);
182  // make sure to account for under-/overflow bins
183  size_t new_area = area * (thisAxis.getNBins() + 2);
184  grid_helper_impl<N - 1>::getLocalBinIndices(bin, axes, new_area, indices);
185  indices.at(N) = bin / area;
186  bin %= area;
187  }
188 
189  template <class... Axes>
190  static void getLowerLeftBinEdge(
191  std::array<double, sizeof...(Axes)>& llEdge,
192  const std::array<size_t, sizeof...(Axes)>& localIndices,
193  const std::tuple<Axes...>& axes) {
194  llEdge.at(N) = std::get<N>(axes).getBinLowerBound(localIndices.at(N));
195  grid_helper_impl<N - 1>::getLowerLeftBinEdge(llEdge, localIndices, axes);
196  }
197 
198  template <class... Axes>
199  static void getLowerLeftBinIndices(
200  std::array<size_t, sizeof...(Axes)>& localIndices,
201  const std::tuple<Axes...>& axes) {
202  localIndices.at(N) = std::get<N>(axes).wrapBin(localIndices.at(N) - 1);
203  grid_helper_impl<N - 1>::getLowerLeftBinIndices(localIndices, axes);
204  }
205 
206  template <class... Axes>
207  static void getNBins(const std::tuple<Axes...>& axes,
208  std::array<size_t, sizeof...(Axes)>& nBinsArray) {
209  // by convention getNBins does not include under-/overflow bins
210  nBinsArray[N] = std::get<N>(axes).getNBins();
211  grid_helper_impl<N - 1>::getNBins(axes, nBinsArray);
212  }
213 
214  template <class... Axes>
215  static void getAxes(const std::tuple<Axes...>& axes,
216  std::array<const IAxis*, sizeof...(Axes)>& axesArr) {
217  axesArr[N] = static_cast<const IAxis*>(&std::get<N>(axes));
218  grid_helper_impl<N - 1>::getAxes(axes, axesArr);
219  }
220 
221  template <class... Axes>
222  static void getUpperRightBinEdge(
223  std::array<double, sizeof...(Axes)>& urEdge,
224  const std::array<size_t, sizeof...(Axes)>& localIndices,
225  const std::tuple<Axes...>& axes) {
226  urEdge.at(N) = std::get<N>(axes).getBinUpperBound(localIndices.at(N));
227  grid_helper_impl<N - 1>::getUpperRightBinEdge(urEdge, localIndices, axes);
228  }
229 
230  template <class... Axes>
231  static void getUpperRightBinIndices(
232  std::array<size_t, sizeof...(Axes)>& localIndices,
233  const std::tuple<Axes...>& axes) {
234  localIndices.at(N) = std::get<N>(axes).wrapBin(localIndices.at(N) + 1);
235  grid_helper_impl<N - 1>::getUpperRightBinIndices(localIndices, axes);
236  }
237 
238  template <class... Axes>
239  static void getMin(const std::tuple<Axes...>& axes,
240  std::array<double, sizeof...(Axes)>& minArray) {
241  minArray[N] = std::get<N>(axes).getMin();
242  grid_helper_impl<N - 1>::getMin(axes, minArray);
243  }
244 
245  template <class... Axes>
246  static void getMax(const std::tuple<Axes...>& axes,
247  std::array<double, sizeof...(Axes)>& maxArray) {
248  maxArray[N] = std::get<N>(axes).getMax();
249  grid_helper_impl<N - 1>::getMax(axes, maxArray);
250  }
251 
252  template <class... Axes>
253  static void getWidth(const std::tuple<Axes...>& axes,
254  std::array<double, sizeof...(Axes)>& widthArray) {
255  widthArray[N] = std::get<N>(axes).getBinWidth();
256  grid_helper_impl<N - 1>::getWidth(axes, widthArray);
257  }
258 
259  template <class Point, class... Axes>
260  static bool isInside(const Point& position, const std::tuple<Axes...>& axes) {
261  bool insideThisAxis = std::get<N>(axes).isInside(position[N]);
262  return insideThisAxis && grid_helper_impl<N - 1>::isInside(position, axes);
263  }
264 
265  template <class... Axes>
266  static void neighborHoodIndices(
267  const std::array<size_t, sizeof...(Axes)>& localIndices,
268  std::pair<size_t, size_t> sizes, const std::tuple<Axes...>& axes,
269  std::array<NeighborHoodIndices, sizeof...(Axes)>& neighborIndices) {
270  // ask n-th axis
271  size_t locIdx = localIndices.at(N);
272  NeighborHoodIndices locNeighbors =
273  std::get<N>(axes).neighborHoodIndices(locIdx, sizes);
274  neighborIndices.at(N) = locNeighbors;
275 
276  grid_helper_impl<N - 1>::neighborHoodIndices(localIndices, sizes, axes,
277  neighborIndices);
278  }
279 
280  template <class... Axes>
281  static void exteriorBinIndices(std::array<size_t, sizeof...(Axes)>& idx,
282  std::array<bool, sizeof...(Axes)> isExterior,
283  std::set<size_t>& combinations,
284  const std::tuple<Axes...>& axes) {
285  // iterate over this axis' bins, remembering which bins are exterior
286  for (size_t i = 0; i < std::get<N>(axes).getNBins() + 2; ++i) {
287  idx.at(N) = i;
288  isExterior.at(N) = (i == 0) || (i == std::get<N>(axes).getNBins() + 1);
289  // vary other axes recursively
290  grid_helper_impl<N - 1>::exteriorBinIndices(idx, isExterior, combinations,
291  axes);
292  }
293  }
294 };
295 
296 template <>
297 struct grid_helper_impl<0u> {
298  template <class... Axes>
299  static void getBinCenter(
300  std::array<double, sizeof...(Axes)>& center,
301  const std::array<size_t, sizeof...(Axes)>& localIndices,
302  const std::tuple<Axes...>& axes) {
303  center.at(0u) = std::get<0u>(axes).getBinCenter(localIndices.at(0u));
304  }
305 
306  template <class... Axes>
307  static void getGlobalBin(const std::array<size_t, sizeof...(Axes)>& localBins,
308  const std::tuple<Axes...>& /*axes*/, size_t& bin,
309  size_t& area) {
310  bin += area * localBins.at(0u);
311  }
312 
313  template <class Point, class... Axes>
314  static void getLocalBinIndices(const Point& point,
315  const std::tuple<Axes...>& axes,
316  std::array<size_t, sizeof...(Axes)>& indices) {
317  const auto& thisAxis = std::get<0u>(axes);
318  indices.at(0u) = thisAxis.getBin(point[0u]);
319  }
320 
321  template <class... Axes>
322  static void getLocalBinIndices(size_t& bin,
323  const std::tuple<Axes...>& /*axes*/,
324  size_t& area,
325  std::array<size_t, sizeof...(Axes)>& indices) {
326  // make sure to account for under-/overflow bins
327  indices.at(0u) = bin / area;
328  bin %= area;
329  }
330 
331  template <class... Axes>
332  static void getLowerLeftBinEdge(
333  std::array<double, sizeof...(Axes)>& llEdge,
334  const std::array<size_t, sizeof...(Axes)>& localIndices,
335  const std::tuple<Axes...>& axes) {
336  llEdge.at(0u) = std::get<0u>(axes).getBinLowerBound(localIndices.at(0u));
337  }
338 
339  template <class... Axes>
340  static void getLowerLeftBinIndices(
341  std::array<size_t, sizeof...(Axes)>& localIndices,
342  const std::tuple<Axes...>& axes) {
343  localIndices.at(0u) = std::get<0u>(axes).wrapBin(localIndices.at(0u) - 1);
344  }
345 
346  template <class... Axes>
347  static void getNBins(const std::tuple<Axes...>& axes,
348  std::array<size_t, sizeof...(Axes)>& nBinsArray) {
349  // by convention getNBins does not include under-/overflow bins
350  nBinsArray[0u] = std::get<0u>(axes).getNBins();
351  }
352 
353  template <class... Axes>
354  static void getAxes(const std::tuple<Axes...>& axes,
355  std::array<const IAxis*, sizeof...(Axes)>& axesArr) {
356  axesArr[0u] = static_cast<const IAxis*>(&std::get<0u>(axes));
357  }
358 
359  template <class... Axes>
360  static void getUpperRightBinEdge(
361  std::array<double, sizeof...(Axes)>& urEdge,
362  const std::array<size_t, sizeof...(Axes)>& localIndices,
363  const std::tuple<Axes...>& axes) {
364  urEdge.at(0u) = std::get<0u>(axes).getBinUpperBound(localIndices.at(0u));
365  }
366 
367  template <class... Axes>
368  static void getUpperRightBinIndices(
369  std::array<size_t, sizeof...(Axes)>& localIndices,
370  const std::tuple<Axes...>& axes) {
371  localIndices.at(0u) = std::get<0u>(axes).wrapBin(localIndices.at(0u) + 1);
372  }
373 
374  template <class... Axes>
375  static void getMin(const std::tuple<Axes...>& axes,
376  std::array<double, sizeof...(Axes)>& minArray) {
377  minArray[0u] = std::get<0u>(axes).getMin();
378  }
379 
380  template <class... Axes>
381  static void getMax(const std::tuple<Axes...>& axes,
382  std::array<double, sizeof...(Axes)>& maxArray) {
383  maxArray[0u] = std::get<0u>(axes).getMax();
384  }
385 
386  template <class... Axes>
387  static void getWidth(const std::tuple<Axes...>& axes,
388  std::array<double, sizeof...(Axes)>& widthArray) {
389  widthArray[0u] = std::get<0u>(axes).getBinWidth();
390  }
391 
392  template <class Point, class... Axes>
393  static bool isInside(const Point& position, const std::tuple<Axes...>& axes) {
394  return std::get<0u>(axes).isInside(position[0u]);
395  }
396 
397  template <class... Axes>
398  static void neighborHoodIndices(
399  const std::array<size_t, sizeof...(Axes)>& localIndices,
400  std::pair<size_t, size_t> sizes, const std::tuple<Axes...>& axes,
401  std::array<NeighborHoodIndices, sizeof...(Axes)>& neighborIndices) {
402  // ask 0-th axis
403  size_t locIdx = localIndices.at(0u);
404  NeighborHoodIndices locNeighbors =
405  std::get<0u>(axes).neighborHoodIndices(locIdx, sizes);
406  neighborIndices.at(0u) = locNeighbors;
407  }
408 
409  template <class... Axes>
410  static void exteriorBinIndices(std::array<size_t, sizeof...(Axes)>& idx,
411  std::array<bool, sizeof...(Axes)> isExterior,
412  std::set<size_t>& combinations,
413  const std::tuple<Axes...>& axes) {
414  // For each exterior bin on this axis, we will do this
415  auto recordExteriorBin = [&](size_t i) {
416  idx.at(0u) = i;
417  // at this point, combinations are complete: save the global bin
418  size_t bin = 0, area = 1;
419  grid_helper_impl<sizeof...(Axes) - 1>::getGlobalBin(idx, axes, bin, area);
420  combinations.insert(bin);
421  };
422 
423  // The first and last bins on this axis are exterior by definition
424  for (size_t i :
425  {static_cast<size_t>(0), std::get<0u>(axes).getNBins() + 1}) {
426  recordExteriorBin(i);
427  }
428 
429  // If no other axis is on an exterior index, stop here
430  bool otherAxisExterior = false;
431  for (size_t N = 1; N < sizeof...(Axes); ++N) {
432  otherAxisExterior = otherAxisExterior | isExterior[N];
433  }
434  if (!otherAxisExterior) {
435  return;
436  }
437 
438  // Otherwise, we're on a grid border: iterate over all the other indices
439  for (size_t i = 1; i <= std::get<0u>(axes).getNBins(); ++i) {
440  recordExteriorBin(i);
441  }
442  }
443 };
445 
447 struct grid_helper {
459  template <class... Axes>
460  static GlobalNeighborHoodIndices<sizeof...(Axes)> closestPointsIndices(
461  const std::array<size_t, sizeof...(Axes)>& localIndices,
462  const std::tuple<Axes...>& axes) {
463  // get neighboring bins, but only increment.
464  return neighborHoodIndices(localIndices, std::make_pair(0, 1), axes);
465  }
466 
476  template <class... Axes>
477  static std::array<double, sizeof...(Axes)> getBinCenter(
478  const std::array<size_t, sizeof...(Axes)>& localIndices,
479  const std::tuple<Axes...>& axes) {
480  std::array<double, sizeof...(Axes)> center;
481  constexpr size_t MAX = sizeof...(Axes) - 1;
482  grid_helper_impl<MAX>::getBinCenter(center, localIndices, axes);
483 
484  return center;
485  }
486 
497  template <class... Axes>
498  static size_t getGlobalBin(
499  const std::array<size_t, sizeof...(Axes)>& localBins,
500  const std::tuple<Axes...>& axes) {
501  constexpr size_t MAX = sizeof...(Axes) - 1;
502  size_t area = 1;
503  size_t bin = 0;
504 
505  grid_helper_impl<MAX>::getGlobalBin(localBins, axes, bin, area);
506 
507  return bin;
508  }
509 
524  template <class Point, class... Axes>
525  static std::array<size_t, sizeof...(Axes)> getLocalBinIndices(
526  const Point& point, const std::tuple<Axes...>& axes) {
527  constexpr size_t MAX = sizeof...(Axes) - 1;
528  std::array<size_t, sizeof...(Axes)> indices;
529 
530  grid_helper_impl<MAX>::getLocalBinIndices(point, axes, indices);
531 
532  return indices;
533  }
534 
545  template <class... Axes>
546  static std::array<size_t, sizeof...(Axes)> getLocalBinIndices(
547  size_t bin, const std::tuple<Axes...>& axes) {
548  constexpr size_t MAX = sizeof...(Axes) - 1;
549  size_t area = 1;
550  std::array<size_t, sizeof...(Axes)> indices;
551 
552  grid_helper_impl<MAX>::getLocalBinIndices(bin, axes, area, indices);
553 
554  return indices;
555  }
556 
566  template <class... Axes>
567  static std::array<double, sizeof...(Axes)> getLowerLeftBinEdge(
568  const std::array<size_t, sizeof...(Axes)>& localIndices,
569  const std::tuple<Axes...>& axes) {
570  std::array<double, sizeof...(Axes)> llEdge;
571  constexpr size_t MAX = sizeof...(Axes) - 1;
572  grid_helper_impl<MAX>::getLowerLeftBinEdge(llEdge, localIndices, axes);
573 
574  return llEdge;
575  }
576 
590  template <class... Axes>
591  static std::array<size_t, sizeof...(Axes)> getLowerLeftBinIndices(
592  const std::array<size_t, sizeof...(Axes)>& localIndices,
593  const std::tuple<Axes...>& axes) {
594  constexpr size_t MAX = sizeof...(Axes) - 1;
595  auto llIndices = localIndices;
596  grid_helper_impl<MAX>::getLowerLeftBinIndices(llIndices, axes);
597 
598  return llIndices;
599  }
600 
609  template <class... Axes>
610  static std::array<size_t, sizeof...(Axes)> getNBins(
611  const std::tuple<Axes...>& axes) {
612  std::array<size_t, sizeof...(Axes)> nBinsArray;
613  grid_helper_impl<sizeof...(Axes) - 1>::getNBins(axes, nBinsArray);
614  return nBinsArray;
615  }
616 
623  template <class... Axes>
624  static std::array<const IAxis*, sizeof...(Axes)> getAxes(
625  const std::tuple<Axes...>& axes) {
626  std::array<const IAxis*, sizeof...(Axes)> arr;
627  grid_helper_impl<sizeof...(Axes) - 1>::getAxes(axes, arr);
628  return arr;
629  }
630 
640  template <class... Axes>
641  static std::array<double, sizeof...(Axes)> getUpperRightBinEdge(
642  const std::array<size_t, sizeof...(Axes)>& localIndices,
643  const std::tuple<Axes...>& axes) {
644  std::array<double, sizeof...(Axes)> urEdge;
645  constexpr size_t MAX = sizeof...(Axes) - 1;
646  grid_helper_impl<MAX>::getUpperRightBinEdge(urEdge, localIndices, axes);
647 
648  return urEdge;
649  }
650 
664  template <class... Axes>
665  static std::array<size_t, sizeof...(Axes)> getUpperRightBinIndices(
666  const std::array<size_t, sizeof...(Axes)>& localIndices,
667  const std::tuple<Axes...>& axes) {
668  constexpr size_t MAX = sizeof...(Axes) - 1;
669  auto urIndices = localIndices;
670  grid_helper_impl<MAX>::getUpperRightBinIndices(urIndices, axes);
671 
672  return urIndices;
673  }
674 
680  template <class... Axes>
681  static std::array<double, sizeof...(Axes)> getMin(
682  const std::tuple<Axes...>& axes) {
683  std::array<double, sizeof...(Axes)> minArray;
684  grid_helper_impl<sizeof...(Axes) - 1>::getMin(axes, minArray);
685  return minArray;
686  }
687 
693  template <class... Axes>
694  static std::array<double, sizeof...(Axes)> getMax(
695  const std::tuple<Axes...>& axes) {
696  std::array<double, sizeof...(Axes)> maxArray;
697  grid_helper_impl<sizeof...(Axes) - 1>::getMax(axes, maxArray);
698  return maxArray;
699  }
700 
706  template <class... Axes>
707  static std::array<double, sizeof...(Axes)> getWidth(
708  const std::tuple<Axes...>& axes) {
709  std::array<double, sizeof...(Axes)> widthArray;
710  grid_helper_impl<sizeof...(Axes) - 1>::getWidth(axes, widthArray);
711  return widthArray;
712  }
713 
734  template <class... Axes>
735  static GlobalNeighborHoodIndices<sizeof...(Axes)> neighborHoodIndices(
736  const std::array<size_t, sizeof...(Axes)>& localIndices,
737  std::pair<size_t, size_t> sizes, const std::tuple<Axes...>& axes) {
738  constexpr size_t MAX = sizeof...(Axes) - 1;
739 
740  // length N array which contains local neighbors based on size par
741  std::array<NeighborHoodIndices, sizeof...(Axes)> neighborIndices;
742  // get local bin indices for neighboring bins
743  grid_helper_impl<MAX>::neighborHoodIndices(localIndices, sizes, axes,
744  neighborIndices);
745 
746  // Query the number of bins
747  std::array<size_t, sizeof...(Axes)> nBinsArray = getNBins(axes);
748 
749  // Produce iterator of global indices
750  return GlobalNeighborHoodIndices(neighborIndices, nBinsArray);
751  }
752 
753  template <class... Axes>
754  static GlobalNeighborHoodIndices<sizeof...(Axes)> neighborHoodIndices(
755  const std::array<size_t, sizeof...(Axes)>& localIndices, size_t size,
756  const std::tuple<Axes...>& axes) {
757  return neighborHoodIndices(localIndices, std::make_pair(size, size), axes);
758  }
759 
765  template <class... Axes>
766  static std::set<size_t> exteriorBinIndices(const std::tuple<Axes...>& axes) {
767  constexpr size_t MAX = sizeof...(Axes) - 1;
768 
769  std::array<size_t, sizeof...(Axes)> idx;
770  std::array<bool, sizeof...(Axes)> isExterior;
771  std::set<size_t> combinations;
772  grid_helper_impl<MAX>::exteriorBinIndices(idx, isExterior, combinations,
773  axes);
774 
775  return combinations;
776  }
777 
791  template <class Point, class... Axes>
792  static bool isInside(const Point& position, const std::tuple<Axes...>& axes) {
793  constexpr size_t MAX = sizeof...(Axes) - 1;
794  return grid_helper_impl<MAX>::isInside(position, axes);
795  }
796 };
797 
798 } // namespace detail
799 
800 } // namespace Acts