Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooNLLVar.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/**
18\file RooNLLVar.cxx
19\class RooNLLVar
20\ingroup Roofitcore
21
22Class RooNLLVar implements a -log(likelihood) calculation from a dataset
23and a PDF. The NLL is calculated as
24\f[
25 \sum_\mathrm{data} -\log( \mathrm{pdf}(x_\mathrm{data}))
26\f]
27In extended mode, a
28\f$ N_\mathrm{expect} - N_\mathrm{observed}*log(N_\mathrm{expect}) \f$ term is added.
29**/
30
31#include <RooNLLVar.h>
32
33#include <RooAbsData.h>
34#include <RooAbsDataStore.h>
35#include <RooAbsPdf.h>
36#include <RooCmdConfig.h>
37#include <RooDataHist.h>
38#include <RooHistPdf.h>
39#include <RooMsgService.h>
40#include <RooNaNPacker.h>
41#include <RooProdPdf.h>
42#include <RooRealMPFE.h>
43#include <RooRealSumPdf.h>
44#include <RooRealVar.h>
45
46#ifdef ROOFIT_CHECK_CACHED_VALUES
47#include <iomanip>
48#endif
49
50#include "TMath.h"
51#include "Math/Util.h"
52
53#include <algorithm>
54
55namespace {
56 template<class ...Args>
57 RooAbsTestStatistic::Configuration makeRooAbsTestStatisticCfg(Args const& ... args) {
59 cfg.rangeName = RooCmdConfig::decodeStringOnTheFly("RooNLLVar::RooNLLVar","RangeWithName",0,"",args...);
60 cfg.addCoefRangeName = RooCmdConfig::decodeStringOnTheFly("RooNLLVar::RooNLLVar","AddCoefRange",0,"",args...);
61 cfg.nCPU = RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","NumCPU",0,1,args...);
63 cfg.verbose = static_cast<bool>(RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","Verbose",0,1,args...));
64 cfg.splitCutRange = static_cast<bool>(RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","SplitRange",0,0,args...));
65 cfg.cloneInputData = static_cast<bool>(RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","CloneData",0,1,args...));
66 cfg.integrateOverBinsPrecision = RooCmdConfig::decodeDoubleOnTheFly("RooNLLVar::RooNLLVar", "IntegrateBins", 0, -1., {args...});
67 return cfg;
68 }
69}
70
72
74
76
77////////////////////////////////////////////////////////////////////////////////
78/// Construct likelihood from given p.d.f and (binned or unbinned dataset)
79///
80/// Argument | Description
81/// -------------------------|------------
82/// Extended() | Include extended term in calculation
83/// NumCPU() | Activate parallel processing feature
84/// Range() | Fit only selected region
85/// SumCoefRange() | Set the range in which to interpret the coefficients of RooAddPdf components
86/// SplitRange() | Fit range is split by index category of simultaneous PDF
87/// ConditionalObservables() | Define conditional observables
88/// Verbose() | Verbose output of GOF framework classes
89/// CloneData() | Clone input dataset for internal use (default is true)
90/// BatchMode() | Evaluate batches of data events (faster if PDFs support it)
91/// IntegrateBins() | Integrate PDF within each bin. This sets the desired precision. Only useful for binned fits.
92RooNLLVar::RooNLLVar(const char *name, const char* title, RooAbsPdf& pdf, RooAbsData& indata,
93 const RooCmdArg& arg1, const RooCmdArg& arg2,const RooCmdArg& arg3,
94 const RooCmdArg& arg4, const RooCmdArg& arg5,const RooCmdArg& arg6,
95 const RooCmdArg& arg7, const RooCmdArg& arg8,const RooCmdArg& arg9) :
96 RooAbsOptTestStatistic(name,title,pdf,indata,
97 *RooCmdConfig::decodeSetOnTheFly(
98 "RooNLLVar::RooNLLVar","ProjectedObservables",0,&_emptySet,
99 arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9),
100 makeRooAbsTestStatisticCfg(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9))
101{
102 RooCmdConfig pc("RooNLLVar::RooNLLVar") ;
103 pc.allowUndefined() ;
104 pc.defineInt("extended","Extended",0,false) ;
105 pc.defineInt("BatchMode", "BatchMode", 0, false);
106
107 pc.process(arg1) ; pc.process(arg2) ; pc.process(arg3) ;
108 pc.process(arg4) ; pc.process(arg5) ; pc.process(arg6) ;
109 pc.process(arg7) ; pc.process(arg8) ; pc.process(arg9) ;
110
111 _extended = pc.getInt("extended") ;
112 _skipZeroWeights = true;
113}
114
115
116////////////////////////////////////////////////////////////////////////////////
117/// Construct likelihood from given p.d.f and (binned or unbinned dataset)
118/// For internal use.
119
120RooNLLVar::RooNLLVar(const char *name, const char *title, RooAbsPdf& pdf, RooAbsData& indata,
121 bool extended, RooAbsTestStatistic::Configuration const& cfg) :
122 RooNLLVar{name, title, pdf, indata, RooArgSet(), extended, cfg} {}
123
124
125////////////////////////////////////////////////////////////////////////////////
126/// Construct likelihood from given p.d.f and (binned or unbinned dataset)
127/// For internal use.
128
129RooNLLVar::RooNLLVar(const char *name, const char *title, RooAbsPdf& pdf, RooAbsData& indata,
130 const RooArgSet& projDeps,
131 bool extended, RooAbsTestStatistic::Configuration const& cfg) :
132 RooAbsOptTestStatistic(name,title,pdf,indata,projDeps, cfg),
133 _extended(extended)
134{
135 // If binned likelihood flag is set, pdf is a RooRealSumPdf representing a yield vector
136 // for a binned likelihood calculation
137 _binnedPdf = cfg.binnedL ? static_cast<RooRealSumPdf*>(_funcClone) : nullptr ;
138
139 // Retrieve and cache bin widths needed to convert un-normalized binnedPdf values back to yields
140 if (_binnedPdf) {
141
142 // The Active label will disable pdf integral calculations
143 _binnedPdf->setAttribute("BinnedLikelihoodActive") ;
144
145 RooArgSet obs;
147 if (obs.size()!=1) {
148 _binnedPdf = nullptr;
149 } else {
150 auto* var = static_cast<RooRealVar*>(obs.first());
151 std::unique_ptr<std::list<double>> boundaries{_binnedPdf->binBoundaries(*var,var->getMin(),var->getMax())};
152 auto biter = boundaries->begin() ;
153 _binw.reserve(boundaries->size()-1) ;
154 double lastBound = (*biter) ;
155 ++biter ;
156 while (biter!=boundaries->end()) {
157 _binw.push_back((*biter) - lastBound);
158 lastBound = (*biter) ;
159 ++biter ;
160 }
161 }
162
163 _skipZeroWeights = false;
164 } else {
165 _skipZeroWeights = true;
166 }
167}
168
169
170
171////////////////////////////////////////////////////////////////////////////////
172/// Copy constructor
173
174RooNLLVar::RooNLLVar(const RooNLLVar& other, const char* name) :
176 _extended(other._extended),
177 _weightSq(other._weightSq),
178 _offsetSaveW2(other._offsetSaveW2),
179 _binw(other._binw),
180 _binnedPdf{other._binnedPdf}
181{
182}
183
184
185////////////////////////////////////////////////////////////////////////////////
186/// Create a test statistic using several properties of the current instance. This is used to duplicate
187/// the test statistic in multi-processing scenarios.
188RooAbsTestStatistic* RooNLLVar::create(const char *name, const char *title, RooAbsReal& pdf, RooAbsData& adata,
189 const RooArgSet& projDeps, RooAbsTestStatistic::Configuration const& cfg) {
190 RooAbsPdf & thePdf = dynamic_cast<RooAbsPdf&>(pdf);
191 // check if pdf can be extended
192 bool extendedPdf = _extended && thePdf.canBeExtended();
193
194 auto testStat = new RooNLLVar(name, title, thePdf, adata, projDeps, extendedPdf, cfg);
195 return testStat;
196}
197
198
199////////////////////////////////////////////////////////////////////////////////
200
202{
203 if (_gofOpMode==Slave) {
204 if (flag != _weightSq) {
205 _weightSq = flag;
206 std::swap(_offset, _offsetSaveW2);
207 }
209 } else if ( _gofOpMode==MPMaster) {
210 for (int i=0 ; i<_nCPU ; i++)
211 _mpfeArray[i]->applyNLLWeightSquared(flag);
212 } else if ( _gofOpMode==SimMaster) {
213 for(auto& gof : _gofArray)
214 static_cast<RooNLLVar&>(*gof).applyWeightSquared(flag);
215 }
216}
217
218
219////////////////////////////////////////////////////////////////////////////////
220/// Calculate and return likelihood on subset of data.
221/// \param[in] firstEvent First event to be processed.
222/// \param[in] lastEvent First event not to be processed, any more.
223/// \param[in] stepSize Steps between events.
224/// \note For batch computations, the step size **must** be one.
225///
226/// If this an extended likelihood, the extended term is added to the return likelihood
227/// in the batch that encounters the event with index 0.
228
229double RooNLLVar::evaluatePartition(std::size_t firstEvent, std::size_t lastEvent, std::size_t stepSize) const
230{
231 // Throughout the calculation, we use Kahan's algorithm for summing to
232 // prevent loss of precision - this is a factor four more expensive than
233 // straight addition, but since evaluating the PDF is usually much more
234 // expensive than that, we tolerate the additional cost...
236 double sumWeight{0.0};
237
238 auto * pdfClone = static_cast<RooAbsPdf*>(_funcClone);
239
240
241 // If pdf is marked as binned - do a binned likelihood calculation here (sum of log-Poisson for each bin)
242 if (_binnedPdf) {
243 ROOT::Math::KahanSum<double> sumWeightKahanSum{0.0};
244 for (auto i=firstEvent ; i<lastEvent ; i+=stepSize) {
245
246 _dataClone->get(i) ;
247
248 double eventWeight = _dataClone->weight();
249
250
251 // Calculate log(Poisson(N|mu) for this bin
252 double N = eventWeight ;
253 double mu = _binnedPdf->getVal()*_binw[i] ;
254 //cout << "RooNLLVar::binnedL(" << GetName() << ") N=" << N << " mu = " << mu << endl ;
255
256 if (mu<=0 && N>0) {
257
258 // Catch error condition: data present where zero events are predicted
259 logEvalError(Form("Observed %f events in bin %lu with zero event yield",N,(unsigned long)i)) ;
260
261 } else if (std::abs(mu)<1e-10 && std::abs(N)<1e-10) {
262
263 // Special handling of this case since log(Poisson(0,0)=0 but can't be calculated with usual log-formula
264 // since log(mu)=0. No update of result is required since term=0.
265
266 } else {
267
268 result += -1*(-mu + N * std::log(mu) - TMath::LnGamma(N+1));
269 sumWeightKahanSum += eventWeight;
270
271 }
272 }
273
274 sumWeight = sumWeightKahanSum.Sum();
275
276 } else { //unbinned PDF
277
278 std::tie(result, sumWeight) = computeScalar(stepSize, firstEvent, lastEvent);
279
280 // include the extended maximum likelihood term, if requested
281 if(_extended && _setNum==_extSet) {
282 result += pdfClone->extendedTerm(*_dataClone, _weightSq, _doBinOffset);
283 }
284 } //unbinned PDF
285
286
287 // If part of simultaneous PDF normalize probability over
288 // number of simultaneous PDFs: -sum(log(p/n)) = -sum(log(p)) + N*log(n)
289 if (_simCount>1) {
290 result += sumWeight * std::log(static_cast<double>(_simCount));
291 }
292
293
294 // At the end of the first full calculation, wire the caches
295 if (_first) {
296 _first = false ;
298 }
299
300
301 // Check if value offset flag is set.
302 if (_doOffset) {
303
304 // If no offset is stored enable this feature now
305 if (_offset.Sum() == 0 && _offset.Carry() == 0 && (result.Sum() != 0 || result.Carry() != 0)) {
306 coutI(Minimization) << "RooNLLVar::evaluatePartition(" << GetName() << ") first = "<< firstEvent << " last = " << lastEvent << " Likelihood offset now set to " << result.Sum() << std::endl ;
307 _offset = result ;
308 }
309
310 // Subtract offset
311 result -= _offset;
312 }
313
314 _evalCarry = result.Carry();
315 return result.Sum() ;
316}
317
318RooNLLVar::ComputeResult RooNLLVar::computeScalar(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const {
319 auto pdfClone = static_cast<const RooAbsPdf*>(_funcClone);
320 return computeScalarFunc(pdfClone, _dataClone, _normSet, _weightSq, stepSize, firstEvent, lastEvent, _offsetPdf.get());
321}
322
323// static function, also used from TestStatistics::RooUnbinnedL
325 RooArgSet *normSet, bool weightSq, std::size_t stepSize,
326 std::size_t firstEvent, std::size_t lastEvent, RooAbsPdf const* offsetPdf)
327{
330 RooNaNPacker packedNaN(0.f);
331
332 for (auto i=firstEvent; i<lastEvent; i+=stepSize) {
333 dataClone->get(i) ;
334
335 double weight = dataClone->weight(); //FIXME
336
337 if (0. == weight * weight) continue ;
338 if (weightSq) weight = dataClone->weightSquared() ;
339
340 double logProba = pdfClone->getLogVal(normSet);
341
342 if(offsetPdf) {
343 logProba -= offsetPdf->getLogVal(normSet);
344 }
345
346 const double term = -weight * logProba;
347
348 kahanWeight.Add(weight);
349 kahanProb.Add(term);
350 packedNaN.accumulate(term);
351 }
352
353 if (packedNaN.getPayload() != 0.) {
354 // Some events with evaluation errors. Return "badness" of errors.
355 return {ROOT::Math::KahanSum<double>{packedNaN.getNaNWithPayload()}, kahanWeight.Sum()};
356 }
357
358 return {kahanProb, kahanWeight.Sum()};
359}
360
361bool RooNLLVar::setData(RooAbsData &data, bool cloneData)
362{
363 bool ret = RooAbsOptTestStatistic::setData(data, cloneData);
364 // To re-create the data template pdf if necessary
365 _offsetPdf.reset();
367 return ret;
368}
369
371{
372 if (on && !_offsetPdf) {
373 std::string name = std::string{GetName()} + "_offsetPdf";
374 std::unique_ptr<RooDataHist> dataTemplate;
375 if (auto dh = dynamic_cast<RooDataHist *>(_dataClone)) {
376 dataTemplate = std::make_unique<RooDataHist>(*dh);
377 } else {
378 dataTemplate = std::unique_ptr<RooDataHist>(static_cast<RooDataSet const &>(*_dataClone).binnedClone());
379 }
380 _offsetPdf = std::make_unique<RooHistPdf>(name.c_str(), name.c_str(), *_funcObsSet, std::move(dataTemplate));
381 _offsetPdf->setOperMode(ADirty);
382 }
384}
#define e(i)
Definition RSha256.hxx:103
#define coutI(a)
#define ClassImp(name)
Definition Rtypes.h:377
#define N
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void on
char name[80]
Definition TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2468
The Kahan summation is a compensated summation algorithm, which significantly reduces numerical error...
Definition Util.h:122
T Sum() const
Definition Util.h:240
T Carry() const
Definition Util.h:250
void Add(T x)
Single-element accumulation. Will not vectorise.
Definition Util.h:165
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
void wireAllCaches()
void setValueDirty()
Mark the element dirty. This forces a re-evaluation when a value is requested.
Definition RooAbsArg.h:491
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
Storage_t::size_type size() const
RooAbsArg * first() const
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:57
virtual double weight() const =0
virtual const RooArgSet * get() const
Definition RooAbsData.h:101
virtual double weightSquared() const =0
RooAbsOptTestStatistic is the abstract base class for test statistics objects that evaluate a functio...
RooAbsReal * _funcClone
Pointer to internal clone of input function.
bool _skipZeroWeights
! Whether to skip entries with weight zero in the evaluation
RooArgSet * _funcObsSet
List of observables in the pdf expression.
RooArgSet * _normSet
Pointer to set with observables used for normalization.
RooAbsData * _dataClone
Pointer to internal clone if input data.
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition RooAbsPdf.h:218
virtual double getLogVal(const RooArgSet *set=nullptr) const
Return the log of the current value with given normalization An error message is printed if the argum...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:59
virtual std::list< double > * binBoundaries(RooAbsRealLValue &obs, double xlo, double xhi) const
Retrieve bin boundaries if this distribution is binned in obs.
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
virtual bool setData(RooAbsData &, bool=true)
Definition RooAbsReal.h:374
void logEvalError(const char *message, const char *serverValueString=nullptr) const
Log evaluation error message.
RooAbsTestStatistic is the abstract base class for all test statistics.
Int_t _setNum
Partition number of this instance in parallel calculation mode.
double _evalCarry
! carry of Kahan sum in evaluatePartition
Int_t _nCPU
Number of processors to use in parallel calculation mode.
GOFOpMode _gofOpMode
Operation mode of test statistic instance.
Int_t _simCount
Total number of component p.d.f.s in RooSimultaneous (if any)
ROOT::Math::KahanSum< double > _offset
! Offset as KahanSum to avoid loss of precision
Int_t _extSet
! Number of designated set to calculated extended term
std::vector< std::unique_ptr< RooAbsTestStatistic > > _gofArray
! Array of sub-contexts representing part of the combined test statistic
pRooRealMPFE * _mpfeArray
! Array of parallel execution frond ends
bool _doOffset
Apply interval value offset to control numeric precision?
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition RooCmdArg.h:26
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
static std::string decodeStringOnTheFly(const char *callerID, const char *cmdArgName, int intIdx, const char *defVal, Args_t &&...args)
Static decoder function allows to retrieve string property from set of RooCmdArgs For use in base mem...
static double decodeDoubleOnTheFly(const char *callerID, const char *cmdArgName, int idx, double defVal, std::initializer_list< std::reference_wrapper< const RooCmdArg > > args)
Find a given double in a list of RooCmdArg.
bool defineInt(const char *name, const char *argName, int intNum, int defValue=0)
Define integer property name 'name' mapped to integer in slot 'intNum' in RooCmdArg with name argName...
void allowUndefined(bool flag=true)
If flag is true the processing of unrecognized RooCmdArgs is not considered an error.
int getInt(const char *name, int defaultValue=0) const
Return integer property registered with name 'name'.
static int decodeIntOnTheFly(const char *callerID, const char *cmdArgName, int intIdx, int defVal, Args_t &&...args)
Static decoder function allows to retrieve integer property from set of RooCmdArgs For use in base me...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition RooDataHist.h:39
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:57
RooFit::OwningPtr< RooDataHist > binnedClone(const char *newName=nullptr, const char *newTitle=nullptr) const
Return binned clone of this dataset.
Class RooNLLVar implements a -log(likelihood) calculation from a dataset and a PDF.
Definition RooNLLVar.h:25
ComputeResult computeScalar(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
bool setData(RooAbsData &data, bool cloneData=true) override
Change dataset that is used to given one.
std::unique_ptr< RooAbsPdf > _offsetPdf
! An optional per-bin likelihood offset
Definition RooNLLVar.h:82
static RooNLLVar::ComputeResult computeScalarFunc(const RooAbsPdf *pdfClone, RooAbsData *dataClone, RooArgSet *normSet, bool weightSq, std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent, RooAbsPdf const *offsetPdf=nullptr)
bool _doBinOffset
Definition RooNLLVar.h:75
ROOT::Math::KahanSum< double > _offsetSaveW2
!
Definition RooNLLVar.h:78
static RooArgSet _emptySet
Definition RooNLLVar.h:69
RooAbsPdf * _binnedPdf
!
Definition RooNLLVar.h:81
void applyWeightSquared(bool flag) override
Disables or enables the usage of squared weights.
std::vector< double > _binw
!
Definition RooNLLVar.h:80
std::pair< ROOT::Math::KahanSum< double >, double > ComputeResult
Definition RooNLLVar.h:56
bool _extended
Definition RooNLLVar.h:74
RooNLLVar(const char *name, const char *title, RooAbsPdf &pdf, RooAbsData &data, const RooCmdArg &arg1={}, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}, const RooCmdArg &arg9={})
Construct likelihood from given p.d.f and (binned or unbinned dataset)
Definition RooNLLVar.cxx:92
RooAbsTestStatistic * create(const char *name, const char *title, RooAbsReal &pdf, RooAbsData &adata, const RooArgSet &projDeps, RooAbsTestStatistic::Configuration const &cfg) override
Create a test statistic using several properties of the current instance.
void enableBinOffsetting(bool on=true)
bool _first
!
Definition RooNLLVar.h:77
bool _weightSq
Apply weights squared?
Definition RooNLLVar.h:76
double evaluatePartition(std::size_t firstEvent, std::size_t lastEvent, std::size_t stepSize) const override
Calculate and return likelihood on subset of data.
The class RooRealSumPdf implements a PDF constructed from a sum of functions:
RooRealVar represents a variable that can be changed from the outside.
Definition RooRealVar.h:37
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
@ BulkPartition
Double_t LnGamma(Double_t z)
Computation of ln[gamma(z)] for all z.
Definition TMath.cxx:509
std::string rangeName
Stores the configuration parameters for RooAbsTestStatistic.
Little struct that can pack a float into the unused bits of the mantissa of a NaN double.
float getPayload() const
Retrieve packed float.
double getNaNWithPayload() const
Retrieve a NaN with the current float payload packed into the mantissa.
void accumulate(double val)
Accumulate a packed float from another NaN into this.