Logo ROOT  
Reference Guide
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 "RooAbsPdf.h"
35#include "RooCmdConfig.h"
36#include "RooMsgService.h"
37#include "RooAbsDataStore.h"
38#include "RooRealMPFE.h"
39#include "RooRealSumPdf.h"
40#include "RooRealVar.h"
41#include "RooProdPdf.h"
42#include "RooNaNPacker.h"
43#include "RunContext.h"
44
45#ifdef ROOFIT_CHECK_CACHED_VALUES
46#include <iomanip>
47#endif
48
49#include "TMath.h"
50#include "Math/Util.h"
51
52#include <algorithm>
53
54namespace {
55 template<class ...Args>
56 RooAbsTestStatistic::Configuration makeRooAbsTestStatisticCfg(Args const& ... args) {
58 cfg.rangeName = RooCmdConfig::decodeStringOnTheFly("RooNLLVar::RooNLLVar","RangeWithName",0,"",args...);
59 cfg.addCoefRangeName = RooCmdConfig::decodeStringOnTheFly("RooNLLVar::RooNLLVar","AddCoefRange",0,"",args...);
60 cfg.nCPU = RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","NumCPU",0,1,args...);
62 cfg.verbose = static_cast<bool>(RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","Verbose",0,1,args...));
63 cfg.splitCutRange = static_cast<bool>(RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","SplitRange",0,0,args...));
64 cfg.cloneInputData = static_cast<bool>(RooCmdConfig::decodeIntOnTheFly("RooNLLVar::RooNLLVar","CloneData",0,1,args...));
65 cfg.integrateOverBinsPrecision = RooCmdConfig::decodeDoubleOnTheFly("RooNLLVar::RooNLLVar", "IntegrateBins", 0, -1., {args...});
66 return cfg;
67 }
68}
69
71
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 *static_cast<const RooArgSet*>(RooCmdConfig::decodeObjOnTheFly(
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 _batchEvaluations = pc.getInt("BatchMode");
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}
164
165
166
167////////////////////////////////////////////////////////////////////////////////
168/// Copy constructor
169
170RooNLLVar::RooNLLVar(const RooNLLVar& other, const char* name) :
172 _extended(other._extended),
173 _batchEvaluations(other._batchEvaluations),
174 _weightSq(other._weightSq),
175 _offsetSaveW2(other._offsetSaveW2),
176 _binw(other._binw),
177 _binnedPdf{other._binnedPdf}
178{
179}
180
181
182////////////////////////////////////////////////////////////////////////////////
183/// Create a test statistic using several properties of the current instance. This is used to duplicate
184/// the test statistic in multi-processing scenarios.
185RooAbsTestStatistic* RooNLLVar::create(const char *name, const char *title, RooAbsReal& pdf, RooAbsData& adata,
186 const RooArgSet& projDeps, RooAbsTestStatistic::Configuration const& cfg) {
187 RooAbsPdf & thePdf = dynamic_cast<RooAbsPdf&>(pdf);
188 // check if pdf can be extended
189 bool extendedPdf = _extended && thePdf.canBeExtended();
190
191 auto testStat = new RooNLLVar(name, title, thePdf, adata, projDeps, extendedPdf, cfg);
192 testStat->batchMode(_batchEvaluations);
193 return testStat;
194}
195
196
197////////////////////////////////////////////////////////////////////////////////
198
200{
201 if (_gofOpMode==Slave) {
202 if (flag != _weightSq) {
203 _weightSq = flag;
205 }
207 } else if ( _gofOpMode==MPMaster) {
208 for (int i=0 ; i<_nCPU ; i++)
209 _mpfeArray[i]->applyNLLWeightSquared(flag);
210 } else if ( _gofOpMode==SimMaster) {
211 for (int i=0 ; i<_nGof ; i++)
212 static_cast<RooNLLVar*>(_gofArray[i])->applyWeightSquared(flag);
213 }
214}
215
216
217////////////////////////////////////////////////////////////////////////////////
218/// Calculate and return likelihood on subset of data.
219/// \param[in] firstEvent First event to be processed.
220/// \param[in] lastEvent First event not to be processed, any more.
221/// \param[in] stepSize Steps between events.
222/// \note For batch computations, the step size **must** be one.
223///
224/// If this an extended likelihood, the extended term is added to the return likelihood
225/// in the batch that encounters the event with index 0.
226
227double RooNLLVar::evaluatePartition(std::size_t firstEvent, std::size_t lastEvent, std::size_t stepSize) const
228{
229 // Throughout the calculation, we use Kahan's algorithm for summing to
230 // prevent loss of precision - this is a factor four more expensive than
231 // straight addition, but since evaluating the PDF is usually much more
232 // expensive than that, we tolerate the additional cost...
234 double sumWeight{0.0};
235
236 auto * pdfClone = static_cast<RooAbsPdf*>(_funcClone);
237
238 // cout << "RooNLLVar::evaluatePartition(" << GetName() << ") projDeps = " << (_projDeps?*_projDeps:RooArgSet()) << endl ;
239
240 _dataClone->store()->recalculateCache( _projDeps, firstEvent, lastEvent, stepSize, (_binnedPdf?false:true) ) ;
241
242
243
244 // If pdf is marked as binned - do a binned likelihood calculation here (sum of log-Poisson for each bin)
245 if (_binnedPdf) {
246 ROOT::Math::KahanSum<double> sumWeightKahanSum{0.0};
247 for (auto i=firstEvent ; i<lastEvent ; i+=stepSize) {
248
249 _dataClone->get(i) ;
250
251 double eventWeight = _dataClone->weight();
252
253
254 // Calculate log(Poisson(N|mu) for this bin
255 double N = eventWeight ;
256 double mu = _binnedPdf->getVal()*_binw[i] ;
257 //cout << "RooNLLVar::binnedL(" << GetName() << ") N=" << N << " mu = " << mu << endl ;
258
259 if (mu<=0 && N>0) {
260
261 // Catch error condition: data present where zero events are predicted
262 logEvalError(Form("Observed %f events in bin %lu with zero event yield",N,(unsigned long)i)) ;
263
264 } else if (std::abs(mu)<1e-10 && std::abs(N)<1e-10) {
265
266 // Special handling of this case since log(Poisson(0,0)=0 but can't be calculated with usual log-formula
267 // since log(mu)=0. No update of result is required since term=0.
268
269 } else {
270
271 result += -1*(-mu + N * std::log(mu) - TMath::LnGamma(N+1));
272 sumWeightKahanSum += eventWeight;
273
274 }
275 }
276
277 sumWeight = sumWeightKahanSum.Sum();
278
279 } else { //unbinned PDF
280
281 if (_batchEvaluations) {
282 std::tie(result, sumWeight) = computeBatched(stepSize, firstEvent, lastEvent);
283#ifdef ROOFIT_CHECK_CACHED_VALUES
284
285 ROOT::Math::KahanSum<double> resultScalar, sumWeightScalar;
286 std::tie(resultScalar, sumWeightScalar) = computeScalar(stepSize, firstEvent, lastEvent);
287 double carryScalar = resultScalar.Carry();
288
289 constexpr bool alwaysPrint = false;
290
291 if (alwaysPrint || std::abs(result - resultScalar)/resultScalar > 5.E-15) {
292 std::cerr << "RooNLLVar: result is off\n\t" << std::setprecision(15) << result
293 << "\n\t" << resultScalar << std::endl;
294 }
295
296 if (alwaysPrint || std::abs(carry - carryScalar)/carryScalar > 500.) {
297 std::cerr << "RooNLLVar: carry is far off\n\t" << std::setprecision(15) << carry
298 << "\n\t" << carryScalar << std::endl;
299 }
300
301 if (alwaysPrint || std::abs(sumWeight - sumWeightScalar)/sumWeightScalar > 1.E-15) {
302 std::cerr << "RooNLLVar: sumWeight is off\n\t" << std::setprecision(15) << sumWeight
303 << "\n\t" << sumWeightScalar << std::endl;
304 }
305
306#endif
307 } else { //scalar mode
308 std::tie(result, sumWeight) = computeScalar(stepSize, firstEvent, lastEvent);
309 }
310
311 // include the extended maximum likelihood term, if requested
312 if(_extended && _setNum==_extSet) {
313 result += pdfClone->extendedTerm(*_dataClone, _weightSq);
314 }
315 } //unbinned PDF
316
317
318 // If part of simultaneous PDF normalize probability over
319 // number of simultaneous PDFs: -sum(log(p/n)) = -sum(log(p)) + N*log(n)
320 if (_simCount>1) {
321 result += sumWeight * std::log(static_cast<double>(_simCount));
322 }
323
324
325 // At the end of the first full calculation, wire the caches
326 if (_first) {
327 _first = false ;
329 }
330
331
332 // Check if value offset flag is set.
333 if (_doOffset) {
334
335 // If no offset is stored enable this feature now
336 if (_offset==0 && result !=0 ) {
337 coutI(Minimization) << "RooNLLVar::evaluatePartition(" << GetName() << ") first = "<< firstEvent << " last = " << lastEvent << " Likelihood offset now set to " << result << std::endl ;
338 _offset = result ;
339 }
340
341 // Subtract offset
342 result -= _offset;
343 }
344
345 _evalCarry = result.Carry();
346 return result.Sum() ;
347}
348
349
350////////////////////////////////////////////////////////////////////////////////
351/// Compute probabilites of all data events. Use faster batch interface.
352/// \param[in] stepSize Stride when moving through the dataset.
353/// \note For batch computations, the step size **must** be one.
354/// \param[in] firstEvent First event to be processed.
355/// \param[in] lastEvent First event not to be processed.
356/// \return Tuple with (Kahan sum of probabilities, carry of kahan sum, sum of weights)
357RooNLLVar::ComputeResult RooNLLVar::computeBatched(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
358{
359 auto pdfClone = static_cast<const RooAbsPdf*>(_funcClone);
360 return computeBatchedFunc(pdfClone, _dataClone, _evalData, _normSet, _weightSq, stepSize, firstEvent, lastEvent);
361}
362
363// static function, also used from TestStatistics::RooUnbinnedL
365 std::unique_ptr<RooBatchCompute::RunContext> &evalData,
366 RooArgSet *normSet, bool weightSq, std::size_t stepSize,
367 std::size_t firstEvent, std::size_t lastEvent)
368{
369 const auto nEvents = lastEvent - firstEvent;
370
371 if (stepSize != 1) {
372 throw std::invalid_argument(std::string("Error in ") + __FILE__ + ": Step size for batch computations can only be 1.");
373 }
374
375 // Create a RunContext that will own the memory where computation results are stored.
376 // Holding on to this struct in between function calls will make sure that the memory
377 // is only allocated once.
378 if (!evalData) {
379 evalData.reset(new RooBatchCompute::RunContext);
380 }
381 evalData->clear();
382 dataClone->getBatches(*evalData, firstEvent, nEvents);
383
384 auto results = pdfClone->getLogProbabilities(*evalData, normSet);
385
386#ifdef ROOFIT_CHECK_CACHED_VALUES
387
388 for (std::size_t evtNo = firstEvent; evtNo < std::min(lastEvent, firstEvent + 10); ++evtNo) {
389 dataClone->get(evtNo);
390 if (dataClone->weight() == 0.) // 0-weight events are not cached, so cannot compare against them.
391 continue;
392
393 try {
394 // Cross check results with strict tolerance and complain
395 BatchInterfaceAccessor::checkBatchComputation(*pdfClone, *evalData, evtNo-firstEvent, normSet, 1.E-13);
396 } catch (std::exception& e) {
397 std::cerr << __FILE__ << ":" << __LINE__ << " ERROR when checking batch computation for event " << evtNo << ":\n"
398 << e.what() << std::endl;
399
400 // It becomes a real problem if it's very wrong. We fail in this case:
401 try {
402 BatchInterfaceAccessor::checkBatchComputation(*pdfClone, *evalData, evtNo-firstEvent, normSet, 1.E-9);
403 } catch (std::exception& e2) {
404 assert(false);
405 }
406 }
407 }
408
409#endif
410
411
412 // Compute sum of event weights. First check if we need squared weights
413 const RooSpan<const double> eventWeights = dataClone->getWeightBatch(firstEvent, nEvents, weightSq);
414
415 //Sum the event weights and probabilities
417 double uniformSingleEventWeight{0.0};
418 double sumOfWeights;
419 if (eventWeights.empty()) {
420 uniformSingleEventWeight = weightSq ? dataClone->weightSquared() : dataClone->weight();
421 sumOfWeights = nEvents * uniformSingleEventWeight;
422 for (std::size_t i = 0; i < results.size(); ++i) { //CHECK_VECTORISE
423 kahanProb.AddIndexed(-uniformSingleEventWeight * results[i], i);
424 }
425 } else {
426 assert(results.size() == eventWeights.size());
428 for (std::size_t i = 0; i < results.size(); ++i) { //CHECK_VECTORISE
429 const double weight = eventWeights[i];
430 kahanProb.AddIndexed(-weight * results[i], i);
431 kahanWeight.AddIndexed(weight, i);
432 }
433 sumOfWeights = kahanWeight.Sum();
434 }
435
436 if (std::isnan(kahanProb.Sum())) {
437 // Special handling of evaluation errors.
438 // We can recover if the bin/event that results in NaN has a weight of zero:
440 RooNaNPacker nanPacker;
441 for (std::size_t i = 0; i < results.size(); ++i) {
442 double weight = eventWeights.empty() ? uniformSingleEventWeight : eventWeights[i];
443
444 if (weight == 0.)
445 continue;
446
447 if (std::isnan(results[i])) {
448 nanPacker.accumulate(results[i]);
449 } else {
450 kahanSanitised += -weight * results[i];
451 }
452 }
453
454 // Some events with evaluation errors. Return "badness" of errors.
455 if (nanPacker.getPayload() > 0.) {
456 return {{nanPacker.getNaNWithPayload()}, sumOfWeights};
457 } else {
458 return {kahanSanitised, sumOfWeights};
459 }
460 }
461
462 return {kahanProb, sumOfWeights};
463}
464
465
466RooNLLVar::ComputeResult RooNLLVar::computeScalar(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const {
467 auto pdfClone = static_cast<const RooAbsPdf*>(_funcClone);
468 return computeScalarFunc(pdfClone, _dataClone, _normSet, _weightSq, stepSize, firstEvent, lastEvent);
469}
470
471// static function, also used from TestStatistics::RooUnbinnedL
473 RooArgSet *normSet, bool weightSq, std::size_t stepSize,
474 std::size_t firstEvent, std::size_t lastEvent)
475{
478 RooNaNPacker packedNaN(0.f);
479
480 for (auto i=firstEvent; i<lastEvent; i+=stepSize) {
481 dataClone->get(i) ;
482
483 double eventWeight = dataClone->weight(); //FIXME
484 if (0. == eventWeight * eventWeight) continue ;
485 if (weightSq) eventWeight = dataClone->weightSquared() ;
486
487 const double term = -eventWeight * pdfClone->getLogVal(normSet);
488
489 kahanWeight.Add(eventWeight);
490 kahanProb.Add(term);
491 packedNaN.accumulate(term);
492 }
493
494 if (packedNaN.getPayload() != 0.) {
495 // Some events with evaluation errors. Return "badness" of errors.
496 return {{packedNaN.getNaNWithPayload()}, kahanWeight.Sum()};
497 }
498
499 return {kahanProb, kahanWeight.Sum()};
500}
#define e(i)
Definition: RSha256.hxx:103
#define coutI(a)
Definition: RooMsgService.h:34
#define ClassImp(name)
Definition: Rtypes.h:375
#define N
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
char name[80]
Definition: TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2447
static void checkBatchComputation(const RooAbsReal &theReal, const RooBatchCompute::RunContext &evalData, std::size_t evtNo, const RooArgSet *normSet=nullptr, double relAccuracy=1.E-13)
Definition: RooAbsReal.h:595
T Sum() const
Definition: Util.h:240
T Carry() const
Definition: Util.h:255
void AddIndexed(T input, std::size_t index)
Add input to the sum.
Definition: Util.h:231
void Add(T x)
Single-element accumulation. Will not vectorise.
Definition: Util.h:165
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:317
void wireAllCaches()
Definition: RooAbsArg.cxx:2384
void setValueDirty()
Mark the element dirty. This forces a re-evaluation when a value is requested.
Definition: RooAbsArg.h:513
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:282
Storage_t::size_type size() const
RooAbsArg * first() const
virtual void recalculateCache(const RooArgSet *, Int_t, Int_t, Int_t, bool)
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:61
virtual double weight() const =0
virtual const RooArgSet * get() const
Definition: RooAbsData.h:105
RooAbsDataStore * store()
Definition: RooAbsData.h:81
virtual RooSpan< const double > getWeightBatch(std::size_t first, std::size_t len, bool sumW2=false) const =0
Return event weights of all events in range [first, first+len).
virtual void getBatches(RooBatchCompute::RunContext &evalData, std::size_t first=0, std::size_t len=std::numeric_limits< std::size_t >::max()) const
Write information to retrieve data columns into evalData.spans.
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.
RooArgSet * _normSet
Pointer to set with observables used for normalization.
RooAbsData * _dataClone
Pointer to internal clone if input data.
RooArgSet * _projDeps
Set of projected observable.
virtual double getLogVal(const RooArgSet *set=0) const
Return the log of the current value with given normalization An error message is printed if the argum...
Definition: RooAbsPdf.cxx:642
RooSpan< const double > getLogProbabilities(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet=nullptr) const
Compute the log-likelihoods for all events in the requested batch.
Definition: RooAbsPdf.cxx:722
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:263
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
void logEvalError(const char *message, const char *serverValueString=0) 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 _nGof
Number of sub-contexts.
Int_t _nCPU
Number of processors to use in parallel calculation mode.
pRooAbsTestStatistic * _gofArray
! Array of sub-contexts representing part of the combined test statistic
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
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:57
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.
Definition: RooCmdConfig.h:31
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.
static std::string decodeStringOnTheFly(const char *callerID, const char *cmdArgName, Int_t intIdx, const char *defVal, Args_t &&...args)
Static decoder function allows to retrieve string property from set of RooCmdArgs For use in base mem...
Definition: RooCmdConfig.h:200
static Int_t decodeIntOnTheFly(const char *callerID, const char *cmdArgName, Int_t intIdx, Int_t defVal, Args_t &&...args)
Static decoder function allows to retrieve integer property from set of RooCmdArgs For use in base me...
Definition: RooCmdConfig.h:184
Class RooNLLVar implements a -log(likelihood) calculation from a dataset and a PDF.
Definition: RooNLLVar.h:30
ComputeResult computeScalar(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
Definition: RooNLLVar.cxx:466
RooRealSumPdf * _binnedPdf
!
Definition: RooNLLVar.h:92
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)
Definition: RooNLLVar.cxx:472
ROOT::Math::KahanSum< double > _offsetSaveW2
!
Definition: RooNLLVar.h:89
static RooNLLVar::ComputeResult computeBatchedFunc(const RooAbsPdf *pdfClone, RooAbsData *dataClone, std::unique_ptr< RooBatchCompute::RunContext > &evalData, RooArgSet *normSet, bool weightSq, std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent)
Definition: RooNLLVar.cxx:364
bool _batchEvaluations
Definition: RooNLLVar.h:86
static RooArgSet _emptySet
Definition: RooNLLVar.h:79
void applyWeightSquared(bool flag) override
Disables or enables the usage of squared weights.
Definition: RooNLLVar.cxx:199
std::vector< double > _binw
!
Definition: RooNLLVar.h:91
std::pair< ROOT::Math::KahanSum< double >, double > ComputeResult
Definition: RooNLLVar.h:64
bool _extended
Definition: RooNLLVar.h:85
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.
Definition: RooNLLVar.cxx:185
~RooNLLVar() override
Definition: RooNLLVar.cxx:73
bool _first
!
Definition: RooNLLVar.h:88
std::unique_ptr< RooBatchCompute::RunContext > _evalData
! Struct to store function evaluation workspaces.
Definition: RooNLLVar.h:93
ComputeResult computeBatched(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
Compute probabilites of all data events.
Definition: RooNLLVar.cxx:357
bool _weightSq
Apply weights squared?
Definition: RooNLLVar.h:87
double evaluatePartition(std::size_t firstEvent, std::size_t lastEvent, std::size_t stepSize) const override
Calculate and return likelihood on subset of data.
Definition: RooNLLVar.cxx:227
The class RooRealSumPdf implements a PDF constructed from a sum of functions:
Definition: RooRealSumPdf.h:24
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Retrieve bin boundaries if this distribution is binned in obs.
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
constexpr std::span< T >::index_type size() const noexcept
Definition: RooSpan.h:121
constexpr bool empty() const noexcept
Definition: RooSpan.h:125
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
RVec< PromoteType< T > > abs(const RVec< T > &v)
Definition: RVec.hxx:1739
RVec< PromoteType< T > > log(const RVec< T > &v)
Definition: RVec.hxx:1748
void swap(RDirectoryEntry &e1, RDirectoryEntry &e2) noexcept
@ BulkPartition
Definition: RooGlobalFunc.h:66
@ Minimization
Definition: RooGlobalFunc.h:63
static constexpr double pc
Double_t LnGamma(Double_t z)
Computation of ln[gamma(z)] for all z.
Definition: TMath.cxx:486
std::string rangeName
Stores the configuration parameters for RooAbsTestStatistic.
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:31
Little struct that can pack a float into the unused bits of the mantissa of a NaN double.
Definition: RooNaNPacker.h:28
float getPayload() const
Retrieve packed float.
Definition: RooNaNPacker.h:85
double getNaNWithPayload() const
Retrieve a NaN with the current float payload packed into the mantissa.
Definition: RooNaNPacker.h:90
void accumulate(double val)
Accumulate a packed float from another NaN into this.
Definition: RooNaNPacker.h:57