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 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 "RooMsgService.h"
66#include "Riostream.h"
67#include "RooResolutionModel.h"
68#include "RooRealVar.h"
69#include "RooFormulaVar.h"
70#include "RooConvGenContext.h"
71#include "RooGenContext.h"
72#include "RooTruthModel.h"
73#include "RooConvCoefVar.h"
74#include "RooNameReg.h"
75
76using namespace std;
77
79
80
81////////////////////////////////////////////////////////////////////////////////
82/// Default constructor, required for persistence
83
85 _isCopy(false),
86 _coefNormMgr(this,10)
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
96RooAbsAnaConvPdf::RooAbsAnaConvPdf(const char *name, const char *title,
97 const RooResolutionModel& model, RooRealVar& cVar) :
98 RooAbsPdf(name,title), _isCopy(false),
99 _model("!model","Original resolution model",this,(RooResolutionModel&)model,false,false),
100 _convVar("!convVar","Convolution variable",this,cVar,false,false),
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(true),
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
158Int_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 allOK(true) ;
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 = false ;
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 true ;
229 }
230
231 // Replace old convolutions with new set
233 _convSet.addOwned(newConvSet) ;
234
235 // Update server link by hand, since _model.setArg() below will not do this
236 replaceServer((RooAbsArg&)_model.arg(),(RooAbsArg&)newModel,false,false) ;
237
238 _model.setArg((RooResolutionModel&)newModel) ;
239 return false ;
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 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(),true,true) ;
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 pdfCanDir = (getGenerator(*convVar(),dummy) != 0) ;
267 bool 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 true ;
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 result(0) ;
330
331 Int_t index(0) ;
332 for (auto convArg : _convSet) {
333 auto conv = static_cast<RooAbsPdf*>(convArg);
334 double coef = coefficient(index++) ;
335 if (coef!=0.) {
336 double c = conv->getVal(0) ;
337 double 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.empty()) return 0 ;
369
370 if (_forceNumInt) return 0 ;
371
372 // Select subset of allVars that are actual dependents
373 RooArgSet allDeps;
374 getObservables(&allVars, allDeps);
375 std::unique_ptr<RooArgSet> normSet{normSet2 ? getObservables(normSet2) : nullptr};
376
377 RooArgSet intSetAll{allDeps,"intSetAll"};
378
379 // Split intSetAll in coef/conv parts
380 auto intCoefSet = std::make_unique<RooArgSet>("intCoefSet");
381 auto intConvSet = std::make_unique<RooArgSet>("intConvSet");
382
383 for (RooAbsArg * arg : intSetAll) {
384 bool ok(true) ;
385 for (RooAbsArg * conv : _convSet) {
386 if (conv->dependsOn(*arg)) ok=false ;
387 }
388
389 if (ok) {
390 intCoefSet->add(*arg) ;
391 } else {
392 intConvSet->add(*arg) ;
393 }
394
395 }
396
397 // Split normSetAll in coef/conv parts
398 auto normCoefSet = std::make_unique<RooArgSet>("normCoefSet");
399 auto normConvSet = std::make_unique<RooArgSet>("normConvSet");
400 if (normSet) {
401 for (RooAbsArg * arg : *normSet) {
402 bool ok(true) ;
403 for (RooAbsArg * conv : _convSet) {
404 if (conv->dependsOn(*arg)) ok=false ;
405 }
406
407 if (ok) {
408 normCoefSet->add(*arg) ;
409 } else {
410 normConvSet->add(*arg) ;
411 }
412
413 }
414 }
415
416 if (intCoefSet->empty()) intCoefSet.reset();
417 if (intConvSet->empty()) intConvSet.reset();
418 if (normCoefSet->empty()) normCoefSet.reset();
419 if (normConvSet->empty()) normConvSet.reset();
420
421
422 // Store integration configuration in registry
423 Int_t masterCode(0) ;
424 std::vector<Int_t> tmp(1, 0) ;
425
426 // takes ownership of all sets
427 masterCode = _codeReg.store(tmp,
428 intCoefSet.release(),
429 intConvSet.release(),
430 normCoefSet.release(),
431 normConvSet.release()) + 1;
432
433 analVars.add(allDeps) ;
434
435 return masterCode ;
436}
437
438
439
440
441////////////////////////////////////////////////////////////////////////////////
442/// Return analytical integral defined by given code, which is returned
443/// by getAnalyticalIntegralWN()
444///
445/// For unnormalized integrals the returned value is
446/// \f[
447/// \mathrm{PDF} = \sum_k \int \mathrm{coef}_k \; \mathrm{d}\bar{x}
448/// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}\bar{y},
449/// \f]
450/// where \f$ \bar{x} \f$ is the set of coefficient dependents to be integrated,
451/// and \f$ \bar{y} \f$ the set of basis function dependents to be integrated.
452///
453/// For normalized integrals this becomes
454/// \f[
455/// \mathrm{PDF} = \frac{\sum_k \int \mathrm{coef}_k \; \mathrm{d}x
456/// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}y}
457/// {\sum_k \int \mathrm{coef}_k \; \mathrm{d}v
458/// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}w},
459/// \f]
460/// where
461/// * \f$ x \f$ is the set of coefficient dependents to be integrated,
462/// * \f$ y \f$ the set of basis function dependents to be integrated,
463/// * \f$ v \f$ is the set of coefficient dependents over which is normalized and
464/// * \f$ w \f$ is the set of basis function dependents over which is normalized.
465///
466/// Set \f$ x \f$ must be contained in \f$ v \f$ and set \f$ y \f$ must be contained in \f$ w \f$.
467///
468
469double RooAbsAnaConvPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
470{
471 // WVE needs adaptation to handle new rangeName feature
472
473 // Handle trivial passthrough scenario
474 if (code==0) return getVal(normSet) ;
475
476 // Unpack master code
477 RooArgSet *intCoefSet, *intConvSet, *normCoefSet, *normConvSet ;
478 _codeReg.retrieve(code-1,intCoefSet,intConvSet,normCoefSet,normConvSet) ;
479
480 Int_t index(0) ;
481 double answer(0) ;
482
483 if (normCoefSet==0&&normConvSet==0) {
484
485 // Integral over unnormalized function
486 double integral(0) ;
487 const TNamed *_rangeName = RooNameReg::ptr(rangeName);
488 for (auto convArg : _convSet) {
489 auto conv = static_cast<RooResolutionModel*>(convArg);
490 double coef = getCoefNorm(index++,intCoefSet,_rangeName) ;
491 //cout << "coefInt[" << index << "] = " << coef << " " ; intCoefSet->Print("1") ;
492 if (coef!=0) {
493 integral += coef* conv->getNormObj(0,intConvSet,_rangeName)->getVal();
494 cxcoutD(Eval) << "RooAbsAnaConv::aiWN(" << GetName() << ") [" << index-1 << "] integral += " << conv->getNorm(intConvSet) << endl ;
495 }
496
497 }
498 answer = integral ;
499
500 } else {
501
502 // Integral over normalized function
503 double integral(0) ;
504 double norm(0) ;
505 const TNamed *_rangeName = RooNameReg::ptr(rangeName);
506 for (auto convArg : _convSet) {
507 auto conv = static_cast<RooResolutionModel*>(convArg);
508
509 double coefInt = getCoefNorm(index,intCoefSet,_rangeName) ;
510 //cout << "coefInt[" << index << "] = " << coefInt << "*" << term << " " << (intCoefSet?*intCoefSet:RooArgSet()) << endl ;
511 if (coefInt!=0) {
512 double term = conv->getNormObj(0,intConvSet,_rangeName)->getVal();
513 integral += coefInt*term ;
514 }
515
516 double coefNorm = getCoefNorm(index,normCoefSet) ;
517 //cout << "coefNorm[" << index << "] = " << coefNorm << "*" << term << " " << (normCoefSet?*normCoefSet:RooArgSet()) << endl ;
518 if (coefNorm!=0) {
519 double term = conv->getNormObj(0,normConvSet)->getVal();
520 norm += coefNorm*term ;
521 }
522
523 index++ ;
524 }
525 answer = integral/norm ;
526 }
527
528 return answer ;
529}
530
531
532
533////////////////////////////////////////////////////////////////////////////////
534/// Default implementation of function advertising integration capabilities. The interface is
535/// similar to that of getAnalyticalIntegral except that an integer code is added that
536/// designates the coefficient number for which the integration capabilities are requested
537///
538/// This default implementation advertises that no internal integrals are supported.
539
540Int_t RooAbsAnaConvPdf::getCoefAnalyticalIntegral(Int_t /* coef*/, RooArgSet& /*allVars*/, RooArgSet& /*analVars*/, const char* /*rangeName*/) const
541{
542 return 0 ;
543}
544
545
546
547////////////////////////////////////////////////////////////////////////////////
548/// Default implementation of function implementing advertised integrals. Only
549/// the pass-through scenario (no integration) is implemented.
550
551double RooAbsAnaConvPdf::coefAnalyticalIntegral(Int_t coef, Int_t code, const char* /*rangeName*/) const
552{
553 if (code==0) return coefficient(coef) ;
554 coutE(InputArguments) << "RooAbsAnaConvPdf::coefAnalyticalIntegral(" << GetName() << ") ERROR: unrecognized integration code: " << code << endl ;
555 assert(0) ;
556 return 1 ;
557}
558
559
560
561////////////////////////////////////////////////////////////////////////////////
562/// This function forces RooRealIntegral to offer all integration dependents
563/// to RooAbsAnaConvPdf::getAnalyticalIntegralWN() for consideration for
564/// internal integration, if RooRealIntegral considers this to be unsafe (e.g. due
565/// to hidden Jacobian terms).
566///
567/// RooAbsAnaConvPdf will not attempt to actually integrate all these dependents
568/// but feed them to the resolution models integration interface, which will
569/// make the final determination on how to integrate these dependents.
570
572{
573 return true ;
574}
575
576
577
578////////////////////////////////////////////////////////////////////////////////
579/// Returns the normalization integral value of the coefficient with number coefIdx over normalization
580/// set nset in range rangeName
581
582double RooAbsAnaConvPdf::getCoefNorm(Int_t coefIdx, const RooArgSet* nset, const TNamed* rangeName) const
583{
584 if (nset==0) return coefficient(coefIdx) ;
585
586 CacheElem* cache = (CacheElem*) _coefNormMgr.getObj(nset,0,0,rangeName) ;
587 if (!cache) {
588
589 cache = new CacheElem ;
590
591 // Make list of coefficient normalizations
592 Int_t i ;
594
595 for (i=0 ; i<cache->_coefVarList.getSize() ; i++) {
596 RooAbsReal* coefInt = static_cast<RooAbsReal&>(*cache->_coefVarList.at(i)).createIntegral(*nset,RooNameReg::str(rangeName)) ;
597 cache->_normList.addOwned(*coefInt) ;
598 }
599
600 _coefNormMgr.setObj(nset,0,cache,rangeName) ;
601 }
602
603 return ((RooAbsReal*)cache->_normList.at(coefIdx))->getVal() ;
604}
605
606
607
608////////////////////////////////////////////////////////////////////////////////
609/// Build complete list of coefficient variables
610
612{
613 // Instantate a coefficient variables
614 for (Int_t i=0 ; i<_convSet.getSize() ; i++) {
615 RooArgSet* cvars = coefVars(i) ;
616 RooAbsReal* coefVar = new RooConvCoefVar(Form("%s_coefVar_%d",GetName(),i),"coefVar",*this,i,cvars) ;
617 varList.addOwned(*coefVar) ;
618 delete cvars ;
619 }
620
621}
622
623
624////////////////////////////////////////////////////////////////////////////////
625/// Return set of parameters with are used exclusively by the coefficient functions
626
628{
629 RooArgSet* cVars = getParameters((RooArgSet*)0) ;
630 std::vector<RooAbsArg*> tmp;
631 for (auto arg : *cVars) {
632 for (auto convSetArg : _convSet) {
633 if (convSetArg->dependsOn(*arg)) {
634 tmp.push_back(arg);
635 }
636 }
637 }
638
639 cVars->remove(tmp.begin(), tmp.end(), true, true);
640
641 return cVars ;
642}
643
644
645
646
647////////////////////////////////////////////////////////////////////////////////
648/// Print info about this object to the specified stream. In addition to the info
649/// from RooAbsPdf::printStream() we add:
650///
651/// Verbose : detailed information on convolution integrals
652
653void RooAbsAnaConvPdf::printMultiline(ostream& os, Int_t contents, bool verbose, TString indent) const
654{
656
657 os << indent << "--- RooAbsAnaConvPdf ---" << endl;
658 for (auto * conv : static_range_cast<RooResolutionModel*>(_convSet)) {
659 conv->printMultiline(os,contents,verbose,indent) ;
660 }
661}
662
663
664///////////////////////////////////////////////////////////////////////////////
665/// Label OK'ed components with cache-and-track
667{
668 for (auto const* carg : static_range_cast<RooAbsArg*>(_convSet)) {
669 if (carg->canNodeBeCached()==Always) {
670 trackNodes.add(*carg) ;
671 //cout << "tracking node RooAddPdf component " << carg->ClassName() << "::" << carg->GetName() << endl ;
672 }
673 }
674}
#define c(i)
Definition: RSha256.hxx:101
#define coutI(a)
Definition: RooMsgService.h:34
#define cxcoutD(a)
Definition: RooMsgService.h:85
#define coutE(a)
Definition: RooMsgService.h:37
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
char name[80]
Definition: TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2456
const std::vector< Int_t > & retrieve(Int_t masterCode) const
Retrieve the array of integer codes associated with the given master code.
Int_t store(const std::vector< Int_t > &codeList, RooArgSet *set1=nullptr, RooArgSet *set2=nullptr, RooArgSet *set3=nullptr, RooArgSet *set4=nullptr)
Store given arrays of integer codes, and up to four RooArgSets in the registry (each setX pointer may...
RooAbsAnaConvPdf is the base class for PDFs that represent a physics model that can be analytically c...
friend class RooConvGenContext
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Advertise capability to perform (analytical) integrals internally.
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Return analytical integral defined by given code, which is returned by getAnalyticalIntegralWN()
virtual double coefAnalyticalIntegral(Int_t coef, Int_t code, const char *rangeName=nullptr) const
Default implementation of function implementing advertised integrals.
virtual bool changeModel(const RooResolutionModel &newModel)
Change the current resolution model to newModel.
double getCoefNorm(Int_t coefIdx, const RooArgSet &nset, const char *rangeName) const
void setCacheAndTrackHints(RooArgSet &) override
Label OK'ed components with cache-and-track.
bool forceAnalyticalInt(const RooAbsArg &dep) const override
This function forces RooRealIntegral to offer all integration dependents to RooAbsAnaConvPdf::getAnal...
RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false) const override
Create a generator context for this p.d.f.
virtual double coefficient(Int_t basisIndex) const =0
RooArgList _basisList
! List of created basis functions
RooObjCacheManager _coefNormMgr
! Coefficient normalization manager
void makeCoefVarList(RooArgList &) const
Build complete list of coefficient variables.
RooAICRegistry _codeReg
! Registry of analytical integration codes
virtual Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=nullptr) const
Default implementation of function advertising integration capabilities.
~RooAbsAnaConvPdf() override
Destructor.
bool isDirectGenSafe(const RooAbsArg &arg) const override
Return true if it is safe to generate the convolution observable from the internal generator (this is...
RooAbsRealLValue * convVar()
Retrieve the convolution variable.
double evaluate() const override
Calculate the current unnormalized value of the PDF.
RooRealProxy _model
Original model.
void printMultiline(std::ostream &stream, Int_t contents, bool verbose=false, TString indent="") const override
Print info about this object to the specified stream.
virtual RooArgSet * coefVars(Int_t coefIdx) const
Return set of parameters with are used exclusively by the coefficient functions.
RooRealProxy _convVar
Convolution variable.
RooAbsAnaConvPdf()
Default constructor, required for persistence.
Int_t declareBasis(const char *expression, const RooArgList &params)
Declare a basis function for use in this physics model.
RooListProxy _convSet
Set of (resModel (x) basisFunc) convolution objects.
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:71
void replaceServer(RooAbsArg &oldServer, RooAbsArg &newServer, bool valueProp, bool shapeProp)
Replace 'oldServer' with 'newServer'.
Definition: RooAbsArg.cxx:427
void setOperMode(OperMode mode, bool recurseADirty=true)
Set the operation mode of this node.
Definition: RooAbsArg.cxx:1866
RooArgSet * getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
Definition: RooAbsArg.h:293
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
Definition: RooAbsArg.cxx:246
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:541
OperMode operMode() const
Query the operation mode of this node.
Definition: RooAbsArg.h:484
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
bool empty() const
RooAbsCollection * snapshot(bool deepCopy=true) const
Take a snap shot of current collection contents.
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Int_t index(const RooAbsArg *arg) const
Returns index of given arg, or -1 if arg is not in the collection.
const_iterator end() const
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
const_iterator begin() const
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
virtual bool isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
Definition: RooAbsPdf.cxx:2357
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print multi line detailed information of this RooAbsPdf.
Definition: RooAbsPdf.cxx:1898
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool staticInitOK=true) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooAbsPdf.cxx:2322
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:62
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables std::liste...
Definition: RooAbsReal.cxx:522
bool _forceNumInt
Force numerical integration if flag set.
Definition: RooAbsReal.h:483
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:110
RooAbsArg * absArg() const
Return pointer to contained argument.
Definition: RooArgProxy.h:47
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:56
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=nullptr)
Setter function without integration set.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=nullptr, const TNamed *isetRangeName=nullptr)
Getter function without integration set.
void removeAll() override
Remove all argument inset using remove(const RooAbsArg&).
bool addOwned(RooAbsArg &var, bool silent=false) override
Overloaded RooCollection_t::addOwned() method insert object into owning set and registers object as s...
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false) override
Remove object 'var' from set and deregister 'var' as server to owner.
RooConvCoefVar is an auxilary class that represents the coefficient of a RooAbsAnaConvPdf implementat...
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:55
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Definition: RooFormulaVar.h:30
Class RooGenContext implement a universal generator context for all RooAbsPdf classes that do not hav...
Definition: RooGenContext.h:30
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition: RooNameReg.h:37
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:81
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
RooResolutionModel is the base class for PDFs that represent a resolution model that can be convolute...
virtual RooResolutionModel * convolution(RooFormulaVar *basis, RooAbsArg *owner) const
Instantiate a clone of this resolution model representing a convolution with given basis function.
virtual RooAbsGenContext * modelGenContext(const RooAbsAnaConvPdf &, const RooArgSet &, const RooDataSet *, const RooArgSet *, bool) const
RooAbsRealLValue & convVar() const
Return the convolution variable of the resolution model.
bool setArg(T &newRef)
Change object held in proxy into newRef.
const T & arg() const
Return reference to object held in proxy.
RooTruthModel is an implementation of RooResolution model that provides a delta-function resolution m...
Definition: RooTruthModel.h:21
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
Basic string class.
Definition: TString.h:136
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition: TString.cxx:451
TString & Append(const char *cs)
Definition: TString.h:564
@ Generation
Definition: RooGlobalFunc.h:61
@ InputArguments
Definition: RooGlobalFunc.h:62