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