axes.hpp 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385
  1. // Copyright 2015-2018 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_DETAIL_AXES_HPP
  7. #define BOOST_HISTOGRAM_DETAIL_AXES_HPP
  8. #include <array>
  9. #include <boost/assert.hpp>
  10. #include <boost/core/nvp.hpp>
  11. #include <boost/histogram/axis/traits.hpp>
  12. #include <boost/histogram/axis/variant.hpp>
  13. #include <boost/histogram/detail/make_default.hpp>
  14. #include <boost/histogram/detail/optional_index.hpp>
  15. #include <boost/histogram/detail/static_if.hpp>
  16. #include <boost/histogram/fwd.hpp>
  17. #include <boost/mp11/algorithm.hpp>
  18. #include <boost/mp11/list.hpp>
  19. #include <boost/mp11/tuple.hpp>
  20. #include <boost/mp11/utility.hpp>
  21. #include <boost/throw_exception.hpp>
  22. #include <stdexcept>
  23. #include <string>
  24. #include <tuple>
  25. #include <type_traits>
  26. /* Most of the histogram code is generic and works for any number of axes. Buffers with a
  27. * fixed maximum capacity are used in some places, which have a size equal to the rank of
  28. * a histogram. The buffers are statically allocated to improve performance, which means
  29. * that they need a preset maximum capacity. 32 seems like a safe upper limit for the rank
  30. * (you can nevertheless increase it here if necessary): the simplest non-trivial axis has
  31. * 2 bins; even if counters are used which need only a byte of storage per bin, this still
  32. * corresponds to 4 GB of storage.
  33. */
  34. #ifndef BOOST_HISTOGRAM_DETAIL_AXES_LIMIT
  35. #define BOOST_HISTOGRAM_DETAIL_AXES_LIMIT 32
  36. #endif
  37. namespace boost {
  38. namespace histogram {
  39. namespace detail {
  40. template <class T>
  41. unsigned axes_rank(const T& axes) {
  42. using std::begin;
  43. using std::end;
  44. return static_cast<unsigned>(std::distance(begin(axes), end(axes)));
  45. }
  46. template <class... Ts>
  47. constexpr unsigned axes_rank(const std::tuple<Ts...>&) {
  48. return static_cast<unsigned>(sizeof...(Ts));
  49. }
  50. template <class T>
  51. void throw_if_axes_is_too_large(const T& axes) {
  52. if (axes_rank(axes) > BOOST_HISTOGRAM_DETAIL_AXES_LIMIT)
  53. BOOST_THROW_EXCEPTION(
  54. std::invalid_argument("length of axis vector exceeds internal buffers, "
  55. "recompile with "
  56. "-DBOOST_HISTOGRAM_DETAIL_AXES_LIMIT=<new max size> "
  57. "to increase internal buffers"));
  58. }
  59. // tuple is never too large because internal buffers adapt to size of tuple
  60. template <class... Ts>
  61. void throw_if_axes_is_too_large(const std::tuple<Ts...>&) {}
  62. template <unsigned N, class... Ts>
  63. decltype(auto) axis_get(std::tuple<Ts...>& axes) {
  64. return std::get<N>(axes);
  65. }
  66. template <unsigned N, class... Ts>
  67. decltype(auto) axis_get(const std::tuple<Ts...>& axes) {
  68. return std::get<N>(axes);
  69. }
  70. template <unsigned N, class T>
  71. decltype(auto) axis_get(T& axes) {
  72. return axes[N];
  73. }
  74. template <unsigned N, class T>
  75. decltype(auto) axis_get(const T& axes) {
  76. return axes[N];
  77. }
  78. template <class... Ts>
  79. auto axis_get(std::tuple<Ts...>& axes, const unsigned i) {
  80. constexpr auto S = sizeof...(Ts);
  81. using V = mp11::mp_unique<axis::variant<Ts*...>>;
  82. return mp11::mp_with_index<S>(i, [&axes](auto i) { return V(&std::get<i>(axes)); });
  83. }
  84. template <class... Ts>
  85. auto axis_get(const std::tuple<Ts...>& axes, const unsigned i) {
  86. constexpr auto S = sizeof...(Ts);
  87. using V = mp11::mp_unique<axis::variant<const Ts*...>>;
  88. return mp11::mp_with_index<S>(i, [&axes](auto i) { return V(&std::get<i>(axes)); });
  89. }
  90. template <class T>
  91. decltype(auto) axis_get(T& axes, const unsigned i) {
  92. return axes[i];
  93. }
  94. template <class T>
  95. decltype(auto) axis_get(const T& axes, const unsigned i) {
  96. return axes[i];
  97. }
  98. template <class... Ts, class... Us>
  99. bool axes_equal(const std::tuple<Ts...>& ts, const std::tuple<Us...>& us) {
  100. using namespace ::boost::mp11;
  101. return static_if<std::is_same<mp_list<Ts...>, mp_list<Us...>>>(
  102. [](const auto& ts, const auto& us) {
  103. using N = mp_size<std::decay_t<decltype(ts)>>;
  104. bool equal = true;
  105. mp_for_each<mp_iota<N>>(
  106. [&](auto I) { equal &= relaxed_equal(std::get<I>(ts), std::get<I>(us)); });
  107. return equal;
  108. },
  109. [](const auto&, const auto&) { return false; }, ts, us);
  110. }
  111. template <class T, class... Us>
  112. bool axes_equal(const T& t, const std::tuple<Us...>& u) {
  113. using namespace ::boost::mp11;
  114. if (t.size() != sizeof...(Us)) return false;
  115. bool equal = true;
  116. mp_for_each<mp_iota_c<sizeof...(Us)>>([&](auto I) { equal &= t[I] == std::get<I>(u); });
  117. return equal;
  118. }
  119. template <class... Ts, class U>
  120. bool axes_equal(const std::tuple<Ts...>& t, const U& u) {
  121. return axes_equal(u, t);
  122. }
  123. template <class T, class U>
  124. bool axes_equal(const T& t, const U& u) {
  125. if (t.size() != u.size()) return false;
  126. return std::equal(t.begin(), t.end(), u.begin());
  127. }
  128. template <class... Ts, class... Us>
  129. void axes_assign(std::tuple<Ts...>& t, const std::tuple<Us...>& u) {
  130. using namespace ::boost::mp11;
  131. static_if<std::is_same<mp_list<Ts...>, mp_list<Us...>>>(
  132. [](auto& a, const auto& b) { a = b; },
  133. [](auto&, const auto&) {
  134. BOOST_THROW_EXCEPTION(
  135. std::invalid_argument("cannot assign axes, types do not match"));
  136. },
  137. t, u);
  138. }
  139. template <class... Ts, class U>
  140. void axes_assign(std::tuple<Ts...>& t, const U& u) {
  141. using namespace ::boost::mp11;
  142. mp_for_each<mp_iota_c<sizeof...(Ts)>>([&](auto I) {
  143. using T = mp_at_c<std::tuple<Ts...>, I>;
  144. std::get<I>(t) = axis::get<T>(u[I]);
  145. });
  146. }
  147. template <class T, class... Us>
  148. void axes_assign(T& t, const std::tuple<Us...>& u) {
  149. // resize instead of reserve, because t may not be empty and we want exact capacity
  150. t.resize(sizeof...(Us));
  151. using namespace ::boost::mp11;
  152. mp_for_each<mp_iota_c<sizeof...(Us)>>([&](auto I) { t[I] = std::get<I>(u); });
  153. }
  154. template <typename T, typename U>
  155. void axes_assign(T& t, const U& u) {
  156. t.assign(u.begin(), u.end());
  157. }
  158. template <class Archive, class T>
  159. void axes_serialize(Archive& ar, T& axes) {
  160. ar& make_nvp("axes", axes);
  161. }
  162. template <class Archive, class... Ts>
  163. void axes_serialize(Archive& ar, std::tuple<Ts...>& axes) {
  164. // needed to keep serialization format backward compatible
  165. struct proxy {
  166. std::tuple<Ts...>& t;
  167. void serialize(Archive& ar, unsigned /* version */) {
  168. mp11::tuple_for_each(t, [&ar](auto& x) { ar& make_nvp("item", x); });
  169. }
  170. };
  171. proxy p{axes};
  172. ar& make_nvp("axes", p);
  173. }
  174. // create empty dynamic axis which can store any axes types from the argument
  175. template <class T>
  176. auto make_empty_dynamic_axes(const T& axes) {
  177. return make_default(axes);
  178. }
  179. template <class... Ts>
  180. auto make_empty_dynamic_axes(const std::tuple<Ts...>&) {
  181. using namespace ::boost::mp11;
  182. using L = mp_unique<axis::variant<Ts...>>;
  183. // return std::vector<axis::variant<Axis0, Axis1, ...>> or std::vector<Axis0>
  184. return std::vector<mp_if_c<(mp_size<L>::value == 1), mp_first<L>, L>>{};
  185. }
  186. template <class T>
  187. void axis_index_is_valid(const T& axes, const unsigned N) {
  188. BOOST_ASSERT_MSG(N < axes_rank(axes), "index out of range");
  189. }
  190. template <class Axes, class V>
  191. void for_each_axis_impl(std::true_type, Axes&& axes, V&& v) {
  192. for (auto&& a : axes) { axis::visit(std::forward<V>(v), a); }
  193. }
  194. template <class Axes, class V>
  195. void for_each_axis_impl(std::false_type, Axes&& axes, V&& v) {
  196. for (auto&& a : axes) std::forward<V>(v)(a);
  197. }
  198. template <class Axes, class V>
  199. void for_each_axis(Axes&& a, V&& v) {
  200. using namespace ::boost::mp11;
  201. using T = mp_first<std::decay_t<Axes>>;
  202. for_each_axis_impl(is_axis_variant<T>(), std::forward<Axes>(a), std::forward<V>(v));
  203. }
  204. template <class V, class... Axis>
  205. void for_each_axis(const std::tuple<Axis...>& a, V&& v) {
  206. mp11::tuple_for_each(a, std::forward<V>(v));
  207. }
  208. template <class V, class... Axis>
  209. void for_each_axis(std::tuple<Axis...>& a, V&& v) {
  210. mp11::tuple_for_each(a, std::forward<V>(v));
  211. }
  212. // total number of bins including *flow bins
  213. template <class T>
  214. std::size_t bincount(const T& axes) {
  215. std::size_t n = 1;
  216. for_each_axis(axes, [&n](const auto& a) {
  217. const auto old = n;
  218. const auto s = axis::traits::extent(a);
  219. n *= s;
  220. if (s > 0 && n < old) BOOST_THROW_EXCEPTION(std::overflow_error("bincount overflow"));
  221. });
  222. return n;
  223. }
  224. // initial offset for the linear index
  225. template <class T>
  226. std::size_t offset(const T& axes) {
  227. std::size_t n = 0;
  228. for_each_axis(axes, [&n, stride = static_cast<std::size_t>(1)](const auto& a) mutable {
  229. if (axis::traits::options(a) & axis::option::growth)
  230. n = invalid_index;
  231. else if (n != invalid_index && axis::traits::options(a) & axis::option::underflow)
  232. n += stride;
  233. stride *= axis::traits::extent(a);
  234. });
  235. return n;
  236. }
  237. template <class T>
  238. using buffer_size_impl = typename std::tuple_size<T>::type;
  239. template <class T>
  240. using buffer_size = mp11::mp_eval_or<
  241. std::integral_constant<std::size_t, BOOST_HISTOGRAM_DETAIL_AXES_LIMIT>,
  242. buffer_size_impl, T>;
  243. template <class T, std::size_t N>
  244. class sub_array : public std::array<T, N> {
  245. using base_type = std::array<T, N>;
  246. public:
  247. explicit sub_array(std::size_t s) noexcept(
  248. std::is_nothrow_default_constructible<T>::value)
  249. : size_(s) {
  250. BOOST_ASSERT_MSG(size_ <= N, "requested size exceeds size of static buffer");
  251. }
  252. sub_array(std::size_t s,
  253. const T& value) noexcept(std::is_nothrow_copy_constructible<T>::value)
  254. : size_(s) {
  255. BOOST_ASSERT_MSG(size_ <= N, "requested size exceeds size of static buffer");
  256. std::array<T, N>::fill(value);
  257. }
  258. // need to override both versions of std::array
  259. auto end() noexcept { return base_type::begin() + size_; }
  260. auto end() const noexcept { return base_type::begin() + size_; }
  261. auto size() const noexcept { return size_; }
  262. private:
  263. std::size_t size_;
  264. };
  265. template <class U, class T>
  266. using stack_buffer = sub_array<U, buffer_size<T>::value>;
  267. // make default-constructed buffer (no initialization for POD types)
  268. template <class U, class T>
  269. auto make_stack_buffer(const T& t) {
  270. return stack_buffer<U, T>(axes_rank(t));
  271. }
  272. // make buffer with elements initialized to v
  273. template <class U, class T, class V>
  274. auto make_stack_buffer(const T& t, V&& v) {
  275. return stack_buffer<U, T>(axes_rank(t), std::forward<V>(v));
  276. }
  277. template <class T>
  278. using has_underflow =
  279. decltype(axis::traits::static_options<T>::test(axis::option::underflow));
  280. template <class T>
  281. using is_growing = decltype(axis::traits::static_options<T>::test(axis::option::growth));
  282. template <class T>
  283. using is_not_inclusive = mp11::mp_not<axis::traits::static_is_inclusive<T>>;
  284. // for vector<T>
  285. template <class T>
  286. struct axis_types_impl {
  287. using type = mp11::mp_list<std::decay_t<T>>;
  288. };
  289. // for vector<variant<Ts...>>
  290. template <class... Ts>
  291. struct axis_types_impl<axis::variant<Ts...>> {
  292. using type = mp11::mp_list<std::decay_t<Ts>...>;
  293. };
  294. // for tuple<Ts...>
  295. template <class... Ts>
  296. struct axis_types_impl<std::tuple<Ts...>> {
  297. using type = mp11::mp_list<std::decay_t<Ts>...>;
  298. };
  299. template <class T>
  300. using axis_types =
  301. typename axis_types_impl<mp11::mp_if<is_vector_like<T>, mp11::mp_first<T>, T>>::type;
  302. template <template <class> class Trait, class Axes>
  303. using has_special_axis = mp11::mp_any_of<axis_types<Axes>, Trait>;
  304. template <class Axes>
  305. using has_growing_axis = mp11::mp_any_of<axis_types<Axes>, is_growing>;
  306. template <class Axes>
  307. using has_non_inclusive_axis = mp11::mp_any_of<axis_types<Axes>, is_not_inclusive>;
  308. template <class T>
  309. constexpr std::size_t type_score() {
  310. return sizeof(T) *
  311. (std::is_integral<T>::value ? 1 : std::is_floating_point<T>::value ? 10 : 100);
  312. }
  313. // arbitrary ordering of types
  314. template <class T, class U>
  315. using type_less = mp11::mp_bool<(type_score<T>() < type_score<U>())>;
  316. template <class Axes>
  317. using value_types = mp11::mp_sort<
  318. mp11::mp_unique<mp11::mp_transform<axis::traits::value_type, axis_types<Axes>>>,
  319. type_less>;
  320. } // namespace detail
  321. } // namespace histogram
  322. } // namespace boost
  323. #endif