algorithm::sum gains coverage arg and improved documentation

This commit is contained in:
Hans Dembinski 2020-01-03 11:07:21 +01:00 committed by GitHub
parent 87e757f9e1
commit 270f6a1afb
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
2 changed files with 83 additions and 35 deletions

View File

@ -9,6 +9,7 @@
#include <boost/histogram/accumulators/sum.hpp> #include <boost/histogram/accumulators/sum.hpp>
#include <boost/histogram/fwd.hpp> #include <boost/histogram/fwd.hpp>
#include <boost/histogram/indexed.hpp>
#include <boost/mp11/utility.hpp> #include <boost/mp11/utility.hpp>
#include <numeric> #include <numeric>
#include <type_traits> #include <type_traits>
@ -16,24 +17,50 @@
namespace boost { namespace boost {
namespace histogram { namespace histogram {
namespace algorithm { namespace algorithm {
/** Compute the sum over all histogram cells, including underflow/overflow bins.
If the value type of the histogram is an integral or floating point type, /** Compute the sum over all histogram cells (underflow/overflow included by default).
boost::accumulators::sum<double> is used to compute the sum, else the original value
type is used. Compilation fails, if the value type does not support operator+=.
Return type is double if the value type of the histogram is integral or floating point, The implementation favors accuracy and protection against overflow over speed. If the
and the original value type otherwise. value type of the histogram is an integral or floating point type,
*/ accumulators::sum<double> is used to compute the sum, else the original value type is
used. Compilation fails, if the value type does not support operator+=. The return type
is double if the value type of the histogram is integral or floating point, and the
original value type otherwise.
If you need a different trade-off, you can write your own loop or use `std::accumulate`:
```
// iterate over all bins
auto sum_all = std::accumulate(hist.begin(), hist.end(), 0.0);
// skip underflow/overflow bins
double sum = 0;
for (auto&& x : indexed(hist))
sum += *x; // dereference accessor
// or:
// auto ind = boost::histogram::indexed(hist);
// auto sum = std::accumulate(ind.begin(), ind.end(), 0.0);
```
@returns accumulator type or double
@param hist Const reference to the histogram.
@param cov Iterate over all or only inner bins (optional, default: all).
*/
template <class A, class S> template <class A, class S>
auto sum(const histogram<A, S>& h) { auto sum(const histogram<A, S>& hist, const coverage cov = coverage::all) {
using T = typename histogram<A, S>::value_type; using T = typename histogram<A, S>::value_type;
using Sum = mp11::mp_if<std::is_arithmetic<T>, accumulators::sum<double>, T>; using sum_type = mp11::mp_if<std::is_arithmetic<T>, accumulators::sum<double>, T>;
Sum sum; sum_type sum;
for (auto&& x : h) sum += x; if (cov == coverage::all)
for (auto&& x : hist) sum += x;
else
// sum += x also works if sum_type::operator+=(const sum_type&) exists
for (auto&& x : indexed(hist)) sum += *x;
using R = mp11::mp_if<std::is_arithmetic<T>, double, T>; using R = mp11::mp_if<std::is_arithmetic<T>, double, T>;
return static_cast<R>(sum); return static_cast<R>(sum);
} }
} // namespace algorithm } // namespace algorithm
} // namespace histogram } // namespace histogram
} // namespace boost } // namespace boost

View File

@ -9,9 +9,9 @@
#include <boost/histogram/accumulators/weighted_sum.hpp> #include <boost/histogram/accumulators/weighted_sum.hpp>
#include <boost/histogram/algorithm/sum.hpp> #include <boost/histogram/algorithm/sum.hpp>
#include <boost/histogram/axis/integer.hpp> #include <boost/histogram/axis/integer.hpp>
#include "throw_exception.hpp"
#include <unordered_map> #include <unordered_map>
#include <vector> #include <vector>
#include "throw_exception.hpp"
#include "utility_histogram.hpp" #include "utility_histogram.hpp"
using namespace boost::histogram; using namespace boost::histogram;
@ -19,35 +19,56 @@ using boost::histogram::algorithm::sum;
template <typename Tag> template <typename Tag>
void run_tests() { void run_tests() {
auto ax = axis::integer<>(0, 100); auto ax = axis::integer<>(0, 10);
auto h1 = make(Tag(), ax); {
for (unsigned i = 0; i < 100; ++i) h1(i); auto h = make(Tag(), ax);
BOOST_TEST_EQ(sum(h1), 100); std::fill(h.begin(), h.end(), 1);
BOOST_TEST_EQ(sum(h), 12);
BOOST_TEST_EQ(sum(h, coverage::inner), 10);
}
auto h2 = make_s(Tag(), std::vector<double>(), ax, ax); {
for (unsigned i = 0; i < 100; ++i) auto h = make_s(Tag(), std::array<int, 12>(), ax);
for (unsigned j = 0; j < 100; ++j) h2(i, j); std::fill(h.begin(), h.end(), 1);
BOOST_TEST_EQ(sum(h2), 10000); BOOST_TEST_EQ(sum(h), 12);
BOOST_TEST_EQ(sum(h, coverage::inner), 10);
}
auto h3 = make_s(Tag(), std::array<int, 102>(), ax); {
for (unsigned i = 0; i < 100; ++i) h3(i); auto h = make_s(Tag(), std::unordered_map<std::size_t, int>(), ax);
BOOST_TEST_EQ(sum(h3), 100); std::fill(h.begin(), h.end(), 1);
BOOST_TEST_EQ(sum(h), 12);
BOOST_TEST_EQ(sum(h, coverage::inner), 10);
}
auto h4 = make_s(Tag(), std::unordered_map<std::size_t, int>(), ax); {
for (unsigned i = 0; i < 100; ++i) h4(i); auto h = make_s(Tag(), std::vector<double>(), ax, ax);
BOOST_TEST_EQ(sum(h4), 100); std::fill(h.begin(), h.end(), 1);
BOOST_TEST_EQ(sum(h), 12 * 12);
BOOST_TEST_EQ(sum(h, coverage::inner), 10 * 10);
}
auto h5 = {
make_s(Tag(), std::vector<accumulators::weighted_sum<>>(), axis::integer<>(0, 1), using W = accumulators::weighted_sum<>;
axis::integer<int, axis::null_type, axis::option::none_t>(2, 4)); auto h = make_s(Tag(), std::vector<W>(), axis::integer<>(0, 2),
h5(weight(2), 0, 2); axis::integer<int, axis::null_type, axis::option::none_t>(2, 4));
h5(-1, 2); W w(0, 2);
h5(1, 3); for (auto&& x : h) {
x = w;
w = W(w.value() + 1, 2);
}
const auto v = algorithm::sum(h5); // x-axis has 4 bins, y-axis has 2 = 8 bins total with 4 inner bins
BOOST_TEST_EQ(v.value(), 4);
BOOST_TEST_EQ(v.variance(), 6); const auto v1 = algorithm::sum(h);
BOOST_TEST_EQ(v1.value(), 0 + 1 + 2 + 3 + 4 + 5 + 6 + 7);
BOOST_TEST_EQ(v1.variance(), 8 * 2);
const auto v2 = algorithm::sum(h, coverage::inner);
BOOST_TEST_EQ(v2.value(), 1 + 2 + 5 + 6);
BOOST_TEST_EQ(v2.variance(), 4 * 2);
}
} }
int main() { int main() {