Logo ROOT   6.10/09
Reference Guide
RooAddModel.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 //
19 // RooAddModel is an efficient implementation of a sum of PDFs of the form
20 //
21 // c_1*PDF_1 + c_2*PDF_2 + ... c_n*PDF_n
22 //
23 // or
24 //
25 // c_1*PDF_1 + c_2*PDF_2 + ... (1-sum(c_1...c_n-1))*PDF_n
26 //
27 // The first form is for extended likelihood fits, where the
28 // expected number of events is Sum(i) c_i. The coefficients c_i
29 // can either be explicitly provided, or, if all components support
30 // extended likelihood fits, they can be calculated the contribution
31 // of each PDF to the total number of expected events.
32 //
33 // In the second form, the sum of the coefficients is enforced to be one,
34 // and the coefficient of the last PDF is calculated from that condition.
35 //
36 // RooAddPdf relies on each component PDF to be normalized and will perform
37 // no normalization other than calculating the proper last coefficient c_n, if requested.
38 // An (enforced) condition for this assuption is that each PDF_i is independent
39 // of each coefficient_i.
40 //
41 //
42 
43 #include "RooFit.h"
44 #include "RooMsgService.h"
45 
46 #include "TIterator.h"
47 #include "TIterator.h"
48 #include "TList.h"
49 #include "RooAddModel.h"
50 #include "RooDataSet.h"
51 #include "RooRealProxy.h"
52 #include "RooPlot.h"
53 #include "RooRealVar.h"
54 #include "RooAddGenContext.h"
55 #include "RooRealConstant.h"
56 #include "RooNameReg.h"
57 #include "RooMsgService.h"
58 #include "RooRealIntegral.h"
59 
60 #include "Riostream.h"
61 
62 
63 
64 using namespace std;
65 
67 ;
68 
69 
70 ////////////////////////////////////////////////////////////////////////////////
71 
73  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
74  _refCoefRangeName(0),
75  _codeReg(10),
76  _snormList(0)
77 {
80 
81  _coefCache = new Double_t[10] ;
83 }
84 
85 
86 
87 ////////////////////////////////////////////////////////////////////////////////
88 /// Generic constructor from list of PDFs and list of coefficients.
89 /// Each pdf list element (i) is paired with coefficient list element (i).
90 /// The number of coefficients must be either equal to the number of PDFs,
91 /// in which case extended MLL fitting is enabled, or be one less.
92 ///
93 /// All PDFs must inherit from RooAbsPdf. All coefficients must inherit from RooAbsReal
94 
95 RooAddModel::RooAddModel(const char *name, const char *title, const RooArgList& inPdfList, const RooArgList& inCoefList, Bool_t ownPdfList) :
96  RooResolutionModel(name,title,((RooResolutionModel*)inPdfList.at(0))->convVar()),
97  _refCoefNorm("!refCoefNorm","Reference coefficient normalization set",this,kFALSE,kFALSE),
100  _projCacheMgr(this,10),
101  _intCacheMgr(this,10),
102  _codeReg(10),
103  _pdfList("!pdfs","List of PDFs",this),
104  _coefList("!coefficients","List of coefficients",this),
107 {
108  if (inPdfList.getSize()>inCoefList.getSize()+1) {
109  coutE(InputArguments) << "RooAddModel::RooAddModel(" << GetName()
110  << ") number of pdfs and coefficients inconsistent, must have Npdf=Ncoef or Npdf=Ncoef+1" << endl ;
111  assert(0) ;
112  }
113 
116 
117  // Constructor with N PDFs and N or N-1 coefs
118  TIterator* pdfIter = inPdfList.createIterator() ;
119  TIterator* coefIter = inCoefList.createIterator() ;
120  RooAbsPdf* pdf ;
121  RooAbsReal* coef ;
122 
123  while((coef = (RooAbsPdf*)coefIter->Next())) {
124  pdf = (RooAbsPdf*) pdfIter->Next() ;
125  if (!pdf) {
126  coutE(InputArguments) << "RooAddModel::RooAddModel(" << GetName()
127  << ") number of pdfs and coefficients inconsistent, must have Npdf=Ncoef or Npdf=Ncoef+1" << endl ;
128  assert(0) ;
129  }
130  if (!dynamic_cast<RooAbsReal*>(coef)) {
131  coutE(InputArguments) << "RooAddModel::RooAddModel(" << GetName() << ") coefficient " << coef->GetName() << " is not of type RooAbsReal, ignored" << endl ;
132  continue ;
133  }
134  if (!dynamic_cast<RooAbsReal*>(pdf)) {
135  coutE(InputArguments) << "RooAddModel::RooAddModel(" << GetName() << ") pdf " << pdf->GetName() << " is not of type RooAbsPdf, ignored" << endl ;
136  continue ;
137  }
138  _pdfList.add(*pdf) ;
139  _coefList.add(*coef) ;
140  }
141 
142  pdf = (RooAbsPdf*) pdfIter->Next() ;
143  if (pdf) {
144  if (!dynamic_cast<RooAbsReal*>(pdf)) {
145  coutE(InputArguments) << "RooAddModel::RooAddModel(" << GetName() << ") last pdf " << coef->GetName() << " is not of type RooAbsPdf, fatal error" << endl ;
146  assert(0) ;
147  }
148  _pdfList.add(*pdf) ;
149  } else {
151  }
152 
153  delete pdfIter ;
154  delete coefIter ;
155 
158 
159  if (ownPdfList) {
161  }
162 
163 }
164 
165 
166 
167 ////////////////////////////////////////////////////////////////////////////////
168 /// Copy constructor
169 
170 RooAddModel::RooAddModel(const RooAddModel& other, const char* name) :
171  RooResolutionModel(other,name),
172  _refCoefNorm("!refCoefNorm",this,other._refCoefNorm),
175  _projCacheMgr(other._projCacheMgr,this),
176  _intCacheMgr(other._intCacheMgr,this),
177  _codeReg(other._codeReg),
178  _pdfList("!pdfs",this,other._pdfList),
179  _coefList("!coefficients",this,other._coefList),
182 {
187 }
188 
189 
190 
191 ////////////////////////////////////////////////////////////////////////////////
192 /// Destructor
193 
195 {
196  delete _pdfIter ;
197  delete _coefIter ;
198 
199  if (_coefCache) delete[] _coefCache ;
200 }
201 
202 
203 
204 ////////////////////////////////////////////////////////////////////////////////
205 /// By default the interpretation of the fraction coefficients is
206 /// performed in the contextual choice of observables. This makes the
207 /// shape of the p.d.f explicitly dependent on the choice of
208 /// observables. This method instructs RooAddPdf to freeze the
209 /// interpretation of the coefficients to be done in the given set of
210 /// observables. If frozen, fractions are automatically transformed
211 /// from the reference normalization set to the contextual normalization
212 /// set by ratios of integrals
213 
215 {
216  if (refCoefNorm.getSize()==0) {
218  return ;
219  }
220  _projectCoefs = kTRUE ;
221 
223  _refCoefNorm.add(refCoefNorm) ;
224 
225  _projCacheMgr.reset() ;
226 }
227 
228 
229 
230 ////////////////////////////////////////////////////////////////////////////////
231 /// By default the interpretation of the fraction coefficients is
232 /// performed in the default range. This make the shape of a RooAddPdf
233 /// explicitly dependent on the range of the observables. To allow
234 /// a range independent definition of the fraction this function
235 /// instructs RooAddPdf to freeze its interpretation in the given
236 /// named range. If the current normalization range is different
237 /// from the reference range, the appropriate fraction coefficients
238 /// are automically calculation from the reference fractions using
239 /// ratios if integrals
240 
241 void RooAddModel::fixCoefRange(const char* rangeName)
242 {
243  _refCoefRangeName = (TNamed*)RooNameReg::ptr(rangeName) ;
245 }
246 
247 
248 
249 ////////////////////////////////////////////////////////////////////////////////
250 /// Instantiate a clone of this resolution model representing a convolution with given
251 /// basis function. The owners object name is incorporated in the clones name
252 /// to avoid multiple convolution objects with the same name in complex PDF structures.
253 ///
254 /// RooAddModel will clone all the component models to create a composite convolution object
255 
257 {
258  // Check that primary variable of basis functions is our convolution variable
259  if (inBasis->getParameter(0) != x.absArg()) {
260  coutE(InputArguments) << "RooAddModel::convolution(" << GetName()
261  << ") convolution parameter of basis function and PDF don't match" << endl ;
262  ccoutE(InputArguments) << "basis->findServer(0) = " << inBasis->findServer(0) << " " << inBasis->findServer(0)->GetName() << endl ;
263  ccoutE(InputArguments) << "x.absArg() = " << x.absArg() << " " << x.absArg()->GetName() << endl ;
264  inBasis->Print("v") ;
265  return 0 ;
266  }
267 
268  TString newName(GetName()) ;
269  newName.Append("_conv_") ;
270  newName.Append(inBasis->GetName()) ;
271  newName.Append("_[") ;
272  newName.Append(owner->GetName()) ;
273  newName.Append("]") ;
274 
275  TString newTitle(GetTitle()) ;
276  newTitle.Append(" convoluted with basis function ") ;
277  newTitle.Append(inBasis->GetName()) ;
278 
279  _pdfIter->Reset() ;
281  RooArgList modelList ;
282  while((model = (RooResolutionModel*)_pdfIter->Next())) {
283  // Create component convolution
284  RooResolutionModel* conv = model->convolution(inBasis,owner) ;
285  modelList.add(*conv) ;
286  }
287 
288  _coefIter->Reset() ;
289  RooAbsReal* coef ;
290  RooArgList theCoefList ;
291  while((coef = (RooAbsReal*)_coefIter->Next())) {
292  theCoefList.add(*coef) ;
293  }
294 
295  RooAddModel* convSum = new RooAddModel(newName,newTitle,modelList,theCoefList,kTRUE) ;
296  for (std::set<std::string>::const_iterator attrIt = _boolAttrib.begin();
297  attrIt != _boolAttrib.end(); ++attrIt) {
298  convSum->setAttribute((*attrIt).c_str()) ;
299  }
300  for (std::map<std::string,std::string>::const_iterator attrIt = _stringAttrib.begin();
301  attrIt != _stringAttrib.end(); ++attrIt) {
302  convSum->setStringAttribute((attrIt->first).c_str(), (attrIt->second).c_str()) ;
303  }
304  convSum->changeBasis(inBasis) ;
305  return convSum ;
306 }
307 
308 
309 
310 ////////////////////////////////////////////////////////////////////////////////
311 /// Return code for basis function representing by 'name' string.
312 /// The basis code of the first component model will be returned,
313 /// if the basis is supported by all components. Otherwise 0
314 /// is returned
315 
316 Int_t RooAddModel::basisCode(const char* name) const
317 {
318  TIterator* mIter = _pdfList.createIterator() ;
320  Bool_t first(kTRUE), code(0) ;
321  while((model = (RooResolutionModel*)mIter->Next())) {
322  Int_t subCode = model->basisCode(name) ;
323  if (first) {
324  code = subCode ;
325  first = kFALSE ;
326  } else if (subCode==0) {
327  code = 0 ;
328  }
329  }
330  delete mIter ;
331 
332  return code ;
333 }
334 
335 
336 
337 ////////////////////////////////////////////////////////////////////////////////
338 /// Retrieve cache element with for calculation of p.d.f value with normalization set nset and integrated over iset
339 /// in range 'rangeName'. If cache element does not exist, create and fill it on the fly. The cache contains
340 /// suplemental normalization terms (in case not all added p.d.f.s have the same observables), projection
341 /// integrals to calculated transformed fraction coefficients when a frozen reference frame is provided
342 /// and projection integrals for similar transformations when a frozen reference range is provided.
343 
344 RooAddModel::CacheElem* RooAddModel::getProjCache(const RooArgSet* nset, const RooArgSet* iset, const char* rangeName) const
345 {
346  // Check if cache already exists
347  CacheElem* cache = (CacheElem*) _projCacheMgr.getObj(nset,iset,0,RooNameReg::ptr(rangeName)) ;
348  if (cache) {
349  return cache ;
350  }
351 
352  //Create new cache
353  cache = new CacheElem ;
354 
355  // *** PART 1 : Create supplemental normalization list ***
356 
357  // Retrieve the combined set of dependents of this PDF ;
358  RooArgSet *fullDepList = getObservables(nset) ;
359  if (iset) {
360  fullDepList->remove(*iset,kTRUE,kTRUE) ;
361  }
362 
363  // Fill with dummy unit RRVs for now
364  _pdfIter->Reset() ;
365  _coefIter->Reset() ;
366  RooAbsPdf* pdf ;
367  RooAbsReal* coef ;
368  while((pdf=(RooAbsPdf*)_pdfIter->Next())) {
369  coef=(RooAbsPdf*)_coefIter->Next() ;
370 
371  // Start with full list of dependents
372  RooArgSet supNSet(*fullDepList) ;
373 
374  // Remove PDF dependents
375  RooArgSet* pdfDeps = pdf->getObservables(nset) ;
376  if (pdfDeps) {
377  supNSet.remove(*pdfDeps,kTRUE,kTRUE) ;
378  delete pdfDeps ;
379  }
380 
381  // Remove coef dependents
382  RooArgSet* coefDeps = coef ? coef->getObservables(nset) : 0 ;
383  if (coefDeps) {
384  supNSet.remove(*coefDeps,kTRUE,kTRUE) ;
385  delete coefDeps ;
386  }
387 
388  RooAbsReal* snorm ;
389  TString name(GetName()) ;
390  name.Append("_") ;
391  name.Append(pdf->GetName()) ;
392  name.Append("_SupNorm") ;
393  if (supNSet.getSize()>0) {
394  snorm = new RooRealIntegral(name,"Supplemental normalization integral",RooRealConstant::value(1.0),supNSet) ;
395  } else {
396  snorm = new RooRealVar(name,"Unit Supplemental normalization integral",1.0) ;
397  }
398  cache->_suppNormList.addOwned(*snorm) ;
399  }
400 
401  delete fullDepList ;
402 
403  if (_verboseEval>1) {
404  cxcoutD(Caching) << "RooAddModel::syncSuppNormList(" << GetName() << ") synching supplemental normalization list for norm" << (nset?*nset:RooArgSet()) << endl ;
405  if (dologD(Caching)) {
406  cache->_suppNormList.Print("v") ;
407  }
408  }
409 
410 
411  // *** PART 2 : Create projection coefficients ***
412 
413  // If no projections required stop here
414  if (!_projectCoefs || _basis!=0 ) {
415  _projCacheMgr.setObj(nset,iset,cache,RooNameReg::ptr(rangeName)) ;
416  return cache ;
417  }
418 
419 
420  // Reduce iset/nset to actual dependents of this PDF
421  RooArgSet* nset2 = nset ? getObservables(nset) : new RooArgSet() ;
422 
423  // Check if requested transformation is not identity
424  if (!nset2->equals(_refCoefNorm) || _refCoefRangeName !=0 || rangeName !=0 ) {
425 
426  coutI(Caching) << "RooAddModel::syncCoefProjList(" << GetName() << ") creating coefficient projection integrals" << endl ;
427  ccoutI(Caching) << " from current normalization: " ; nset2->Print("1") ;
428  ccoutI(Caching) << " with current range: " << (rangeName?rangeName:"<none>") << endl ;
429  ccoutI(Caching) << " to reference normalization: " ; _refCoefNorm.Print("1") ;
430  ccoutI(Caching) << " with reference range: " << (_refCoefRangeName?RooNameReg::str(_refCoefRangeName):"<none>") << endl ;
431 
432  // Recalculate projection integrals of PDFs
433  _pdfIter->Reset() ;
434  RooAbsPdf* thePdf ;
435 
436  while((thePdf=(RooAbsPdf*)_pdfIter->Next())) {
437 
438  // Calculate projection integral
439  RooAbsReal* pdfProj ;
440  if (!nset2->equals(_refCoefNorm)) {
441  pdfProj = thePdf->createIntegral(*nset2,_refCoefNorm) ;
442  pdfProj->setOperMode(operMode()) ;
443  } else {
444  TString name(GetName()) ;
445  name.Append("_") ;
446  name.Append(thePdf->GetName()) ;
447  name.Append("_ProjectNorm") ;
448  pdfProj = new RooRealVar(name,"Unit Projection normalization integral",1.0) ;
449  }
450 
451  cache->_projList.addOwned(*pdfProj) ;
452 
453  // Calculation optional supplemental normalization term
454  RooArgSet supNormSet(_refCoefNorm) ;
455  RooArgSet* deps = thePdf->getParameters(RooArgSet()) ;
456  supNormSet.remove(*deps,kTRUE,kTRUE) ;
457  delete deps ;
458 
459  RooAbsReal* snorm ;
460  TString name(GetName()) ;
461  name.Append("_") ;
462  name.Append(thePdf->GetName()) ;
463  name.Append("_ProjSupNorm") ;
464  if (supNormSet.getSize()>0) {
465  snorm = new RooRealIntegral(name,"Projection Supplemental normalization integral",
466  RooRealConstant::value(1.0),supNormSet) ;
467  } else {
468  snorm = new RooRealVar(name,"Unit Projection Supplemental normalization integral",1.0) ;
469  }
470  cache->_suppProjList.addOwned(*snorm) ;
471 
472  // Calculate reference range adjusted projection integral
473  RooAbsReal* rangeProj1 ;
475  rangeProj1 = thePdf->createIntegral(_refCoefNorm,_refCoefNorm,RooNameReg::str(_refCoefRangeName)) ;
476  rangeProj1->setOperMode(operMode()) ;
477  } else {
478  TString theName(GetName()) ;
479  theName.Append("_") ;
480  theName.Append(thePdf->GetName()) ;
481  theName.Append("_RangeNorm1") ;
482  rangeProj1 = new RooRealVar(theName,"Unit range normalization integral",1.0) ;
483  }
484  cache->_refRangeProjList.addOwned(*rangeProj1) ;
485 
486 
487  // Calculate range adjusted projection integral
488  RooAbsReal* rangeProj2 ;
489  if (rangeName && _refCoefNorm.getSize()>0) {
490  rangeProj2 = thePdf->createIntegral(_refCoefNorm,_refCoefNorm,rangeName) ;
491  rangeProj2->setOperMode(operMode()) ;
492  } else {
493  TString theName(GetName()) ;
494  theName.Append("_") ;
495  theName.Append(thePdf->GetName()) ;
496  theName.Append("_RangeNorm2") ;
497  rangeProj2 = new RooRealVar(theName,"Unit range normalization integral",1.0) ;
498  }
499  cache->_rangeProjList.addOwned(*rangeProj2) ;
500 
501  }
502 
503  }
504 
505  delete nset2 ;
506 
507  _projCacheMgr.setObj(nset,iset,cache,RooNameReg::ptr(rangeName)) ;
508 
509  return cache ;
510 }
511 
512 
513 
514 ////////////////////////////////////////////////////////////////////////////////
515 /// Update the coefficient values in the given cache element: calculate new remainder
516 /// fraction, normalize fractions obtained from extended ML terms to unity and
517 /// multiply these the various range and dimensional corrections needed in the
518 /// current use context
519 
520 void RooAddModel::updateCoefficients(CacheElem& cache, const RooArgSet* nset) const
521 {
522  // cxcoutD(ChangeTracking) << "RooAddModel::updateCoefficients(" << GetName() << ") update coefficients" << endl ;
523 
524  Int_t i ;
525 
526  // Straight coefficients
527  if (_allExtendable) {
528 
529  // coef[i] = expectedEvents[i] / SUM(expectedEvents)
530  Double_t coefSum(0) ;
531  for (i=0 ; i<_pdfList.getSize() ; i++) {
533  coefSum += _coefCache[i] ;
534  }
535  if (coefSum==0.) {
536  coutW(Eval) << "RooAddModel::updateCoefCache(" << GetName() << ") WARNING: total number of expected events is 0" << endl ;
537  } else {
538  for (i=0 ; i<_pdfList.getSize() ; i++) {
539  _coefCache[i] /= coefSum ;
540  }
541  }
542 
543  } else {
544  if (_haveLastCoef) {
545 
546  // coef[i] = coef[i] / SUM(coef)
547  Double_t coefSum(0) ;
548  for (i=0 ; i<_coefList.getSize() ; i++) {
549  _coefCache[i] = ((RooAbsPdf*)_coefList.at(i))->getVal(nset) ;
550  coefSum += _coefCache[i] ;
551  }
552  for (i=0 ; i<_coefList.getSize() ; i++) {
553  _coefCache[i] /= coefSum ;
554  }
555  } else {
556 
557  // coef[i] = coef[i] ; coef[n] = 1-SUM(coef[0...n-1])
558  Double_t lastCoef(1) ;
559  for (i=0 ; i<_coefList.getSize() ; i++) {
560  _coefCache[i] = ((RooAbsPdf*)_coefList.at(i))->getVal(nset) ;
561  cxcoutD(Caching) << "SYNC: orig coef[" << i << "] = " << _coefCache[i] << endl ;
562  lastCoef -= _coefCache[i] ;
563  }
564  _coefCache[_coefList.getSize()] = lastCoef ;
565  cxcoutD(Caching) << "SYNC: orig coef[" << _coefList.getSize() << "] = " << _coefCache[_coefList.getSize()] << endl ;
566 
567 
568  // Warn about coefficient degeneration
569  if ((lastCoef<-1e-05 || (lastCoef-1)>1e-5) && _coefErrCount-->0) {
570  coutW(Eval) << "RooAddModel::updateCoefCache(" << GetName()
571  << " WARNING: sum of PDF coefficients not in range [0-1], value="
572  << 1-lastCoef << endl ;
573  if (_coefErrCount==0) {
574  coutW(Eval) << " (no more will be printed)" << endl ;
575  }
576  }
577  }
578  }
579 
580 
581 
582  // Stop here if not projection is required or needed
583  if ((!_projectCoefs) || cache._projList.getSize()==0) {
584  // cout << "SYNC no projection required rangeName = " << (rangeName?rangeName:"<none>") << endl ;
585  return ;
586  }
587 
588  // Adjust coefficients for given projection
589  Double_t coefSum(0) ;
590  for (i=0 ; i<_pdfList.getSize() ; i++) {
592 
593  RooAbsReal* pp = ((RooAbsReal*)cache._projList.at(i)) ;
594  RooAbsReal* sn = ((RooAbsReal*)cache._suppProjList.at(i)) ;
595  RooAbsReal* r1 = ((RooAbsReal*)cache._refRangeProjList.at(i)) ;
596  RooAbsReal* r2 = ((RooAbsReal*)cache._rangeProjList.at(i)) ;
597 
598  if (dologD(Eval)) {
599  cxcoutD(Eval) << "pp = " << pp->GetName() << endl
600  << "sn = " << sn->GetName() << endl
601  << "r1 = " << r1->GetName() << endl
602  << "r2 = " << r2->GetName() << endl ;
605  }
606 
607  Double_t proj = pp->getVal()/sn->getVal()*(r2->getVal()/r1->getVal()) ;
608 
610 
611  _coefCache[i] *= proj ;
612  coefSum += _coefCache[i] ;
613  }
614  for (i=0 ; i<_pdfList.getSize() ; i++) {
615  _coefCache[i] /= coefSum ;
616 // cout << "POST-SYNC coef[" << i << "] = " << _coefCache[i] << endl ;
617  }
618 
619 
620 
621 }
622 
623 
624 
625 ////////////////////////////////////////////////////////////////////////////////
626 /// Calculate the current value
627 
629 {
630  const RooArgSet* nset = _normSet ;
631  CacheElem* cache = getProjCache(nset) ;
632 
633  updateCoefficients(*cache,nset) ;
634 
635 
636  // Do running sum of coef/pdf pairs, calculate lastCoef.
637  _pdfIter->Reset() ;
638  _coefIter->Reset() ;
639  RooAbsPdf* pdf ;
640 
641  Double_t snormVal ;
642  Double_t value(0) ;
643  Int_t i(0) ;
644  while((pdf = (RooAbsPdf*)_pdfIter->Next())) {
645  if (_coefCache[i]!=0.) {
646  snormVal = nset ? ((RooAbsReal*)cache->_suppNormList.at(i))->getVal() : 1.0 ;
647  Double_t pdfVal = pdf->getVal(nset) ;
648  // Double_t pdfNorm = pdf->getNorm(nset) ;
649  if (pdf->isSelectedComp()) {
650  value += pdfVal*_coefCache[i]/snormVal ;
651  cxcoutD(Eval) << "RooAddModel::evaluate(" << GetName() << ") value += ["
652  << pdf->GetName() << "] " << pdfVal << " * " << _coefCache[i] << " / " << snormVal << endl ;
653  }
654  }
655  i++ ;
656  }
657 
658  return value ;
659 }
660 
661 
662 
663 ////////////////////////////////////////////////////////////////////////////////
664 /// Reset error counter to given value, limiting the number
665 /// of future error messages for this pdf to 'resetValue'
666 
668 {
669  RooAbsPdf::resetErrorCounters(resetValue) ;
670  _coefErrCount = resetValue ;
671 }
672 
673 
674 
675 ////////////////////////////////////////////////////////////////////////////////
676 /// Check if PDF is valid for given normalization set.
677 /// Coeffient and PDF must be non-overlapping, but pdf-coefficient
678 /// pairs may overlap each other
679 
681 {
682  Bool_t ret(kFALSE) ;
683 
684  _pdfIter->Reset() ;
685  _coefIter->Reset() ;
686  RooAbsReal* coef ;
687  RooAbsReal* pdf ;
688  while((coef=(RooAbsReal*)_coefIter->Next())) {
689  pdf = (RooAbsReal*)_pdfIter->Next() ;
690  if (pdf->observableOverlaps(nset,*coef)) {
691  coutE(InputArguments) << "RooAddModel::checkObservables(" << GetName() << "): ERROR: coefficient " << coef->GetName()
692  << " and PDF " << pdf->GetName() << " have one or more dependents in common" << endl ;
693  ret = kTRUE ;
694  }
695  }
696 
697  return ret ;
698 }
699 
700 
701 
702 ////////////////////////////////////////////////////////////////////////////////
703 
705  const RooArgSet* normSet, const char* rangeName) const
706 {
707  if (_forceNumInt) return 0 ;
708 
709  // Declare that we can analytically integrate all requested observables
710  analVars.add(allVars) ;
711 
712  // Retrieve (or create) the required component integral list
713  Int_t code ;
714  RooArgList *cilist ;
715  getCompIntList(normSet,&allVars,cilist,code,rangeName) ;
716 
717  return code+1 ;
718 
719 }
720 
721 
722 
723 ////////////////////////////////////////////////////////////////////////////////
724 /// Check if this configuration was created before
725 
726 void RooAddModel::getCompIntList(const RooArgSet* nset, const RooArgSet* iset, pRooArgList& compIntList, Int_t& code, const char* isetRangeName) const
727 {
728  Int_t sterileIdx(-1) ;
729 
730  IntCacheElem* cache = (IntCacheElem*) _intCacheMgr.getObj(nset,iset,&sterileIdx,RooNameReg::ptr(isetRangeName)) ;
731  if (cache) {
732  code = _intCacheMgr.lastIndex() ;
733  compIntList = &cache->_intList ;
734 
735  return ;
736  }
737 
738  // Create containers for partial integral components to be generated
739  cache = new IntCacheElem ;
740 
741  // Fill Cache
742  _pdfIter->Reset() ;
744  while ((model=(RooResolutionModel*)_pdfIter->Next())) {
745  RooAbsReal* intPdf = model->createIntegral(*iset,nset,0,isetRangeName) ;
746  cache->_intList.addOwned(*intPdf) ;
747  }
748 
749  // Store the partial integral list and return the assigned code ;
750  code = _intCacheMgr.setObj(nset,iset,(RooAbsCacheElement*)cache,RooNameReg::ptr(isetRangeName)) ;
751 
752  // Fill references to be returned
753  compIntList = &cache->_intList ;
754 }
755 
756 
757 
758 ////////////////////////////////////////////////////////////////////////////////
759 /// Return analytical integral defined by given scenario code
760 
761 Double_t RooAddModel::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
762 {
763  // No integration scenario
764  if (code==0) {
765  return getVal(normSet) ;
766  }
767 
768  // Partial integration scenarios
770 
771  RooArgList* compIntList ;
772 
773  // If cache has been sterilized, revive this slot
774  if (cache==0) {
775  RooArgSet* vars = getParameters(RooArgSet()) ;
776  RooArgSet* nset = _intCacheMgr.nameSet1ByIndex(code-1)->select(*vars) ;
777  RooArgSet* iset = _intCacheMgr.nameSet2ByIndex(code-1)->select(*vars) ;
778 
779  Int_t code2(-1) ;
780  getCompIntList(nset,iset,compIntList,code2,rangeName) ;
781 
782  delete vars ;
783  delete nset ;
784  delete iset ;
785  } else {
786 
787  compIntList = &cache->_intList ;
788 
789  }
790 
791  // Calculate the current value
792  const RooArgSet* nset = _normSet ;
793  CacheElem* pcache = getProjCache(nset) ;
794 
795  updateCoefficients(*pcache,nset) ;
796 
797  // Do running sum of coef/pdf pairs, calculate lastCoef.
798  TIterator* compIntIter = compIntList->createIterator() ;
799  _coefIter->Reset() ;
800  RooAbsReal* pdfInt ;
801 
802  Double_t snormVal ;
803  Double_t value(0) ;
804  Int_t i(0) ;
805  while((pdfInt = (RooAbsReal*)compIntIter->Next())) {
806  if (_coefCache[i]!=0.) {
807  snormVal = nset ? ((RooAbsReal*)pcache->_suppNormList.at(i))->getVal() : 1.0 ;
808  Double_t intVal = pdfInt->getVal(nset) ;
809  value += intVal*_coefCache[i]/snormVal ;
810  cxcoutD(Eval) << "RooAddModel::evaluate(" << GetName() << ") value += ["
811  << pdfInt->GetName() << "] " << intVal << " * " << _coefCache[i] << " / " << snormVal << endl ;
812  }
813  i++ ;
814  }
815 
816  delete compIntIter ;
817 
818  return value ;
819 
820 }
821 
822 
823 
824 ////////////////////////////////////////////////////////////////////////////////
825 /// Return the number of expected events, which is either the sum of all coefficients
826 /// or the sum of the components extended terms
827 
829 {
830  Double_t expectedTotal(0.0);
831  RooAbsPdf* pdf ;
832 
833  if (_allExtendable) {
834 
835  // Sum of the extended terms
836  _pdfIter->Reset() ;
837  while((pdf = (RooAbsPdf*)_pdfIter->Next())) {
838  expectedTotal += pdf->expectedEvents(nset) ;
839  }
840 
841  } else {
842 
843  // Sum the coefficients
844  _coefIter->Reset() ;
845  RooAbsReal* coef ;
846  while((coef=(RooAbsReal*)_coefIter->Next())) {
847  expectedTotal += coef->getVal() ;
848  }
849  }
850 
851  return expectedTotal;
852 }
853 
854 
855 
856 ////////////////////////////////////////////////////////////////////////////////
857 /// Interface function used by test statistics to freeze choice of observables
858 /// for interpretation of fraction coefficients
859 
861 {
862  if (!force && _refCoefNorm.getSize()!=0) {
863  return ;
864  }
865 
866  if (!depSet) {
868  return ;
869  }
870 
871  RooArgSet* myDepSet = getObservables(depSet) ;
872  fixCoefNormalization(*myDepSet) ;
873  delete myDepSet ;
874 }
875 
876 
877 
878 ////////////////////////////////////////////////////////////////////////////////
879 /// Interface function used by test statistics to freeze choice of range
880 /// for interpretation of fraction coefficients
881 
882 void RooAddModel::selectNormalizationRange(const char* rangeName, Bool_t force)
883 {
884  if (!force && _refCoefRangeName) {
885  return ;
886  }
887 
888  fixCoefRange(rangeName) ;
889 }
890 
891 
892 
893 ////////////////////////////////////////////////////////////////////////////////
894 /// Return specialized context to efficiently generate toy events from RooAddPdfs
895 
897  const RooArgSet* auxProto, Bool_t verbose) const
898 {
899  return new RooAddGenContext(*this,vars,prototype,auxProto,verbose) ;
900 }
901 
902 
903 
904 ////////////////////////////////////////////////////////////////////////////////
905 /// Direct generation is safe if all components say so
906 
908 {
909  _pdfIter->Reset() ;
910  RooAbsPdf* pdf ;
911  while((pdf=(RooAbsPdf*)_pdfIter->Next())) {
912  if (!pdf->isDirectGenSafe(arg)) {
913  return kFALSE ;
914  }
915  }
916  return kTRUE ;
917 }
918 
919 
920 
921 ////////////////////////////////////////////////////////////////////////////////
922 /// Return pseud-code that indicates if all components can do internal generation (1) or not (0)
923 
924 Int_t RooAddModel::getGenerator(const RooArgSet& directVars, RooArgSet &/*generateVars*/, Bool_t /*staticInitOK*/) const
925 {
926  _pdfIter->Reset() ;
927  RooAbsPdf* pdf ;
928  while((pdf=(RooAbsPdf*)_pdfIter->Next())) {
929  RooArgSet tmp ;
930  if (pdf->getGenerator(directVars,tmp)==0) {
931  return 0 ;
932  }
933  }
934  return 1 ;
935 }
936 
937 
938 
939 
940 ////////////////////////////////////////////////////////////////////////////////
941 /// This function should never be called as RooAddModel implements a custom generator context
942 
944 {
945  assert(0) ;
946 }
947 
948 
949 
950 
951 ////////////////////////////////////////////////////////////////////////////////
952 /// List all RooAbsArg derived contents in this cache element
953 
955 {
956  RooArgList allNodes;
957  allNodes.add(_projList) ;
958  allNodes.add(_suppProjList) ;
959  allNodes.add(_refRangeProjList) ;
960  allNodes.add(_rangeProjList) ;
961 
962  return allNodes ;
963 }
964 
965 
966 
967 ////////////////////////////////////////////////////////////////////////////////
968 /// List all RooAbsArg derived contents in this cache element
969 
971 {
972  RooArgList allNodes(_intList) ;
973  return allNodes ;
974 }
975 
976 
977 ////////////////////////////////////////////////////////////////////////////////
978 /// Customized printing of arguments of a RooAddModel to more intuitively reflect the contents of the
979 /// product operator construction
980 
981 void RooAddModel::printMetaArgs(ostream& os) const
982 {
983  _pdfIter->Reset() ;
984  _coefIter->Reset() ;
985 
986  Bool_t first(kTRUE) ;
987 
988  os << "(" ;
989  RooAbsArg* coef, *pdf ;
990  while((coef=(RooAbsArg*)_coefIter->Next())) {
991  if (!first) {
992  os << " + " ;
993  } else {
994  first = kFALSE ;
995  }
996  pdf=(RooAbsArg*)_pdfIter->Next() ;
997  os << coef->GetName() << " * " << pdf->GetName() ;
998  }
999  pdf = (RooAbsArg*) _pdfIter->Next() ;
1000  if (pdf) {
1001  os << " + [%] * " << pdf->GetName() ;
1002  }
1003  os << ") " ;
1004 }
1005 
void setAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:266
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
TIterator * createIterator(Bool_t dir=kIterForward) const
Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Direct generation is safe if all components say so.
#define coutE(a)
Definition: RooMsgService.h:34
static void globalSelectComp(Bool_t flag)
Global switch controlling the activation of the selectComp() functionality.
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer, which is interpreted as an OR of &#39;enum ContentsOptions&#39; values and in the style given by &#39;enum StyleOption&#39;.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
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
TIterator * _pdfIter
List of supplemental normalization factors.
Definition: RooAddModel.h:136
virtual void Reset()=0
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...
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Definition: RooAbsArg.h:194
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
Variant of getAnalyticalIntegral that is also passed the normalization set that should be applied to ...
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...
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
#define coutI(a)
Definition: RooMsgService.h:31
virtual void removeAll()
Remove all argument inset using remove(const RooAbsArg&).
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
#define cxcoutD(a)
Definition: RooMsgService.h:79
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.cxx:135
Int_t _errorCount
Definition: RooAbsPdf.h:324
void fixCoefNormalization(const RooArgSet &refCoefNorm)
By default the interpretation of the fraction coefficients is performed in the contextual choice of o...
RooRealVar & convVar() const
Return the convolution variable of the resolution model.
T * getObjByIndex(Int_t index) const
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
Bool_t equals(const RooAbsCollection &otherColl) const
Check if this and other collection have identically named contents.
std::set< std::string > _boolAttrib
Definition: RooAbsArg.h:523
virtual RooResolutionModel * convolution(RooFormulaVar *basis, RooAbsArg *owner) const
Instantiate a clone of this resolution model representing a convolution with given basis function...
STL namespace.
#define coutW(a)
Definition: RooMsgService.h:33
RooArgList _rangeProjList
Definition: RooAddModel.h:111
#define ccoutI(a)
Definition: RooMsgService.h:38
Iterator abstract base class.
Definition: TIterator.h:30
void updateCoefficients(CacheElem &cache, const RooArgSet *nset) const
Update the coefficient values in the given cache element: calculate new remainder fraction...
void getCompIntList(const RooArgSet *nset, const RooArgSet *iset, pRooArgList &compIntList, Int_t &code, const char *isetRangeName) const
Check if this configuration was created before.
RooArgList _suppProjList
Definition: RooAddModel.h:109
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...
RooAbsArg * findServer(const char *name) const
Definition: RooAbsArg.h:122
const RooNameSet * nameSet1ByIndex(Int_t index) const
TNamed * _refCoefRangeName
Reference observable set for coefficient interpretation.
Definition: RooAddModel.h:96
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string &#39;value&#39; to this object under key &#39;key&#39;.
Definition: RooAbsArg.cxx:298
RooListProxy _coefList
Definition: RooAddModel.h:134
friend class CacheElem
Definition: RooAbsPdf.h:314
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:817
void printCompactTree(const char *indent="", const char *fileName=0, const char *namePat=0, RooAbsArg *client=0)
Print tree structure of expression tree on stdout, or to file if filename is specified.
Definition: RooAbsArg.cxx:1783
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
virtual Int_t basisCode(const char *name) const
Return code for basis function representing by &#39;name&#39; string.
RooArgList _suppNormList
Definition: RooAddModel.h:104
RooArgSet _ownedComps
Coefficient error counter.
Definition: RooAddModel.h:144
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
virtual ~RooAddModel()
Destructor.
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooAbsArg.h:227
friend class RooArgSet
Definition: RooAbsArg.h:469
virtual void Print(Option_t *options=0) const
This method must be overridden when a class wants to print itself.
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
RooListProxy _pdfList
Registry of component analytical integration codes.
Definition: RooAddModel.h:133
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events from this p.d.f for use in extended likelihood calculations.
Definition: RooAbsPdf.cxx:2930
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
RooFormulaVar * _basis
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
CacheElem * getProjCache(const RooArgSet *nset, const RooArgSet *iset=0, const char *rangeName=0) const
Retrieve cache element with for calculation of p.d.f value with normalization set nset and integrated...
virtual void changeBasis(RooFormulaVar *basis)
Change the basis function we convolute with.
virtual RooArgList containedArgs(Action)
List all RooAbsArg derived contents in this cache element.
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
#define ccoutE(a)
Definition: RooMsgService.h:41
Bool_t _forceNumInt
Definition: RooAbsReal.h:392
static Int_t _verboseEval
Definition: RooAbsPdf.h:295
RooAbsArg * at(Int_t idx) const
Definition: RooArgList.h:84
RooArgList _refRangeProjList
Definition: RooAddModel.h:110
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:125
unsigned int r1[N_CITIES]
Definition: simanTSP.cxx:321
TIterator * _coefIter
Iterator over PDF list.
Definition: RooAddModel.h:137
friend class RooAddGenContext
Definition: RooAddModel.h:88
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
RooObjCacheManager _projCacheMgr
Definition: RooAddModel.h:116
bool verbose
Int_t lastIndex() const
RooAbsArg * absArg() const
Definition: RooArgProxy.h:37
virtual Bool_t addOwned(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling addOwned() for each element in the source...
Definition: RooArgSet.h:90
OperMode operMode() const
Definition: RooAbsArg.h:399
RooObjCacheManager _intCacheMgr
Definition: RooAddModel.h:129
Int_t _coefErrCount
Definition: RooAddModel.h:142
virtual void resetErrorCounters(Int_t resetValue=10)
Reset error counter to given value, limiting the number of future error messages for this pdf to &#39;res...
Definition: RooAbsPdf.cxx:568
Bool_t _allExtendable
Definition: RooAddModel.h:140
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:29
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Return pseud-code that indicates if all components can do internal generation (1) or not (0) ...
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooAddModel to more intuitively reflect the contents of the pro...
void generateEvent(Int_t code)
This function should never be called as RooAddModel implements a custom generator context...
virtual Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
Definition: RooAbsPdf.cxx:2136
const Bool_t kFALSE
Definition: RtypesCore.h:92
virtual RooArgList containedArgs(Action)
List all RooAbsArg derived contents in this cache element.
const RooNameSet * nameSet2ByIndex(Int_t index) const
#define ClassImp(name)
Definition: Rtypes.h:336
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
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooAbsPdf.cxx:2101
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:560
std::map< std::string, std::string > _stringAttrib
Definition: RooAbsArg.h:524
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.
Double_t evaluate() const
Calculate the current value.
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
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:189
RooAICRegistry _codeReg
Definition: RooAddModel.h:131
Transiet cache with transformed values of coefficients.
Definition: RooAddModel.h:102
Bool_t _haveLastCoef
Iterator over coefficient list.
Definition: RooAddModel.h:139
Double_t * _coefCache
Definition: RooAddModel.h:99
#define dologD(a)
Definition: RooMsgService.h:63
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
friend class RooRealIntegral
Definition: RooAbsPdf.h:294
virtual TObject * Next()=0
static RooConstVar & value(Double_t value)
Return a constant value object with given value.
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Change cache operation mode to given mode.
Definition: RooAbsArg.cxx:1754
Bool_t _projectCoefs
Reference range name for coefficient interpreation.
Definition: RooAddModel.h:98
Bool_t isSelectedComp() const
If true, the current pdf is a selected component (for use in plotting)
RooArgSet * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:302
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Definition: first.py:1
virtual Int_t basisCode(const char *name) const =0
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
void fixCoefRange(const char *rangeName)
By default the interpretation of the fraction coefficients is performed in the default range...
#define ccoutD(a)
Definition: RooMsgService.h:37
virtual RooResolutionModel * convolution(RooFormulaVar *basis, RooAbsArg *owner) const
Instantiate a clone of this resolution model representing a convolution with given basis function...
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooAbsArg * getParameter(const char *name) const
Definition: RooFormulaVar.h:39
virtual void resetErrorCounters(Int_t resetValue=10)
Reset error counter to given value, limiting the number of future error messages for this pdf to &#39;res...
RooSetProxy _refCoefNorm
Definition: RooAddModel.h:95
unsigned int r2[N_CITIES]
Definition: simanTSP.cxx:322
virtual Bool_t checkObservables(const RooArgSet *nset) const
Check if PDF is valid for given normalization set.
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:48
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Overloaded RooArgSet::add() method inserts &#39;var&#39; into set and registers &#39;var&#39; as server to owner with...
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return analytical integral defined by given scenario code.