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