indexed.hpp 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372
  1. // Copyright 2015-2016 Hans Dembinski
  2. //
  3. // Distributed under the Boost Software License, Version 1.0.
  4. // (See accompanying file LICENSE_1_0.txt
  5. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  6. #ifndef BOOST_HISTOGRAM_INDEXED_HPP
  7. #define BOOST_HISTOGRAM_INDEXED_HPP
  8. #include <array>
  9. #include <boost/config.hpp> // BOOST_ATTRIBUTE_NODISCARD
  10. #include <boost/histogram/axis/traits.hpp>
  11. #include <boost/histogram/detail/axes.hpp>
  12. #include <boost/histogram/detail/iterator_adaptor.hpp>
  13. #include <boost/histogram/detail/operators.hpp>
  14. #include <boost/histogram/fwd.hpp>
  15. #include <iterator>
  16. #include <type_traits>
  17. #include <utility>
  18. namespace boost {
  19. namespace histogram {
  20. /** Coverage mode of the indexed range generator.
  21. Defines options for the iteration strategy.
  22. */
  23. enum class coverage {
  24. inner, /*!< iterate over inner bins, exclude underflow and overflow */
  25. all, /*!< iterate over all bins, including underflow and overflow */
  26. };
  27. /** Input iterator range over histogram bins with multi-dimensional index.
  28. The iterator returned by begin() can only be incremented. begin() may only be called
  29. once, calling it a second time returns the end() iterator. If several copies of the
  30. input iterators exist, the other copies become invalid if one of them is incremented.
  31. */
  32. template <class Histogram>
  33. class BOOST_ATTRIBUTE_NODISCARD indexed_range {
  34. private:
  35. using histogram_type = Histogram;
  36. static constexpr unsigned buffer_size =
  37. detail::buffer_size<typename std::decay_t<histogram_type>::axes_type>::value;
  38. public:
  39. using value_iterator = std::conditional_t<std::is_const<histogram_type>::value,
  40. typename histogram_type::const_iterator,
  41. typename histogram_type::iterator>;
  42. using value_reference = typename std::iterator_traits<value_iterator>::reference;
  43. using value_type = typename std::iterator_traits<value_iterator>::value_type;
  44. class iterator;
  45. using range_iterator [[deprecated("use iterator instead")]] = iterator; ///< deprecated
  46. /** Lightweight view to access value and index of current cell.
  47. The methods provide access to the current cell indices and bins. It acts like a
  48. pointer to the cell value, and in a limited way also like a reference. To interoperate
  49. with the algorithms of the standard library, the accessor is implicitly convertible to
  50. a cell value. Assignments and comparisons are passed through to the cell. An accessor
  51. is coupled to its parent indexed_range::iterator. Moving the parent iterator
  52. forward also updates the linked accessor. Accessors are not copyable. They cannot be
  53. stored in containers, but indexed_range::iterator can be stored.
  54. */
  55. class BOOST_ATTRIBUTE_NODISCARD accessor : detail::mirrored<accessor, void> {
  56. public:
  57. /// Array-like view into the current multi-dimensional index.
  58. class index_view {
  59. using index_pointer = const typename iterator::index_data*;
  60. public:
  61. using const_reference = const axis::index_type&;
  62. using reference [[deprecated("use const_reference instead")]] =
  63. const_reference; ///< deprecated
  64. /// implementation detail
  65. class const_iterator
  66. : public detail::iterator_adaptor<const_iterator, index_pointer,
  67. const_reference> {
  68. public:
  69. const_reference operator*() const noexcept { return const_iterator::base()->idx; }
  70. private:
  71. explicit const_iterator(index_pointer i) noexcept
  72. : const_iterator::iterator_adaptor_(i) {}
  73. friend class index_view;
  74. };
  75. const_iterator begin() const noexcept { return const_iterator{begin_}; }
  76. const_iterator end() const noexcept { return const_iterator{end_}; }
  77. std::size_t size() const noexcept {
  78. return static_cast<std::size_t>(end_ - begin_);
  79. }
  80. const_reference operator[](unsigned d) const noexcept { return begin_[d].idx; }
  81. const_reference at(unsigned d) const { return begin_[d].idx; }
  82. private:
  83. /// implementation detail
  84. index_view(index_pointer b, index_pointer e) : begin_(b), end_(e) {}
  85. index_pointer begin_, end_;
  86. friend class accessor;
  87. };
  88. // assignment is pass-through
  89. accessor& operator=(const accessor& o) {
  90. get() = o.get();
  91. return *this;
  92. }
  93. // assignment is pass-through
  94. template <class T>
  95. accessor& operator=(const T& x) {
  96. get() = x;
  97. return *this;
  98. }
  99. /// Returns the cell reference.
  100. value_reference get() const noexcept { return *(iter_.iter_); }
  101. /// @copydoc get()
  102. value_reference operator*() const noexcept { return get(); }
  103. /// Access fields and methods of the cell object.
  104. value_iterator operator->() const noexcept { return iter_.iter_; }
  105. /// Access current index.
  106. /// @param d axis dimension.
  107. axis::index_type index(unsigned d = 0) const noexcept {
  108. return iter_.indices_[d].idx;
  109. }
  110. /// Access indices as an iterable range.
  111. index_view indices() const noexcept {
  112. assert(iter_.indices_.hist_);
  113. return {iter_.indices_.begin(), iter_.indices_.end()};
  114. }
  115. /// Access current bin.
  116. /// @tparam N axis dimension.
  117. template <unsigned N = 0>
  118. decltype(auto) bin(std::integral_constant<unsigned, N> = {}) const {
  119. assert(iter_.indices_.hist_);
  120. return iter_.indices_.hist_->axis(std::integral_constant<unsigned, N>())
  121. .bin(index(N));
  122. }
  123. /// Access current bin.
  124. /// @param d axis dimension.
  125. decltype(auto) bin(unsigned d) const {
  126. assert(iter_.indices_.hist_);
  127. return iter_.indices_.hist_->axis(d).bin(index(d));
  128. }
  129. /** Computes density in current cell.
  130. The density is computed as the cell value divided by the product of bin widths. Axes
  131. without bin widths, like axis::category, are treated as having unit bin with.
  132. */
  133. double density() const {
  134. assert(iter_.indices_.hist_);
  135. double x = 1;
  136. unsigned d = 0;
  137. iter_.indices_.hist_->for_each_axis([&](const auto& a) {
  138. const auto w = axis::traits::width_as<double>(a, this->index(d++));
  139. x *= w != 0 ? w : 1;
  140. });
  141. return get() / x;
  142. }
  143. // forward all comparison operators to the value
  144. bool operator<(const accessor& o) noexcept { return get() < o.get(); }
  145. bool operator>(const accessor& o) noexcept { return get() > o.get(); }
  146. bool operator==(const accessor& o) noexcept { return get() == o.get(); }
  147. bool operator!=(const accessor& o) noexcept { return get() != o.get(); }
  148. bool operator<=(const accessor& o) noexcept { return get() <= o.get(); }
  149. bool operator>=(const accessor& o) noexcept { return get() >= o.get(); }
  150. template <class U>
  151. bool operator<(const U& o) const noexcept {
  152. return get() < o;
  153. }
  154. template <class U>
  155. bool operator>(const U& o) const noexcept {
  156. return get() > o;
  157. }
  158. template <class U>
  159. bool operator==(const U& o) const noexcept {
  160. return get() == o;
  161. }
  162. template <class U>
  163. bool operator!=(const U& o) const noexcept {
  164. return get() != o;
  165. }
  166. template <class U>
  167. bool operator<=(const U& o) const noexcept {
  168. return get() <= o;
  169. }
  170. template <class U>
  171. bool operator>=(const U& o) const noexcept {
  172. return get() >= o;
  173. }
  174. operator value_type() const noexcept { return get(); }
  175. private:
  176. accessor(iterator& i) noexcept : iter_(i) {}
  177. accessor(const accessor&) = default; // only callable by indexed_range::iterator
  178. iterator& iter_;
  179. friend class iterator;
  180. };
  181. /// implementation detail
  182. class iterator {
  183. public:
  184. using value_type = typename indexed_range::value_type;
  185. using reference = accessor;
  186. private:
  187. struct pointer_proxy {
  188. reference* operator->() noexcept { return std::addressof(ref_); }
  189. reference ref_;
  190. };
  191. public:
  192. using pointer = pointer_proxy;
  193. using difference_type = std::ptrdiff_t;
  194. using iterator_category = std::forward_iterator_tag;
  195. reference operator*() noexcept { return *this; }
  196. pointer operator->() noexcept { return pointer_proxy{operator*()}; }
  197. iterator& operator++() {
  198. assert(iter_ < indices_.hist_->end());
  199. const auto cbeg = indices_.begin();
  200. auto c = cbeg;
  201. ++iter_;
  202. ++c->idx;
  203. if (c->idx < c->end) return *this;
  204. while (c->idx == c->end) {
  205. iter_ += c->end_skip;
  206. if (++c == indices_.end()) return *this;
  207. ++c->idx;
  208. }
  209. while (c-- != cbeg) {
  210. c->idx = c->begin;
  211. iter_ += c->begin_skip;
  212. }
  213. return *this;
  214. }
  215. iterator operator++(int) {
  216. auto prev = *this;
  217. operator++();
  218. return prev;
  219. }
  220. bool operator==(const iterator& x) const noexcept { return iter_ == x.iter_; }
  221. bool operator!=(const iterator& x) const noexcept { return !operator==(x); }
  222. // make iterator ready for C++17 sentinels
  223. bool operator==(const value_iterator& x) const noexcept { return iter_ == x; }
  224. bool operator!=(const value_iterator& x) const noexcept { return !operator==(x); }
  225. // useful for iterator debugging
  226. std::size_t offset() const noexcept { return iter_ - indices_.hist_->begin(); }
  227. private:
  228. iterator(value_iterator i, histogram_type& h) : iter_(i), indices_(&h) {}
  229. value_iterator iter_;
  230. struct index_data {
  231. axis::index_type idx, begin, end;
  232. std::size_t begin_skip, end_skip;
  233. };
  234. struct indices_t : private std::array<index_data, buffer_size> {
  235. using base_type = std::array<index_data, buffer_size>;
  236. using pointer = index_data*;
  237. using const_pointer = const index_data*;
  238. indices_t(histogram_type* h) noexcept : hist_{h} {}
  239. using base_type::operator[];
  240. unsigned size() const noexcept { return hist_->rank(); }
  241. pointer begin() noexcept { return base_type::data(); }
  242. const_pointer begin() const noexcept { return base_type::data(); }
  243. pointer end() noexcept { return begin() + size(); }
  244. const_pointer end() const noexcept { return begin() + size(); }
  245. histogram_type* hist_;
  246. } indices_;
  247. friend class indexed_range;
  248. };
  249. indexed_range(histogram_type& hist, coverage cov)
  250. : begin_(hist.begin(), hist), end_(hist.end(), hist) {
  251. begin_.indices_.hist_->for_each_axis([ca = begin_.indices_.begin(), cov,
  252. stride = std::size_t{1},
  253. this](const auto& a) mutable {
  254. using opt = axis::traits::get_options<std::decay_t<decltype(a)>>;
  255. constexpr axis::index_type under = opt::test(axis::option::underflow);
  256. constexpr axis::index_type over = opt::test(axis::option::overflow);
  257. const auto size = a.size();
  258. // -1 if underflow and cover all, else 0
  259. ca->begin = cov == coverage::all ? -under : 0;
  260. // size + 1 if overflow and cover all, else size
  261. ca->end = cov == coverage::all ? size + over : size;
  262. ca->idx = ca->begin;
  263. // if axis has *flow and coverage::all OR axis has no *flow:
  264. // begin + under == 0, size + over - end == 0
  265. // if axis has *flow and coverage::inner:
  266. // begin + under == 1, size + over - end == 1
  267. ca->begin_skip = static_cast<std::size_t>(ca->begin + under) * stride;
  268. ca->end_skip = static_cast<std::size_t>(size + over - ca->end) * stride;
  269. begin_.iter_ += ca->begin_skip;
  270. stride *= size + under + over;
  271. ++ca;
  272. });
  273. }
  274. iterator begin() noexcept { return begin_; }
  275. iterator end() noexcept { return end_; }
  276. private:
  277. iterator begin_, end_;
  278. };
  279. /** Generates an indexed range of <a
  280. href="https://en.cppreference.com/w/cpp/named_req/ForwardIterator">forward iterators</a>
  281. over the histogram cells.
  282. Use this in a range-based for loop:
  283. ```
  284. for (auto&& x : indexed(hist)) { ... }
  285. ```
  286. This generates an optimized loop which is nearly always faster than a hand-written loop
  287. over the histogram cells. The iterators dereference to an indexed_range::accessor, which
  288. has methods to query the current indices and bins and acts like a pointer to the cell
  289. value. The returned iterators are forward iterators. They can be stored in a container,
  290. but may not be used after the life-time of the histogram ends.
  291. @returns indexed_range
  292. @param hist Reference to the histogram.
  293. @param cov Iterate over all or only inner bins (optional, default: inner).
  294. */
  295. template <class Histogram>
  296. auto indexed(Histogram&& hist, coverage cov = coverage::inner) {
  297. return indexed_range<std::remove_reference_t<Histogram>>{std::forward<Histogram>(hist),
  298. cov};
  299. }
  300. } // namespace histogram
  301. } // namespace boost
  302. #endif