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