Acts/Utilities/GridIterator.ipp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2016-2024 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 http://mozilla.org/MPL/2.0/.
0008
0009 namespace Acts {
0010 // Global Iterator
0011 template <typename T, class... Axes>
0012 GridGlobalIterator<T, Axes...>::GridGlobalIterator(
0013     const Acts::Grid<T, Axes...>& grid, std::size_t idx)
0014     : m_grid(&grid), m_idx(idx) {}
0015
0016 template <typename T, class... Axes>
0017 GridGlobalIterator<T, Axes...>::GridGlobalIterator(
0018     GridGlobalIterator<T, Axes...>&& other) noexcept
0019     : m_grid(std::exchange(other.m_grid.ptr, nullptr)), m_idx(other.m_idx) {}
0020
0021 template <typename T, class... Axes>
0022 GridGlobalIterator<T, Axes...>& GridGlobalIterator<T, Axes...>::operator=(
0023     GridGlobalIterator<T, Axes...>&& other) noexcept {
0024   m_grid.ptr = std::exchange(other.m_grid.ptr, nullptr);
0025   m_idx = other.m_idx;
0026   return *this;
0027 }
0028
0029 template <typename T, class... Axes>
0030 bool GridGlobalIterator<T, Axes...>::operator==(
0031     const GridGlobalIterator<T, Axes...>& other) const {
0032   // This will always return false if we are comparing two iterators from
0033   // different grids.
0034   // As such a loop from itrStart (from grid A) to itrStop (from grid B) will
0035   // never complete since itrStop will not be reachable from itrStart
0036   return (m_grid.ptr == other.m_grid.ptr) && m_idx == other.m_idx;
0037 }
0038
0039 template <typename T, class... Axes>
0040 bool GridGlobalIterator<T, Axes...>::operator!=(
0041     const GridGlobalIterator<T, Axes...>& other) const {
0042   return !(*this == other);
0043 }
0044
0045 template <typename T, class... Axes>
0046 bool GridGlobalIterator<T, Axes...>::operator<(
0047     const GridGlobalIterator<T, Axes...>& other) const {
0048   // This operator only makes sense if the two iterators we are comparing
0049   // are using the same grid
0050   assert(m_grid.ptr == other.m_grid.ptr);
0051   return m_idx < other.m_idx;
0052 }
0053
0054 template <typename T, class... Axes>
0055 bool GridGlobalIterator<T, Axes...>::operator>(
0056     const GridGlobalIterator<T, Axes...>& other) const {
0057   // This operator only makes sense if the two iterators we are comparing
0058   // are using the same grid
0059   assert(m_grid.ptr == other.m_grid.ptr);
0060   return m_idx > other.m_idx;
0061 }
0062
0063 template <typename T, class... Axes>
0064 bool GridGlobalIterator<T, Axes...>::operator<=(
0065     const GridGlobalIterator<T, Axes...>& other) const {
0066   return !(*this > other);
0067 }
0068
0069 template <typename T, class... Axes>
0070 bool GridGlobalIterator<T, Axes...>::operator>=(
0071     const GridGlobalIterator<T, Axes...>& other) const {
0072   return !(*this < other);
0073 }
0074
0075 template <typename T, class... Axes>
0076 GridGlobalIterator<T, Axes...>& GridGlobalIterator<T, Axes...>::operator+=(
0077     const std::size_t offset) {
0078   m_idx += offset;
0079   return *this;
0080 }
0081
0082 template <typename T, class... Axes>
0083 GridGlobalIterator<T, Axes...>& GridGlobalIterator<T, Axes...>::operator-=(
0084     const std::size_t offset) {
0085   m_idx -= offset;
0086   return *this;
0087 }
0088
0089 template <typename T, class... Axes>
0090 GridGlobalIterator<T, Axes...> GridGlobalIterator<T, Axes...>::operator+(
0091     const std::size_t offset) const {
0092   return {*m_grid, m_idx + offset};
0093 }
0094
0095 template <typename T, class... Axes>
0096 GridGlobalIterator<T, Axes...> GridGlobalIterator<T, Axes...>::operator-(
0097     const std::size_t offset) const {
0098   return {*m_grid, m_idx - offset};
0099 }
0100
0101 template <typename T, class... Axes>
0102 typename GridGlobalIterator<T, Axes...>::difference_type
0103 GridGlobalIterator<T, Axes...>::operator-(
0104     const GridGlobalIterator<T, Axes...>& other) const {
0105   assert(m_grid.ptr == other.m_grid.ptr);
0106   assert(other <= *this);
0107   return m_idx - other.m_idx;
0108 }
0109
0110 template <typename T, class... Axes>
0111 const typename GridGlobalIterator<T, Axes...>::value_type&
0112 GridGlobalIterator<T, Axes...>::operator*() const {
0113   return m_grid->at(m_idx);
0114 }
0115
0116 template <typename T, class... Axes>
0117 GridGlobalIterator<T, Axes...>& GridGlobalIterator<T, Axes...>::operator++() {
0118   ++m_idx;
0119   return *this;
0120 }
0121
0122 template <typename T, class... Axes>
0123 GridGlobalIterator<T, Axes...> GridGlobalIterator<T, Axes...>::operator++(int) {
0124   GridGlobalIterator<T, Axes...> output(*m_grid, m_idx++);
0125   return output;
0126 }
0127
0128 template <typename T, class... Axes>
0129 std::size_t GridGlobalIterator<T, Axes...>::globalBinIndex() const {
0130   return m_idx;
0131 }
0132
0133 template <typename T, class... Axes>
0134 std::array<std::size_t, GridGlobalIterator<T, Axes...>::DIM>
0135 GridGlobalIterator<T, Axes...>::localBinsIndices() const {
0136   return m_grid->localBinsFromGlobalBin(m_idx);
0137 }
0138
0139 // Local Iterator
0140 template <typename T, class... Axes>
0141 Acts::GridLocalIterator<T, Axes...>::GridLocalIterator(
0142     const Acts::Grid<T, Axes...>& grid,
0143     const std::array<std::size_t, DIM>& indices)
0144     : m_grid(&grid),
0145       m_numLocalBins(grid.numLocalBins()),
0146       m_currentIndex(indices) {
0147   // Since the user has not defined a custom navigation pattern, we tell the
0148   // iterator we want to iterate on all the local bins in ascending order from
0149   // 1ul to numLocalBin for that specific axis.
0150   for (std::size_t i(0); i < DIM; ++i) {
0151     m_navigationIndex[i].resize(m_numLocalBins[i]);
0152     std::iota(m_navigationIndex[i].begin(), m_navigationIndex[i].end(), 1ul);
0153   }
0154 }
0155
0156 template <typename T, class... Axes>
0157 Acts::GridLocalIterator<T, Axes...>::GridLocalIterator(
0158     const Acts::Grid<T, Axes...>& grid,
0159     const std::array<std::size_t, DIM>& indices,
0160     std::array<std::vector<std::size_t>, DIM> navigation)
0161     : m_grid(&grid),
0162       m_numLocalBins(grid.numLocalBins()),
0163       m_currentIndex(indices),
0164       m_navigationIndex(std::move(navigation)) {
0165   /// We can allow navigation on only a subset of bins.
0166   /// If the number of specified bins in the navigation for one axis is not
0167   /// zero then override the maximum number of navigation bins instead of using
0168   /// the total number of available bins in the axis
0169   for (std::size_t i(0ul); i < DIM; ++i) {
0170     /// We do not allow empty bin sequences
0171     if (m_navigationIndex[i].size() == 0) {
0172       throw std::invalid_argument(
0173           "Invalid navigation sequence in local grid iterator. No bins "
0174           "specified.");
0175     }
0176     /// Too many bins
0177     if (m_navigationIndex[i].size() > m_numLocalBins[i]) {
0178       throw std::invalid_argument(
0179           "Invalid navigation sequence in local grid iterator. Too many bins "
0180           "specified.");
0181     }
0182     m_numLocalBins[i] = m_navigationIndex[i].size();
0183   }
0184 }
0185
0186 template <typename T, class... Axes>
0187 Acts::GridLocalIterator<T, Axes...>::GridLocalIterator(
0188     Acts::GridLocalIterator<T, Axes...>&& other) noexcept
0189     : m_grid(std::exchange(other.m_grid.ptr, nullptr)),
0190       m_numLocalBins(other.m_numLocalBins),
0191       m_currentIndex(other.m_currentIndex),
0192       m_navigationIndex(std::move(other.m_navigationIndex)) {}
0193
0194 template <typename T, class... Axes>
0195 Acts::GridLocalIterator<T, Axes...>&
0196 Acts::GridLocalIterator<T, Axes...>::operator=(
0197     Acts::GridLocalIterator<T, Axes...>&& other) noexcept {
0198   m_grid.ptr = std::exchange(other.m_grid.ptr, nullptr);
0199   m_numLocalBins = other.m_numLocalBins;
0200   m_currentIndex = other.m_currentIndex;
0201   m_navigationIndex = std::move(other.m_navigationIndex);
0202   return *this;
0203 }
0204
0205 template <typename T, class... Axes>
0206 bool Acts::GridLocalIterator<T, Axes...>::operator==(
0207     const Acts::GridLocalIterator<T, Axes...>& other) const {
0208   // This will always return false if we are comparing two iterators from
0209   // different grids.
0210   // As such a loop from itrStart (from grid A) to itrStop (from grid B) will
0211   // never complete since itrStop will not be reachable from itrStart
0212   if (m_grid.ptr != other.m_grid.ptr) {
0213     return false;
0214   }
0215
0216   for (std::size_t i(0); i < DIM; ++i) {
0217     if (m_currentIndex[i] != other.m_currentIndex[i]) {
0218       return false;
0219     }
0220   }
0221
0222   return true;
0223 }
0224
0225 template <typename T, class... Axes>
0226 bool Acts::GridLocalIterator<T, Axes...>::operator!=(
0227     const Acts::GridLocalIterator<T, Axes...>& other) const {
0228   return !(*this == other);
0229 }
0230
0231 template <typename T, class... Axes>
0232 const typename Acts::GridLocalIterator<T, Axes...>::value_type&
0233 Acts::GridLocalIterator<T, Axes...>::operator*() const {
0234   std::array<std::size_t, DIM> localPositionBin{};
0235   for (std::size_t i(0); i < DIM; ++i) {
0236     localPositionBin[i] = m_navigationIndex[i][m_currentIndex[i]];
0237   }
0238   return m_grid->atLocalBins(localPositionBin);
0239 }
0240
0241 template <typename T, class... Axes>
0242 GridLocalIterator<T, Axes...>& GridLocalIterator<T, Axes...>::operator++() {
0243   increment<DIM - 1>();
0244   return *this;
0245 }
0246
0247 template <typename T, class... Axes>
0248 GridLocalIterator<T, Axes...> GridLocalIterator<T, Axes...>::operator++(int) {
0249   GridLocalIterator<T, Axes...> output(*this);
0250   this->operator++();
0251   return output;
0252 }
0253
0254 template <typename T, class... Axes>
0255 template <std::size_t N>
0256 void GridLocalIterator<T, Axes...>::increment() {
0257   // Check if the current local bin can be incremented, or we reached the end
0258   // of bins in the axis
0259   if (++m_currentIndex[N] < m_numLocalBins[N]) {
0260     return;
0261   }
0262   // We have reached the last bin in the axis, we set the position to 0ul and
0263   // try to increment another axis
0264   if constexpr (N != 0) {
0265     m_currentIndex[N] = 0;
0266     increment<N - 1>();
0267   } else {
0268     m_currentIndex = m_numLocalBins;
0269   }
0270 }
0271
0272 template <typename T, class... Axes>
0273 std::size_t GridLocalIterator<T, Axes...>::globalBinIndex() const {
0274   return m_grid->globalBinFromLocalBins(localBinsIndices());
0275 }
0276
0277 template <typename T, class... Axes>
0278 std::array<std::size_t, GridLocalIterator<T, Axes...>::DIM>
0279 GridLocalIterator<T, Axes...>::localBinsIndices() const {
0280   std::array<std::size_t, DIM> output{};
0281   for (std::size_t i(0); i < DIM; ++i) {
0282     output[i] = m_navigationIndex[i][m_currentIndex[i]];
0283   }
0284   return output;
0285 }
0286
0287 }  // namespace Acts