15 #ifndef ROOT7_THistImpl 16 #define ROOT7_THistImpl 29 namespace Experimental {
32 template <
int D_,
class P_,
template <
class P__>
class STORAGE>
class... STAT>
38 template<
int NDIM>
using AxisIter_t = std::array<TAxisBase::const_iterator, NDIM>;
51 return static_cast<int>(
a) & static_cast<int>(b);
66 template <
int DIMENSIONS>
67 class THistImplPrecisionAgnosticBase {
81 static constexpr
int GetNDim() {
return DIMENSIONS; }
84 virtual int GetNBins()
const noexcept = 0;
87 const std::string&
GetTitle()
const {
return fTitle; }
93 virtual int GetBinIndexAndGrow(
const CoordArray_t& x) = 0;
104 virtual double GetBinUncertainty(
int binidx)
const = 0;
108 virtual bool HasBinUncertainty()
const = 0;
111 virtual double GetBinContentAsDouble(
int binidx)
const = 0;
116 virtual TAxisView GetAxis(
int iAxis)
const = 0;
124 GetRange(
const std::array<Hist::EOverflow, DIMENSIONS>& withOverUnder)
const = 0;
167 virtual void FillN(
const std::array_view<CoordArray_t> xN,
168 const std::array_view<Weight_t> weightN) = 0;
171 virtual void FillN(
const std::array_view<CoordArray_t> xN) = 0;
194 virtual double GetBinUncertainty(
const CoordArray_t& x)
const = 0;
198 int GetNBins() const noexcept final {
return fStatistics.size(); }
215 return (
double) GetBinContent(binidx);
237 return std::get<0>(axes).GetNBins();
242 template <
int I,
class AXES>
245 return std::get<I>(axes).GetNBins() *
TGetBinCount<
I - 1, AXES>()(axes);
250 template<
class... AXISCONFIG>
252 using axesTuple = std::tuple<AXISCONFIG...>;
253 return TGetBinCount<
sizeof...(AXISCONFIG) - 1, axesTuple>()(axesTuple{axisArgs...});
257 template <
int IDX,
class HISTIMPL,
class AXES,
bool GROW>
261 template <
class HISTIMPL,
class AXES,
bool GROW>
270 template <
int I,
class HISTIMPL,
class AXES,
bool GROW>
274 constexpr
const int thisAxis = HISTIMPL::GetNDim() -
I - 1;
275 int bin = std::get<thisAxis>(axes).FindBin(x[thisAxis]);
276 if (GROW && std::get<thisAxis>(axes).CanGrow()
277 && (bin < 0 || bin > std::get<thisAxis>(axes).GetNBinsNoOver())) {
278 hist->GrowAxis(
I, x[thisAxis]);
284 return bin +
TGetBinIndex<
I - 1, HISTIMPL, AXES, GROW>()(hist, axes, x, status)
285 * std::get<thisAxis>(axes).GetNBins();
296 const std::array<
Hist::EOverflow, std::tuple_size<AXES>::value>& )
const {}
302 template<
int I,
class AXES>
306 const std::array<
Hist::EOverflow, std::tuple_size<AXES>::value> &over)
const {
308 range[0][
I] = std::get<I>(axes).begin_with_underflow();
310 range[0][
I] = std::get<I>(axes).begin();
312 range[1][
I] = std::get<I>(axes).end_with_overflow();
314 range[1][
I] = std::get<I>(axes).end();
336 template<
int I,
class COORD,
class AXES>
339 int axisbin = binidx % std::get<I>(axes).GetNBins();
340 size_t coordidx = std::tuple_size<AXES>::value -
I - 1;
342 case EBinCoord::kBinFrom:
343 coord[coordidx] = std::get<I>(axes).GetBinFrom(axisbin);
345 case EBinCoord::kBinCenter:
346 coord[coordidx] = std::get<I>(axes).GetBinCenter(axisbin);
348 case EBinCoord::kBinTo:
349 coord[coordidx] = std::get<I>(axes).GetBinTo(axisbin);
353 binidx / std::get<I>(axes).GetNBins());
359 template <
class... AXISCONFIG>
360 static std::array<
TAxisView,
sizeof...(AXISCONFIG)>
362 std::array<TAxisView,
sizeof...(AXISCONFIG)> axisViews = {
374 template <
class DATA,
class... AXISCONFIG>
376 static_assert(
sizeof...(AXISCONFIG) == DATA::GetNDim(),
377 "Number of axes must equal histogram dimension");
379 friend typename DATA::Hist_t;
411 for (
auto&& binref: *
this)
412 op(binref.GetCenter(), binref.GetContent());
418 for (
auto&& binref: *
this)
419 op(binref.GetCenter(), binref.GetContent(), binref.GetUncertainty());
424 const std::tuple<AXISCONFIG...>&
GetAxes()
const {
return fAxes; }
428 return std::apply(Internal::GetAxisView<AXISCONFIG...>, fAxes)[iAxis];
437 decltype(fAxes),
false>()(
nullptr, fAxes,
x, status);
451 (
this, fAxes,
x, status);
485 void FillN(
const std::array_view<CoordArray_t> xN,
486 const std::array_view<Weight_t> weightN)
final {
488 if (xN.size() != weightN.size()) {
489 R__ERROR_HERE(
"HIST") <<
"Not the same number of points and weights!";
494 for (
size_t i = 0; i < xN.size(); ++i) {
495 Fill(xN[i], weightN[i]);
502 void FillN(
const std::array_view<CoordArray_t> xN)
final {
510 int bin = GetBinIndexAndGrow(x);
511 this->GetStat().Fill(x, bin, w);
516 int bin = GetBinIndex(
x);
518 return ImplBase_t::GetBinContent(bin);
524 return this->GetStat().GetBinUncertainty(binidx);
529 const int bin = GetBinIndex(
x);
530 return this->GetBinUncertainty(bin);
536 return this->GetStat().HasBinUncertainty();
545 std::array<std::array<TAxisBase::const_iterator, DATA::GetNDim()>, 2> ret;
570 template <
class DATA,
class... AXISCONFIG>
575 template <
class DATA,
class... AXISCONFIG>
582 template <
class DATA,
class... AXISCONFIG>
592 template <
class DATA>
595 THistImplRuntime(std::array<
TAxisConfig, DATA::GetNDim()>&& axisCfg);
Weight_t GetBinContent(const CoordArray_t &x) const final
Get the content of the bin at position x.
Exclude under- and overflows.
DATA Stat_t
Type of the statistics (bin content, uncertainties etc).
CoordArray_t GetBinTo(int binidx) const final
Get the coordinate of the high limit of the bin.
void AddBinContent(int binidx, Weight_t w)
Add w to the bin at index bin.
std::array< double, DIMENSIONS > CoordArray_t
This namespace contains pre-defined functions to be used in conjuction with TExecutor::Map and TExecu...
Weight_t GetBinContent(int binidx) const
Get the bin content (sum of weights) for bin index binidx.
const Stat_t & GetStat() const noexcept
Const access to statistics.
Interface class for THistImpl.
int GetBinIndexAndGrow(const CoordArray_t &x) final
Gets the bin index for coordinate x, growing the axes as needed and possible.
typename Hist::AxisIterRange_t< NDIM > AxisIterRange_t
virtual void ApplyXCE(std::function< void(const CoordArray_t &, Weight_t, double)> op) const final
Apply a function (lambda) to all bins of the histogram.
void operator()(COORD &, const AXES &, EBinCoord, int) const
int operator()(HISTIMPL *, const AXES &, const typename HISTIMPL::CoordArray_t &, TAxisBase::EFindStatus &status) const
bool operator&(EOverflow a, EOverflow b)
Weight_t & GetBinContent(int binidx)
Get the bin content (sum of weights) for bin index binidx (non-const).
Iterates over the bins of a THist or THistImpl.
int GetNBins() const noexcept final
Get the number of bins in this histogram, including possible under- and overflow bins.
std::tuple< AXISCONFIG... > fAxes
The histogram's axes.
void GrowAxis(int, double)
Grow the axis number iAxis to fit the coordinate x.
Common view on a TAxis, no matter what its kind.
CoordArray_t GetBinFrom(int binidx) const final
Get the coordinate of the low limit of the bin.
void operator()(Hist::AxisIterRange_t< std::tuple_size< AXES >::value > &, const AXES &, const std::array< Hist::EOverflow, std::tuple_size< AXES >::value > &) const
const std::tuple< AXISCONFIG... > & GetAxes() const
Get the axes of this histogram.
EFindStatus
Status of FindBin(x)
CoordArray_t GetBinCenter(int binidx) const final
Get the center coordinate of the bin.
AxisIterRange_t< DATA::GetNDim()> GetRange(const std::array< Hist::EOverflow, DATA::GetNDim()> &withOverUnder) const final
Get the begin() and end() for each axis.
const std::string & GetTitle() const
Get the histogram title.
FillFunc_t GetFillFunc() const final
Retrieve the fill function for this histogram implementation, to prevent the virtual function call fo...
Stat_t & GetStat() noexcept
Non-const access to statistics.
Hist::CoordArray_t< DIMENSIONS > CoordArray_t
Type of the coordinate: a DIMENSIONS-dimensional array of doubles.
int GetBinIndex(const CoordArray_t &x) const final
Gets the bin index for coordinate x; returns -1 if there is no such bin, e.g.
int GetNBinsFromAxes(AXISCONFIG...axisArgs)
std::string fTitle
Histogram title.
int operator()(const AXES &axes) const
void FillN(const std::array_view< CoordArray_t > xN) final
Fill an array of weightN to the bins specified by coordinates xN.
static std::array< TAxisView, sizeof...(AXISCONFIG)> GetAxisView(const AXISCONFIG &...axes) noexcept
void ApplyXC(std::function< void(const CoordArray_t &, Weight_t)> op) const final
Apply a function (lambda) to all bins of the histogram.
void Apply(std::function< void(THistBinRef< const ImplBase_t >)> op) const final
Apply a function (lambda) to all bins of the histogram.
static void GetRange(const char *comments, Double_t &xmin, Double_t &xmax, Double_t &factor)
Parse comments to search for a range specifier of the style: [xmin,xmax] or [xmin,xmax,nbits] [0,1] [-10,100]; [-pi,pi], [-pi/2,pi/4],[-2pi,2*pi] [-10,100,16] [0,0,8] if nbits is not specified, or nbits <2 or nbits>32 it is set to 32 if (xmin==0 and xmax==0 and nbits <=16) the double word will be converted to a float and its mantissa truncated to nbits significative bits.
iterator begin() noexcept
bool HasBinUncertainty() const final
Whether this histogram's statistics provide storage for uncertainties, or whether uncertainties are d...
THistImplPrecisionAgnosticBase(std::string_view title)
Include both under- and overflows.
int operator()(HISTIMPL *hist, const AXES &axes, const typename HISTIMPL::CoordArray_t &x, TAxisBase::EFindStatus &status) const
int operator()(const AXES &axes) const
void Fill(const CoordArray_t &x, Weight_t w=1.)
Add a single weight w to the bin at coordinate x.
Objects used to configure the different axis types.
typename DATA::Weight_t Weight_t
Type of the bin content (and thus weights).
Represents a bin reference.
Hist::CoordArray_t< DATA::GetNDim()> CoordArray_t
Type of the coordinate: a DIMENSIONS-dimensional array of doubles.
static constexpr int GetNDim()
Number of dimensions of the coordinates.
const_iterator end() const noexcept
THistImplBase(std::string_view title, size_t numBins)
void(THistImplBase::*)(const CoordArray_t &x, Weight_t w) FillFunc_t
Type of the Fill(x, w) function.
The returned bin index is valid.
double GetBinUncertainty(int binidx) const final
Return the uncertainties for the given bin.
virtual ~THistImplPrecisionAgnosticBase()
double GetBinContentAsDouble(int binidx) const final
Get the bin content (sum of weights) for bin index binidx, cast to double.
Base class for THistImplBase that abstracts out the histogram's PRECISION.
Hist::AxisIterRange_t< DIMENSIONS > AxisIterRange_t
Range type.
Stat_t fStatistics
The histogram's bin content, uncertainties etc.
void operator()(Hist::AxisIterRange_t< std::tuple_size< AXES >::value > &range, const AXES &axes, const std::array< Hist::EOverflow, std::tuple_size< AXES >::value > &over) const
void FillN(const std::array_view< CoordArray_t > xN, const std::array_view< Weight_t > weightN) final
Fill an array of weightN to the bins specified by coordinates xN.
typedef void((*Func_t)())
Fill range with begin() and end() of all axes, including under/overflow as specified by over...
void operator()(COORD &coord, const AXES &axes, EBinCoord kind, int binidx) const
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
const_iterator begin() const noexcept
Fill coord with low bin edge or center or high bin edge of all axes.
Coordinate could fit after growing the axis.
std::array< AxisIter_t< NDIM >, 2 > AxisIterRange_t
Range over n dimensional axes - a pair of arrays of n axis iterators.
EOverflow
Kinds of under- and overflow handling.
TAxisView GetAxis(int iAxis) const final
Normalized axes access, converting the actual axis to TAxisConfig.
decltype(auto) constexpr apply(F &&f, Tuple &&t)
#define R__ERROR_HERE(GROUP)
double GetBinUncertainty(const CoordArray_t &x) const final
Get the bin uncertainty for the bin at coordinate x.
THistImplBase(size_t numBins)
std::array< TAxisBase::const_iterator, NDIM > AxisIter_t
Iterator over n dimensional axes - an array of n axis iterators.