Logo ROOT  
Reference Guide
RooRealSumPdf.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 RooRealSumPdf
19 \ingroup Roofitcore
20
21
22The class RooRealSumPdf implements a PDF constructed from a sum of functions:
23\f[
24 \mathrm{PDF}(x) = \frac{ \sum_{i=1}^{n-1} \mathrm{coef}_i * \mathrm{func}_i(x) + \left[ 1 - \sum_{i=1}^{n-1} \mathrm{coef}_i \right] * \mathrm{func}_n(x) }
25 {\sum_{i=1}^{n-1} \mathrm{coef}_i * \int \mathrm{func}_i(x)dx + \left[ 1 - \sum_{i=1}^{n-1} \mathrm{coef}_i \right] * \int \mathrm{func}_n(x) dx }
26\f]
27
28where \f$\mathrm{coef}_i\f$ and \f$\mathrm{func}_i\f$ are RooAbsReal objects, and \f$ x \f$ is the collection of dependents.
29In the present version \f$\mathrm{coef}_i\f$ may not depend on \f$ x \f$, but this limitation could be removed should the need arise.
30
31If the number of coefficients is one less than the number of functions, the PDF is assumed to be normalised. Due to this additional constraint,
32\f$\mathrm{coef}_n\f$ is computed from the other coefficients.
33
34### Extending the PDF
35If an \f$ n^\mathrm{th} \f$ coefficient is provided, the PDF **can** be used as an extended PDF, *i.e.* the total number of events will be measured in addition
36to the fractions of the various functions. **This requires setting the last argument of the constructor to `true`.**
37\note For the RooAddPdf, the extension happens automatically.
38
39### Difference to RooAddPdf / RooRealSumFunc
40- RooAddPdf is a PDF of PDFs, *i.e.* its components need to be normalised and non-negative.
41- RooRealSumPdf is a PDF of functions, *i.e.*, its components can be negative, but their sum cannot be. The normalisation
42 is computed automatically, unless the PDF is extended (see above).
43- RooRealSumFunc is a sum of functions. It is neither normalised, nor need it be positive.
44
45*/
46
47#include "RooFit.h"
48#include "Riostream.h"
49
50#include "TError.h"
51#include "TIterator.h"
52#include "TList.h"
53#include "RooRealSumPdf.h"
54#include "RooRealProxy.h"
55#include "RooPlot.h"
56#include "RooRealVar.h"
57#include "RooAddGenContext.h"
58#include "RooRealConstant.h"
59#include "RooRealIntegral.h"
60#include "RooMsgService.h"
61#include "RooNameReg.h"
62
63#include <algorithm>
64#include <memory>
65
66using namespace std;
67
69;
70
72
73////////////////////////////////////////////////////////////////////////////////
74/// Default constructor
75/// coverity[UNINIT_CTOR]
76
78{
81}
82
83
84
85////////////////////////////////////////////////////////////////////////////////
86/// Constructor with name and title
87
88RooRealSumPdf::RooRealSumPdf(const char *name, const char *title) :
89 RooAbsPdf(name,title),
90 _normIntMgr(this,10),
91 _funcList("!funcList","List of functions",this),
92 _coefList("!coefList","List of coefficients",this),
93 _extended(kFALSE),
94 _doFloor(kFALSE)
95{
96
97}
98
99
100
101////////////////////////////////////////////////////////////////////////////////
102/// Construct p.d.f consisting of \f$ \mathrm{coef}_1 * \mathrm{func}_1 + (1-\mathrm{coef}_1) * \mathrm{func}_2 \f$.
103/// The input coefficients and functions are allowed to be negative
104/// but the resulting sum is not, which is enforced at runtime.
105
106RooRealSumPdf::RooRealSumPdf(const char *name, const char *title,
107 RooAbsReal& func1, RooAbsReal& func2, RooAbsReal& coef1) :
108 RooAbsPdf(name,title),
109 _normIntMgr(this,10),
110 _funcList("!funcList","List of functions",this),
111 _coefList("!coefList","List of coefficients",this),
112 _extended(kFALSE),
113 _doFloor(kFALSE)
114{
115 // Special constructor with two functions and one coefficient
116
117 _funcList.add(func1) ;
118 _funcList.add(func2) ;
119 _coefList.add(coef1) ;
120
121}
122
123
124////////////////////////////////////////////////////////////////////////////////
125/// Constructor for a PDF from a list of functions and coefficients.
126/// It implements
127/// \f[
128/// \sum_i \mathrm{coef}_i \cdot \mathrm{func}_i,
129/// \f]
130/// if \f$ N_\mathrm{coef} = N_\mathrm{func} \f$. With `extended=true`, the coefficients can take any values. With `extended=false`,
131/// there is the danger of getting a degenerate minimisation problem because a PDF has to be normalised, which needs one degree
132/// of freedom less.
133///
134/// A plain (normalised) PDF can therefore be implemented with one less coefficient. RooFit then computes
135/// \f[
136/// \sum_i^{N-1} \mathrm{coef}_i \cdot \mathrm{func}_i + (1 - \sum_i \mathrm{coef}_i ) \cdot \mathrm{func}_N,
137/// \f]
138/// if \f$ N_\mathrm{coef} = N_\mathrm{func} - 1 \f$.
139///
140/// All coefficients and functions are allowed to be negative
141/// but the sum (*i.e.* the PDF) is not, which is enforced at runtime.
142///
143/// \param name Name of the PDF
144/// \param title Title (for plotting)
145/// \param inFuncList List of functions to sum
146/// \param inCoefList List of coefficients
147/// \param extended Interpret as extended PDF (requires equal number of functions and coefficients)
148
149RooRealSumPdf::RooRealSumPdf(const char *name, const char *title,
150 const RooArgList& inFuncList, const RooArgList& inCoefList, Bool_t extended) :
151 RooAbsPdf(name,title),
152 _normIntMgr(this,10),
153 _funcList("!funcList","List of functions",this),
154 _coefList("!coefList","List of coefficients",this),
155 _extended(extended),
156 _doFloor(kFALSE)
157{
158 if (!(inFuncList.getSize()==inCoefList.getSize()+1 || inFuncList.getSize()==inCoefList.getSize())) {
159 coutE(InputArguments) << "RooRealSumPdf::RooRealSumPdf(" << GetName()
160 << ") number of pdfs and coefficients inconsistent, must have Nfunc=Ncoef or Nfunc=Ncoef+1" << endl ;
161 assert(0) ;
162 }
163
164 // Constructor with N functions and N or N-1 coefs
165 for (unsigned int i = 0; i < inCoefList.size(); ++i) {
166 const auto& func = inFuncList[i];
167 const auto& coef = inCoefList[i];
168
169 if (!dynamic_cast<const RooAbsReal*>(&coef)) {
170 coutW(InputArguments) << "RooRealSumPdf::RooRealSumPdf(" << GetName() << ") coefficient " << coef.GetName() << " is not of type RooAbsReal, ignored" << endl ;
171 continue ;
172 }
173 if (!dynamic_cast<const RooAbsReal*>(&func)) {
174 coutW(InputArguments) << "RooRealSumPdf::RooRealSumPdf(" << GetName() << ") func " << func.GetName() << " is not of type RooAbsReal, ignored" << endl ;
175 continue ;
176 }
177 _funcList.add(func) ;
178 _coefList.add(coef) ;
179 }
180
181 if (inFuncList.size() == inCoefList.size() + 1) {
182 const auto& func = inFuncList[inFuncList.size()-1];
183 if (!dynamic_cast<const RooAbsReal*>(&func)) {
184 coutE(InputArguments) << "RooRealSumPdf::RooRealSumPdf(" << GetName() << ") last func " << func.GetName() << " is not of type RooAbsReal, fatal error" << endl ;
185 assert(0) ;
186 }
187 _funcList.add(func);
188 }
189
190}
191
192
193
194
195////////////////////////////////////////////////////////////////////////////////
196/// Copy constructor
197
199 RooAbsPdf(other,name),
200 _normIntMgr(other._normIntMgr,this),
201 _funcList("!funcList",this,other._funcList),
202 _coefList("!coefList",this,other._coefList),
203 _extended(other._extended),
204 _doFloor(other._doFloor)
205{
206
207}
208
209
210
211////////////////////////////////////////////////////////////////////////////////
212/// Destructor
213
215{
216
217}
218
219
220
221
222
223////////////////////////////////////////////////////////////////////////////////
224
226{
228}
229
230
231
232
233////////////////////////////////////////////////////////////////////////////////
234/// Calculate the current value
235
237{
238 Double_t value(0) ;
239
240 // Do running sum of coef/func pairs, calculate lastCoef.
241
242 // N funcs, N-1 coefficients
243 Double_t lastCoef(1) ;
244 auto funcIt = _funcList.begin();
245 for (const auto coefArg : _coefList) {
246 assert(funcIt != _funcList.end());
247 auto func = static_cast<const RooAbsReal*>(*funcIt++);
248 auto coef = static_cast<const RooAbsReal*>(coefArg);
249
250 Double_t coefVal = coef->getVal() ;
251 if (coefVal) {
252 cxcoutD(Eval) << "RooRealSumPdf::eval(" << GetName() << ") coefVal = " << coefVal << " funcVal = " << func->IsA()->GetName() << "::" << func->GetName() << " = " << func->getVal() << endl ;
253 if (func->isSelectedComp()) {
254 value += func->getVal()*coefVal ;
255 }
256 lastCoef -= coef->getVal() ;
257 }
258 }
259
260 if (!haveLastCoef()) {
261 assert(funcIt != _funcList.end());
262 // Add last func with correct coefficient
263 auto func = static_cast<const RooAbsReal*>(*funcIt);
264 if (func->isSelectedComp()) {
265 value += func->getVal()*lastCoef ;
266 }
267
268 cxcoutD(Eval) << "RooRealSumPdf::eval(" << GetName() << ") lastCoef = " << lastCoef << " funcVal = " << func->getVal() << endl ;
269
270 // Warn about coefficient degeneration
271 if (lastCoef<0 || lastCoef>1) {
272 coutW(Eval) << "RooRealSumPdf::evaluate(" << GetName()
273 << ") WARNING: sum of FUNC coefficients not in range [0-1], value="
274 << 1-lastCoef << ". This means that the PDF is not properly normalised. If the PDF was meant to be extended, provide as many coefficients as functions." << endl ;
275 }
276 }
277
278 // Introduce floor if so requested
279 if (value<0 && (_doFloor || _doFloorGlobal)) {
280 value = 0 ;
281 }
282
283 return value ;
284}
285
286
287
288
289////////////////////////////////////////////////////////////////////////////////
290/// Check if FUNC is valid for given normalization set.
291/// Coefficient and FUNC must be non-overlapping, but func-coefficient
292/// pairs may overlap each other.
293///
294/// In the present implementation, coefficients may not be observables or derive
295/// from observables.
296
298{
299 Bool_t ret(kFALSE) ;
300
301 for (unsigned int i=0; i < _coefList.size(); ++i) {
302 const auto& coef = _coefList[i];
303 const auto& func = _funcList[i];
304
305 if (func.observableOverlaps(nset, coef)) {
306 coutE(InputArguments) << "RooRealSumPdf::checkObservables(" << GetName() << "): ERROR: coefficient " << coef.GetName()
307 << " and FUNC " << func.GetName() << " have one or more observables in common" << endl ;
308 ret = kTRUE ;
309 }
310 if (coef.dependsOn(*nset)) {
311 coutE(InputArguments) << "RooRealPdf::checkObservables(" << GetName() << "): ERROR coefficient " << coef.GetName()
312 << " depends on one or more of the following observables" ; nset->Print("1") ;
313 ret = kTRUE ;
314 }
315 }
316
317 return ret ;
318}
319
320
321
322
323////////////////////////////////////////////////////////////////////////////////
324/// Advertise that all integrals can be handled internally.
325
327 const RooArgSet* normSet2, const char* rangeName) const
328{
329 // Handle trivial no-integration scenario
330 if (allVars.getSize()==0) return 0 ;
331 if (_forceNumInt) return 0 ;
332
333 // Select subset of allVars that are actual dependents
334 analVars.add(allVars) ;
335 RooArgSet* normSet = normSet2 ? getObservables(normSet2) : 0 ;
336
337
338 // Check if this configuration was created before
339 Int_t sterileIdx(-1) ;
340 CacheElem* cache = (CacheElem*) _normIntMgr.getObj(normSet,&analVars,&sterileIdx,RooNameReg::ptr(rangeName)) ;
341 if (cache) {
342 //cout << "RooRealSumPdf("<<this<<")::getAnalyticalIntegralWN:"<<GetName()<<"("<<allVars<<","<<analVars<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << " -> " << _normIntMgr.lastIndex()+1 << " (cached)" << endl;
343 return _normIntMgr.lastIndex()+1 ;
344 }
345
346 // Create new cache element
347 cache = new CacheElem ;
348
349 // Make list of function projection and normalization integrals
350 for (const auto elm : _funcList) {
351 const auto func = static_cast<RooAbsReal*>(elm);
352
353 RooAbsReal* funcInt = func->createIntegral(analVars,rangeName) ;
354 if(funcInt->InheritsFrom(RooRealIntegral::Class())) ((RooRealIntegral*)funcInt)->setAllowComponentSelection(true);
355 cache->_funcIntList.addOwned(*funcInt) ;
356 if (normSet && normSet->getSize()>0) {
357 RooAbsReal* funcNorm = func->createIntegral(*normSet) ;
358 cache->_funcNormList.addOwned(*funcNorm) ;
359 }
360 }
361
362 // Store cache element
363 Int_t code = _normIntMgr.setObj(normSet,&analVars,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
364
365 if (normSet) {
366 delete normSet ;
367 }
368
369 //cout << "RooRealSumPdf("<<this<<")::getAnalyticalIntegralWN:"<<GetName()<<"("<<allVars<<","<<analVars<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << " -> " << code+1 << endl;
370 return code+1 ;
371}
372
373
374
375
376////////////////////////////////////////////////////////////////////////////////
377/// Implement analytical integrations by deferring integration of component
378/// functions to integrators of components.
379
380Double_t RooRealSumPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet2, const char* rangeName) const
381{
382 // Handle trivial passthrough scenario
383 if (code==0) return getVal(normSet2) ;
384
385
386 // WVE needs adaptation for rangeName feature
387 CacheElem* cache = (CacheElem*) _normIntMgr.getObjByIndex(code-1) ;
388 if (cache==0) { // revive the (sterilized) cache
389 //cout << "RooRealSumPdf("<<this<<")::analyticalIntegralWN:"<<GetName()<<"("<<code<<","<<(normSet2?*normSet2:RooArgSet())<<","<<(rangeName?rangeName:"<none>") << ": reviving cache "<< endl;
390 std::unique_ptr<RooArgSet> vars( getParameters(RooArgSet()) );
391 std::unique_ptr<RooArgSet> iset( _normIntMgr.nameSet2ByIndex(code-1)->select(*vars) );
392 std::unique_ptr<RooArgSet> nset( _normIntMgr.nameSet1ByIndex(code-1)->select(*vars) );
394 Int_t code2 = getAnalyticalIntegralWN(*iset,dummy,nset.get(),rangeName);
395 R__ASSERT(code==code2); // must have revived the right (sterilized) slot...
396 cache = (CacheElem*) _normIntMgr.getObjByIndex(code-1) ;
397 R__ASSERT(cache!=0);
398 }
399
400 Double_t value(0) ;
401
402 // N funcs, N-1 coefficients
403 Double_t lastCoef(1) ;
404 auto funcIt = _funcList.begin();
405 auto funcIntIt = cache->_funcIntList.begin();
406 for (const auto coefArg : _coefList) {
407 assert(funcIt != _funcList.end());
408 const auto coef = static_cast<const RooAbsReal*>(coefArg);
409 const auto func = static_cast<const RooAbsReal*>(*funcIt++);
410 const auto funcInt = static_cast<RooAbsReal*>(*funcIntIt++);
411
412 Double_t coefVal = coef->getVal(normSet2) ;
413 if (coefVal) {
414 assert(func);
415 if (normSet2 ==0 || func->isSelectedComp()) {
416 assert(funcInt);
417 value += funcInt->getVal()*coefVal ;
418 }
419 lastCoef -= coef->getVal(normSet2) ;
420 }
421 }
422
423 if (!haveLastCoef()) {
424 // Add last func with correct coefficient
425 const auto func = static_cast<const RooAbsReal*>(*funcIt);
426 const auto funcInt = static_cast<RooAbsReal*>(*funcIntIt);
427 assert(func);
428
429 if (normSet2 ==0 || func->isSelectedComp()) {
430 assert(funcInt);
431 value += funcInt->getVal()*lastCoef ;
432 }
433
434 // Warn about coefficient degeneration
435 if (lastCoef<0 || lastCoef>1) {
436 coutW(Eval) << "RooRealSumPdf::evaluate(" << GetName()
437 << " WARNING: sum of FUNC coefficients not in range [0-1], value="
438 << 1-lastCoef << endl ;
439 }
440 }
441
442 Double_t normVal(1) ;
443 if (normSet2 && normSet2->getSize()>0) {
444 normVal = 0 ;
445
446 // N funcs, N-1 coefficients
447 auto funcNormIter = cache->_funcNormList.begin();
448 for (const auto coefAsArg : _coefList) {
449 auto coef = static_cast<RooAbsReal*>(coefAsArg);
450 auto funcNorm = static_cast<RooAbsReal*>(*funcNormIter++);
451
452 Double_t coefVal = coef->getVal(normSet2);
453 if (coefVal) {
454 assert(funcNorm);
455 normVal += funcNorm->getVal()*coefVal ;
456 }
457 }
458
459 // Add last func with correct coefficient
460 if (!haveLastCoef()) {
461 auto funcNorm = static_cast<RooAbsReal*>(*funcNormIter);
462 assert(funcNorm);
463
464 normVal += funcNorm->getVal()*lastCoef;
465 }
466 }
467
468 return value / normVal;
469}
470
471
472////////////////////////////////////////////////////////////////////////////////
473
475{
476 Double_t n = getNorm(nset) ;
477 if (n<0) {
478 logEvalError("Expected number of events is negative") ;
479 }
480 return n ;
481}
482
483
484////////////////////////////////////////////////////////////////////////////////
485
486std::list<Double_t>* RooRealSumPdf::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
487{
488 list<Double_t>* sumBinB = 0 ;
489 Bool_t needClean(kFALSE) ;
490
491 // Loop over components pdf
492 for (const auto elm : _funcList) {
493 auto func = static_cast<RooAbsReal*>(elm);
494
495 list<Double_t>* funcBinB = func->binBoundaries(obs,xlo,xhi) ;
496
497 // Process hint
498 if (funcBinB) {
499 if (!sumBinB) {
500 // If this is the first hint, then just save it
501 sumBinB = funcBinB ;
502 } else {
503
504 list<Double_t>* newSumBinB = new list<Double_t>(sumBinB->size()+funcBinB->size()) ;
505
506 // Merge hints into temporary array
507 merge(funcBinB->begin(),funcBinB->end(),sumBinB->begin(),sumBinB->end(),newSumBinB->begin()) ;
508
509 // Copy merged array without duplicates to new sumBinBArrau
510 delete sumBinB ;
511 delete funcBinB ;
512 sumBinB = newSumBinB ;
513 needClean = kTRUE ;
514 }
515 }
516 }
517
518 // Remove consecutive duplicates
519 if (needClean) {
520 list<Double_t>::iterator new_end = unique(sumBinB->begin(),sumBinB->end()) ;
521 sumBinB->erase(new_end,sumBinB->end()) ;
522 }
523
524 return sumBinB ;
525}
526
527
528
529//_____________________________________________________________________________B
531{
532 // If all components that depend on obs are binned that so is the product
533
534 for (const auto elm : _funcList) {
535 auto func = static_cast<RooAbsReal*>(elm);
536
537 if (func->dependsOn(obs) && !func->isBinnedDistribution(obs)) {
538 return kFALSE ;
539 }
540 }
541
542 return kTRUE ;
543}
544
545
546
547
548
549////////////////////////////////////////////////////////////////////////////////
550
551std::list<Double_t>* RooRealSumPdf::plotSamplingHint(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
552{
553 list<Double_t>* sumHint = 0 ;
554 Bool_t needClean(kFALSE) ;
555
556 // Loop over components pdf
557 for (const auto elm : _funcList) {
558 auto func = static_cast<RooAbsReal*>(elm);
559
560 list<Double_t>* funcHint = func->plotSamplingHint(obs,xlo,xhi) ;
561
562 // Process hint
563 if (funcHint) {
564 if (!sumHint) {
565
566 // If this is the first hint, then just save it
567 sumHint = funcHint ;
568
569 } else {
570
571 list<Double_t>* newSumHint = new list<Double_t>(sumHint->size()+funcHint->size()) ;
572
573 // Merge hints into temporary array
574 merge(funcHint->begin(),funcHint->end(),sumHint->begin(),sumHint->end(),newSumHint->begin()) ;
575
576 // Copy merged array without duplicates to new sumHintArrau
577 delete sumHint ;
578 sumHint = newSumHint ;
579 needClean = kTRUE ;
580 }
581 }
582 }
583
584 // Remove consecutive duplicates
585 if (needClean) {
586 list<Double_t>::iterator new_end = unique(sumHint->begin(),sumHint->end()) ;
587 sumHint->erase(new_end,sumHint->end()) ;
588 }
589
590 return sumHint ;
591}
592
593
594
595
596////////////////////////////////////////////////////////////////////////////////
597/// Label OK'ed components of a RooRealSumPdf with cache-and-track
598
600{
601 for (const auto sarg : _funcList) {
602 if (sarg->canNodeBeCached()==Always) {
603 trackNodes.add(*sarg) ;
604 //cout << "tracking node RealSumPdf component " << sarg->IsA()->GetName() << "::" << sarg->GetName() << endl ;
605 }
606 }
607}
608
609
610
611////////////////////////////////////////////////////////////////////////////////
612/// Customized printing of arguments of a RooRealSumPdf to more intuitively reflect the contents of the
613/// product operator construction
614
615void RooRealSumPdf::printMetaArgs(ostream& os) const
616{
617
619
620 if (_coefList.getSize()!=0) {
621 auto funcIter = _funcList.begin();
622
623 for (const auto coef : _coefList) {
624 if (!first) {
625 os << " + " ;
626 } else {
627 first = kFALSE ;
628 }
629 const auto func = *(funcIter++);
630 os << coef->GetName() << " * " << func->GetName();
631 }
632
633 if (funcIter != _funcList.end()) {
634 os << " + [%] * " << (*funcIter)->GetName() ;
635 }
636 } else {
637
638 for (const auto func : _funcList) {
639 if (!first) {
640 os << " + " ;
641 } else {
642 first = kFALSE ;
643 }
644 os << func->GetName() ;
645 }
646 }
647
648 os << " " ;
649}
void Class()
Definition: Class.C:29
static RooMathCoreReg dummy
#define cxcoutD(a)
Definition: RooMsgService.h:82
#define coutW(a)
Definition: RooMsgService.h:33
#define coutE(a)
Definition: RooMsgService.h:34
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define ClassImp(name)
Definition: Rtypes.h:365
#define R__ASSERT(e)
Definition: TError.h:96
char name[80]
Definition: TGX11.cxx:109
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:240
friend class RooArgSet
Definition: RooAbsArg.h:551
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:548
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
Int_t getSize() const
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Storage_t::size_type size() const
const_iterator begin() const
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
TIterator end() and range-based for loops.")
Double_t getNorm(const RooArgSet &nset) const
Definition: RooAbsPdf.h:208
friend class CacheElem
Definition: RooAbsPdf.h:334
@ CanBeExtended
Definition: RooAbsPdf.h:230
@ CanNotBeExtended
Definition: RooAbsPdf.h:230
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:59
Bool_t _forceNumInt
Definition: RooAbsReal.h:450
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &, Double_t, Double_t) const
Definition: RooAbsReal.h:311
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Definition: RooAbsReal.h:310
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:87
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:562
void logEvalError(const char *message, const char *serverValueString=0) const
Log evaluation error message.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
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
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
const RooNameSet * nameSet2ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Int_t lastIndex() const
const RooNameSet * nameSet1ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooArgSet * select(const RooArgSet &list) const
Construct a RooArgSet of objects in input 'list' whose names match to those in the internal name list...
Definition: RooNameSet.cxx:187
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
The class RooRealSumPdf implements a PDF constructed from a sum of functions:
Definition: RooRealSumPdf.h:24
virtual ExtendMode extendMode() const
RooObjCacheManager _normIntMgr
Definition: RooRealSumPdf.h:82
Double_t evaluate() const
Calculate the current value.
RooListProxy _coefList
Definition: RooRealSumPdf.h:86
RooListProxy _funcList
Definition: RooRealSumPdf.h:85
RooRealSumPdf()
Default constructor coverity[UNINIT_CTOR].
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events from this p.d.f for use in extended likelihood calculations.
Bool_t isBinnedDistribution(const RooArgSet &obs) const
static Bool_t _doFloorGlobal
Definition: RooRealSumPdf.h:90
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &, Double_t, Double_t) const
virtual ~RooRealSumPdf()
Destructor.
virtual Bool_t checkObservables(const RooArgSet *nset) const
Check if FUNC is valid for given normalization set.
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooRealSumPdf to more intuitively reflect the contents of the p...
virtual void setCacheAndTrackHints(RooArgSet &)
Label OK'ed components of a RooRealSumPdf with cache-and-track.
Bool_t _extended
Definition: RooRealSumPdf.h:87
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
Advertise that all integrals can be handled internally.
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Implement analytical integrations by deferring integration of component functions to integrators of c...
bool haveLastCoef() const
Definition: RooRealSumPdf.h:94
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:443
const Int_t n
Definition: legend1.C:16
@ InputArguments
Definition: RooGlobalFunc.h:68
Definition: first.py:1