Logo ROOT  
Reference Guide
RooAddPdf.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 /** \class RooAddPdf
19  \ingroup Roofitcore
20 
21 RooAddPdf is an efficient implementation of a sum of PDFs of the form
22 
23 \f[
24  \sum_{i=1}^{n} c_i \cdot \mathrm{PDF}_i
25 \f]
26 
27 or
28 \f[
29  c_1\cdot\mathrm{PDF}_1 + c_2\cdot\mathrm{PDF}_2 \; + \; ... \; + \; \left( 1-\sum_{i=1}^{n-1}c_i \right) \cdot \mathrm{PDF}_n
30 \f]
31 
32 The first form is for extended likelihood fits, where the
33 expected number of events is \f$ \sum_i c_i \f$. The coefficients \f$ c_i \f$
34 can either be explicitly provided, or, if all components support
35 extended likelihood fits, they can be calculated from the contribution
36 of each PDF to the total expected number of events.
37 
38 In the second form, the sum of the coefficients is required to be 1 or less,
39 and the coefficient of the last PDF is calculated automatically from the condition
40 that the sum of all coefficients has to be 1.
41 
42 ### Recursive coefficients
43 It is also possible to parameterise the coefficients recursively
44 
45 \f[
46  \sum_{i=1}^n c_i \prod_{j=1}^{i-1} \left[ (1-c_j) \right] \cdot \mathrm{PDF}_i \\
47  = c_1 \cdot \mathrm{PDF}_1 + (1-c_1)\, c_2 \cdot \mathrm{PDF}_2 + \ldots + (1-c_1)\ldots(1-c_{n-1}) \cdot 1 \cdot \mathrm{PDF}_n \\
48 \f]
49 
50 In this form the sum of the coefficients is always less than 1.0
51 for all possible values of the individual coefficients between 0 and 1.
52 \note Don't pass the \f$ n^\mathrm{th} \f$ coefficient. It is always 1, since the normalisation condition removes one degree of freedom.
53 
54 RooAddPdf relies on each component PDF to be normalized and will perform
55 no normalization other than calculating the proper last coefficient \f$ c_n \f$, if requested.
56 An (enforced) condition for this assumption is that each \f$ \mathrm{PDF}_i \f$ is independent of each \f$ c_i \f$.
57 
58 ## Difference between RooAddPdf / RooRealSumFunc / RooRealSumPdf
59 - RooAddPdf is a PDF of PDFs, *i.e.* its components need to be normalised and non-negative.
60 - RooRealSumPdf is a PDF of functions, *i.e.*, its components can be negative, but their sum cannot be. The normalisation
61  is computed automatically, unless the PDF is extended (see above).
62 - RooRealSumFunc is a sum of functions. It is neither normalised, nor need it be positive.
63 
64 */
65 
66 #include "RooAddPdf.h"
67 
68 #include "RooDataSet.h"
69 #include "RooRealProxy.h"
70 #include "RooRealVar.h"
71 #include "RooAddGenContext.h"
72 #include "RooRealConstant.h"
73 #include "RooRecursiveFraction.h"
74 #include "RooGlobalFunc.h"
75 #include "RooRealIntegral.h"
76 #include "RooNaNPacker.h"
77 #include "RooBatchCompute.h"
78 
79 #include <algorithm>
80 #include <sstream>
81 
82 using namespace std;
83 
85 
86 
87 ////////////////////////////////////////////////////////////////////////////////
88 /// Default constructor used for persistence
89 
91  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
92  _refCoefRangeName(0),
93  _projectCoefs(false),
94  _codeReg(10),
95  _snormList(0),
96  _haveLastCoef(false),
97  _allExtendable(false),
98  _recursive(false)
99 {
102 }
103 
104 
105 
106 ////////////////////////////////////////////////////////////////////////////////
107 /// Dummy constructor
108 
109 RooAddPdf::RooAddPdf(const char *name, const char *title) :
110  RooAbsPdf(name,title),
111  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
112  _refCoefRangeName(0),
113  _projectCoefs(kFALSE),
114  _projCacheMgr(this,10),
115  _codeReg(10),
116  _pdfList("!pdfs","List of PDFs",this),
117  _coefList("!coefficients","List of coefficients",this),
118  _snormList(0),
119  _haveLastCoef(kFALSE),
120  _allExtendable(kFALSE),
121  _recursive(kFALSE)
122 {
125 }
126 
127 
128 
129 ////////////////////////////////////////////////////////////////////////////////
130 /// Constructor with two PDFs and one coefficient
131 
132 RooAddPdf::RooAddPdf(const char *name, const char *title,
133  RooAbsPdf& pdf1, RooAbsPdf& pdf2, RooAbsReal& coef1) :
134  RooAbsPdf(name,title),
135  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
136  _refCoefRangeName(0),
137  _projectCoefs(kFALSE),
138  _projCacheMgr(this,10),
139  _codeReg(10),
140  _pdfList("!pdfs","List of PDFs",this),
141  _coefList("!coefficients","List of coefficients",this),
142  _haveLastCoef(kFALSE),
143  _allExtendable(kFALSE),
144  _recursive(kFALSE)
145 {
146  _pdfList.add(pdf1) ;
147  _pdfList.add(pdf2) ;
148  _coefList.add(coef1) ;
149 
150  _coefCache.resize(_pdfList.size());
153 }
154 
155 
156 
157 ////////////////////////////////////////////////////////////////////////////////
158 /// Generic constructor from list of PDFs and list of coefficients.
159 /// Each pdf list element (i) is paired with coefficient list element (i).
160 /// The number of coefficients must be either equal to the number of PDFs,
161 /// in which case extended MLL fitting is enabled, or be one less.
162 ///
163 /// All PDFs must inherit from RooAbsPdf. All coefficients must inherit from RooAbsReal
164 ///
165 /// If the recursiveFraction flag is true, the coefficients are interpreted as recursive
166 /// coefficients as explained in the class description.
167 
168 RooAddPdf::RooAddPdf(const char *name, const char *title, const RooArgList& inPdfList, const RooArgList& inCoefList, Bool_t recursiveFractions) :
169  RooAbsPdf(name,title),
170  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
171  _refCoefRangeName(0),
172  _projectCoefs(kFALSE),
173  _projCacheMgr(this,10),
174  _codeReg(10),
175  _pdfList("!pdfs","List of PDFs",this),
176  _coefList("!coefficients","List of coefficients",this),
177  _haveLastCoef(kFALSE),
178  _allExtendable(kFALSE),
179  _recursive(recursiveFractions)
180 {
181  if (inPdfList.getSize()>inCoefList.getSize()+1 || inPdfList.getSize()<inCoefList.getSize()) {
182  std::stringstream errorMsg;
183  errorMsg << "RooAddPdf::RooAddPdf(" << GetName()
184  << ") number of pdfs and coefficients inconsistent, must have Npdf=Ncoef or Npdf=Ncoef+1." << endl ;
185  coutE(InputArguments) << errorMsg.str();
186  throw std::invalid_argument(errorMsg.str().c_str());
187  }
188 
189  if (recursiveFractions && inPdfList.getSize()!=inCoefList.getSize()+1) {
190  std::stringstream errorMsg;
191  errorMsg << "RooAddPdf::RooAddPdf(" << GetName()
192  << "): Recursive fractions option can only be used if Npdf=Ncoef+1." << endl;
193  coutE(InputArguments) << errorMsg.str();
194  throw std::invalid_argument(errorMsg.str());
195  }
196 
197  // Constructor with N PDFs and N or N-1 coefs
198  RooArgList partinCoefList ;
199 
200  Bool_t first(kTRUE) ;
201 
202  for (auto i = 0u; i < inCoefList.size(); ++i) {
203  auto coef = dynamic_cast<RooAbsReal*>(inCoefList.at(i));
204  auto pdf = dynamic_cast<RooAbsPdf*>(inPdfList.at(i));
205  if (inPdfList.at(i) == nullptr) {
206  std::stringstream errorMsg;
207  errorMsg << "RooAddPdf::RooAddPdf(" << GetName()
208  << ") number of pdfs and coefficients inconsistent, must have Npdf=Ncoef or Npdf=Ncoef+1" << endl ;
209  coutE(InputArguments) << errorMsg.str();
210  throw std::invalid_argument(errorMsg.str());
211  }
212  if (!coef) {
213  std::stringstream errorMsg;
214  errorMsg << "RooAddPdf::RooAddPdf(" << GetName() << ") coefficient " << (coef ? coef->GetName() : "") << " is not of type RooAbsReal, ignored" << endl ;
215  coutE(InputArguments) << errorMsg.str();
216  throw std::invalid_argument(errorMsg.str());
217  }
218  if (!pdf) {
219  std::stringstream errorMsg;
220  errorMsg << "RooAddPdf::RooAddPdf(" << GetName() << ") pdf " << (pdf ? pdf->GetName() : "") << " is not of type RooAbsPdf, ignored" << endl ;
221  coutE(InputArguments) << errorMsg.str();
222  throw std::invalid_argument(errorMsg.str());
223  }
224  _pdfList.add(*pdf) ;
225 
226  // Process recursive fraction mode separately
227  if (recursiveFractions) {
228  partinCoefList.add(*coef) ;
229  if (first) {
230 
231  // The first fraction is the first plain fraction
232  first = kFALSE ;
233  _coefList.add(*coef) ;
234 
235  } else {
236 
237  // The i-th recursive fraction = (1-f1)*(1-f2)*...(fi) and is calculated from the list (f1,...,fi) by RooRecursiveFraction)
238  RooAbsReal* rfrac = new RooRecursiveFraction(Form("%s_recursive_fraction_%s",GetName(),pdf->GetName()),"Recursive Fraction",partinCoefList) ;
239  addOwnedComponents(*rfrac) ;
240  _coefList.add(*rfrac) ;
241 
242  }
243 
244  } else {
245  _coefList.add(*coef) ;
246  }
247  }
248 
249  if (inPdfList.size() == inCoefList.size() + 1) {
250  auto pdf = dynamic_cast<RooAbsPdf*>(inPdfList.at(inCoefList.size()));
251 
252  if (!pdf) {
253  coutE(InputArguments) << "RooAddPdf::RooAddPdf(" << GetName() << ") last argument " << inPdfList.at(inCoefList.size())->GetName() << " is not of type RooAbsPdf." << endl ;
254  throw std::invalid_argument("Last argument for RooAddPdf is not a PDF.");
255  }
256  _pdfList.add(*pdf) ;
257 
258  // Process recursive fractions mode. Above, we verified that we don't have a last coefficient
259  if (recursiveFractions) {
260 
261  // The last recursive fraction = (1-f1)*(1-f2)*...(1-fN) and is calculated from the list (f1,...,fN,1) by RooRecursiveFraction
262  partinCoefList.add(RooFit::RooConst(1)) ;
263  RooAbsReal* rfrac = new RooRecursiveFraction(Form("%s_recursive_fraction_%s",GetName(),pdf->GetName()),"Recursive Fraction",partinCoefList) ;
264  addOwnedComponents(*rfrac) ;
265  _coefList.add(*rfrac) ;
266 
267  // In recursive mode we always have Ncoef=Npdf, since we added it just above
269  }
270 
271  } else {
273  }
274 
275 
276  _coefCache.resize(_pdfList.size());
278 
280 }
281 
282 
283 
284 ////////////////////////////////////////////////////////////////////////////////
285 /// Generic constructor from list of extended PDFs. There are no coefficients as the expected
286 /// number of events from each components determine the relative weight of the PDFs.
287 ///
288 /// All PDFs must inherit from RooAbsPdf.
289 
290 RooAddPdf::RooAddPdf(const char *name, const char *title, const RooArgList& inPdfList) :
291  RooAbsPdf(name,title),
292  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
293  _refCoefRangeName(0),
294  _projectCoefs(kFALSE),
295  _projCacheMgr(this,10),
296  _pdfList("!pdfs","List of PDFs",this),
297  _coefList("!coefficients","List of coefficients",this),
298  _haveLastCoef(kFALSE),
299  _allExtendable(kTRUE),
300  _recursive(kFALSE)
301 {
302  // Constructor with N PDFs
303  for (const auto pdfArg : inPdfList) {
304  auto pdf = dynamic_cast<const RooAbsPdf*>(pdfArg);
305 
306  if (!pdf) {
307  coutE(InputArguments) << "RooAddPdf::RooAddPdf(" << GetName() << ") pdf " << (pdf ? pdf->GetName() : "") << " is not of type RooAbsPdf, ignored" << endl ;
308  continue ;
309  }
310  if (!pdf->canBeExtended()) {
311  coutE(InputArguments) << "RooAddPdf::RooAddPdf(" << GetName() << ") pdf " << pdf->GetName() << " is not extendable, ignored" << endl ;
312  continue ;
313  }
314  _pdfList.add(*pdf) ;
315  }
316 
317  _coefCache.resize(_pdfList.size());
320 }
321 
322 
323 
324 
325 ////////////////////////////////////////////////////////////////////////////////
326 /// Copy constructor
327 
328 RooAddPdf::RooAddPdf(const RooAddPdf& other, const char* name) :
329  RooAbsPdf(other,name),
330  _refCoefNorm("!refCoefNorm",this,other._refCoefNorm),
331  _refCoefRangeName((TNamed*)other._refCoefRangeName),
332  _projectCoefs(other._projectCoefs),
333  _projCacheMgr(other._projCacheMgr,this),
334  _codeReg(other._codeReg),
335  _pdfList("!pdfs",this,other._pdfList),
336  _coefList("!coefficients",this,other._coefList),
337  _haveLastCoef(other._haveLastCoef),
338  _allExtendable(other._allExtendable),
339  _recursive(other._recursive)
340 {
341  _coefCache.resize(_pdfList.size());
344 }
345 
346 
347 
348 ////////////////////////////////////////////////////////////////////////////////
349 /// Destructor
350 
352 {
354 }
355 
356 
357 
358 ////////////////////////////////////////////////////////////////////////////////
359 /// By default the interpretation of the fraction coefficients is
360 /// performed in the contextual choice of observables. This makes the
361 /// shape of the p.d.f explicitly dependent on the choice of
362 /// observables. This method instructs RooAddPdf to freeze the
363 /// interpretation of the coefficients to be done in the given set of
364 /// observables. If frozen, fractions are automatically transformed
365 /// from the reference normalization set to the contextual normalization
366 /// set by ratios of integrals.
367 
369 {
370  if (refCoefNorm.getSize()==0) {
372  return ;
373  }
374  _projectCoefs = kTRUE ;
375 
377  _refCoefNorm.add(refCoefNorm) ;
378 
379  _projCacheMgr.reset() ;
380 }
381 
382 
383 
384 ////////////////////////////////////////////////////////////////////////////////
385 /// By default, fraction coefficients are assumed to refer to the default
386 /// fit range. This makes the shape of a RooAddPdf
387 /// explicitly dependent on the range of the observables. Calling this function
388 /// allows for a range-independent definition of the fractions, because it
389 /// ties all coefficients to the given
390 /// named range. If the normalisation range is different
391 /// from this reference range, the appropriate fraction coefficients
392 /// are automatically calculated from the reference fractions by
393 /// integrating over the ranges, and comparing these integrals.
394 
395 void RooAddPdf::fixCoefRange(const char* rangeName)
396 {
397  _refCoefRangeName = (TNamed*)RooNameReg::ptr(rangeName) ;
399 }
400 
401 
402 
403 ////////////////////////////////////////////////////////////////////////////////
404 /// Retrieve cache element for the computation of the PDF normalisation.
405 /// \param[in] nset Current normalisation set (integration over these variables yields 1).
406 /// \param[in] iset Integration set. Variables to be integrated over (if integrations are performed).
407 /// \param[in] rangeName Reference range for the integrals.
408 ///
409 /// If a cache element does not exist, create and fill it on the fly. The cache also contains
410 /// - Supplemental normalization terms (in case not all added p.d.f.s have the same observables)
411 /// - Projection integrals to calculate transformed fraction coefficients when a frozen reference frame is provided
412 /// - Projection integrals for similar transformations when a frozen reference range is provided.
413 
414 RooAddPdf::CacheElem* RooAddPdf::getProjCache(const RooArgSet* nset, const RooArgSet* iset, const char* rangeName) const
415 {
416 
417  // Check if cache already exists
418  CacheElem* cache = (CacheElem*) _projCacheMgr.getObj(nset,iset,0,rangeName) ;
419  if (cache) {
420  return cache ;
421  }
422 
423  //Create new cache
424  cache = new CacheElem ;
425 
426  // *** PART 1 : Create supplemental normalization list ***
427 
428  // Retrieve the combined set of dependents of this PDF ;
429  RooArgSet *fullDepList = getObservables(nset) ;
430  if (iset) {
431  fullDepList->remove(*iset,kTRUE,kTRUE) ;
432  }
433 
434  // Fill with dummy unit RRVs for now
435  for (int i = 0; i < _pdfList.getSize(); ++i) {
436  auto pdf = static_cast<const RooAbsPdf *>(_pdfList.at(i));
437  auto coef = static_cast<const RooAbsReal*>(_coefList.at(i));
438 
439  // Start with full list of dependents
440  RooArgSet supNSet(*fullDepList) ;
441 
442  // Remove PDF dependents
443  RooArgSet* pdfDeps = pdf->getObservables(nset) ;
444  if (pdfDeps) {
445  supNSet.remove(*pdfDeps,kTRUE,kTRUE) ;
446  delete pdfDeps ;
447  }
448 
449  // Remove coef dependents
450  RooArgSet* coefDeps = coef ? coef->getObservables(nset) : 0 ;
451  if (coefDeps) {
452  supNSet.remove(*coefDeps,kTRUE,kTRUE) ;
453  delete coefDeps ;
454  }
455 
456  RooAbsReal* snorm ;
457  TString name(GetName()) ;
458  name.Append("_") ;
459  name.Append(pdf->GetName()) ;
460  name.Append("_SupNorm") ;
461  cache->_needSupNorm = kFALSE ;
462  if (supNSet.getSize()>0) {
463  snorm = new RooRealIntegral(name,"Supplemental normalization integral",RooRealConstant::value(1.0),supNSet) ;
464  cxcoutD(Caching) << "RooAddPdf " << GetName() << " making supplemental normalization set " << supNSet << " for pdf component " << pdf->GetName() << endl ;
465  cache->_needSupNorm = kTRUE ;
466  } else {
467  snorm = new RooRealVar(name,"Unit Supplemental normalization integral",1.0) ;
468  }
469  cache->_suppNormList.addOwned(*snorm) ;
470  }
471 
472  delete fullDepList ;
473 
474  if (_verboseEval>1) {
475  cxcoutD(Caching) << "RooAddPdf::syncSuppNormList(" << GetName() << ") synching supplemental normalization list for norm" << (nset?*nset:RooArgSet()) << endl ;
476  if dologD(Caching) {
477  cache->_suppNormList.Print("v") ;
478  }
479  }
480 
481 
482  // *** PART 2 : Create projection coefficients ***
483 
484 // cout << " this = " << this << " (" << GetName() << ")" << endl ;
485 // cout << "projectCoefs = " << (_projectCoefs?"T":"F") << endl ;
486 // cout << "_normRange.Length() = " << _normRange.Length() << endl ;
487 
488  // If no projections required stop here
489  if (!_projectCoefs && !rangeName) {
490  _projCacheMgr.setObj(nset,iset,cache,RooNameReg::ptr(rangeName)) ;
491 // cout << " no projection required" << endl ;
492  return cache ;
493  }
494 
495 
496 // cout << "calculating projection" << endl ;
497 
498  // Reduce iset/nset to actual dependents of this PDF
499  RooArgSet* nset2 = nset ? getObservables(nset) : new RooArgSet() ;
500  cxcoutD(Caching) << "RooAddPdf(" << GetName() << ")::getPC nset = " << (nset?*nset:RooArgSet()) << " nset2 = " << *nset2 << endl ;
501 
502  if (nset2->getSize()==0 && _refCoefNorm.getSize()!=0) {
503  //cout << "WVE: evaluating RooAddPdf without normalization, but have reference normalization for coefficient definition" << endl ;
504 
505  nset2->add(_refCoefNorm) ;
506  if (_refCoefRangeName) {
507  rangeName = RooNameReg::str(_refCoefRangeName) ;
508  }
509  }
510 
511 
512  // Check if requested transformation is not identity
513  if (!nset2->equals(_refCoefNorm) || _refCoefRangeName !=0 || rangeName !=0 || _normRange.Length()>0) {
514 
515  cxcoutD(Caching) << "ALEX: RooAddPdf::syncCoefProjList(" << GetName() << ") projecting coefficients from "
516  << *nset2 << (rangeName?":":"") << (rangeName?rangeName:"")
517  << " to " << ((_refCoefNorm.getSize()>0)?_refCoefNorm:*nset2) << (_refCoefRangeName?":":"") << (_refCoefRangeName?RooNameReg::str(_refCoefRangeName):"") << endl ;
518 
519  // Recalculate projection integrals of PDFs
520  for (auto arg : _pdfList) {
521  auto thePdf = static_cast<const RooAbsPdf*>(arg);
522 
523  // Calculate projection integral
524  RooAbsReal* pdfProj ;
525  if (!nset2->equals(_refCoefNorm)) {
526  pdfProj = thePdf->createIntegral(*nset2,_refCoefNorm,_normRange.Length()>0?_normRange.Data():0) ;
527  pdfProj->setOperMode(operMode()) ;
528  cxcoutD(Caching) << "RooAddPdf(" << GetName() << ")::getPC nset2(" << *nset2 << ")!=_refCoefNorm(" << _refCoefNorm << ") --> pdfProj = " << pdfProj->GetName() << endl ;
529  } else {
530  TString name(GetName()) ;
531  name.Append("_") ;
532  name.Append(thePdf->GetName()) ;
533  name.Append("_ProjectNorm") ;
534  pdfProj = new RooRealVar(name,"Unit Projection normalization integral",1.0) ;
535  cxcoutD(Caching) << "RooAddPdf(" << GetName() << ")::getPC nset2(" << *nset2 << ")==_refCoefNorm(" << _refCoefNorm << ") --> pdfProj = " << pdfProj->GetName() << endl ;
536  }
537 
538  cache->_projList.addOwned(*pdfProj) ;
539  cxcoutD(Caching) << " RooAddPdf::syncCoefProjList(" << GetName() << ") PP = " << pdfProj->GetName() << endl ;
540 
541  // Calculation optional supplemental normalization term
542  RooArgSet supNormSet(_refCoefNorm) ;
543  RooArgSet* deps = thePdf->getParameters(RooArgSet()) ;
544  supNormSet.remove(*deps,kTRUE,kTRUE) ;
545  delete deps ;
546 
547  RooAbsReal* snorm ;
548  TString name(GetName()) ;
549  name.Append("_") ;
550  name.Append(thePdf->GetName()) ;
551  name.Append("_ProjSupNorm") ;
552  if (supNormSet.getSize()>0 && !nset2->equals(_refCoefNorm) ) {
553  snorm = new RooRealIntegral(name,"Projection Supplemental normalization integral",
554  RooRealConstant::value(1.0),supNormSet) ;
555  } else {
556  snorm = new RooRealVar(name,"Unit Projection Supplemental normalization integral",1.0) ;
557  }
558  cxcoutD(Caching) << " RooAddPdf::syncCoefProjList(" << GetName() << ") SN = " << snorm->GetName() << endl ;
559  cache->_suppProjList.addOwned(*snorm) ;
560 
561  // Calculate reference range adjusted projection integral
562  RooAbsReal* rangeProj1 ;
563 
564  // cout << "ALEX >>>> RooAddPdf(" << GetName() << ")::getPC _refCoefRangeName WVE = "
565 // <<(_refCoefRangeName?":":"") << (_refCoefRangeName?RooNameReg::str(_refCoefRangeName):"")
566 // <<" _refCoefRangeName AK = " << (_refCoefRangeName?_refCoefRangeName->GetName():"")
567 // << " && _refCoefNorm" << _refCoefNorm << " with size = _refCoefNorm.getSize() " << _refCoefNorm.getSize() << endl ;
568 
569  // Check if _refCoefRangeName is identical to default range for all observables,
570  // If so, substitute by unit integral
571 
572  // ----------
573  RooArgSet* tmpObs = thePdf->getObservables(_refCoefNorm) ;
574  RooAbsArg* obsArg ;
575  TIterator* iter = tmpObs->createIterator() ;
576  Bool_t allIdent = kTRUE ;
577  while((obsArg=(RooAbsArg*)iter->Next())) {
578  RooRealVar* rvarg = dynamic_cast<RooRealVar*>(obsArg) ;
579  if (rvarg) {
580  if (rvarg->getMin(RooNameReg::str(_refCoefRangeName))!=rvarg->getMin() ||
581  rvarg->getMax(RooNameReg::str(_refCoefRangeName))!=rvarg->getMax()) {
582  allIdent=kFALSE ;
583  }
584  }
585  }
586  delete iter ;
587  delete tmpObs ;
588  // -------------
589 
590  if (_refCoefRangeName && _refCoefNorm.getSize()>0 && !allIdent) {
591 
592 
593  RooArgSet* tmp = thePdf->getObservables(_refCoefNorm) ;
594  rangeProj1 = thePdf->createIntegral(*tmp,*tmp,RooNameReg::str(_refCoefRangeName)) ;
595 
596  //rangeProj1->setOperMode(operMode()) ;
597 
598  delete tmp ;
599  } else {
600 
601  TString theName(GetName()) ;
602  theName.Append("_") ;
603  theName.Append(thePdf->GetName()) ;
604  theName.Append("_RangeNorm1") ;
605  rangeProj1 = new RooRealVar(theName,"Unit range normalization integral",1.0) ;
606 
607  }
608  cxcoutD(Caching) << " RooAddPdf::syncCoefProjList(" << GetName() << ") R1 = " << rangeProj1->GetName() << endl ;
609  cache->_refRangeProjList.addOwned(*rangeProj1) ;
610 
611 
612  // Calculate range adjusted projection integral
613  RooAbsReal* rangeProj2 ;
614  cxcoutD(Caching) << "RooAddPdf::syncCoefProjList(" << GetName() << ") rangename = " << (rangeName?rangeName:"<null>")
615  << " nset = " << (nset?*nset:RooArgSet()) << endl ;
616  if (rangeName && _refCoefNorm.getSize()>0) {
617 
618  rangeProj2 = thePdf->createIntegral(_refCoefNorm,_refCoefNorm,rangeName) ;
619  //rangeProj2->setOperMode(operMode()) ;
620 
621  } else if (_normRange.Length()>0) {
622 
623  RooArgSet* tmp = thePdf->getObservables(_refCoefNorm) ;
624  rangeProj2 = thePdf->createIntegral(*tmp,*tmp,_normRange.Data()) ;
625  delete tmp ;
626 
627  } else {
628 
629  TString theName(GetName()) ;
630  theName.Append("_") ;
631  theName.Append(thePdf->GetName()) ;
632  theName.Append("_RangeNorm2") ;
633  rangeProj2 = new RooRealVar(theName,"Unit range normalization integral",1.0) ;
634 
635  }
636  cxcoutD(Caching) << " RooAddPdf::syncCoefProjList(" << GetName() << ") R2 = " << rangeProj2->GetName() << endl ;
637  cache->_rangeProjList.addOwned(*rangeProj2) ;
638 
639  }
640 
641  }
642 
643  delete nset2 ;
644 
645  _projCacheMgr.setObj(nset,iset,cache,RooNameReg::ptr(rangeName)) ;
646 
647  return cache ;
648 }
649 
650 
651 ////////////////////////////////////////////////////////////////////////////////
652 /// Update the coefficient values in the given cache element: calculate new remainder
653 /// fraction, normalize fractions obtained from extended ML terms to unity, and
654 /// multiply the various range and dimensional corrections needed in the
655 /// current use context.
656 
657 void RooAddPdf::updateCoefficients(CacheElem& cache, const RooArgSet* nset) const
658 {
659  // Since this function updates the cache, it obviously needs write access:
660  auto& myCoefCache = const_cast<std::vector<double>&>(_coefCache);
661  myCoefCache.resize(_haveLastCoef ? _coefList.size() : _pdfList.size(), 0.);
662 
663  // Straight coefficients
664  if (_allExtendable) {
665 
666  // coef[i] = expectedEvents[i] / SUM(expectedEvents)
667  Double_t coefSum(0) ;
668  std::size_t i = 0;
669  for (auto arg : _pdfList) {
670  auto pdf = static_cast<RooAbsPdf*>(arg);
671  myCoefCache[i] = pdf->expectedEvents(_refCoefNorm.getSize()>0?&_refCoefNorm:nset) ;
672  coefSum += myCoefCache[i] ;
673  i++ ;
674  }
675 
676  if (coefSum==0.) {
677  coutW(Eval) << "RooAddPdf::updateCoefCache(" << GetName() << ") WARNING: total number of expected events is 0" << endl ;
678  } else {
679  for (int j=0; j < _pdfList.getSize(); j++) {
680  myCoefCache[j] /= coefSum ;
681  }
682  }
683 
684  } else {
685  if (_haveLastCoef) {
686 
687  // coef[i] = coef[i] / SUM(coef)
688  Double_t coefSum(0) ;
689  std::size_t i=0;
690  for (auto coefArg : _coefList) {
691  auto coef = static_cast<RooAbsReal*>(coefArg);
692  myCoefCache[i] = coef->getVal(nset) ;
693  coefSum += myCoefCache[i++];
694  }
695  if (coefSum==0.) {
696  coutW(Eval) << "RooAddPdf::updateCoefCache(" << GetName() << ") WARNING: sum of coefficients is zero 0" << endl ;
697  } else {
698  const double invCoefSum = 1./coefSum;
699  for (int j=0; j < _coefList.getSize(); j++) {
700  myCoefCache[j] *= invCoefSum;
701  }
702  }
703  } else {
704 
705  // coef[i] = coef[i] ; coef[n] = 1-SUM(coef[0...n-1])
706  Double_t lastCoef(1) ;
707  std::size_t i=0;
708  for (auto coefArg : _coefList) {
709  auto coef = static_cast<RooAbsReal*>(coefArg);
710  myCoefCache[i] = coef->getVal(nset) ;
711  lastCoef -= myCoefCache[i++];
712  }
713  myCoefCache[_coefList.getSize()] = lastCoef ;
714 
715  // Treat coefficient degeneration
716  const float coefDegen = lastCoef < 0. ? -lastCoef : (lastCoef > 1. ? lastCoef - 1. : 0.);
717  if (coefDegen > 1.E-5) {
718  myCoefCache[_coefList.getSize()] = RooNaNPacker::packFloatIntoNaN(100.f*coefDegen);
719 
720  std::stringstream msg;
721  if (_coefErrCount-->0) {
722  msg << "RooAddPdf::updateCoefCache(" << GetName()
723  << " WARNING: sum of PDF coefficients not in range [0-1], value="
724  << 1-lastCoef ;
725  if (_coefErrCount==0) {
726  msg << " (no more will be printed)" ;
727  }
728  coutW(Eval) << msg.str() << std::endl;
729  }
730  }
731  }
732  }
733 
734 
735  // Stop here if not projection is required or needed
736  if ((!_projectCoefs && _normRange.Length()==0) || cache._projList.getSize()==0) {
737  return ;
738  }
739 
740  // Adjust coefficients for given projection
741  Double_t coefSum(0) ;
742  {
744 
745  for (int i = 0; i < _pdfList.getSize(); i++) {
746 
747  RooAbsReal* pp = ((RooAbsReal*)cache._projList.at(i)) ;
748  RooAbsReal* sn = ((RooAbsReal*)cache._suppProjList.at(i)) ;
749  RooAbsReal* r1 = ((RooAbsReal*)cache._refRangeProjList.at(i)) ;
750  RooAbsReal* r2 = ((RooAbsReal*)cache._rangeProjList.at(i)) ;
751 
752  Double_t proj = pp->getVal()/sn->getVal()*(r2->getVal()/r1->getVal()) ;
753 
754  myCoefCache[i] *= proj ;
755  coefSum += myCoefCache[i] ;
756  }
757  }
758 
759 
761  for (int i=0; i < _pdfList.getSize(); ++i) {
762  ccoutD(Caching) << " ALEX: POST-SYNC coef[" << i << "] = " << myCoefCache[i]
763  << " ( _coefCache[i]/coefSum = " << myCoefCache[i]*coefSum << "/" << coefSum << " ) "<< endl ;
764  }
765  }
766 
767  if (coefSum==0.) {
768  coutE(Eval) << "RooAddPdf::updateCoefCache(" << GetName() << ") sum of coefficients is zero." << endl ;
769  }
770 
771  for (int i=0; i < _pdfList.getSize(); i++) {
772  myCoefCache[i] /= coefSum ;
773  }
774 
775 }
776 
777 ////////////////////////////////////////////////////////////////////////////////
778 /// Look up projection cache and per-PDF norm sets. If a PDF doesn't have a special
779 /// norm set, use the `defaultNorm`. If `defaultNorm == nullptr`, use the member
780 /// _normSet.
781 std::pair<const RooArgSet*, RooAddPdf::CacheElem*> RooAddPdf::getNormAndCache(const RooArgSet* defaultNorm) const {
782  const RooArgSet* nset = defaultNorm ? defaultNorm : _normSet;
783 
784  if (nset==0 || nset->getSize()==0) {
785  if (_refCoefNorm.getSize()!=0) {
786  nset = &_refCoefNorm ;
787  }
788  }
789 
790  CacheElem* cache = getProjCache(nset) ;
791  updateCoefficients(*cache,nset) ;
792 
793  return {nset, cache};
794 }
795 
796 
797 ////////////////////////////////////////////////////////////////////////////////
798 /// Calculate and return the current value
799 
801 {
802  auto normAndCache = getNormAndCache();
803  const RooArgSet* nset = normAndCache.first;
804  CacheElem* cache = normAndCache.second;
805  // check that nset is not a null pointer
806  // in this case interpretation of coefficient is arbitrary
807  if (!nset)
808  oocoutW(this, Eval) << "Evaluating RooAddPdf without a defined normalization set. This can lead to ambiguos coefficients definition and incorrect results."
809  << " Use RooAddPdf::fixCoefNormalization(nset) to provide a normalization set for defining uniquely RooAddPdf coefficients!" << std::endl;
810 
811  // Do running sum of coef/pdf pairs, calculate lastCoef.
812  Double_t value(0);
813 
814  for (unsigned int i=0; i < _pdfList.size(); ++i) {
815  const auto& pdf = static_cast<RooAbsPdf&>(_pdfList[i]);
816  double snormVal = 1.;
817  if (cache->_needSupNorm) {
818  snormVal = ((RooAbsReal*)cache->_suppNormList.at(i))->getVal();
819  }
820 
821  Double_t pdfVal = pdf.getVal(nset);
822  if (pdf.isSelectedComp()) {
823  value += pdfVal*_coefCache[i]/snormVal;
824  }
825  }
826 
827  return value;
828 }
829 
830 
831 ////////////////////////////////////////////////////////////////////////////////
832 /// Compute addition of PDFs in batches.
834  auto normAndCache = getNormAndCache(normSet);
835  const RooArgSet* nset = normAndCache.first;
836  CacheElem* cache = normAndCache.second;
837 
838 
840 
841  for (unsigned int pdfNo = 0; pdfNo < _pdfList.size(); ++pdfNo) {
842  const auto& pdf = static_cast<RooAbsPdf&>(_pdfList[pdfNo]);
843  auto pdfOutputs = pdf.getValues(evalData, nset);
844  if (output.empty()) {
845  output = evalData.makeBatch(this, pdfOutputs.size());
846  for (double& val : output) { //CHECK_VECTORISE
847  val = 0.;
848  }
849  }
850  assert(output.size() == pdfOutputs.size());
851 
852  const double coef = _coefCache[pdfNo] / (cache->_needSupNorm ?
853  static_cast<RooAbsReal*>(cache->_suppNormList.at(pdfNo))->getVal() :
854  1.);
855 
856  if (pdf.isSelectedComp()) {
857  for (std::size_t i = 0; i < output.size(); ++i) { //CHECK_VECTORISE
858  output[i] += pdfOutputs[i] * coef;
859  }
860  }
861  }
862 
863  return output;
864 }
865 
866 
867 ////////////////////////////////////////////////////////////////////////////////
868 /// Reset error counter to given value, limiting the number
869 /// of future error messages for this pdf to 'resetValue'
870 
872 {
873  RooAbsPdf::resetErrorCounters(resetValue) ;
874  _coefErrCount = resetValue ;
875 }
876 
877 
878 
879 ////////////////////////////////////////////////////////////////////////////////
880 /// Check if PDF is valid for given normalization set.
881 /// Coeffient and PDF must be non-overlapping, but pdf-coefficient
882 /// pairs may overlap each other
883 
885 {
886  Bool_t ret(kFALSE) ;
887 
888  for (int i = 0; i < _pdfList.getSize(); ++i) {
889  auto pdf = static_cast<const RooAbsPdf *>(_pdfList.at(i));
890  auto coef = static_cast<const RooAbsReal*>(_coefList.at(i));
891  if (pdf->observableOverlaps(nset,*coef)) {
892  coutE(InputArguments) << "RooAddPdf::checkObservables(" << GetName() << "): ERROR: coefficient " << coef->GetName()
893  << " and PDF " << pdf->GetName() << " have one or more dependents in common" << endl ;
894  ret = kTRUE ;
895  }
896  }
897 
898  return ret ;
899 }
900 
901 
902 ////////////////////////////////////////////////////////////////////////////////
903 /// Determine which part (if any) of given integral can be performed analytically.
904 /// If any analytical integration is possible, return integration scenario code
905 ///
906 /// RooAddPdf queries each component PDF for its analytical integration capability of the requested
907 /// set ('allVars'). It finds the largest common set of variables that can be integrated
908 /// by all components. If such a set exists, it reconfirms that each component is capable of
909 /// analytically integrating the common set, and combines the components individual integration
910 /// codes into a single integration code valid for RooAddPdf.
911 
913  const RooArgSet* normSet, const char* rangeName) const
914 {
915 
916  RooArgSet* allDepVars = getObservables(allVars) ;
917  RooArgSet allAnalVars(*allDepVars) ;
918  delete allDepVars ;
919 
920  Int_t n(0) ;
921 
922  // First iteration, determine what each component can integrate analytically
923  for (const auto pdfArg : _pdfList) {
924  auto pdf = static_cast<const RooAbsPdf *>(pdfArg);
925  RooArgSet subAnalVars ;
926  pdf->getAnalyticalIntegralWN(allVars,subAnalVars,normSet,rangeName) ;
927 
928  // Observables that cannot be integrated analytically by this component are dropped from the common list
929  for (const auto arg : allVars) {
930  if (!subAnalVars.find(arg->GetName()) && pdf->dependsOn(*arg)) {
931  allAnalVars.remove(*arg,kTRUE,kTRUE) ;
932  }
933  }
934  n++ ;
935  }
936 
937  // If no observables can be integrated analytically, return code 0 here
938  if (allAnalVars.getSize()==0) {
939  return 0 ;
940  }
941 
942 
943  // Now retrieve codes for integration over common set of analytically integrable observables for each component
944  n=0 ;
945  std::vector<Int_t> subCode(_pdfList.getSize());
946  Bool_t allOK(kTRUE) ;
947  for (const auto arg : _pdfList) {
948  auto pdf = static_cast<const RooAbsPdf *>(arg);
949  RooArgSet subAnalVars ;
950  RooArgSet* allAnalVars2 = pdf->getObservables(allAnalVars) ;
951  subCode[n] = pdf->getAnalyticalIntegralWN(*allAnalVars2,subAnalVars,normSet,rangeName) ;
952  if (subCode[n]==0 && allAnalVars2->getSize()>0) {
953  coutE(InputArguments) << "RooAddPdf::getAnalyticalIntegral(" << GetName() << ") WARNING: component PDF " << pdf->GetName()
954  << " advertises inconsistent set of integrals (e.g. (X,Y) but not X or Y individually."
955  << " Distributed analytical integration disabled. Please fix PDF" << endl ;
956  allOK = kFALSE ;
957  }
958  delete allAnalVars2 ;
959  n++ ;
960  }
961  if (!allOK) {
962  return 0 ;
963  }
964 
965  // Mare all analytically integrated observables as such
966  analVars.add(allAnalVars) ;
967 
968  // Store set of variables analytically integrated
969  RooArgSet* intSet = new RooArgSet(allAnalVars) ;
970  Int_t masterCode = _codeReg.store(subCode,intSet)+1 ;
971 
972  return masterCode ;
973 }
974 
975 
976 
977 ////////////////////////////////////////////////////////////////////////////////
978 /// Return analytical integral defined by given scenario code
979 
980 Double_t RooAddPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
981 {
982  // WVE needs adaptation to handle new rangeName feature
983  if (code==0) {
984  return getVal(normSet) ;
985  }
986 
987  // Retrieve analytical integration subCodes and set of observabels integrated over
988  RooArgSet* intSet ;
989  const std::vector<Int_t>& subCode = _codeReg.retrieve(code-1,intSet) ;
990  if (subCode.empty()) {
991  coutE(InputArguments) << "RooAddPdf::analyticalIntegral(" << GetName() << "): ERROR unrecognized integration code, " << code << endl ;
992  assert(0) ;
993  }
994 
995  cxcoutD(Caching) << "RooAddPdf::aiWN(" << GetName() << ") calling getProjCache with nset = " << (normSet?*normSet:RooArgSet()) << endl ;
996 
997  if ((normSet==0 || normSet->getSize()==0) && _refCoefNorm.getSize()>0) {
998 // cout << "WVE integration of RooAddPdf without normalization, but have reference set, using ref set for normalization" << endl ;
999  normSet = &_refCoefNorm ;
1000  }
1001 
1002  CacheElem* cache = getProjCache(normSet,intSet,0) ; // WVE rangename here?
1003  updateCoefficients(*cache,normSet) ;
1004 
1005  // Calculate the current value of this object
1006  Double_t value(0) ;
1007 
1008  // Do running sum of coef/pdf pairs, calculate lastCoef.
1009  Double_t snormVal ;
1010 
1011  //cout << "ROP::aIWN updateCoefCache with rangeName = " << (rangeName?rangeName:"<null>") << endl ;
1012  RooArgList* snormSet = (cache->_suppNormList.getSize()>0) ? &cache->_suppNormList : 0 ;
1013  for (int i = 0; i < _pdfList.getSize(); ++i ) {
1014  auto pdf = static_cast<const RooAbsPdf*>(_pdfList.at(i));
1015 
1016  if (_coefCache[i]) {
1017  snormVal = snormSet ? ((RooAbsReal*) cache->_suppNormList.at(i))->getVal() : 1.0 ;
1018 
1019  // WVE swap this?
1020  Double_t val = pdf->analyticalIntegralWN(subCode[i],normSet,rangeName) ;
1021  if (pdf->isSelectedComp()) {
1022  value += val*_coefCache[i]/snormVal ;
1023  }
1024  }
1025  }
1026 
1027  return value ;
1028 }
1029 
1030 
1031 
1032 ////////////////////////////////////////////////////////////////////////////////
1033 /// Return the number of expected events, which is either the sum of all coefficients
1034 /// or the sum of the components extended terms, multiplied with the fraction that
1035 /// is in the current range w.r.t the reference range
1036 
1038 {
1039  Double_t expectedTotal(0.0);
1040 
1041  cxcoutD(Caching) << "RooAddPdf::expectedEvents(" << GetName() << ") calling getProjCache with nset = " << (nset?*nset:RooArgSet()) << endl ;
1042  CacheElem* cache = getProjCache(nset) ;
1043  updateCoefficients(*cache,nset) ;
1044 
1045  if (cache->_rangeProjList.getSize()>0) {
1046 
1047  RooFIter iter1 = cache->_refRangeProjList.fwdIterator() ;
1048  RooFIter iter2 = cache->_rangeProjList.fwdIterator() ;
1049  RooFIter iter3 = _pdfList.fwdIterator() ;
1050 
1051  if (_allExtendable) {
1052 
1053  RooAbsPdf* pdf ;
1054  while ((pdf=(RooAbsPdf*)iter3.next())) {
1055  RooAbsReal* r1 = (RooAbsReal*)iter1.next() ;
1056  RooAbsReal* r2 = (RooAbsReal*)iter2.next() ;
1057  expectedTotal += (r2->getVal()/r1->getVal()) * pdf->expectedEvents(nset) ;
1058  }
1059 
1060  } else {
1061 
1062  RooFIter citer = _coefList.fwdIterator() ;
1063  RooAbsReal* coef ;
1064  while((coef=(RooAbsReal*)citer.next())) {
1065  Double_t ncomp = coef->getVal(nset) ;
1066  RooAbsReal* r1 = (RooAbsReal*)iter1.next() ;
1067  RooAbsReal* r2 = (RooAbsReal*)iter2.next() ;
1068  expectedTotal += (r2->getVal()/r1->getVal()) * ncomp ;
1069  }
1070 
1071  }
1072 
1073 
1074 
1075  } else {
1076 
1077  if (_allExtendable) {
1078 
1079  RooFIter iter = _pdfList.fwdIterator() ;
1080  RooAbsPdf* pdf ;
1081  while((pdf=(RooAbsPdf*)iter.next())) {
1082  expectedTotal += pdf->expectedEvents(nset) ;
1083  }
1084 
1085  } else {
1086 
1087  RooFIter citer = _coefList.fwdIterator() ;
1088  RooAbsReal* coef ;
1089  while((coef=(RooAbsReal*)citer.next())) {
1090  Double_t ncomp = coef->getVal(nset) ;
1091  expectedTotal += ncomp ;
1092  }
1093 
1094  }
1095 
1096  }
1097  return expectedTotal ;
1098 }
1099 
1100 
1101 
1102 ////////////////////////////////////////////////////////////////////////////////
1103 /// Interface function used by test statistics to freeze choice of observables
1104 /// for interpretation of fraction coefficients
1105 
1107 {
1108 
1109  if (!force && _refCoefNorm.getSize()!=0) {
1110  return ;
1111  }
1112 
1113  if (!depSet) {
1115  return ;
1116  }
1117 
1118  RooArgSet* myDepSet = getObservables(depSet) ;
1119  fixCoefNormalization(*myDepSet) ;
1120  delete myDepSet ;
1121 }
1122 
1123 
1124 
1125 ////////////////////////////////////////////////////////////////////////////////
1126 /// Interface function used by test statistics to freeze choice of range
1127 /// for interpretation of fraction coefficients
1128 
1129 void RooAddPdf::selectNormalizationRange(const char* rangeName, Bool_t force)
1130 {
1131  if (!force && _refCoefRangeName) {
1132  return ;
1133  }
1134 
1135  fixCoefRange(rangeName) ;
1136 }
1137 
1138 
1139 
1140 ////////////////////////////////////////////////////////////////////////////////
1141 /// Return specialized context to efficiently generate toy events from RooAddPdfs
1142 /// return RooAbsPdf::genContext(vars,prototype,auxProto,verbose) ; // WVE DEBUG
1143 
1145  const RooArgSet* auxProto, Bool_t verbose) const
1146 {
1147  return new RooAddGenContext(*this,vars,prototype,auxProto,verbose) ;
1148 }
1149 
1150 
1151 
1152 ////////////////////////////////////////////////////////////////////////////////
1153 /// List all RooAbsArg derived contents in this cache element
1154 
1156 {
1157  RooArgList allNodes;
1158  allNodes.add(_projList) ;
1159  allNodes.add(_suppProjList) ;
1160  allNodes.add(_refRangeProjList) ;
1161  allNodes.add(_rangeProjList) ;
1162 
1163  return allNodes ;
1164 }
1165 
1166 
1167 
1168 ////////////////////////////////////////////////////////////////////////////////
1169 /// Loop over components for plot sampling hints and merge them if there are multiple
1170 
1171 std::list<Double_t>* RooAddPdf::plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
1172 {
1173  list<Double_t>* sumHint = 0 ;
1174  Bool_t needClean(kFALSE) ;
1175 
1176  // Loop over components pdf
1177  for (const auto arg : _pdfList) {
1178  auto pdf = static_cast<const RooAbsPdf*>(arg);
1179 
1180  list<Double_t>* pdfHint = pdf->plotSamplingHint(obs,xlo,xhi) ;
1181 
1182  // Process hint
1183  if (pdfHint) {
1184  if (!sumHint) {
1185 
1186  // If this is the first hint, then just save it
1187  sumHint = pdfHint ;
1188 
1189  } else {
1190 
1191  list<Double_t>* newSumHint = new list<Double_t>(sumHint->size()+pdfHint->size()) ;
1192 
1193  // Merge hints into temporary array
1194  merge(pdfHint->begin(),pdfHint->end(),sumHint->begin(),sumHint->end(),newSumHint->begin()) ;
1195 
1196  // Copy merged array without duplicates to new sumHintArrau
1197  delete sumHint ;
1198  sumHint = newSumHint ;
1199  needClean = kTRUE ;
1200 
1201  }
1202  }
1203  }
1204  if (needClean) {
1205  list<Double_t>::iterator new_end = unique(sumHint->begin(),sumHint->end()) ;
1206  sumHint->erase(new_end,sumHint->end()) ;
1207  }
1208 
1209  return sumHint ;
1210 }
1211 
1212 
1213 ////////////////////////////////////////////////////////////////////////////////
1214 /// Loop over components for plot sampling hints and merge them if there are multiple
1215 
1216 std::list<Double_t>* RooAddPdf::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
1217 {
1218  list<Double_t>* sumBinB = 0 ;
1219  Bool_t needClean(kFALSE) ;
1220 
1221  // Loop over components pdf
1222  for (auto arg : _pdfList) {
1223  auto pdf = static_cast<const RooAbsPdf *>(arg);
1224  list<Double_t>* pdfBinB = pdf->binBoundaries(obs,xlo,xhi) ;
1225 
1226  // Process hint
1227  if (pdfBinB) {
1228  if (!sumBinB) {
1229 
1230  // If this is the first hint, then just save it
1231  sumBinB = pdfBinB ;
1232 
1233  } else {
1234 
1235  list<Double_t>* newSumBinB = new list<Double_t>(sumBinB->size()+pdfBinB->size()) ;
1236 
1237  // Merge hints into temporary array
1238  merge(pdfBinB->begin(),pdfBinB->end(),sumBinB->begin(),sumBinB->end(),newSumBinB->begin()) ;
1239 
1240  // Copy merged array without duplicates to new sumBinBArrau
1241  delete sumBinB ;
1242  delete pdfBinB ;
1243  sumBinB = newSumBinB ;
1244  needClean = kTRUE ;
1245  }
1246  }
1247  }
1248 
1249  // Remove consecutive duplicates
1250  if (needClean) {
1251  list<Double_t>::iterator new_end = unique(sumBinB->begin(),sumBinB->end()) ;
1252  sumBinB->erase(new_end,sumBinB->end()) ;
1253  }
1254 
1255  return sumBinB ;
1256 }
1257 
1258 
1259 ////////////////////////////////////////////////////////////////////////////////
1260 /// If all components that depend on obs are binned, so is their sum.
1262 {
1263  for (const auto arg : _pdfList) {
1264  auto pdf = static_cast<const RooAbsPdf*>(arg);
1265  if (pdf->dependsOn(obs) && !pdf->isBinnedDistribution(obs)) {
1266  return kFALSE ;
1267  }
1268  }
1269 
1270  return kTRUE ;
1271 }
1272 
1273 
1274 ////////////////////////////////////////////////////////////////////////////////
1275 /// Label OK'ed components of a RooAddPdf with cache-and-track
1276 
1278 {
1279  RooFIter aiter = pdfList().fwdIterator() ;
1280  RooAbsArg* aarg ;
1281  while ((aarg=aiter.next())) {
1282  if (aarg->canNodeBeCached()==Always) {
1283  trackNodes.add(*aarg) ;
1284  //cout << "tracking node RooAddPdf component " << aarg->IsA()->GetName() << "::" << aarg->GetName() << endl ;
1285  }
1286  }
1287 }
1288 
1289 
1290 
1291 ////////////////////////////////////////////////////////////////////////////////
1292 /// Customized printing of arguments of a RooAddPdf to more intuitively reflect the contents of the
1293 /// product operator construction
1294 
1295 void RooAddPdf::printMetaArgs(ostream& os) const
1296 {
1297  Bool_t first(kTRUE) ;
1298 
1299  if (_coefList.getSize() != 0) {
1300  for (int i = 0; i < _pdfList.getSize(); ++i ) {
1301  const RooAbsArg * coef = _coefList.at(i);
1302  const RooAbsArg * pdf = _pdfList.at(i);
1303  if (!first) {
1304  os << " + " ;
1305  } else {
1306  first = kFALSE ;
1307  }
1308 
1309  if (i < _coefList.getSize()) {
1310  os << coef->GetName() << " * " << pdf->GetName();
1311  } else {
1312  os << "[%] * " << pdf->GetName();
1313  }
1314  }
1315  } else {
1316 
1317  for (const auto pdf : _pdfList) {
1318  if (!first) {
1319  os << " + " ;
1320  } else {
1321  first = kFALSE ;
1322  }
1323  os << pdf->GetName() ;
1324  }
1325  }
1326 
1327  os << " " ;
1328 }
RooAddPdf::checkObservables
virtual Bool_t checkObservables(const RooArgSet *nset) const
Check if PDF is valid for given normalization set.
Definition: RooAddPdf.cxx:884
RooCacheManager::reset
void reset()
Clear the cache.
Definition: RooCacheManager.h:192
RooAbsPdf::_normSet
RooArgSet * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:322
RooCacheManager::setObj
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:49
RooBatchCompute.h
n
const Int_t n
Definition: legend1.C:16
RooAbsPdf::CacheElem
friend class CacheElem
Definition: RooAbsPdf.h:334
RooAbsPdf::_normRange
TString _normRange
MC generator configuration specific for this object.
Definition: RooAbsPdf.h:355
first
Definition: first.py:1
RooAbsArg::operMode
OperMode operMode() const
Query the operation mode of this node.
Definition: RooAbsArg.h:482
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooBatchCompute::RunContext::makeBatch
RooSpan< double > makeBatch(const RooAbsReal *owner, std::size_t size)
Create a writable batch.
Definition: RunContext.cxx:58
RooAddPdf
RooAddPdf is an efficient implementation of a sum of PDFs of the form.
Definition: RooAddPdf.h:32
ccoutD
#define ccoutD(a)
Definition: RooMsgService.h:37
RooAddPdf::updateCoefficients
void updateCoefficients(CacheElem &cache, const RooArgSet *nset) const
Update the coefficient values in the given cache element: calculate new remainder fraction,...
Definition: RooAddPdf.cxx:657
RooAbsCollection::first
RooAbsArg * first() const
Definition: RooAbsCollection.h:177
RooMsgService::isActive
Bool_t isActive(const RooAbsArg *self, RooFit::MsgTopic facility, RooFit::MsgLevel level)
Check if logging is active for given object/topic/RooFit::MsgLevel combination.
Definition: RooMsgService.cxx:398
RooAbsReal::plotSamplingHint
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Interface for returning an optional hint for initial sampling points when constructing a curve projec...
Definition: RooAbsReal.cxx:3840
RooNameReg::ptr
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooAddPdf::RooAddGenContext
friend class RooAddGenContext
Definition: RooAddPdf.h:128
RooAbsRealLValue::getMax
virtual Double_t getMax(const char *name=0) const
Get maximum of currently defined range.
Definition: RooAbsRealLValue.h:89
RooAddPdf::isBinnedDistribution
Bool_t isBinnedDistribution(const RooArgSet &obs) const
If all components that depend on obs are binned, so is their sum.
Definition: RooAddPdf.cxx:1261
RooAddPdf::CacheElem::_projList
RooArgList _projList
Definition: RooAddPdf.h:115
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:68
RooAbsArg::Always
@ Always
Definition: RooAbsArg.h:387
RooCacheManager::getObj
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:44
RooNaNPacker::packFloatIntoNaN
static double packFloatIntoNaN(float payload)
Pack float into mantissa of a NaN.
Definition: RooNaNPacker.h:109
RooAddPdf::plotSamplingHint
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Loop over components for plot sampling hints and merge them if there are multiple.
Definition: RooAddPdf.cxx:1171
RooAddPdf::fixCoefRange
void fixCoefRange(const char *rangeName)
By default, fraction coefficients are assumed to refer to the default fit range.
Definition: RooAddPdf.cxx:395
RooFit::DEBUG
@ DEBUG
Definition: RooGlobalFunc.h:65
TString::Data
const char * Data() const
Definition: TString.h:369
RooAddPdf::_haveLastCoef
Bool_t _haveLastCoef
List of supplemental normalization factors.
Definition: RooAddPdf.h:143
RooAbsReal::createIntegral
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:573
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
RooAddPdf::selectNormalizationRange
virtual void selectNormalizationRange(const char *rangeName=0, Bool_t force=kFALSE)
Interface function used by test statistics to freeze choice of range for interpretation of fraction c...
Definition: RooAddPdf.cxx:1129
RooAddPdf::_projectCoefs
Bool_t _projectCoefs
Definition: RooAddPdf.h:104
RooAbsCollection::fwdIterator
RooFIter fwdIterator() const
One-time forward iterator.
Definition: RooAbsCollection.h:134
coutE
#define coutE(a)
Definition: RooMsgService.h:33
RooAbsReal::binBoundaries
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Retrieve bin boundaries if this distribution is binned in obs.
Definition: RooAbsReal.cxx:3829
RooAddPdf::evaluate
Double_t evaluate() const
Calculate and return the current value.
Definition: RooAddPdf.cxx:800
coutW
#define coutW(a)
Definition: RooMsgService.h:32
output
static void output(int code)
Definition: gifencode.c:226
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooAddPdf::analyticalIntegralWN
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return analytical integral defined by given scenario code.
Definition: RooAddPdf.cxx:980
RooAddPdf::_projCacheMgr
RooObjCacheManager _projCacheMgr
Definition: RooAddPdf.h:123
RooRealConstant::value
static RooConstVar & value(Double_t value)
Return a constant value object with given value.
Definition: RooRealConstant.cxx:52
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooAddPdf::_refCoefRangeName
TNamed * _refCoefRangeName
Definition: RooAddPdf.h:102
RooAbsCollection::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Definition: RooAbsCollection.cxx:584
RooAddPdf::selectNormalization
virtual void selectNormalization(const RooArgSet *depSet=0, Bool_t force=kFALSE)
Interface function used by test statistics to freeze choice of observables for interpretation of frac...
Definition: RooAddPdf.cxx:1106
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
RooAbsArg::setOperMode
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1780
TString::Length
Ssiz_t Length() const
Definition: TString.h:410
RooAddPdf::pdfList
const RooArgList & pdfList() const
Definition: RooAddPdf.h:70
RooAddPdf.h
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooArgList::at
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:70
RooAddPdf::CacheElem::_suppProjList
RooArgList _suppProjList
Definition: RooAddPdf.h:116
RooAddPdf::CacheElem::_suppNormList
RooArgList _suppNormList
Definition: RooAddPdf.h:110
RooAddPdf::CacheElem
Transient cache with transformed values of coefficients.
Definition: RooAddPdf.h:108
TString
Basic string class.
Definition: TString.h:136
RooAddPdf::CacheElem::_refRangeProjList
RooArgList _refRangeProjList
Definition: RooAddPdf.h:117
RooAbsCollection::GetName
const char * GetName() const
Returns name of object.
Definition: RooAbsCollection.h:227
RooDataSet.h
bool
TIterator
Iterator abstract base class.
Definition: TIterator.h:30
RooAbsRealLValue::getMin
virtual Double_t getMin(const char *name=0) const
Get miniminum of currently defined range.
Definition: RooAbsRealLValue.h:86
RooAddPdf::resetErrorCounters
virtual void resetErrorCounters(Int_t resetValue=10)
Reset error counter to given value, limiting the number of future error messages for this pdf to 'res...
Definition: RooAddPdf.cxx:871
RooRecursiveFraction
Class RooRecursiveFraction is a RooAbsReal implementation that calculates the plain fraction of sum o...
Definition: RooRecursiveFraction.h:25
RooAbsReal::GlobalSelectComponentRAII
Definition: RooAbsReal.h:546
RooAddPdf::binBoundaries
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Loop over components for plot sampling hints and merge them if there are multiple.
Definition: RooAddPdf.cxx:1216
RooAddPdf::getNormAndCache
std::pair< const RooArgSet *, CacheElem * > getNormAndCache(const RooArgSet *defaultNorm=nullptr) const
Coefficient error counter.
Definition: RooAddPdf.cxx:781
RooAddPdf::_allExtendable
Bool_t _allExtendable
Definition: RooAddPdf.h:144
RooAbsPdf::_errorCount
Int_t _errorCount
Definition: RooAbsPdf.h:344
RooAICRegistry::store
Int_t store(const std::vector< Int_t > &codeList, RooArgSet *set1=0, RooArgSet *set2=0, RooArgSet *set3=0, RooArgSet *set4=0)
Store given arrays of integer codes, and up to four RooArgSets in the registry (each setX pointer may...
Definition: RooAICRegistry.cxx:100
RooAddPdf::expectedEvents
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return the number of expected events, which is either the sum of all coefficients or the sum of the c...
Definition: RooAddPdf.cxx:1037
RooAddPdf::_refCoefNorm
RooSetProxy _refCoefNorm
Definition: RooAddPdf.h:101
TRACE_DESTROY
#define TRACE_DESTROY
Definition: RooTrace.h:24
RooAddPdf::evaluateSpan
RooSpan< double > evaluateSpan(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Compute addition of PDFs in batches.
Definition: RooAddPdf.cxx:833
RooNaNPacker.h
RooAddGenContext.h
TNamed
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
RooFIter
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
Definition: RooLinkedListIter.h:40
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
RooAddPdf::_coefList
RooListProxy _coefList
Definition: RooAddPdf.h:140
RooArgSet::add
Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add element to non-owning set.
Definition: RooArgSet.cxx:261
RooAbsCollection::createIterator
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:119
RooFIter::next
RooAbsArg * next()
Return next element or nullptr if at end.
Definition: RooLinkedListIter.h:49
RooRecursiveFraction.h
RooRealProxy.h
RooAbsArg::RooArgSet
friend class RooArgSet
Definition: RooAbsArg.h:586
RooAbsGenContext
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
Definition: RooAbsGenContext.h:26
TRACE_CREATE
#define TRACE_CREATE
Definition: RooTrace.h:23
RooAbsCollection::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:437
RooAbsCollection::size
Storage_t::size_type size() const
Definition: RooAbsCollection.h:155
oocoutW
#define oocoutW(o, a)
Definition: RooMsgService.h:47
RooRealVar.h
RooAbsCollection::addOwned
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:389
TIterator::Next
virtual TObject * Next()=0
RooAddPdf::getAnalyticalIntegralWN
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
Determine which part (if any) of given integral can be performed analytically.
Definition: RooAddPdf.cxx:912
RooGlobalFunc.h
RooSetProxy::removeAll
virtual void removeAll() override
Remove all argument inset using remove(const RooAbsArg&).
Definition: RooSetProxy.cxx:229
RooAddPdf::CacheElem::_rangeProjList
RooArgList _rangeProjList
Definition: RooAddPdf.h:118
RooAddPdf::CacheElem::containedArgs
virtual RooArgList containedArgs(Action)
List all RooAbsArg derived contents in this cache element.
Definition: RooAddPdf.cxx:1155
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:294
RooAbsArg::addOwnedComponents
Bool_t addOwnedComponents(const RooArgSet &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2168
RooAbsCollection::equals
Bool_t equals(const RooAbsCollection &otherColl) const
Check if this and other collection have identically-named contents.
Definition: RooAbsCollection.cxx:763
RooAddPdf::printMetaArgs
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooAddPdf to more intuitively reflect the contents of the produ...
Definition: RooAddPdf.cxx:1295
RooAddPdf::fixCoefNormalization
void fixCoefNormalization(const RooArgSet &refCoefNorm)
By default the interpretation of the fraction coefficients is performed in the contextual choice of o...
Definition: RooAddPdf.cxx:368
RooAbsCacheElement::Action
Action
Definition: RooAbsCacheElement.h:39
RooAddPdf::~RooAddPdf
virtual ~RooAddPdf()
Destructor.
Definition: RooAddPdf.cxx:351
RooAddPdf::genContext
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Return specialized context to efficiently generate toy events from RooAddPdfs return RooAbsPdf::genCo...
Definition: RooAddPdf.cxx:1144
RooAddPdf::_pdfList
RooListProxy _pdfList
Registry of component analytical integration codes.
Definition: RooAddPdf.h:139
RooAbsPdf::getValues
RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Compute batch of values for given input data, and normalise by integrating over the observables in no...
Definition: RooAbsPdf.cxx:351
RooMsgService::_debugCount
static Int_t _debugCount
Definition: RooMsgService.h:172
RooNameReg::str
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.cxx:103
RooAbsPdf::_verboseEval
static Int_t _verboseEval
Definition: RooAbsPdf.h:315
RooAddPdf::RooAddPdf
RooAddPdf()
Default constructor used for persistence.
Definition: RooAddPdf.cxx:90
name
char name[80]
Definition: TGX11.cxx:110
RooAbsArg::dependsOn
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:764
RooAddPdf::setCacheAndTrackHints
virtual void setCacheAndTrackHints(RooArgSet &)
Label OK'ed components of a RooAddPdf with cache-and-track.
Definition: RooAddPdf.cxx:1277
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
RooAddPdf::_coefCache
std::vector< double > _coefCache
Definition: RooAddPdf.h:105
RooFit::Caching
@ Caching
Definition: RooGlobalFunc.h:68
RooDataSet
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
RooAddPdf::CacheElem::_needSupNorm
Bool_t _needSupNorm
Definition: RooAddPdf.h:113
dologD
#define dologD(a)
Definition: RooMsgService.h:65
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooListProxy::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Reimplementation of standard RooArgList::add()
Definition: RooListProxy.cxx:104
RooAbsArg::canNodeBeCached
virtual CacheMode canNodeBeCached() const
Definition: RooAbsArg.h:410
RooRealIntegral.h
RooAbsCollection::Print
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
Definition: RooAbsCollection.h:200
RooAbsPdf
Definition: RooAbsPdf.h:43
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooAddPdf::getProjCache
CacheElem * getProjCache(const RooArgSet *nset, const RooArgSet *iset=0, const char *rangeName=0) const
Retrieve cache element for the computation of the PDF normalisation.
Definition: RooAddPdf.cxx:414
RooAICRegistry::retrieve
const std::vector< Int_t > & retrieve(Int_t masterCode) const
Retrieve the array of integer codes associated with the given master code.
Definition: RooAICRegistry.cxx:153
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooFit::Eval
@ Eval
Definition: RooGlobalFunc.h:68
RooAbsRealLValue
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Definition: RooAbsRealLValue.h:31
for
for(Int_t i=0;i< n;i++)
Definition: legend1.C:18
RooAddPdf::_coefErrCount
Int_t _coefErrCount
Definition: RooAddPdf.h:147
RooMsgService::instance
static RooMsgService & instance()
Return reference to singleton instance.
Definition: RooMsgService.cxx:366
RooAbsPdf::resetErrorCounters
virtual void resetErrorCounters(Int_t resetValue=10)
Reset error counter to given value, limiting the number of future error messages for this pdf to 'res...
Definition: RooAbsPdf.cxx:604
RooAbsPdf::RooRealIntegral
friend class RooRealIntegral
Definition: RooAbsPdf.h:314
RooBatchCompute::RunContext
Data that has to be passed around when evaluating functions / PDFs.
Definition: RunContext.h:32
RooAddPdf::_codeReg
RooAICRegistry _codeReg
Definition: RooAddPdf.h:137
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:172
RooSpan
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
RooRealConstant.h
RooSetProxy::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Overloaded RooArgSet::add() method inserts 'var' into set and registers 'var' as server to owner with...
Definition: RooSetProxy.cxx:165
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
cxcoutD
#define cxcoutD(a)
Definition: RooMsgService.h:81
int
RooAbsPdf::expectedEvents
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events from this p.d.f for use in extended likelihood calculations.
Definition: RooAbsPdf.cxx:3364
RooFit::RooConst
RooConstVar & RooConst(Double_t val)
Definition: RooGlobalFunc.cxx:347