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 
22 Class RooNLLVar implements a -log(likelihood) calculation from a dataset
23 and a PDF. The NLL is calculated as
24 \f[
25  \sum_\mathrm{data} -\log( \mathrm{pdf}(x_\mathrm{data}))
26 \f]
27 In 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 
54 namespace {
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 
73 
75 { }
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 kTRUE)
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.
92 RooNLLVar::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,kFALSE) ;
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  _weightSq = kFALSE ;
114  _first = kTRUE ;
115  _offsetSaveW2 = 0.;
116 
117  _binnedPdf = 0 ;
118 }
119 
120 
121 
122 ////////////////////////////////////////////////////////////////////////////////
123 /// Construct likelihood from given p.d.f and (binned or unbinned dataset)
124 /// For internal use.
125 
126 RooNLLVar::RooNLLVar(const char *name, const char *title, RooAbsPdf& pdf, RooAbsData& indata,
127  RooAbsTestStatistic::Configuration const& cfg, bool extended) :
128  RooAbsOptTestStatistic(name,title,pdf,indata,RooArgSet(),cfg),
129  _extended(extended),
130  _weightSq(kFALSE),
131  _first(kTRUE)
132 {
133  // If binned likelihood flag is set, pdf is a RooRealSumPdf representing a yield vector
134  // for a binned likelihood calculation
136 
137  // Retrieve and cache bin widths needed to convert un-normalized binnedPdf values back to yields
138  if (_binnedPdf) {
139 
140  // The Active label will disable pdf integral calculations
141  _binnedPdf->setAttribute("BinnedLikelihoodActive") ;
142 
144  if (obs->getSize()!=1) {
145  _binnedPdf = 0 ;
146  } else {
147  RooRealVar* var = (RooRealVar*) obs->first() ;
148  std::list<Double_t>* boundaries = _binnedPdf->binBoundaries(*var,var->getMin(),var->getMax()) ;
149  std::list<Double_t>::iterator biter = boundaries->begin() ;
150  _binw.resize(boundaries->size()-1) ;
151  Double_t lastBound = (*biter) ;
152  ++biter ;
153  int ibin=0 ;
154  while (biter!=boundaries->end()) {
155  _binw[ibin] = (*biter) - lastBound ;
156  lastBound = (*biter) ;
157  ibin++ ;
158  ++biter ;
159  }
160  }
161  }
162 }
163 
164 
165 
166 ////////////////////////////////////////////////////////////////////////////////
167 /// Construct likelihood from given p.d.f and (binned or unbinned dataset)
168 /// For internal use.
169 
170 RooNLLVar::RooNLLVar(const char *name, const char *title, RooAbsPdf& pdf, RooAbsData& indata,
171  const RooArgSet& projDeps,
172  RooAbsTestStatistic::Configuration const& cfg, bool extended) :
173  RooAbsOptTestStatistic(name,title,pdf,indata,projDeps, cfg),
174  _extended(extended),
175  _weightSq(kFALSE),
176  _first(kTRUE)
177 {
178  // If binned likelihood flag is set, pdf is a RooRealSumPdf representing a yield vector
179  // for a binned likelihood calculation
181 
182  // Retrieve and cache bin widths needed to convert un-normalized binnedPdf values back to yields
183  if (_binnedPdf) {
184 
186  if (obs->getSize()!=1) {
187  _binnedPdf = 0 ;
188  } else {
189  RooRealVar* var = (RooRealVar*) obs->first() ;
190  std::list<Double_t>* boundaries = _binnedPdf->binBoundaries(*var,var->getMin(),var->getMax()) ;
191  std::list<Double_t>::iterator biter = boundaries->begin() ;
192  _binw.resize(boundaries->size()-1) ;
193  Double_t lastBound = (*biter) ;
194  ++biter ;
195  int ibin=0 ;
196  while (biter!=boundaries->end()) {
197  _binw[ibin] = (*biter) - lastBound ;
198  lastBound = (*biter) ;
199  ibin++ ;
200  ++biter ;
201  }
202  }
203  }
204 }
205 
206 
207 
208 ////////////////////////////////////////////////////////////////////////////////
209 /// Copy constructor
210 
211 RooNLLVar::RooNLLVar(const RooNLLVar& other, const char* name) :
213  _extended(other._extended),
214  _batchEvaluations(other._batchEvaluations),
215  _weightSq(other._weightSq),
216  _first(kTRUE),
217  _offsetSaveW2(other._offsetSaveW2),
218  _binw(other._binw) {
220 }
221 
222 
223 ////////////////////////////////////////////////////////////////////////////////
224 /// Create a test statistic using several properties of the current instance. This is used to duplicate
225 /// the test statistic in multi-processing scenarios.
226 RooAbsTestStatistic* RooNLLVar::create(const char *name, const char *title, RooAbsReal& pdf, RooAbsData& adata,
227  const RooArgSet& projDeps, RooAbsTestStatistic::Configuration const& cfg) {
228  RooAbsPdf & thePdf = dynamic_cast<RooAbsPdf&>(pdf);
229  // check if pdf can be extended
230  bool extendedPdf = _extended && thePdf.canBeExtended();
231 
232  auto testStat = new RooNLLVar(name, title, thePdf, adata, projDeps, cfg, extendedPdf);
233  testStat->batchMode(_batchEvaluations);
234  return testStat;
235 }
236 
237 
238 ////////////////////////////////////////////////////////////////////////////////
239 /// Destructor
240 
242 {
243 }
244 
245 
246 
247 
248 ////////////////////////////////////////////////////////////////////////////////
249 
251 {
252  if (_gofOpMode==Slave) {
253  if (flag != _weightSq) {
254  _weightSq = flag;
256  }
257  setValueDirty();
258  } else if ( _gofOpMode==MPMaster) {
259  for (Int_t i=0 ; i<_nCPU ; i++)
260  _mpfeArray[i]->applyNLLWeightSquared(flag);
261  } else if ( _gofOpMode==SimMaster) {
262  for (Int_t i=0 ; i<_nGof ; i++)
263  ((RooNLLVar*)_gofArray[i])->applyWeightSquared(flag);
264  }
265 }
266 
267 
268 ////////////////////////////////////////////////////////////////////////////////
269 /// Calculate and return likelihood on subset of data.
270 /// \param[in] firstEvent First event to be processed.
271 /// \param[in] lastEvent First event not to be processed, any more.
272 /// \param[in] stepSize Steps between events.
273 /// \note For batch computations, the step size **must** be one.
274 ///
275 /// If this an extended likelihood, the extended term is added to the return likelihood
276 /// in the batch that encounters the event with index 0.
277 
278 Double_t RooNLLVar::evaluatePartition(std::size_t firstEvent, std::size_t lastEvent, std::size_t stepSize) const
279 {
280  // Throughout the calculation, we use Kahan's algorithm for summing to
281  // prevent loss of precision - this is a factor four more expensive than
282  // straight addition, but since evaluating the PDF is usually much more
283  // expensive than that, we tolerate the additional cost...
284  ROOT::Math::KahanSum<double> result{0.0};
285  double sumWeight{0.0};
286 
287  RooAbsPdf* pdfClone = (RooAbsPdf*) _funcClone ;
288 
289  // cout << "RooNLLVar::evaluatePartition(" << GetName() << ") projDeps = " << (_projDeps?*_projDeps:RooArgSet()) << endl ;
290 
291  _dataClone->store()->recalculateCache( _projDeps, firstEvent, lastEvent, stepSize, (_binnedPdf?kFALSE:kTRUE) ) ;
292 
293 
294 
295  // If pdf is marked as binned - do a binned likelihood calculation here (sum of log-Poisson for each bin)
296  if (_binnedPdf) {
297  ROOT::Math::KahanSum<double> sumWeightKahanSum{0.0};
298  for (auto i=firstEvent ; i<lastEvent ; i+=stepSize) {
299 
300  _dataClone->get(i) ;
301 
302  if (!_dataClone->valid()) continue;
303 
304  Double_t eventWeight = _dataClone->weight();
305 
306 
307  // Calculate log(Poisson(N|mu) for this bin
308  Double_t N = eventWeight ;
309  Double_t mu = _binnedPdf->getVal()*_binw[i] ;
310  //cout << "RooNLLVar::binnedL(" << GetName() << ") N=" << N << " mu = " << mu << endl ;
311 
312  if (mu<=0 && N>0) {
313 
314  // Catch error condition: data present where zero events are predicted
315  logEvalError(Form("Observed %f events in bin %lu with zero event yield",N,(unsigned long)i)) ;
316 
317  } else if (fabs(mu)<1e-10 && fabs(N)<1e-10) {
318 
319  // Special handling of this case since log(Poisson(0,0)=0 but can't be calculated with usual log-formula
320  // since log(mu)=0. No update of result is required since term=0.
321 
322  } else {
323 
324  result += -1*(-mu + N*log(mu) - TMath::LnGamma(N+1));
325  sumWeightKahanSum += eventWeight;
326 
327  }
328  }
329 
330  sumWeight = sumWeightKahanSum.Sum();
331 
332  } else { //unbinned PDF
333 
334  if (_batchEvaluations) {
335  std::tie(result, sumWeight) = computeBatched(stepSize, firstEvent, lastEvent);
336 #ifdef ROOFIT_CHECK_CACHED_VALUES
337 
338  ROOT::Math::KahanSum<double> resultScalar, sumWeightScalar;
339  std::tie(resultScalar, sumWeightScalar) = computeScalar(stepSize, firstEvent, lastEvent);
340  double carryScalar = resultScalar.Carry();
341 
342  constexpr bool alwaysPrint = false;
343 
344  if (alwaysPrint || fabs(result - resultScalar)/resultScalar > 5.E-15) {
345  std::cerr << "RooNLLVar: result is off\n\t" << std::setprecision(15) << result
346  << "\n\t" << resultScalar << std::endl;
347  }
348 
349  if (alwaysPrint || fabs(carry - carryScalar)/carryScalar > 500.) {
350  std::cerr << "RooNLLVar: carry is far off\n\t" << std::setprecision(15) << carry
351  << "\n\t" << carryScalar << std::endl;
352  }
353 
354  if (alwaysPrint || fabs(sumWeight - sumWeightScalar)/sumWeightScalar > 1.E-15) {
355  std::cerr << "RooNLLVar: sumWeight is off\n\t" << std::setprecision(15) << sumWeight
356  << "\n\t" << sumWeightScalar << std::endl;
357  }
358 
359 #endif
360  } else { //scalar mode
361  std::tie(result, sumWeight) = computeScalar(stepSize, firstEvent, lastEvent);
362  }
363 
364  // include the extended maximum likelihood term, if requested
365  if(_extended && _setNum==_extSet) {
366  if (_weightSq) {
367 
368 
369  // Calculate sum of weights-squared here for extended term
370  Double_t sumW2;
371  if (_batchEvaluations) {
372  const RooSpan<const double> eventWeights = _dataClone->getWeightBatch(0, _nEvents);
373  if (eventWeights.empty()) {
374  sumW2 = (lastEvent - firstEvent) * _dataClone->weightSquared();
375  } else {
377  for (std::size_t i = 0; i < eventWeights.size(); ++i) {
378  kahanWeight.AddIndexed(eventWeights[i] * eventWeights[i], i);
379  }
380  sumW2 = kahanWeight.Sum();
381  }
382  } else { // scalar mode
383  ROOT::Math::KahanSum<double> sumW2KahanSum;
384  for (decltype(_dataClone->numEntries()) i = 0; i < _dataClone->numEntries() ; i++) {
385  _dataClone->get(i);
386  sumW2KahanSum += _dataClone->weightSquared();
387  }
388  sumW2 = sumW2KahanSum.Sum();
389  }
390 
391  Double_t expected= pdfClone->expectedEvents(_dataClone->get());
392 
393  // Adjust calculation of extended term with W^2 weighting: adjust poisson such that
394  // estimate of Nexpected stays at the same value, but has a different variance, rescale
395  // both the observed and expected count of the Poisson with a factor sum[w] / sum[w^2] which is
396  // the effective weight of the Poisson term.
397  // i.e. change Poisson(Nobs = sum[w]| Nexp ) --> Poisson( sum[w] * sum[w] / sum[w^2] | Nexp * sum[w] / sum[w^2] )
398  // weighted by the effective weight sum[w^2]/ sum[w] in the likelihood.
399  // Since here we compute the likelihood with the weight square we need to multiply by the
400  // square of the effective weight
401  // expectedW = expected * sum[w] / sum[w^2] : effective expected entries
402  // observedW = sum[w] * sum[w] / sum[w^2] : effective observed entries
403  // The extended term for the likelihood weighted by the square of the weight will be then:
404  // (sum[w^2]/ sum[w] )^2 * expectedW - (sum[w^2]/ sum[w] )^2 * observedW * log (expectedW) and this is
405  // using the previous expressions for expectedW and observedW
406  // sum[w^2] / sum[w] * expected - sum[w^2] * log (expectedW)
407  // and since the weights are constants in the likelihood we can use log(expected) instead of log(expectedW)
408 
409  Double_t expectedW2 = expected * sumW2 / _dataClone->sumEntries() ;
410  Double_t extra= expectedW2 - sumW2*log(expected );
411 
412  // Double_t extra = pdfClone->extendedTerm(sumW2, _dataClone->get());
413 
414  result += extra;
415  } else {
416  result += pdfClone->extendedTerm(_dataClone->sumEntries(), _dataClone->get());
417  }
418  }
419  } //unbinned PDF
420 
421 
422  // If part of simultaneous PDF normalize probability over
423  // number of simultaneous PDFs: -sum(log(p/n)) = -sum(log(p)) + N*log(n)
424  if (_simCount>1) {
425  result += sumWeight * log(1.0*_simCount);
426  }
427 
428 
429  // At the end of the first full calculation, wire the caches
430  if (_first) {
431  _first = kFALSE ;
433  }
434 
435 
436  // Check if value offset flag is set.
437  if (_doOffset) {
438 
439  // If no offset is stored enable this feature now
440  if (_offset==0 && result !=0 ) {
441  coutI(Minimization) << "RooNLLVar::evaluatePartition(" << GetName() << ") first = "<< firstEvent << " last = " << lastEvent << " Likelihood offset now set to " << result << std::endl ;
442  _offset = result ;
443  }
444 
445  // Subtract offset
446  result -= _offset;
447  }
448 
449  _evalCarry = result.Carry();
450  return result.Sum() ;
451 }
452 
453 
454 ////////////////////////////////////////////////////////////////////////////////
455 /// Compute probabilites of all data events. Use faster batch interface.
456 /// \param[in] stepSize Stride when moving through the dataset.
457 /// \note For batch computations, the step size **must** be one.
458 /// \param[in] firstEvent First event to be processed.
459 /// \param[in] lastEvent First event not to be processed.
460 /// \return Tuple with (Kahan sum of probabilities, carry of kahan sum, sum of weights)
461 RooNLLVar::ComputeResult RooNLLVar::computeBatched(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
462 {
463  const auto nEvents = lastEvent - firstEvent;
464 
465  if (stepSize != 1) {
466  throw std::invalid_argument(std::string("Error in ") + __FILE__ + ": Step size for batch computations can only be 1.");
467  }
468 
469  auto pdfClone = static_cast<const RooAbsPdf*>(_funcClone);
470 
471  // Create a RunContext that will own the memory where computation results are stored.
472  // Holding on to this struct in between function calls will make sure that the memory
473  // is only allocated once.
474  if (!_evalData) {
476  }
477  _evalData->clear();
478  _dataClone->getBatches(*_evalData, firstEvent, nEvents);
479 
480  auto results = pdfClone->getLogProbabilities(*_evalData, _normSet);
481 
482 #ifdef ROOFIT_CHECK_CACHED_VALUES
483 
484  for (std::size_t evtNo = firstEvent; evtNo < std::min(lastEvent, firstEvent + 10); ++evtNo) {
485  _dataClone->get(evtNo);
486  if (_dataClone->weight() == 0.) // 0-weight events are not cached, so cannot compare against them.
487  continue;
488 
489  assert(_dataClone->valid());
490  try {
491  // Cross check results with strict tolerance and complain
492  BatchInterfaceAccessor::checkBatchComputation(*pdfClone, *_evalData, evtNo-firstEvent, _normSet, 1.E-13);
493  } catch (std::exception& e) {
494  std::cerr << __FILE__ << ":" << __LINE__ << " ERROR when checking batch computation for event " << evtNo << ":\n"
495  << e.what() << std::endl;
496 
497  // It becomes a real problem if it's very wrong. We fail in this case:
498  try {
499  BatchInterfaceAccessor::checkBatchComputation(*pdfClone, *_evalData, evtNo-firstEvent, _normSet, 1.E-9);
500  } catch (std::exception& e2) {
501  assert(false);
502  }
503  }
504  }
505 
506 #endif
507 
508 
509  // Compute sum of event weights. First check if we need squared weights
510  const RooSpan<const double> eventWeights = _dataClone->getWeightBatch(firstEvent, nEvents);
511  //Capture member for lambda:
512  const bool retrieveSquaredWeights = _weightSq;
513  auto retrieveWeight = [&eventWeights, retrieveSquaredWeights](std::size_t i) {
514  return retrieveSquaredWeights ? eventWeights[i] * eventWeights[i] : eventWeights[i];
515  };
516 
517  //Sum the event weights and probabilities
519  double uniformSingleEventWeight{0.0};
520  double sumOfWeights;
521  if (eventWeights.empty()) {
522  uniformSingleEventWeight = retrieveSquaredWeights ? _dataClone->weightSquared() : _dataClone->weight();
523  sumOfWeights = nEvents * uniformSingleEventWeight;
524  for (std::size_t i = 0; i < results.size(); ++i) { //CHECK_VECTORISE
525  kahanProb.AddIndexed(-uniformSingleEventWeight * results[i], i);
526  }
527  } else {
528  assert(results.size() == eventWeights.size());
530  for (std::size_t i = 0; i < results.size(); ++i) { //CHECK_VECTORISE
531  const double weight = retrieveWeight(i);
532  kahanProb.AddIndexed(-weight * results[i], i);
533  kahanWeight.AddIndexed(weight, i);
534  }
535  sumOfWeights = kahanWeight.Sum();
536  }
537 
538  if (std::isnan(kahanProb.Sum())) {
539  // Special handling of evaluation errors.
540  // We can recover if the bin/event that results in NaN has a weight of zero:
541  ROOT::Math::KahanSum<double, 4u> kahanSanitised;
542  RooNaNPacker nanPacker;
543  for (std::size_t i = 0; i < results.size(); ++i) {
544  double weight = eventWeights.empty() ? uniformSingleEventWeight : retrieveWeight(i);
545 
546  if (weight == 0.)
547  continue;
548 
549  if (std::isnan(results[i])) {
550  nanPacker.accumulate(results[i]);
551  } else {
552  kahanSanitised += -weight * results[i];
553  }
554  }
555 
556  // Some events with evaluation errors. Return "badness" of errors.
557  if (nanPacker.getPayload() > 0.) {
558  return {{nanPacker.getNaNWithPayload()}, sumOfWeights};
559  } else {
560  return {kahanSanitised, sumOfWeights};
561  }
562  }
563 
564  return {kahanProb, sumOfWeights};
565 }
566 
567 
568 RooNLLVar::ComputeResult RooNLLVar::computeScalar(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const {
569  auto pdfClone = static_cast<const RooAbsPdf*>(_funcClone);
570 
571  ROOT::Math::KahanSum<double> kahanWeight;
573  RooNaNPacker packedNaN(0.f);
574 
575  for (auto i=firstEvent; i<lastEvent; i+=stepSize) {
576  _dataClone->get(i) ;
577 
578  if (!_dataClone->valid()) continue;
579 
580  Double_t eventWeight = _dataClone->weight(); //FIXME
581  if (0. == eventWeight * eventWeight) continue ;
582  if (_weightSq) eventWeight = _dataClone->weightSquared() ;
583 
584  const double term = -eventWeight * pdfClone->getLogVal(_normSet);
585 
586  kahanWeight.Add(eventWeight);
587  kahanProb.Add(term);
588  packedNaN.accumulate(term);
589  }
590 
591  if (packedNaN.getPayload() != 0.) {
592  // Some events with evaluation errors. Return "badness" of errors.
593  return {{packedNaN.getNaNWithPayload()}, kahanWeight.Sum()};
594  }
595 
596  return {kahanProb, kahanWeight.Sum()};
597 }
RooAbsData::getBatches
virtual void getBatches(RooBatchCompute::RunContext &evalData, std::size_t first=0, std::size_t len=std::numeric_limits< std::size_t >::max()) const =0
Retrieve batches of data for each real-valued variable in this dataset.
Util.h
RooCmdArg
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:27
RooNaNPacker::getNaNWithPayload
double getNaNWithPayload() const
Retrieve a NaN with the current float payload packed into the mantissa.
Definition: RooNaNPacker.h:90
RooFit::Minimization
@ Minimization
Definition: RooGlobalFunc.h:60
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
RooCmdConfig.h
e
#define e(i)
Definition: RSha256.hxx:103
RooMsgService.h
RooNLLVar::evaluatePartition
virtual Double_t evaluatePartition(std::size_t firstEvent, std::size_t lastEvent, std::size_t stepSize) const
Calculate and return likelihood on subset of data.
Definition: RooNLLVar.cxx:278
RooAbsCollection::first
RooAbsArg * first() const
Definition: RooAbsCollection.h:236
RooAbsTestStatistic::_mpfeArray
pRooRealMPFE * _mpfeArray
Definition: RooAbsTestStatistic.h:154
RooAbsReal::logEvalError
void logEvalError(const char *message, const char *serverValueString=0) const
Log evaluation error message.
Definition: RooAbsReal.cxx:3710
RooNLLVar::computeBatched
ComputeResult computeBatched(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
Compute probabilites of all data events.
Definition: RooNLLVar.cxx:461
RooAbsData
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:49
RooAbsRealLValue::getMax
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Definition: RooAbsRealLValue.h:89
RooAbsTestStatistic::_gofOpMode
GOFOpMode _gofOpMode
Is object initialized
Definition: RooAbsTestStatistic.h:140
RooNLLVar::_binnedPdf
RooRealSumPdf * _binnedPdf
Definition: RooNLLVar.h:83
RooNaNPacker
Little struct that can pack a float into the unused bits of the mantissa of a NaN double.
Definition: RooNaNPacker.h:28
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
RooAbsTestStatistic::Configuration::cloneInputData
bool cloneInputData
Definition: RooAbsTestStatistic.h:51
RooAbsOptTestStatistic::_funcClone
RooAbsReal * _funcClone
Definition: RooAbsOptTestStatistic.h:81
RooRealMPFE.h
RooNaNPacker::accumulate
void accumulate(double val)
Accumulate a packed float from another NaN into this.
Definition: RooNaNPacker.h:57
RooAbsTestStatistic::_nEvents
Int_t _nEvents
Definition: RooAbsTestStatistic.h:142
RooAbsData::weight
virtual Double_t weight() const =0
BatchInterfaceAccessor::checkBatchComputation
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:588
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooAbsTestStatistic::_nGof
Int_t _nGof
Number of designated set to calculated extended term.
Definition: RooAbsTestStatistic.h:148
log
double log(double)
RooAbsTestStatistic::Configuration::rangeName
std::string rangeName
Stores the configuration parameters for RooAbsTestStatistic.
Definition: RooAbsTestStatistic.h:45
N
#define N
RooAbsTestStatistic::Configuration::verbose
bool verbose
Definition: RooAbsTestStatistic.h:49
TGeant4Unit::pc
static constexpr double pc
Definition: TGeant4SystemOfUnits.h:130
RooAbsData::valid
virtual Bool_t valid() const
Definition: RooAbsData.h:98
RooAbsData::store
RooAbsDataStore * store()
Definition: RooAbsData.h:68
coutI
#define coutI(a)
Definition: RooMsgService.h:30
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooAbsOptTestStatistic::_dataClone
RooAbsData * _dataClone
Definition: RooAbsOptTestStatistic.h:80
ROOT::Math::KahanSum::AddIndexed
void AddIndexed(T input, std::size_t index)
Add input to the sum.
Definition: Util.h:208
RooAbsTestStatistic::_evalCarry
Double_t _evalCarry
Offset as KahanSum to avoid loss of precision.
Definition: RooAbsTestStatistic.h:159
TMath::LnGamma
Double_t LnGamma(Double_t z)
Computation of ln[gamma(z)] for all z.
Definition: TMath.cxx:486
RooNLLVar.h
RooNLLVar::_extended
Bool_t _extended
Definition: RooNLLVar.h:76
RooAbsData::sumEntries
virtual Double_t sumEntries() const =0
Return effective number of entries in dataset, i.e., sum all weights.
bool
RooAbsRealLValue::getMin
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
Definition: RooAbsRealLValue.h:86
RooAbsTestStatistic::Configuration::addCoefRangeName
std::string addCoefRangeName
Definition: RooAbsTestStatistic.h:46
RooNLLVar::create
virtual RooAbsTestStatistic * create(const char *name, const char *title, RooAbsReal &pdf, RooAbsData &adata, const RooArgSet &projDeps, RooAbsTestStatistic::Configuration const &cfg)
Create a test statistic using several properties of the current instance.
Definition: RooNLLVar.cxx:226
RooSpan::size
constexpr std::span< T >::index_type size() const noexcept
Definition: RooSpan.h:121
RooCmdConfig
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
Definition: RooCmdConfig.h:27
RooNLLVar::_evalData
std::unique_ptr< RooBatchCompute::RunContext > _evalData
Definition: RooNLLVar.h:84
RooAbsPdf::extendedTerm
virtual Double_t extendedTerm(Double_t observedEvents, const RooArgSet *nset=0) const
Return the extended likelihood term ( ) of this PDF for the given number of observed events.
Definition: RooAbsPdf.cxx:784
RooAbsTestStatistic
RooAbsTestStatistic is the abstract base class for all test statistics.
Definition: RooAbsTestStatistic.h:39
RooAbsTestStatistic::Configuration::interleave
RooFit::MPSplit interleave
Definition: RooAbsTestStatistic.h:48
RooNLLVar::ComputeResult
std::pair< ROOT::Math::KahanSum< double >, double > ComputeResult
Definition: RooNLLVar.h:63
ROOT::Math::fabs
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
Definition: UnaryOperators.h:131
RooCmdConfig::decodeDoubleOnTheFly
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.
Definition: RooCmdConfig.cxx:897
RooProdPdf.h
RooNaNPacker::getPayload
float getPayload() const
Retrieve packed float.
Definition: RooNaNPacker.h:85
RooAbsTestStatistic::_gofArray
pRooAbsTestStatistic * _gofArray
Definition: RooAbsTestStatistic.h:149
RooNaNPacker.h
RooAbsPdf.h
RooAbsTestStatistic::_offset
ROOT::Math::KahanSum< double > _offset
Definition: RooAbsTestStatistic.h:158
RooSpan::empty
constexpr bool empty() const noexcept
Definition: RooSpan.h:125
RooAbsOptTestStatistic::_projDeps
RooArgSet * _projDeps
Definition: RooAbsOptTestStatistic.h:82
RooAbsDataStore::recalculateCache
virtual void recalculateCache(const RooArgSet *, Int_t, Int_t, Int_t, Bool_t)
Definition: RooAbsDataStore.h:111
RooNLLVar::_offsetSaveW2
ROOT::Math::KahanSum< double > _offsetSaveW2
Definition: RooNLLVar.h:80
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
RooNLLVar::RooNLLVar
RooNLLVar()
Definition: RooNLLVar.cxx:74
RooAbsData::get
virtual const RooArgSet * get() const
Definition: RooAbsData.h:92
RooNLLVar
Class RooNLLVar implements a -log(likelihood) calculation from a dataset and a PDF.
Definition: RooNLLVar.h:30
RooAbsTestStatistic::Configuration
Definition: RooAbsTestStatistic.h:43
isnan
int isnan(double)
RooAbsData::numEntries
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
Definition: RooAbsData.cxx:304
RooAbsTestStatistic::Slave
@ Slave
Definition: RooAbsTestStatistic.h:112
RooAbsTestStatistic::_doOffset
Bool_t _doOffset
Definition: RooAbsTestStatistic.h:157
RooNLLVar::_binw
std::vector< Double_t > _binw
Definition: RooNLLVar.h:82
RooNLLVar::_batchEvaluations
bool _batchEvaluations
Definition: RooNLLVar.h:77
RooRealVar.h
RooFit::BulkPartition
@ BulkPartition
Definition: RooGlobalFunc.h:63
RooNLLVar::_first
Bool_t _first
Definition: RooNLLVar.h:79
RooAbsDataStore.h
RooAbsTestStatistic::_nCPU
Int_t _nCPU
GOF MP Split mode specified by component (when Auto is active)
Definition: RooAbsTestStatistic.h:153
ROOT::Math::KahanSum::Add
void Add(T x)
Single-element accumulation. Will not vectorise.
Definition: Util.h:142
RooCmdConfig::decodeStringOnTheFly
static std::string decodeStringOnTheFly(const char *callerID, const char *cmdArgName, Int_t intIdx, const char *defVal, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg())
Static decoder function allows to retrieve string property from set of RooCmdArgs For use in base mem...
Definition: RooCmdConfig.cxx:854
RooAbsData.h
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:295
RooNLLVar::applyWeightSquared
void applyWeightSquared(Bool_t flag)
Definition: RooNLLVar.cxx:250
RooAbsTestStatistic::SimMaster
@ SimMaster
Definition: RooAbsTestStatistic.h:112
RooAbsTestStatistic::Configuration::nCPU
int nCPU
Definition: RooAbsTestStatistic.h:47
RooAbsData::weightSquared
virtual Double_t weightSquared() const =0
RooAbsOptTestStatistic
RooAbsOptTestStatistic is the abstract base class for test statistics objects that evaluate a functio...
Definition: RooAbsOptTestStatistic.h:28
RooAbsData::getWeightBatch
virtual RooSpan< const double > getWeightBatch(std::size_t first, std::size_t len) const =0
Return event weights of all events in range [first, first+len).
RooAbsArg::setAttribute
void setAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:291
RooCmdConfig::decodeIntOnTheFly
static Int_t decodeIntOnTheFly(const char *callerID, const char *cmdArgName, Int_t intIdx, Int_t defVal, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg())
Static decoder function allows to retrieve integer property from set of RooCmdArgs For use in base me...
Definition: RooCmdConfig.cxx:834
ROOT::Math::KahanSum::Sum
T Sum() const
Definition: Util.h:217
name
char name[80]
Definition: TGX11.cxx:110
ROOT::Experimental::Internal::swap
void swap(RDirectoryEntry &e1, RDirectoryEntry &e2) noexcept
Definition: RDirectoryEntry.hxx:94
RooAbsArg::wireAllCaches
void wireAllCaches()
Definition: RooAbsArg.cxx:2307
RunContext.h
RooNLLVar::_weightSq
Bool_t _weightSq
Definition: RooNLLVar.h:78
RooRealSumPdf::binBoundaries
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Retrieve bin boundaries if this distribution is binned in obs.
Definition: RooRealSumPdf.cxx:519
RooAbsTestStatistic::Configuration::integrateOverBinsPrecision
double integrateOverBinsPrecision
Definition: RooAbsTestStatistic.h:52
RooAbsArg::setValueDirty
void setValueDirty()
Mark the element dirty. This forces a re-evaluation when a value is requested.
Definition: RooAbsArg.h:490
RooAbsPdf
Definition: RooAbsPdf.h:41
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooAbsTestStatistic::MPMaster
@ MPMaster
Definition: RooAbsTestStatistic.h:112
RooNLLVar::_emptySet
static RooArgSet _emptySet
Definition: RooNLLVar.h:70
RooAbsOptTestStatistic::_normSet
RooArgSet * _normSet
Definition: RooAbsOptTestStatistic.h:78
RooRealSumPdf.h
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooAbsTestStatistic::Configuration::binnedL
bool binnedL
Definition: RooAbsTestStatistic.h:53
ROOT::Math::KahanSum< double >
RooAbsTestStatistic::Configuration::splitCutRange
bool splitCutRange
Definition: RooAbsTestStatistic.h:50
ROOT::Math::KahanSum::Carry
T Carry() const
Definition: Util.h:232
RooAbsTestStatistic::_extSet
Int_t _extSet
Definition: RooAbsTestStatistic.h:145
RooRealSumPdf
The class RooRealSumPdf implements a PDF constructed from a sum of functions:
Definition: RooRealSumPdf.h:24
RooAbsTestStatistic::_simCount
Int_t _simCount
Definition: RooAbsTestStatistic.h:125
RooBatchCompute::RunContext
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:31
RooAbsTestStatistic::_setNum
Int_t _setNum
Definition: RooAbsTestStatistic.h:143
RooAbsPdf::canBeExtended
Bool_t canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:236
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:231
RooSpan
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
TMath.h
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:33
int
RooAbsPdf::expectedEvents
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.cxx:3393
RooNLLVar::computeScalar
ComputeResult computeScalar(std::size_t stepSize, std::size_t firstEvent, std::size_t lastEvent) const
Definition: RooNLLVar.cxx:568
RooNLLVar::~RooNLLVar
virtual ~RooNLLVar()
Destructor.
Definition: RooNLLVar.cxx:241