Logo ROOT  
Reference Guide
RooAbsAnaConvPdf.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 RooAbsAnaConvPdf
19 /// \ingroup Roofitcore
20 ///
21 /// RooAbsAnaConvPdf is the base class for PDFs that represent a
22 /// physics model that can be analytically convolved with a resolution model.
23 ///
24 /// To achieve factorization between the physics model and the resolution
25 /// model, each physics model must be able to be written in the form
26 /// \f[
27 /// \mathrm{Phys}(x, \bar{a}, \bar{b}) = \sum_k \mathrm{coef}_k(\bar{a}) * \mathrm{basis}_k(x,\bar{b})
28 /// \f]
29 ///
30 /// where \f$ \mathrm{basis}_k \f$ are a limited number of functions in terms of the variable
31 /// to be convoluted, and \f$ \mathrm{coef}_k \f$ are coefficients independent of the convolution
32 /// variable.
33 ///
34 /// Classes derived from RooResolutionModel implement
35 /// \f[
36 /// R_k(x,\bar{b},\bar{c}) = \int \mathrm{basis}_k(x', \bar{b}) \cdot \mathrm{resModel}(x-x',\bar{c}) \; \mathrm{d}x',
37 /// \f]
38 ///
39 /// which RooAbsAnaConvPdf uses to construct the pdf for [ Phys (x) R ] :
40 /// \f[
41 /// \mathrm{PDF}(x,\bar{a},\bar{b},\bar{c}) = \sum_k \mathrm{coef}_k(\bar{a}) * R_k(x,\bar{b},\bar{c})
42 /// \f]
43 ///
44 /// A minimal implementation of a RooAbsAnaConvPdf physics model consists of
45 ///
46 /// - A constructor that declares the required basis functions using the declareBasis() method.
47 /// The declareBasis() function assigns a unique identifier code to each declare basis
48 ///
49 /// - An implementation of `coefficient(Int_t code)` returning the coefficient value for each
50 /// declared basis function
51 ///
52 /// Optionally, analytical integrals can be provided for the coefficient functions. The
53 /// interface for this is quite similar to that for integrals of regular PDFs. Two functions,
54 /// \code{.cpp}
55 /// Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
56 /// Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char* rangeName) const
57 /// \endcode
58 ///
59 /// advertise the coefficient integration capabilities and implement them respectively.
60 /// Please see RooAbsPdf for additional details. Advertised analytical integrals must be
61 /// valid for all coefficients.
62 
63 #include "RooAbsAnaConvPdf.h"
64 
65 #include "RooFit.h"
66 #include "RooMsgService.h"
67 #include "Riostream.h"
68 #include "RooResolutionModel.h"
69 #include "RooRealVar.h"
70 #include "RooFormulaVar.h"
71 #include "RooConvGenContext.h"
72 #include "RooGenContext.h"
73 #include "RooTruthModel.h"
74 #include "RooConvCoefVar.h"
75 #include "RooNameReg.h"
76 
77 using namespace std;
78 
80 
81 
82 ////////////////////////////////////////////////////////////////////////////////
83 /// Default constructor, required for persistence
84 
86  _isCopy(kFALSE)
87 {
88 }
89 
90 
91 
92 ////////////////////////////////////////////////////////////////////////////////
93 /// Constructor. The supplied resolution model must be constructed with the same
94 /// convoluted variable as this physics model ('convVar')
95 
96 RooAbsAnaConvPdf::RooAbsAnaConvPdf(const char *name, const char *title,
97  const RooResolutionModel& model, RooRealVar& cVar) :
98  RooAbsPdf(name,title), _isCopy(kFALSE),
99  _model("!model","Original resolution model",this,(RooResolutionModel&)model,kFALSE,kFALSE),
100  _convVar("!convVar","Convolution variable",this,cVar,kFALSE,kFALSE),
101  _convSet("!convSet","Set of resModel X basisFunc convolutions",this),
102  _coefNormMgr(this,10),
103  _codeReg(10)
104 {
105  _model.absArg()->setAttribute("NOCacheAndTrack") ;
106 }
107 
108 
109 
110 ////////////////////////////////////////////////////////////////////////////////
111 
113  RooAbsPdf(other,name), _isCopy(kTRUE),
114  _model("!model",this,other._model),
115  _convVar("!convVar",this,other._convVar),
116  _convSet("!convSet",this,other._convSet),
117  _coefNormMgr(other._coefNormMgr,this),
118  _codeReg(other._codeReg)
119 {
120  // Copy constructor
121  if (_model.absArg()) {
122  _model.absArg()->setAttribute("NOCacheAndTrack") ;
123  }
125 }
126 
127 
128 
129 ////////////////////////////////////////////////////////////////////////////////
130 /// Destructor
131 
133 {
134  if (!_isCopy) {
135  std::vector<RooAbsArg*> tmp(_convSet.begin(), _convSet.end());
136 
137  for (auto arg : tmp) {
138  _convSet.remove(*arg) ;
139  delete arg ;
140  }
141  }
142 
143 }
144 
145 
146 ////////////////////////////////////////////////////////////////////////////////
147 /// Declare a basis function for use in this physics model. The string expression
148 /// must be a valid RooFormulVar expression representing the basis function, referring
149 /// to the convolution variable as '@0', and any additional parameters (supplied in
150 /// 'params' as '@1','@2' etc.
151 ///
152 /// The return value is a unique identifier code, that will be passed to coefficient()
153 /// to identify the basis function for which the coefficient is requested. If the
154 /// resolution model used does not support the declared basis function, code -1 is
155 /// returned.
156 ///
157 
158 Int_t RooAbsAnaConvPdf::declareBasis(const char* expression, const RooArgList& params)
159 {
160  // Sanity check
161  if (_isCopy) {
162  coutE(InputArguments) << "RooAbsAnaConvPdf::declareBasis(" << GetName() << "): ERROR attempt to "
163  << " declare basis functions in a copied RooAbsAnaConvPdf" << endl ;
164  return -1 ;
165  }
166 
167  // Resolution model must support declared basis
168  if (!((RooResolutionModel*)_model.absArg())->isBasisSupported(expression)) {
169  coutE(InputArguments) << "RooAbsAnaConvPdf::declareBasis(" << GetName() << "): resolution model "
170  << _model.absArg()->GetName()
171  << " doesn't support basis function " << expression << endl ;
172  return -1 ;
173  }
174 
175  // Instantiate basis function
176  RooArgList basisArgs(_convVar.arg()) ;
177  basisArgs.add(params) ;
178 
179  TString basisName(expression) ;
180  for (const auto arg : basisArgs) {
181  basisName.Append("_") ;
182  basisName.Append(arg->GetName()) ;
183  }
184 
185  RooFormulaVar* basisFunc = new RooFormulaVar(basisName, expression, basisArgs);
186  basisFunc->setAttribute("RooWorkspace::Recycle") ;
187  basisFunc->setAttribute("NOCacheAndTrack") ;
188  basisFunc->setOperMode(operMode()) ;
189  _basisList.addOwned(*basisFunc) ;
190 
191  // Instantiate resModel x basisFunc convolution
192  RooAbsReal* conv = ((RooResolutionModel*)_model.absArg())->convolution(basisFunc,this) ;
193  if (!conv) {
194  coutE(InputArguments) << "RooAbsAnaConvPdf::declareBasis(" << GetName() << "): unable to construct convolution with basis function '"
195  << expression << "'" << endl ;
196  return -1 ;
197  }
198  _convSet.add(*conv) ;
199 
200  return _convSet.index(conv) ;
201 }
202 
203 
204 
205 ////////////////////////////////////////////////////////////////////////////////
206 /// Change the current resolution model to newModel
207 
209 {
210  RooArgList newConvSet ;
211  Bool_t allOK(kTRUE) ;
212  for (auto convArg : _convSet) {
213  auto conv = static_cast<RooResolutionModel*>(convArg);
214 
215  // Build new resolution model
216  RooResolutionModel* newConv = newModel.convolution((RooFormulaVar*)&conv->basis(),this) ;
217  if (!newConvSet.add(*newConv)) {
218  allOK = kFALSE ;
219  break ;
220  }
221  }
222 
223  // Check if all convolutions were successfully built
224  if (!allOK) {
225  // Delete new basis functions created sofar
226  std::for_each(newConvSet.begin(), newConvSet.end(), [](RooAbsArg* arg){delete arg;});
227 
228  return kTRUE ;
229  }
230 
231  // Replace old convolutions with new set
232  _convSet.removeAll() ;
233  _convSet.addOwned(newConvSet) ;
234 
235  // Update server link by hand, since _model.setArg() below will not do this
237 
238  _model.setArg((RooResolutionModel&)newModel) ;
239  return kFALSE ;
240 }
241 
242 
243 
244 
245 ////////////////////////////////////////////////////////////////////////////////
246 /// Create a generator context for this p.d.f. If both the p.d.f and the resolution model
247 /// support internal generation of the convolution observable on an infinite domain,
248 /// deploy a specialized convolution generator context, which generates the physics distribution
249 /// and the smearing separately, adding them a posteriori. If this is not possible return
250 /// a (slower) generic generation context that uses accept/reject sampling
251 
253  const RooArgSet* auxProto, Bool_t verbose) const
254 {
255  // Check if the resolution model specifies a special context to be used.
256  RooResolutionModel* conv = dynamic_cast<RooResolutionModel*>(_model.absArg());
257  assert(conv);
258 
259  RooArgSet* modelDep = _model.absArg()->getObservables(&vars) ;
260  modelDep->remove(*convVar(),kTRUE,kTRUE) ;
261  Int_t numAddDep = modelDep->getSize() ;
262  delete modelDep ;
263 
264  // Check if physics PDF and resolution model can both directly generate the convolution variable
265  RooArgSet dummy ;
266  Bool_t pdfCanDir = (getGenerator(*convVar(),dummy) != 0) ;
267  Bool_t resCanDir = conv && (conv->getGenerator(*convVar(),dummy)!=0) && conv->isDirectGenSafe(*convVar()) ;
268 
269  if (numAddDep>0 || !pdfCanDir || !resCanDir) {
270  // Any resolution model with more dependents than the convolution variable
271  // or pdf or resmodel do not support direct generation
272  string reason ;
273  if (numAddDep>0) reason += "Resolution model has more observables than the convolution variable. " ;
274  if (!pdfCanDir) reason += "PDF does not support internal generation of convolution observable. " ;
275  if (!resCanDir) reason += "Resolution model does not support internal generation of convolution observable. " ;
276 
277  coutI(Generation) << "RooAbsAnaConvPdf::genContext(" << GetName() << ") Using regular accept/reject generator for convolution p.d.f because: " << reason.c_str() << endl ;
278  return new RooGenContext(*this,vars,prototype,auxProto,verbose) ;
279  }
280 
281  RooAbsGenContext* context = conv->modelGenContext(*this, vars, prototype, auxProto, verbose);
282  if (context) return context;
283 
284  // Any other resolution model: use specialized generator context
285  return new RooConvGenContext(*this,vars,prototype,auxProto,verbose) ;
286 }
287 
288 
289 
290 ////////////////////////////////////////////////////////////////////////////////
291 /// Return true if it is safe to generate the convolution observable
292 /// from the internal generator (this is the case if the chosen resolution
293 /// model is the truth model)
294 
296 {
297 
298  // All direct generation of convolution arg if model is truth model
299  if (!TString(_convVar.absArg()->GetName()).CompareTo(arg.GetName()) &&
300  dynamic_cast<RooTruthModel*>(_model.absArg())) {
301  return kTRUE ;
302  }
303 
304  return RooAbsPdf::isDirectGenSafe(arg) ;
305 }
306 
307 
308 
309 ////////////////////////////////////////////////////////////////////////////////
310 /// Return a pointer to the convolution variable instance used in the resolution model
311 
313 {
315  if (!conv) return 0 ;
316  return &conv->convVar() ;
317 }
318 
319 
320 
321 ////////////////////////////////////////////////////////////////////////////////
322 /// Calculate the current unnormalized value of the PDF
323 ///
324 /// PDF = sum_k coef_k * [ basis_k (x) ResModel ]
325 ///
326 
328 {
329  Double_t result(0) ;
330 
331  Int_t index(0) ;
332  for (auto convArg : _convSet) {
333  auto conv = static_cast<RooAbsPdf*>(convArg);
334  Double_t coef = coefficient(index++) ;
335  if (coef!=0.) {
336  Double_t c = conv->getVal(0) ;
337  Double_t r = coef ;
338  cxcoutD(Eval) << "RooAbsAnaConvPdf::evaluate(" << GetName() << ") val += coef*conv [" << index-1 << "/"
339  << _convSet.getSize() << "] coef = " << r << " conv = " << c << endl ;
340  result += conv->getVal(0)*coef ;
341  } else {
342  cxcoutD(Eval) << "RooAbsAnaConvPdf::evaluate(" << GetName() << ") [" << index-1 << "/" << _convSet.getSize() << "] coef = 0" << endl ;
343  }
344  }
345 
346  return result ;
347 }
348 
349 
350 
351 ////////////////////////////////////////////////////////////////////////////////
352 /// Advertise capability to perform (analytical) integrals
353 /// internally. For a given integration request over allVars while
354 /// normalized over normSet2 and in range 'rangeName', returns
355 /// largest subset that can be performed internally in analVars
356 /// Return code is unique integer code identifying integration scenario
357 /// to be passed to analyticalIntegralWN() to calculate requeste integral
358 ///
359 /// Class RooAbsAnaConv defers analytical integration request to
360 /// resolution model and/or coefficient implementations and
361 /// aggregates results into composite configuration with a unique
362 /// code assigned by RooAICRegistry
363 
365  RooArgSet& analVars, const RooArgSet* normSet2, const char* /*rangeName*/) const
366 {
367  // Handle trivial no-integration scenario
368  if (allVars.getSize()==0) return 0 ;
369 
370  if (_forceNumInt) return 0 ;
371 
372  // Select subset of allVars that are actual dependents
373  RooArgSet* allDeps = getObservables(allVars) ;
374  RooArgSet* normSet = normSet2 ? getObservables(normSet2) : 0 ;
375 
376  RooAbsArg *arg ;
377  RooResolutionModel *conv ;
378 
379  RooArgSet* intSetAll = new RooArgSet(*allDeps,"intSetAll") ;
380 
381  // Split intSetAll in coef/conv parts
382  RooArgSet* intCoefSet = new RooArgSet("intCoefSet") ;
383  RooArgSet* intConvSet = new RooArgSet("intConvSet") ;
384  TIterator* varIter = intSetAll->createIterator() ;
385  TIterator* convIter = _convSet.createIterator() ;
386 
387  while(((arg=(RooAbsArg*) varIter->Next()))) {
388  Bool_t ok(kTRUE) ;
389  convIter->Reset() ;
390  while(((conv=(RooResolutionModel*) convIter->Next()))) {
391  if (conv->dependsOn(*arg)) ok=kFALSE ;
392  }
393 
394  if (ok) {
395  intCoefSet->add(*arg) ;
396  } else {
397  intConvSet->add(*arg) ;
398  }
399 
400  }
401  delete varIter ;
402 
403 
404  // Split normSetAll in coef/conv parts
405  RooArgSet* normCoefSet = new RooArgSet("normCoefSet") ;
406  RooArgSet* normConvSet = new RooArgSet("normConvSet") ;
407  RooArgSet* normSetAll = normSet ? (new RooArgSet(*normSet,"normSetAll")) : 0 ;
408  if (normSetAll) {
409  varIter = normSetAll->createIterator() ;
410  while(((arg=(RooAbsArg*) varIter->Next()))) {
411  Bool_t ok(kTRUE) ;
412  convIter->Reset() ;
413  while(((conv=(RooResolutionModel*) convIter->Next()))) {
414  if (conv->dependsOn(*arg)) ok=kFALSE ;
415  }
416 
417  if (ok) {
418  normCoefSet->add(*arg) ;
419  } else {
420  normConvSet->add(*arg) ;
421  }
422 
423  }
424  delete varIter ;
425  }
426  delete convIter ;
427 
428  if (intCoefSet->getSize()==0) {
429  delete intCoefSet ; intCoefSet=0 ;
430  }
431  if (intConvSet->getSize()==0) {
432  delete intConvSet ; intConvSet=0 ;
433  }
434  if (normCoefSet->getSize()==0) {
435  delete normCoefSet ; normCoefSet=0 ;
436  }
437  if (normConvSet->getSize()==0) {
438  delete normConvSet ; normConvSet=0 ;
439  }
440 
441 
442 
443  // Store integration configuration in registry
444  Int_t masterCode(0) ;
445  std::vector<Int_t> tmp(1, 0) ;
446 
447  masterCode = _codeReg.store(tmp, intCoefSet, intConvSet, normCoefSet, normConvSet) + 1 ; // takes ownership of all sets
448 
449  analVars.add(*allDeps) ;
450  delete allDeps ;
451  if (normSet) delete normSet ;
452  if (normSetAll) delete normSetAll ;
453  delete intSetAll ;
454 
455 // cout << this << "---> masterCode = " << masterCode << endl ;
456 
457  return masterCode ;
458 }
459 
460 
461 
462 
463 ////////////////////////////////////////////////////////////////////////////////
464 /// Return analytical integral defined by given code, which is returned
465 /// by getAnalyticalIntegralWN()
466 ///
467 /// For unnormalized integrals the returned value is
468 /// \f[
469 /// \mathrm{PDF} = \sum_k \int \mathrm{coef}_k \; \mathrm{d}\bar{x}
470 /// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}\bar{y},
471 /// \f]
472 /// where \f$ \bar{x} \f$ is the set of coefficient dependents to be integrated,
473 /// and \f$ \bar{y} \f$ the set of basis function dependents to be integrated.
474 ///
475 /// For normalized integrals this becomes
476 /// \f[
477 /// \mathrm{PDF} = \frac{\sum_k \int \mathrm{coef}_k \; \mathrm{d}x
478 /// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}y}
479 /// {\sum_k \int \mathrm{coef}_k \; \mathrm{d}v
480 /// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}w},
481 /// \f]
482 /// where
483 /// * \f$ x \f$ is the set of coefficient dependents to be integrated,
484 /// * \f$ y \f$ the set of basis function dependents to be integrated,
485 /// * \f$ v \f$ is the set of coefficient dependents over which is normalized and
486 /// * \f$ w \f$ is the set of basis function dependents over which is normalized.
487 ///
488 /// Set \f$ x \f$ must be contained in \f$ v \f$ and set \f$ y \f$ must be contained in \f$ w \f$.
489 ///
490 
491 Double_t RooAbsAnaConvPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
492 {
493  // WVE needs adaptation to handle new rangeName feature
494 
495  // Handle trivial passthrough scenario
496  if (code==0) return getVal(normSet) ;
497 
498  // Unpack master code
499  RooArgSet *intCoefSet, *intConvSet, *normCoefSet, *normConvSet ;
500  _codeReg.retrieve(code-1,intCoefSet,intConvSet,normCoefSet,normConvSet) ;
501 
502  Int_t index(0) ;
503  Double_t answer(0) ;
504 
505  if (normCoefSet==0&&normConvSet==0) {
506 
507  // Integral over unnormalized function
508  Double_t integral(0) ;
509  const TNamed *_rangeName = RooNameReg::ptr(rangeName);
510  for (auto convArg : _convSet) {
511  auto conv = static_cast<RooResolutionModel*>(convArg);
512  Double_t coef = getCoefNorm(index++,intCoefSet,_rangeName) ;
513  //cout << "coefInt[" << index << "] = " << coef << " " ; intCoefSet->Print("1") ;
514  if (coef!=0) {
515  integral += coef*(_rangeName ? conv->getNormObj(0,intConvSet,_rangeName)->getVal() : conv->getNorm(intConvSet) ) ;
516  cxcoutD(Eval) << "RooAbsAnaConv::aiWN(" << GetName() << ") [" << index-1 << "] integral += " << conv->getNorm(intConvSet) << endl ;
517  }
518 
519  }
520  answer = integral ;
521 
522  } else {
523 
524  // Integral over normalized function
525  Double_t integral(0) ;
526  Double_t norm(0) ;
527  const TNamed *_rangeName = RooNameReg::ptr(rangeName);
528  for (auto convArg : _convSet) {
529  auto conv = static_cast<RooResolutionModel*>(convArg);
530 
531  Double_t coefInt = getCoefNorm(index,intCoefSet,_rangeName) ;
532  //cout << "coefInt[" << index << "] = " << coefInt << "*" << term << " " << (intCoefSet?*intCoefSet:RooArgSet()) << endl ;
533  if (coefInt!=0) {
534  Double_t term = (_rangeName ? conv->getNormObj(0,intConvSet,_rangeName)->getVal() : conv->getNorm(intConvSet) ) ;
535  integral += coefInt*term ;
536  }
537 
538  Double_t coefNorm = getCoefNorm(index,normCoefSet) ;
539  //cout << "coefNorm[" << index << "] = " << coefNorm << "*" << term << " " << (normCoefSet?*normCoefSet:RooArgSet()) << endl ;
540  if (coefNorm!=0) {
541  Double_t term = conv->getNorm(normConvSet) ;
542  norm += coefNorm*term ;
543  }
544 
545  index++ ;
546  }
547  answer = integral/norm ;
548  }
549 
550  return answer ;
551 }
552 
553 
554 
555 ////////////////////////////////////////////////////////////////////////////////
556 /// Default implementation of function advertising integration capabilities. The interface is
557 /// similar to that of getAnalyticalIntegral except that an integer code is added that
558 /// designates the coefficient number for which the integration capabilities are requested
559 ///
560 /// This default implementation advertises that no internal integrals are supported.
561 
562 Int_t RooAbsAnaConvPdf::getCoefAnalyticalIntegral(Int_t /* coef*/, RooArgSet& /*allVars*/, RooArgSet& /*analVars*/, const char* /*rangeName*/) const
563 {
564  return 0 ;
565 }
566 
567 
568 
569 ////////////////////////////////////////////////////////////////////////////////
570 /// Default implementation of function implementing advertised integrals. Only
571 /// the pass-through scenario (no integration) is implemented.
572 
573 Double_t RooAbsAnaConvPdf::coefAnalyticalIntegral(Int_t coef, Int_t code, const char* /*rangeName*/) const
574 {
575  if (code==0) return coefficient(coef) ;
576  coutE(InputArguments) << "RooAbsAnaConvPdf::coefAnalyticalIntegral(" << GetName() << ") ERROR: unrecognized integration code: " << code << endl ;
577  assert(0) ;
578  return 1 ;
579 }
580 
581 
582 
583 ////////////////////////////////////////////////////////////////////////////////
584 /// This function forces RooRealIntegral to offer all integration dependents
585 /// to RooAbsAnaConvPdf::getAnalyticalIntegralWN() for consideration for
586 /// internal integration, if RooRealIntegral considers this to be unsafe (e.g. due
587 /// to hidden Jacobian terms).
588 ///
589 /// RooAbsAnaConvPdf will not attempt to actually integrate all these dependents
590 /// but feed them to the resolution models integration interface, which will
591 /// make the final determination on how to integrate these dependents.
592 
594 {
595  return kTRUE ;
596 }
597 
598 
599 
600 ////////////////////////////////////////////////////////////////////////////////
601 /// Returns the normalization integral value of the coefficient with number coefIdx over normalization
602 /// set nset in range rangeName
603 
604 Double_t RooAbsAnaConvPdf::getCoefNorm(Int_t coefIdx, const RooArgSet* nset, const TNamed* rangeName) const
605 {
606  if (nset==0) return coefficient(coefIdx) ;
607 
608  CacheElem* cache = (CacheElem*) _coefNormMgr.getObj(nset,0,0,rangeName) ;
609  if (!cache) {
610 
611  cache = new CacheElem ;
612 
613  // Make list of coefficient normalizations
614  Int_t i ;
615  makeCoefVarList(cache->_coefVarList) ;
616 
617  for (i=0 ; i<cache->_coefVarList.getSize() ; i++) {
618  RooAbsReal* coefInt = static_cast<RooAbsReal&>(*cache->_coefVarList.at(i)).createIntegral(*nset,RooNameReg::str(rangeName)) ;
619  cache->_normList.addOwned(*coefInt) ;
620  }
621 
622  _coefNormMgr.setObj(nset,0,cache,rangeName) ;
623  }
624 
625  return ((RooAbsReal*)cache->_normList.at(coefIdx))->getVal() ;
626 }
627 
628 
629 
630 ////////////////////////////////////////////////////////////////////////////////
631 /// Build complete list of coefficient variables
632 
634 {
635  // Instantate a coefficient variables
636  for (Int_t i=0 ; i<_convSet.getSize() ; i++) {
637  RooArgSet* cvars = coefVars(i) ;
638  RooAbsReal* coefVar = new RooConvCoefVar(Form("%s_coefVar_%d",GetName(),i),"coefVar",*this,i,cvars) ;
639  varList.addOwned(*coefVar) ;
640  delete cvars ;
641  }
642 
643 }
644 
645 
646 ////////////////////////////////////////////////////////////////////////////////
647 /// Return set of parameters with are used exclusively by the coefficient functions
648 
650 {
651  RooArgSet* cVars = getParameters((RooArgSet*)0) ;
652  std::vector<RooAbsArg*> tmp;
653  for (auto arg : *cVars) {
654  for (auto convSetArg : _convSet) {
655  if (convSetArg->dependsOn(*arg)) {
656  tmp.push_back(arg);
657  }
658  }
659  }
660 
661  cVars->remove(tmp.begin(), tmp.end(), true, true);
662 
663  return cVars ;
664 }
665 
666 
667 
668 
669 ////////////////////////////////////////////////////////////////////////////////
670 /// Print info about this object to the specified stream. In addition to the info
671 /// from RooAbsPdf::printStream() we add:
672 ///
673 /// Verbose : detailed information on convolution integrals
674 
676 {
678 
679  os << indent << "--- RooAbsAnaConvPdf ---" << endl;
680  TIter iter = _convSet.createIterator() ;
681  RooResolutionModel* conv ;
682  while (((conv=(RooResolutionModel*)iter.Next()))) {
683  conv->printMultiline(os,contents,verbose,indent) ;
684  }
685 }
686 
687 
688 ///////////////////////////////////////////////////////////////////////////////
689 /// Label OK'ed components with cache-and-track
691 {
692  RooFIter citer = _convSet.fwdIterator() ;
693  RooAbsArg* carg ;
694  while ((carg=citer.next())) {
695  if (carg->canNodeBeCached()==Always) {
696  trackNodes.add(*carg) ;
697  //cout << "tracking node RooAddPdf component " << carg->IsA()->GetName() << "::" << carg->GetName() << endl ;
698  }
699  }
700 }
RooFormulaVar.h
c
#define c(i)
Definition: RSha256.hxx:101
RooCacheManager::setObj
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:50
RooTruthModel.h
RooAbsPdf::CacheElem
friend class CacheElem
Definition: RooAbsPdf.h:336
RooAbsAnaConvPdf::convVar
RooAbsRealLValue * convVar()
Retrieve the convolution variable.
Definition: RooAbsAnaConvPdf.cxx:312
RooAbsArg::operMode
OperMode operMode() const
Query the operation mode of this node.
Definition: RooAbsArg.h:484
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
RooAbsAnaConvPdf::RooConvGenContext
friend class RooConvGenContext
Definition: RooAbsAnaConvPdf.h:94
RooAbsAnaConvPdf::getCoefNorm
Double_t getCoefNorm(Int_t coefIdx, const RooArgSet &nset, const char *rangeName) const
Definition: RooAbsAnaConvPdf.h:50
RooAbsAnaConvPdf::getAnalyticalIntegralWN
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Advertise capability to perform (analytical) integrals internally.
Definition: RooAbsAnaConvPdf.cxx:364
RooMsgService.h
RooNameReg::ptr
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooConvCoefVar
RooConvCoefVar is an auxilary class that represents the coefficient of a RooAbsAnaConvPdf implementat...
Definition: RooConvCoefVar.h:28
RooAbsAnaConvPdf::coefficient
virtual Double_t coefficient(Int_t basisIndex) const =0
RooFit.h
RooAbsAnaConvPdf::declareBasis
Int_t declareBasis(const char *expression, const RooArgList &params)
Declare a basis function for use in this physics model.
Definition: RooAbsAnaConvPdf.cxx:158
RooTemplateProxy::arg
const T & arg() const
Return reference to object held in proxy.
Definition: RooTemplateProxy.h:290
RooFit::InputArguments
@ InputArguments
Definition: RooGlobalFunc.h:61
RooAbsArg::Always
@ Always
Definition: RooAbsArg.h:389
RooAbsArg::getParameters
RooArgSet * getParameters(const RooAbsData *data, bool stripDisconnected=true) 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:578
RooCacheManager::getObj
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:45
RooResolutionModel::convVar
RooAbsRealLValue & convVar() const
Return the convolution variable of the resolution model.
Definition: RooResolutionModel.h:43
RooAbsAnaConvPdf::_basisList
RooArgList _basisList
Definition: RooAbsAnaConvPdf.h:102
RooAbsAnaConvPdf::~RooAbsAnaConvPdf
virtual ~RooAbsAnaConvPdf()
Destructor.
Definition: RooAbsAnaConvPdf.cxx:132
RooAbsAnaConvPdf::coefAnalyticalIntegral
virtual Double_t coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=0) const
Default implementation of function implementing advertised integrals.
Definition: RooAbsAnaConvPdf.cxx:573
RooAbsReal::createIntegral
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:548
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
RooAbsAnaConvPdf::changeModel
virtual Bool_t changeModel(const RooResolutionModel &newModel)
Change the current resolution model to newModel.
Definition: RooAbsAnaConvPdf.cxx:208
RooAbsAnaConvPdf::makeCoefVarList
void makeCoefVarList(RooArgList &) const
Build complete list of coefficient variables.
Definition: RooAbsAnaConvPdf.cxx:633
r
ROOT::R::TRInterface & r
Definition: Object.C:4
RooAbsCollection::fwdIterator
RooFIter fwdIterator() const
One-time forward iterator.
Definition: RooAbsCollection.h:193
coutE
#define coutE(a)
Definition: RooMsgService.h:33
RooAbsAnaConvPdf::CacheElem
List of created basis functions.
Definition: RooAbsAnaConvPdf.h:105
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooGenContext
Class RooGenContext implement a universal generator context for all RooAbsPdf classes that do not hav...
Definition: RooGenContext.h:30
RooAbsAnaConvPdf::CacheElem::_coefVarList
RooArgList _coefVarList
Definition: RooAbsAnaConvPdf.h:115
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooResolutionModel.h
RooAbsPdf::printMultiline
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Print multi line detailed information of this RooAbsPdf.
Definition: RooAbsPdf.cxx:2039
RooAbsCollection::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Definition: RooAbsCollection.cxx:585
RooAbsArg::setOperMode
void setOperMode(OperMode mode, Bool_t recurseADirty=kTRUE)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1817
RooAbsAnaConvPdf::getCoefAnalyticalIntegral
virtual Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Default implementation of function advertising integration capabilities.
Definition: RooAbsAnaConvPdf.cxx:562
RooAbsReal::_forceNumInt
Bool_t _forceNumInt
Definition: RooAbsReal.h:480
RooAbsPdf::getGenerator
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:2474
coutI
#define coutI(a)
Definition: RooMsgService.h:30
indent
static void indent(ostringstream &buf, int indent_level)
Definition: TClingCallFunc.cxx:87
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooArgList::at
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:72
TString
Basic string class.
Definition: TString.h:136
RooAbsPdf::isDirectGenSafe
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:2509
RooListProxy::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE) override
Reimplementation of standard RooArgList::remove()
Definition: RooListProxy.cxx:143
bool
TIterator
Iterator abstract base class.
Definition: TIterator.h:30
RooAbsAnaConvPdf::_isCopy
Bool_t _isCopy
Definition: RooAbsAnaConvPdf.h:88
RooConvCoefVar.h
RooResolutionModel::convolution
virtual RooResolutionModel * convolution(RooFormulaVar *basis, RooAbsArg *owner) const
Instantiate a clone of this resolution model representing a convolution with given basis function.
Definition: RooResolutionModel.cxx:154
RooResolutionModel::printMultiline
virtual void printMultiline(std::ostream &os, Int_t content, Bool_t verbose=kFALSE, TString indent="") const
Print info about this object to the specified stream.
Definition: RooResolutionModel.cxx:345
RooAICRegistry::store
Int_t store(const std::vector< Int_t > &codeList, RooArgSet *set1=0, RooArgSet *set2=0, RooArgSet *set3=0, RooArgSet *set4=0)
Store given arrays of integer codes, and up to four RooArgSets in the registry (each setX pointer may...
Definition: RooAICRegistry.cxx:100
RooFormulaVar
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
RooResolutionModel
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
Definition: RooResolutionModel.h:26
RooAbsArg::replaceServer
void replaceServer(RooAbsArg &oldServer, RooAbsArg &newServer, Bool_t valueProp, Bool_t shapeProp)
Replace 'oldServer' with 'newServer'.
Definition: RooAbsArg.cxx:464
RooAbsAnaConvPdf::_coefNormMgr
RooObjCacheManager _coefNormMgr
Definition: RooAbsAnaConvPdf.h:118
RooAbsAnaConvPdf::RooAbsAnaConvPdf
RooAbsAnaConvPdf()
Default constructor, required for persistence.
Definition: RooAbsAnaConvPdf.cxx:85
RooAbsAnaConvPdf::_convSet
RooListProxy _convSet
Definition: RooAbsAnaConvPdf.h:101
Double_t
RooAbsAnaConvPdf::printMultiline
virtual void printMultiline(std::ostream &stream, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Print info about this object to the specified stream.
Definition: RooAbsAnaConvPdf.cxx:675
RooAbsAnaConvPdf.h
RooAbsAnaConvPdf::evaluate
virtual Double_t evaluate() const
Calculate the current unnormalized value of the PDF.
Definition: RooAbsAnaConvPdf.cxx:327
TNamed
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
RooFIter
A one-time forward iterator working on RooLinkedList or RooAbsCollection.
Definition: RooLinkedListIter.h:40
RooFit::Generation
@ Generation
Definition: RooGlobalFunc.h:60
RooArgProxy::absArg
RooAbsArg * absArg() const
Definition: RooArgProxy.h:37
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
RooAbsAnaConvPdf::isDirectGenSafe
virtual Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Return true if it is safe to generate the convolution observable from the internal generator (this is...
Definition: RooAbsAnaConvPdf.cxx:295
RooAbsAnaConvPdf::setCacheAndTrackHints
virtual void setCacheAndTrackHints(RooArgSet &)
Label OK'ed components with cache-and-track.
Definition: RooAbsAnaConvPdf.cxx:690
RooAbsCollection::createIterator
TIterator * createIterator(Bool_t dir=kIterForward) const
TIterator-style iteration over contained elements.
Definition: RooAbsCollection.h:178
RooFIter::next
RooAbsArg * next()
Return next element or nullptr if at end.
Definition: RooLinkedListIter.h:49
RooAbsArg::RooArgSet
friend class RooArgSet
Definition: RooAbsArg.h:599
RooAbsGenContext
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
Definition: RooAbsGenContext.h:26
RooAbsAnaConvPdf::CacheElem::_normList
RooArgList _normList
Definition: RooAbsAnaConvPdf.h:116
RooAbsCollection::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:455
RooResolutionModel::modelGenContext
virtual RooAbsGenContext * modelGenContext(const RooAbsAnaConvPdf &, const RooArgSet &, const RooDataSet *, const RooArgSet *, Bool_t) const
Definition: RooResolutionModel.h:36
RooGenContext.h
RooAbsAnaConvPdf::_codeReg
RooAICRegistry _codeReg
Definition: RooAbsAnaConvPdf.h:120
RooTemplateProxy::setArg
bool setArg(T &newRef)
Change object held in proxy into newRef.
Definition: RooTemplateProxy.h:227
RooRealVar.h
RooAbsCollection::end
const_iterator end() const
Definition: RooAbsCollection.h:202
RooAbsCollection::addOwned
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
Definition: RooAbsCollection.cxx:403
TIterator::Next
virtual TObject * Next()=0
TIterator::Reset
virtual void Reset()=0
TIter::Next
TObject * Next()
Definition: TCollection.h:249
RooAbsCollection::begin
const_iterator begin() const
Definition: RooAbsCollection.h:198
TString::CompareTo
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition: TString.cxx:442
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:295
RooAbsArg::setAttribute
void setAttribute(const Text_t *name, Bool_t value=kTRUE)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:291
RooNameReg::str
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.cxx:103
RooListProxy::addOwned
virtual Bool_t addOwned(RooAbsArg &var, Bool_t silent=kFALSE) override
Reimplementation of standard RooArgList::addOwned()
Definition: RooListProxy.cxx:114
name
char name[80]
Definition: TGX11.cxx:110
RooAbsAnaConvPdf::analyticalIntegralWN
virtual Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return analytical integral defined by given code, which is returned by getAnalyticalIntegralWN()
Definition: RooAbsAnaConvPdf.cxx:491
RooAbsArg::dependsOn
Bool_t dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=0, Bool_t valueOnly=kFALSE) const
Test whether we depend on (ie, are served by) any object in the specified collection.
Definition: RooAbsArg.cxx:799
genreflex::verbose
bool verbose
Definition: rootcling_impl.cxx:133
TIter
Definition: TCollection.h:233
RooDataSet
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:33
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
RooListProxy::add
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Reimplementation of standard RooArgList::add()
Definition: RooListProxy.cxx:104
RooAbsArg::canNodeBeCached
virtual CacheMode canNodeBeCached() const
Definition: RooAbsArg.h:412
RooAbsAnaConvPdf::forceAnalyticalInt
virtual Bool_t forceAnalyticalInt(const RooAbsArg &dep) const
This function forces RooRealIntegral to offer all integration dependents to RooAbsAnaConvPdf::getAnal...
Definition: RooAbsAnaConvPdf.cxx:593
RooAbsPdf
Definition: RooAbsPdf.h:41
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooAICRegistry::retrieve
const std::vector< Int_t > & retrieve(Int_t masterCode) const
Retrieve the array of integer codes associated with the given master code.
Definition: RooAICRegistry.cxx:153
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:39
RooFit::Eval
@ Eval
Definition: RooGlobalFunc.h:61
Riostream.h
RooAbsAnaConvPdf::_model
RooRealProxy _model
Definition: RooAbsAnaConvPdf.h:96
RooAbsRealLValue
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
Definition: RooAbsRealLValue.h:31
RooNameReg.h
RooConvGenContext.h
RooAbsAnaConvPdf::coefVars
virtual RooArgSet * coefVars(Int_t coefIdx) const
Return set of parameters with are used exclusively by the coefficient functions.
Definition: RooAbsAnaConvPdf.cxx:649
RooAbsAnaConvPdf
RooAbsAnaConvPdf is the base class for PDFs that represent a physics model that can be analytically c...
Definition: RooAbsAnaConvPdf.h:34
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:231
RooAbsAnaConvPdf::_convVar
RooRealProxy _convVar
Definition: RooAbsAnaConvPdf.h:97
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:33
RooAbsCollection::index
Int_t index(const RooAbsArg *arg) const
Returns index of given arg, or -1 if arg is not in the collection.
Definition: RooAbsCollection.h:247
RooAbsCollection::snapshot
RooAbsCollection * snapshot(Bool_t deepCopy=kTRUE) const
Take a snap shot of current collection contents.
Definition: RooAbsCollection.cxx:215
cxcoutD
#define cxcoutD(a)
Definition: RooMsgService.h:81
RooTruthModel
RooTruthModel is an implementation of RooResolution model that provides a delta-function resolution m...
Definition: RooTruthModel.h:21
int
RooAbsAnaConvPdf::genContext
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Create a generator context for this p.d.f.
Definition: RooAbsAnaConvPdf.cxx:252
RooListProxy::removeAll
virtual void removeAll() override
Reimplementation of standard RooArgList::removeAll()
Definition: RooListProxy.cxx:157