20 #include <boost/container/small_vector.hpp>
27 constexpr
T ipow(
T num,
unsigned int pow) {
28 return (pow >=
sizeof(
unsigned int) * 8) ? 0
30 : num *
ipow(num, pow - 1);
48 std::array<NeighborHoodIndices, DIM>& neighborIndices,
49 const std::array<size_t, DIM>& nBinsArray)
54 size_t globalStride = 1;
55 for (
long i = DIM - 2;
i >= 0; --
i) {
56 globalStride *= (nBinsArray[
i + 1] + 2);
66 std::array<NeighborHoodIndices::iterator, DIM>&& localIndicesIter)
74 for (
size_t i = 0;
i < DIM - 1; ++
i) {
86 for (
long i = DIM - 1;
i > 0; --
i) {
120 std::array<NeighborHoodIndices::iterator, DIM> localIndicesIter{};
121 for (
size_t i = 0;
i < DIM; ++
i) {
132 for (
size_t i = 1;
i < DIM; ++
i) {
140 boost::container::small_vector<size_t, ipow(3, DIM)> result;
141 result.reserve(this->
size());
142 for (
size_t idx : *
this) {
143 result.push_back(idx);
150 return {result.begin(), result.end()};
167 template <
class... Axes>
171 const std::tuple<Axes...>& axes) {
176 template <
class... Axes>
178 const std::tuple<Axes...>& axes,
size_t& bin,
180 const auto& thisAxis = std::get<N>(axes);
181 bin += area * localBins.at(
N);
183 area *= (thisAxis.getNBins() + 2);
187 template <
class Point,
class... Axes>
189 const std::tuple<Axes...>& axes,
190 std::array<
size_t,
sizeof...(Axes)>& indices) {
191 const auto& thisAxis = std::get<N>(axes);
192 indices.at(
N) =
static_cast<size_t>(thisAxis.getBin(point[
N]));
196 template <
class... Axes>
199 std::array<
size_t,
sizeof...(Axes)>& indices) {
200 const auto& thisAxis = std::get<N>(axes);
202 size_t new_area = area * (thisAxis.getNBins() + 2);
204 indices.at(
N) = bin / area;
208 template <
class... Axes>
212 const std::tuple<Axes...>& axes) {
213 llEdge.at(
N) = std::get<N>(axes).getBinLowerBound(
localIndices.at(
N));
217 template <
class... Axes>
220 const std::tuple<Axes...>& axes) {
225 template <
class... Axes>
226 static void getNBins(
const std::tuple<Axes...>& axes,
227 std::array<
size_t,
sizeof...(Axes)>& nBinsArray) {
229 nBinsArray[
N] = std::get<N>(axes).
getNBins();
233 template <
class... Axes>
234 static void getAxes(
const std::tuple<Axes...>& axes,
236 axesArr[
N] =
static_cast<const IAxis*
>(&std::get<N>(axes));
240 template <
class... Axes>
244 const std::tuple<Axes...>& axes) {
245 urEdge.at(
N) = std::get<N>(axes).getBinUpperBound(
localIndices.at(
N));
249 template <
class... Axes>
252 const std::tuple<Axes...>& axes) {
257 template <
class... Axes>
258 static void getMin(
const std::tuple<Axes...>& axes,
260 minArray[
N] = std::get<N>(axes).
getMin();
264 template <
class... Axes>
265 static void getMax(
const std::tuple<Axes...>& axes,
267 maxArray[
N] = std::get<N>(axes).
getMax();
271 template <
class... Axes>
272 static void getWidth(
const std::tuple<Axes...>& axes,
274 widthArray[
N] = std::get<N>(axes).getBinWidth();
278 template <
class Point,
class... Axes>
280 bool insideThisAxis = std::get<N>(axes).
isInside(position[
N]);
284 template <
class... Axes>
287 std::pair<int, int> sizes,
const std::tuple<Axes...>& axes,
293 neighborIndices.at(
N) = locNeighbors;
299 template <
class... Axes>
302 std::array<std::pair<int, int>,
sizeof...(Axes)> sizes,
303 const std::tuple<Axes...>& axes,
309 neighborIndices.at(
N) = locNeighbors;
315 template <
class... Axes>
318 std::set<size_t>& combinations,
319 const std::tuple<Axes...>& axes) {
321 for (
size_t i = 0; i < std::get<N>(axes).
getNBins() + 2; ++
i) {
323 isExterior.at(
N) = (
i == 0) || (
i == std::get<N>(axes).getNBins() + 1);
333 template <
class... Axes>
337 const std::tuple<Axes...>& axes) {
341 template <
class... Axes>
343 const std::tuple<Axes...>& ,
size_t& bin,
345 bin += area * localBins.at(0
u);
348 template <
class Point,
class... Axes>
350 const std::tuple<Axes...>& axes,
351 std::array<
size_t,
sizeof...(Axes)>& indices) {
352 const auto& thisAxis = std::get<0u>(axes);
353 indices.at(0
u) = thisAxis.getBin(point[0
u]);
356 template <
class... Axes>
358 const std::tuple<Axes...>& ,
360 std::array<
size_t,
sizeof...(Axes)>& indices) {
362 indices.at(0
u) = bin / area;
366 template <
class... Axes>
370 const std::tuple<Axes...>& axes) {
371 llEdge.at(0
u) = std::get<0u>(axes).getBinLowerBound(
localIndices.at(0
u));
374 template <
class... Axes>
377 const std::tuple<Axes...>& axes) {
381 template <
class... Axes>
382 static void getNBins(
const std::tuple<Axes...>& axes,
383 std::array<
size_t,
sizeof...(Axes)>& nBinsArray) {
385 nBinsArray[0
u] = std::get<0u>(axes).
getNBins();
388 template <
class... Axes>
389 static void getAxes(
const std::tuple<Axes...>& axes,
391 axesArr[0
u] =
static_cast<const IAxis*
>(&std::get<0u>(axes));
394 template <
class... Axes>
398 const std::tuple<Axes...>& axes) {
399 urEdge.at(0
u) = std::get<0u>(axes).getBinUpperBound(
localIndices.at(0
u));
402 template <
class... Axes>
405 const std::tuple<Axes...>& axes) {
409 template <
class... Axes>
410 static void getMin(
const std::tuple<Axes...>& axes,
412 minArray[0
u] = std::get<0u>(axes).
getMin();
415 template <
class... Axes>
416 static void getMax(
const std::tuple<Axes...>& axes,
418 maxArray[0
u] = std::get<0u>(axes).
getMax();
421 template <
class... Axes>
422 static void getWidth(
const std::tuple<Axes...>& axes,
424 widthArray[0
u] = std::get<0u>(axes).getBinWidth();
427 template <
class Point,
class... Axes>
429 return std::get<0u>(axes).
isInside(position[0
u]);
432 template <
class... Axes>
435 std::pair<int, int> sizes,
const std::tuple<Axes...>& axes,
441 neighborIndices.at(0
u) = locNeighbors;
444 template <
class... Axes>
447 std::array<std::pair<int, int>,
sizeof...(Axes)> sizes,
448 const std::tuple<Axes...>& axes,
454 neighborIndices.at(0
u) = locNeighbors;
457 template <
class... Axes>
460 std::set<size_t>& combinations,
461 const std::tuple<Axes...>& axes) {
463 auto recordExteriorBin = [&](
size_t i) {
466 size_t bin = 0, area = 1;
468 combinations.insert(bin);
473 {
static_cast<size_t>(0), std::get<0u>(axes).getNBins() + 1}) {
474 recordExteriorBin(
i);
478 bool otherAxisExterior =
false;
479 for (
size_t N = 1;
N <
sizeof...(Axes); ++
N) {
480 otherAxisExterior = otherAxisExterior | isExterior[
N];
482 if (!otherAxisExterior) {
487 for (
size_t i = 1; i <= std::get<0u>(axes).
getNBins(); ++
i) {
488 recordExteriorBin(
i);
507 template <
class... Axes>
510 const std::tuple<Axes...>& axes) {
524 template <
class... Axes>
527 const std::tuple<Axes...>& axes) {
529 constexpr
size_t MAX =
sizeof...(Axes) - 1;
545 template <
class... Axes>
547 const std::array<
size_t,
sizeof...(Axes)>& localBins,
548 const std::tuple<Axes...>& axes) {
549 constexpr
size_t MAX =
sizeof...(Axes) - 1;
572 template <
class Point,
class... Axes>
574 const Point& point,
const std::tuple<Axes...>& axes) {
575 constexpr
size_t MAX =
sizeof...(Axes) - 1;
576 std::array<size_t,
sizeof...(Axes)> indices{};
593 template <
class... Axes>
595 size_t bin,
const std::tuple<Axes...>& axes) {
596 constexpr
size_t MAX =
sizeof...(Axes) - 1;
598 std::array<size_t,
sizeof...(Axes)> indices{};
614 template <
class... Axes>
617 const std::tuple<Axes...>& axes) {
619 constexpr
size_t MAX =
sizeof...(Axes) - 1;
638 template <
class... Axes>
641 const std::tuple<Axes...>& axes) {
642 constexpr
size_t MAX =
sizeof...(Axes) - 1;
657 template <
class... Axes>
659 const std::tuple<Axes...>& axes) {
660 std::array<size_t,
sizeof...(Axes)> nBinsArray{};
671 template <
class... Axes>
673 const std::tuple<Axes...>& axes) {
688 template <
class... Axes>
691 const std::tuple<Axes...>& axes) {
693 constexpr
size_t MAX =
sizeof...(Axes) - 1;
712 template <
class... Axes>
715 const std::tuple<Axes...>& axes) {
716 constexpr
size_t MAX =
sizeof...(Axes) - 1;
728 template <
class... Axes>
730 const std::tuple<Axes...>& axes) {
741 template <
class... Axes>
743 const std::tuple<Axes...>& axes) {
754 template <
class... Axes>
756 const std::tuple<Axes...>& axes) {
782 template <
class... Axes>
785 std::pair<size_t, size_t> sizes,
const std::tuple<Axes...>& axes) {
786 constexpr
size_t MAX =
sizeof...(Axes) - 1;
801 template <
class... Axes>
804 const std::tuple<Axes...>& axes) {
830 template <
class... Axes>
833 std::array<std::pair<int, int>,
sizeof...(Axes)>& sizes,
834 const std::tuple<Axes...>& axes) {
835 constexpr
size_t MAX =
sizeof...(Axes) - 1;
855 template <
class... Axes>
857 constexpr
size_t MAX =
sizeof...(Axes) - 1;
860 std::array<bool,
sizeof...(Axes)> isExterior{};
861 std::set<size_t> combinations;
881 template <
class Point,
class... Axes>
883 constexpr
size_t MAX =
sizeof...(Axes) - 1;