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 7 #ifndef BOOST_HISTOGRAM_AXIS_VARIABLE_HPP 8 #define BOOST_HISTOGRAM_AXIS_VARIABLE_HPP 9 10 #include <algorithm> 11 #include <boost/core/nvp.hpp> 12 #include <boost/histogram/axis/interval_view.hpp> 13 #include <boost/histogram/axis/iterator.hpp> 14 #include <boost/histogram/axis/metadata_base.hpp> 15 #include <boost/histogram/axis/option.hpp> 16 #include <boost/histogram/detail/convert_integer.hpp> 17 #include <boost/histogram/detail/detect.hpp> 18 #include <boost/histogram/detail/limits.hpp> 19 #include <boost/histogram/detail/relaxed_equal.hpp> 20 #include <boost/histogram/detail/replace_type.hpp> 21 #include <boost/histogram/fwd.hpp> 22 #include <boost/throw_exception.hpp> 23 #include <cassert> 24 #include <cmath> 25 #include <limits> 26 #include <memory> 27 #include <stdexcept> 28 #include <string> 29 #include <type_traits> 30 #include <utility> 31 #include <vector> 32 33 namespace boost { 34 namespace histogram { 35 namespace axis { 36 37 /** 38 Axis for non-equidistant bins on the real line. 39 40 Binning is a O(log(N)) operation. If speed matters and the problem domain 41 allows it, prefer a regular axis, possibly with a transform. 42 43 @tparam Value input value type, must be floating point. 44 @tparam MetaData type to store meta data. 45 @tparam Options see boost::histogram::axis::option (all values allowed). 46 @tparam Allocator allocator to use for dynamic memory management. 47 */ 48 template <class Value, class MetaData, class Options, class Allocator> 49 class variable : public iterator_mixin<variable<Value, MetaData, Options, Allocator>>, 50 public metadata_base_t<MetaData> { 51 // these must be private, so that they are not automatically inherited 52 using value_type = Value; 53 using metadata_base = metadata_base_t<MetaData>; 54 using metadata_type = typename metadata_base::metadata_type; 55 using options_type = 56 detail::replace_default<Options, decltype(option::underflow | option::overflow)>; 57 using allocator_type = Allocator; 58 using vector_type = std::vector<Value, allocator_type>; 59 60 static_assert( 61 std::is_floating_point<value_type>::value, 62 "current version of variable axis requires floating point type; " 63 "if you need a variable axis with an integral type, please submit an issue"); 64 65 static_assert( 66 (!options_type::test(option::circular) && !options_type::test(option::growth)) || 67 (options_type::test(option::circular) ^ options_type::test(option::growth)), 68 "circular and growth options are mutually exclusive"); 69 70 public: 71 constexpr variable() = default; variable(allocator_type alloc)72 explicit variable(allocator_type alloc) : vec_(alloc) {} 73 74 /** Construct from iterator range of bin edges. 75 * 76 * \param begin begin of edge sequence. 77 * \param end end of edge sequence. 78 * \param meta description of the axis. 79 * \param alloc allocator instance to use. 80 */ 81 template <class It, class = detail::requires_iterator<It>> variable(It begin,It end,metadata_type meta={},allocator_type alloc={})82 variable(It begin, It end, metadata_type meta = {}, allocator_type alloc = {}) 83 : metadata_base(std::move(meta)), vec_(std::move(alloc)) { 84 if (std::distance(begin, end) < 2) 85 BOOST_THROW_EXCEPTION(std::invalid_argument("bins > 0 required")); 86 87 vec_.reserve(std::distance(begin, end)); 88 vec_.emplace_back(*begin++); 89 bool strictly_ascending = true; 90 while (begin != end) { 91 if (*begin <= vec_.back()) strictly_ascending = false; 92 vec_.emplace_back(*begin++); 93 } 94 if (!strictly_ascending) 95 BOOST_THROW_EXCEPTION( 96 std::invalid_argument("input sequence must be strictly ascending")); 97 } 98 99 /** Construct variable axis from iterable range of bin edges. 100 * 101 * \param iterable iterable range of bin edges. 102 * \param meta description of the axis. 103 * \param alloc allocator instance to use. 104 */ 105 template <class U, class = detail::requires_iterable<U>> variable(const U & iterable,metadata_type meta={},allocator_type alloc={})106 variable(const U& iterable, metadata_type meta = {}, allocator_type alloc = {}) 107 : variable(std::begin(iterable), std::end(iterable), std::move(meta), 108 std::move(alloc)) {} 109 110 /** Construct variable axis from initializer list of bin edges. 111 * 112 * @param list `std::initializer_list` of bin edges. 113 * @param meta description of the axis. 114 * @param alloc allocator instance to use. 115 */ 116 template <class U> variable(std::initializer_list<U> list,metadata_type meta={},allocator_type alloc={})117 variable(std::initializer_list<U> list, metadata_type meta = {}, 118 allocator_type alloc = {}) 119 : variable(list.begin(), list.end(), std::move(meta), std::move(alloc)) {} 120 121 /// Constructor used by algorithm::reduce to shrink and rebin (not for users). variable(const variable & src,index_type begin,index_type end,unsigned merge)122 variable(const variable& src, index_type begin, index_type end, unsigned merge) 123 : metadata_base(src), vec_(src.get_allocator()) { 124 assert((end - begin) % merge == 0); 125 if (options_type::test(option::circular) && !(begin == 0 && end == src.size())) 126 BOOST_THROW_EXCEPTION(std::invalid_argument("cannot shrink circular axis")); 127 vec_.reserve((end - begin) / merge); 128 const auto beg = src.vec_.begin(); 129 for (index_type i = begin; i <= end; i += merge) vec_.emplace_back(*(beg + i)); 130 } 131 132 /// Return index for value argument. index(value_type x) const133 index_type index(value_type x) const noexcept { 134 if (options_type::test(option::circular)) { 135 const auto a = vec_[0]; 136 const auto b = vec_[size()]; 137 x -= std::floor((x - a) / (b - a)) * (b - a); 138 } 139 return static_cast<index_type>(std::upper_bound(vec_.begin(), vec_.end(), x) - 140 vec_.begin() - 1); 141 } 142 update(value_type x)143 std::pair<index_type, index_type> update(value_type x) noexcept { 144 const auto i = index(x); 145 if (std::isfinite(x)) { 146 if (0 <= i) { 147 if (i < size()) return std::make_pair(i, 0); 148 const auto d = value(size()) - value(size() - 0.5); 149 x = std::nextafter(x, (std::numeric_limits<value_type>::max)()); 150 x = (std::max)(x, vec_.back() + d); 151 vec_.push_back(x); 152 return {i, -1}; 153 } 154 const auto d = value(0.5) - value(0); 155 x = (std::min)(x, value(0) - d); 156 vec_.insert(vec_.begin(), x); 157 return {0, -i}; 158 } 159 return {x < 0 ? -1 : size(), 0}; 160 } 161 162 /// Return value for fractional index argument. value(real_index_type i) const163 value_type value(real_index_type i) const noexcept { 164 if (options_type::test(option::circular)) { 165 auto shift = std::floor(i / size()); 166 i -= shift * size(); 167 double z; 168 const auto k = static_cast<index_type>(std::modf(i, &z)); 169 const auto a = vec_[0]; 170 const auto b = vec_[size()]; 171 return (1.0 - z) * vec_[k] + z * vec_[k + 1] + shift * (b - a); 172 } 173 if (i < 0) return detail::lowest<value_type>(); 174 if (i == size()) return vec_.back(); 175 if (i > size()) return detail::highest<value_type>(); 176 const auto k = static_cast<index_type>(i); // precond: i >= 0 177 const real_index_type z = i - k; 178 return (1.0 - z) * vec_[k] + z * vec_[k + 1]; 179 } 180 181 /// Return bin for index argument. bin(index_type idx) const182 auto bin(index_type idx) const noexcept { return interval_view<variable>(*this, idx); } 183 184 /// Returns the number of bins, without over- or underflow. size() const185 index_type size() const noexcept { return static_cast<index_type>(vec_.size()) - 1; } 186 187 /// Returns the options. options()188 static constexpr unsigned options() noexcept { return options_type::value; } 189 190 template <class V, class M, class O, class A> operator ==(const variable<V,M,O,A> & o) const191 bool operator==(const variable<V, M, O, A>& o) const noexcept { 192 const auto& a = vec_; 193 const auto& b = o.vec_; 194 return std::equal(a.begin(), a.end(), b.begin(), b.end()) && 195 detail::relaxed_equal{}(this->metadata(), o.metadata()); 196 } 197 198 template <class V, class M, class O, class A> operator !=(const variable<V,M,O,A> & o) const199 bool operator!=(const variable<V, M, O, A>& o) const noexcept { 200 return !operator==(o); 201 } 202 203 /// Return allocator instance. get_allocator() const204 auto get_allocator() const { return vec_.get_allocator(); } 205 206 template <class Archive> serialize(Archive & ar,unsigned)207 void serialize(Archive& ar, unsigned /* version */) { 208 ar& make_nvp("seq", vec_); 209 ar& make_nvp("meta", this->metadata()); 210 } 211 212 private: 213 vector_type vec_; 214 215 template <class V, class M, class O, class A> 216 friend class variable; 217 }; 218 219 #if __cpp_deduction_guides >= 201606 220 221 template <class T> 222 variable(std::initializer_list<T>) 223 ->variable<detail::convert_integer<T, double>, null_type>; 224 225 template <class T, class M> 226 variable(std::initializer_list<T>, M) 227 ->variable<detail::convert_integer<T, double>, 228 detail::replace_type<std::decay_t<M>, const char*, std::string>>; 229 230 template <class Iterable, class = detail::requires_iterable<Iterable>> 231 variable(Iterable) 232 ->variable< 233 detail::convert_integer< 234 std::decay_t<decltype(*std::begin(std::declval<Iterable&>()))>, double>, 235 null_type>; 236 237 template <class Iterable, class M> 238 variable(Iterable, M) 239 ->variable< 240 detail::convert_integer< 241 std::decay_t<decltype(*std::begin(std::declval<Iterable&>()))>, double>, 242 detail::replace_type<std::decay_t<M>, const char*, std::string>>; 243 244 #endif 245 246 } // namespace axis 247 } // namespace histogram 248 } // namespace boost 249 250 #endif 251