Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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/// 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}) \;
37/// \mathrm{d}x',
38/// \f]
39///
40/// which RooAbsAnaConvPdf uses to construct the pdf for [ Phys (x) R ] :
41/// \f[
42/// \mathrm{PDF}(x,\bar{a},\bar{b},\bar{c}) = \sum_k \mathrm{coef}_k(\bar{a}) * R_k(x,\bar{b},\bar{c})
43/// \f]
44///
45/// A minimal implementation of a RooAbsAnaConvPdf physics model consists of
46///
47/// - A constructor that declares the required basis functions using the declareBasis() method.
48/// The declareBasis() function assigns a unique identifier code to each declare basis
49///
50/// - An implementation of `coefficient(Int_t code)` returning the coefficient value for each
51/// declared basis function
52///
53/// Optionally, analytical integrals can be provided for the coefficient functions. The
54/// interface for this is quite similar to that for integrals of regular PDFs. Two functions,
55/// \code{.cpp}
56/// Int_t getCoefAnalyticalIntegral(Int_t coef, RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
57/// double coefAnalyticalIntegral(Int_t coef, Int_t code, const char* rangeName) const
58/// \endcode
59///
60/// advertise the coefficient integration capabilities and implement them respectively.
61/// Please see RooAbsPdf for additional details. Advertised analytical integrals must be
62/// valid for all coefficients.
63///
64/// ### The resolution model is a configuration object, not a graph node
65///
66/// The resolution model passed to the constructor is **not** a node of the
67/// computation graph of the RooAbsAnaConvPdf. It is never evaluated directly;
68/// it only serves as a *configuration* object that specifies which resolution
69/// model should be convolved with the basis functions. From it, the
70/// RooAbsAnaConvPdf builds its own internal \f$ \mathrm{basis}_k \otimes
71/// \mathrm{resModel} \f$ convolution objects (one per declared basis function),
72/// and it is *those* convolutions that are the actual value servers of the pdf
73/// and that get evaluated.
74///
75/// Consequently, the resolution model itself is not a server of the
76/// RooAbsAnaConvPdf. It remains accessible via getModel() (for example for
77/// serialization), but it does not appear in the pdf's `servers()` list, in
78/// `getParameters()` / `getVariables()`, or in the printed computation graph.
79///
80/// \note **Behavior change in ROOT 6.42:** in earlier releases the resolution
81/// model was kept as a (non-value, non-shape) server of the RooAbsAnaConvPdf.
82/// As a side effect, importing a RooAbsAnaConvPdf into a RooWorkspace also
83/// dragged the original resolution model into the workspace (and into HS3/JSON
84/// exports), even though it played no role in the computation. As of ROOT 6.42
85/// this is no longer the case: a resolution model that is only used as the
86/// configuration of a RooAbsAnaConvPdf is not imported into the workspace on
87/// its own anymore.
88///
89/// Objects written with older ROOT versions are read back correctly via schema
90/// evolution: the resolution model is dropped as a server, so the *pdf's*
91/// computation graph is the same as for a freshly constructed one. Note,
92/// however, that if such an old file already stored the resolution model as a
93/// standalone RooWorkspace member (which used to happen on import), that member
94/// is *not* retroactively removed from the workspace on read-back -- it simply
95/// is no longer wired into the RooAbsAnaConvPdf. Only newly created workspaces
96/// are guaranteed to be free of the standalone resolution model.
97
98#include "RooAbsAnaConvPdf.h"
99
101#include "RooMsgService.h"
102#include "Riostream.h"
103#include "RooResolutionModel.h"
104#include "RooRealVar.h"
105#include "RooFormulaVar.h"
106#include "RooConvGenContext.h"
107#include "RooGenContext.h"
108#include "RooTruthModel.h"
109#include "RooConvCoefVar.h"
110#include "RooNameReg.h"
111
112using std::endl, std::string, std::ostream;
113
114
115
116////////////////////////////////////////////////////////////////////////////////
117/// Default constructor, required for persistence
118
120 _isCopy(false),
121 _coefNormMgr(this,10)
122{
123}
124
125
126
127////////////////////////////////////////////////////////////////////////////////
128/// Constructor. The supplied resolution model must be constructed with the same
129/// convoluted variable as this physics model ('convVar')
130
131RooAbsAnaConvPdf::RooAbsAnaConvPdf(const char *name, const char *title, const RooResolutionModel &model,
133 : RooAbsPdf(name, title),
134 _isCopy(false),
135 _model{static_cast<RooResolutionModel *>(model.clone(model.GetName()))},
136 _ownModel{true},
137 _convVar("!convVar", "Convolution variable", this, cVar, false, false),
138 _convSet("!convSet", "Set of resModel X basisFunc convolutions", this),
139 _coefNormMgr(this, 10),
140 _codeReg(10)
141{
142 _model->setAttribute("NOCacheAndTrack");
143}
144
145
146
147////////////////////////////////////////////////////////////////////////////////
148
150 : RooAbsPdf(other, name),
151 _isCopy(true),
152 _model{other._model ? static_cast<RooResolutionModel *>(other._model->clone(other._model->GetName())) : nullptr},
153 _ownModel{true},
154 _convVar("!convVar", this, other._convVar),
155 _convSet("!convSet", this, other._convSet),
156 _coefNormMgr(other._coefNormMgr, this),
157 _codeReg(other._codeReg)
158{
159 // Copy constructor
160 if (_model) {
161 _model->setAttribute("NOCacheAndTrack");
162 }
163 other._basisList.snapshot(_basisList);
164}
165
166
167
168////////////////////////////////////////////////////////////////////////////////
169/// Destructor
170
172{
173 if (!_isCopy) {
174 std::vector<RooAbsArg*> tmp(_convSet.begin(), _convSet.end());
175
176 for (auto arg : tmp) {
177 _convSet.remove(*arg) ;
178 delete arg ;
179 }
180 }
181
182 if (_ownModel) {
183 delete _model;
184 }
185}
186
187////////////////////////////////////////////////////////////////////////////////
188/// Forward server redirection to the original resolution model. The resolution
189/// model is not a server of this pdf (it is only used to build the convolutions
190/// and for generation), so it is not redirected by the standard machinery. We
191/// keep its servers in sync here, analogous to RooResolutionModel forwarding
192/// the redirection to its basis function.
193
195 bool isRecursive)
196{
197 if (_model) {
198 // Pass mustReplaceAll=false: the model may legitimately reference servers
199 // that are not part of this particular redirection, and it is never
200 // evaluated as part of the computation graph anyway.
202 }
203
205}
206
207////////////////////////////////////////////////////////////////////////////////
208/// Second-pass schema evolution, called by the RooWorkspace after reading.
209///
210/// In class version <= 3, the original resolution model was held in a
211/// RooRealProxy and was therefore a (non-value, non-shape) server of this pdf.
212/// The schema evolution read rule (see LinkDef.h) already recovered the model
213/// pointer into _model, but the stale server link is also restored from the
214/// file via the RooAbsArg server list. We remove it here, once the full graph
215/// is live, so that an object read from an old file has the same clean server
216/// structure as a freshly constructed one. This is safe because there is no
217/// longer a proxy that could resurrect the server link on copy.
218
220{
222
223 // Force removal (the link may have been added with a reference count > 1):
224 // the model must be completely severed from the server list.
225 if (_model && findServer(*_model)) {
226 removeServer(*_model, true);
227 }
228}
229
230////////////////////////////////////////////////////////////////////////////////
231/// Declare a basis function for use in this physics model. The string expression
232/// must be a valid RooFormulVar expression representing the basis function, referring
233/// to the convolution variable as '@0', and any additional parameters (supplied in
234/// 'params' as '@1','@2' etc.
235///
236/// The return value is a unique identifier code, that will be passed to coefficient()
237/// to identify the basis function for which the coefficient is requested. If the
238/// resolution model used does not support the declared basis function, code -1 is
239/// returned.
240///
241
242Int_t RooAbsAnaConvPdf::declareBasis(const char* expression, const RooArgList& params)
243{
244 // Sanity check
245 if (_isCopy) {
246 coutE(InputArguments) << "RooAbsAnaConvPdf::declareBasis(" << GetName() << "): ERROR attempt to "
247 << " declare basis functions in a copied RooAbsAnaConvPdf" << std::endl ;
248 return -1 ;
249 }
250
251 // Resolution model must support declared basis
252 if (!_model->isBasisSupported(expression)) {
253 coutE(InputArguments) << "RooAbsAnaConvPdf::declareBasis(" << GetName() << "): resolution model "
254 << _model->GetName() << " doesn't support basis function " << expression << std::endl;
255 return -1;
256 }
257
258 // Instantiate basis function
260 basisArgs.add(params) ;
261
262 TString basisName(expression) ;
263 for (const auto arg : basisArgs) {
264 basisName.Append("_") ;
265 basisName.Append(arg->GetName()) ;
266 }
267
268 auto basisFunc = std::make_unique<RooFormulaVar>(basisName, expression, basisArgs);
269 basisFunc->setAttribute("RooWorkspace::Recycle") ;
270 basisFunc->setAttribute("NOCacheAndTrack") ;
271 basisFunc->setOperMode(operMode()) ;
272
273 // Instantiate resModel x basisFunc convolution
274 RooAbsReal *conv = _model->convolution(basisFunc.get(), this);
275 _basisList.addOwned(std::move(basisFunc));
276 if (!conv) {
277 coutE(InputArguments) << "RooAbsAnaConvPdf::declareBasis(" << GetName() << "): unable to construct convolution with basis function '"
278 << expression << "'" << std::endl ;
279 return -1 ;
280 }
281 _convSet.add(*conv) ;
282
283 return _convSet.index(conv) ;
284}
285
286
287
288////////////////////////////////////////////////////////////////////////////////
289/// Change the current resolution model to newModel
290
292{
294 bool allOK(true) ;
296
297 // Build new resolution model
298 std::unique_ptr<RooResolutionModel> newConv{newModel.convolution(const_cast<RooFormulaVar*>(&conv->basis()),this)};
299 if (!newConvSet.addOwned(std::move(newConv))) {
300 allOK = false ;
301 break ;
302 }
303 }
304
305 // Check if all convolutions were successfully built
306 if (!allOK) {
307 return true ;
308 }
309
310 // Replace old convolutions with new set
312 _convSet.addOwned(std::move(newConvSet));
313
314 // Replace the stored original resolution model. Since it is not a server of
315 // this pdf, it cannot (and need not) be redirected via redirectServers(): we
316 // simply own a fresh clone of the new model.
317 if (_ownModel) {
318 delete _model;
319 }
320 _model = static_cast<RooResolutionModel *>(newModel.clone(newModel.GetName()));
321 _ownModel = true;
322 _model->setAttribute("NOCacheAndTrack");
323
324 return false ;
325}
326
327
328
329
330////////////////////////////////////////////////////////////////////////////////
331/// Create a generator context for this p.d.f. If both the p.d.f and the resolution model
332/// support internal generation of the convolution observable on an infinite domain,
333/// deploy a specialized convolution generator context, which generates the physics distribution
334/// and the smearing separately, adding them a posteriori. If this is not possible return
335/// a (slower) generic generation context that uses accept/reject sampling
336
338 const RooArgSet* auxProto, bool verbose) const
339{
340 // Check if the resolution model specifies a special context to be used.
342 assert(conv);
343
344 std::unique_ptr<RooArgSet> modelDep {_model->getObservables(&vars)};
345 modelDep->remove(*convVar(),true,true) ;
346 Int_t numAddDep = modelDep->size() ;
347
348 // Check if physics PDF and resolution model can both directly generate the convolution variable
349 RooArgSet dummy ;
350 bool pdfCanDir = (getGenerator(*convVar(),dummy) != 0) ;
351 bool resCanDir = conv && (conv->getGenerator(*convVar(),dummy)!=0) && conv->isDirectGenSafe(*convVar()) ;
352
353 if (numAddDep>0 || !pdfCanDir || !resCanDir) {
354 // Any resolution model with more dependents than the convolution variable
355 // or pdf or resmodel do not support direct generation
356 string reason ;
357 if (numAddDep>0) reason += "Resolution model has more observables than the convolution variable. " ;
358 if (!pdfCanDir) reason += "PDF does not support internal generation of convolution observable. " ;
359 if (!resCanDir) reason += "Resolution model does not support internal generation of convolution observable. " ;
360
361 coutI(Generation) << "RooAbsAnaConvPdf::genContext(" << GetName() << ") Using regular accept/reject generator for convolution p.d.f because: " << reason.c_str() << std::endl ;
362 return new RooGenContext(*this,vars,prototype,auxProto,verbose) ;
363 }
364
365 RooAbsGenContext* context = conv->modelGenContext(*this, vars, prototype, auxProto, verbose);
366 if (context) return context;
367
368 // Any other resolution model: use specialized generator context
369 return new RooConvGenContext(*this,vars,prototype,auxProto,verbose) ;
370}
371
372
373
374////////////////////////////////////////////////////////////////////////////////
375/// Return true if it is safe to generate the convolution observable
376/// from the internal generator (this is the case if the chosen resolution
377/// model is the truth model)
378
380{
381
382 // All direct generation of convolution arg if model is truth model
383 if (!TString(_convVar.absArg()->GetName()).CompareTo(arg.GetName()) && dynamic_cast<RooTruthModel *>(_model)) {
384 return true;
385 }
386
387 return RooAbsPdf::isDirectGenSafe(arg) ;
388}
389
390
391
392////////////////////////////////////////////////////////////////////////////////
393/// Return a pointer to the convolution variable instance used in the resolution model
394
396{
397 auto* conv = static_cast<RooResolutionModel*>(_convSet.at(0));
398 if (!conv) return nullptr;
399 return &conv->convVar() ;
400}
401
402
403
404////////////////////////////////////////////////////////////////////////////////
405/// Calculate the current unnormalized value of the PDF
406///
407/// PDF = sum_k coef_k * [ basis_k (x) ResModel ]
408///
409
411{
412 double result(0) ;
413
414 Int_t index(0) ;
415 for (auto *conv : static_range_cast<RooAbsPdf*>(_convSet)) {
416 double coef = coefficient(index++) ;
417 if (coef!=0.) {
418 const double c = conv->getVal(nullptr);
419 cxcoutD(Eval) << "RooAbsAnaConvPdf::evaluate(" << GetName() << ") val += coef*conv [" << index-1 << "/"
420 << _convSet.size() << "] coef = " << coef << " conv = " << c << std::endl ;
421 result += c * coef;
422 } else {
423 cxcoutD(Eval) << "RooAbsAnaConvPdf::evaluate(" << GetName() << ") [" << index-1 << "/" << _convSet.size() << "] coef = 0" << std::endl ;
424 }
425 }
426
427 return result ;
428}
429
430
431
432////////////////////////////////////////////////////////////////////////////////
433/// Advertise capability to perform (analytical) integrals
434/// internally. For a given integration request over allVars while
435/// normalized over normSet2 and in range 'rangeName', returns
436/// largest subset that can be performed internally in analVars
437/// Return code is unique integer code identifying integration scenario
438/// to be passed to analyticalIntegralWN() to calculate requeste integral
439///
440/// Class RooAbsAnaConv defers analytical integration request to
441/// resolution model and/or coefficient implementations and
442/// aggregates results into composite configuration with a unique
443/// code assigned by RooAICRegistry
444
446 RooArgSet& analVars, const RooArgSet* normSet2, const char* /*rangeName*/) const
447{
448 // Handle trivial no-integration scenario
449 if (allVars.empty()) return 0 ;
450
451 if (_forceNumInt) return 0 ;
452
453 // Select subset of allVars that are actual dependents
455 getObservables(&allVars, allDeps);
456 std::unique_ptr<RooArgSet> normSet{normSet2 ? getObservables(normSet2) : nullptr};
457
458 RooArgSet intSetAll{allDeps,"intSetAll"};
459
460 // Split intSetAll in coef/conv parts
461 auto intCoefSet = std::make_unique<RooArgSet>("intCoefSet");
462 auto intConvSet = std::make_unique<RooArgSet>("intConvSet");
463
464 for (RooAbsArg * arg : intSetAll) {
465 bool ok(true) ;
466 for (RooAbsArg * conv : _convSet) {
467 if (conv->dependsOn(*arg)) ok=false ;
468 }
469
470 if (ok) {
471 intCoefSet->add(*arg) ;
472 } else {
473 intConvSet->add(*arg) ;
474 }
475
476 }
477
478 // Split normSetAll in coef/conv parts
479 auto normCoefSet = std::make_unique<RooArgSet>("normCoefSet");
480 auto normConvSet = std::make_unique<RooArgSet>("normConvSet");
481 if (normSet) {
482 for (RooAbsArg * arg : *normSet) {
483 bool ok(true) ;
484 for (RooAbsArg * conv : _convSet) {
485 if (conv->dependsOn(*arg)) ok=false ;
486 }
487
488 if (ok) {
489 normCoefSet->add(*arg) ;
490 } else {
491 normConvSet->add(*arg) ;
492 }
493
494 }
495 }
496
497 if (intCoefSet->empty()) intCoefSet.reset();
498 if (intConvSet->empty()) intConvSet.reset();
499 if (normCoefSet->empty()) normCoefSet.reset();
500 if (normConvSet->empty()) normConvSet.reset();
501
502
503 // Store integration configuration in registry
504 Int_t masterCode(0) ;
505 std::vector<Int_t> tmp(1, 0) ;
506
507 // takes ownership of all sets
509 intCoefSet.release(),
510 intConvSet.release(),
511 normCoefSet.release(),
512 normConvSet.release()) + 1;
513
514 analVars.add(allDeps) ;
515
516 return masterCode ;
517}
518
519
520
521
522////////////////////////////////////////////////////////////////////////////////
523/// Return analytical integral defined by given code, which is returned
524/// by getAnalyticalIntegralWN()
525///
526/// For unnormalized integrals the returned value is
527/// \f[
528/// \mathrm{PDF} = \sum_k \int \mathrm{coef}_k \; \mathrm{d}\bar{x}
529/// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}\bar{y},
530/// \f]
531/// where \f$ \bar{x} \f$ is the set of coefficient dependents to be integrated,
532/// and \f$ \bar{y} \f$ the set of basis function dependents to be integrated.
533///
534/// For normalized integrals this becomes
535/// \f[
536/// \mathrm{PDF} = \frac{\sum_k \int \mathrm{coef}_k \; \mathrm{d}x
537/// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}y}
538/// {\sum_k \int \mathrm{coef}_k \; \mathrm{d}v
539/// \cdot \int \mathrm{basis}_k (x) \mathrm{ResModel} \; \mathrm{d}w},
540/// \f]
541/// where
542/// * \f$ x \f$ is the set of coefficient dependents to be integrated,
543/// * \f$ y \f$ the set of basis function dependents to be integrated,
544/// * \f$ v \f$ is the set of coefficient dependents over which is normalized and
545/// * \f$ w \f$ is the set of basis function dependents over which is normalized.
546///
547/// Set \f$ x \f$ must be contained in \f$ v \f$ and set \f$ y \f$ must be contained in \f$ w \f$.
548///
549
551{
552 // WVE needs adaptation to handle new rangeName feature
553
554 // Handle trivial passthrough scenario
555 if (code == 0)
556 return getVal(normSet);
557
558 // Unpack master code
564
565 Int_t index(0);
566
567 if (normCoefSet == nullptr && normConvSet == nullptr) {
568 // Integral over unnormalized function
569 double integral(0);
571 for (auto *conv : static_range_cast<RooAbsPdf *>(_convSet)) {
572 double coef = getCoefNorm(index++, intCoefSet, rangeNamePtr);
573 if (coef != 0) {
574 const double term = coef * conv->getNormObj(nullptr, intConvSet, rangeNamePtr)->getVal();
575 integral += term;
576 cxcoutD(Eval) << "RooAbsAnaConv::aiWN(" << GetName() << ") [" << index - 1 << "] integral += " << term
577 << std::endl;
578 }
579 }
580 return integral;
581 }
582
583 // Integral over normalized function
584 double integral(0);
585 double norm(0);
587 for (auto *conv : static_range_cast<RooAbsPdf *>(_convSet)) {
588
590 if (coefInt != 0) {
591 double term = conv->getNormObj(nullptr, intConvSet, rangeNamePtr)->getVal();
592 integral += coefInt * term;
593 }
594
596 if (coefNorm != 0) {
597 double term = conv->getNormObj(nullptr, normConvSet)->getVal();
598 norm += coefNorm * term;
599 }
600
601 index++;
602 }
603 return integral / norm;
604}
605
606
607
608////////////////////////////////////////////////////////////////////////////////
609/// Default implementation of function advertising integration capabilities. The interface is
610/// similar to that of getAnalyticalIntegral except that an integer code is added that
611/// designates the coefficient number for which the integration capabilities are requested
612///
613/// This default implementation advertises that no internal integrals are supported.
614
615Int_t RooAbsAnaConvPdf::getCoefAnalyticalIntegral(Int_t /* coef*/, RooArgSet& /*allVars*/, RooArgSet& /*analVars*/, const char* /*rangeName*/) const
616{
617 return 0 ;
618}
619
620
621
622////////////////////////////////////////////////////////////////////////////////
623/// Default implementation of function implementing advertised integrals. Only
624/// the pass-through scenario (no integration) is implemented.
625
626double RooAbsAnaConvPdf::coefAnalyticalIntegral(Int_t coef, Int_t code, const char* /*rangeName*/) const
627{
628 if (code==0) return coefficient(coef) ;
629 coutE(InputArguments) << "RooAbsAnaConvPdf::coefAnalyticalIntegral(" << GetName() << ") ERROR: unrecognized integration code: " << code << std::endl ;
630 assert(0) ;
631 return 1 ;
632}
633
634
635
636////////////////////////////////////////////////////////////////////////////////
637/// This function forces RooRealIntegral to offer all integration dependents
638/// to RooAbsAnaConvPdf::getAnalyticalIntegralWN() for consideration for
639/// internal integration, if RooRealIntegral considers this to be unsafe (e.g. due
640/// to hidden Jacobian terms).
641///
642/// RooAbsAnaConvPdf will not attempt to actually integrate all these dependents
643/// but feed them to the resolution models integration interface, which will
644/// make the final determination on how to integrate these dependents.
645
647{
648 return true ;
649}
650
651
652
653////////////////////////////////////////////////////////////////////////////////
654/// Returns the normalization integral value of the coefficient with number coefIdx over normalization
655/// set nset in range rangeName
656
658{
659 if (nset==nullptr) return coefficient(coefIdx) ;
660
661 CacheElem* cache = static_cast<CacheElem*>(_coefNormMgr.getObj(nset,nullptr,nullptr,rangeName)) ;
662 if (!cache) {
663
664 cache = new CacheElem ;
665
666 // Make list of coefficient normalizations
668
669 for (std::size_t i=0 ; i<cache->_coefVarList.size() ; i++) {
670 cache->_normList.addOwned(std::unique_ptr<RooAbsReal>{static_cast<RooAbsReal&>(*cache->_coefVarList.at(i)).createIntegral(*nset,RooNameReg::str(rangeName))});
671 }
672
673 _coefNormMgr.setObj(nset,nullptr,cache,rangeName) ;
674 }
675
676 return (static_cast<RooAbsReal*>(cache->_normList.at(coefIdx)))->getVal() ;
677}
678
679
680
681////////////////////////////////////////////////////////////////////////////////
682/// Build complete list of coefficient variables
683
685{
686 // Instantiate a coefficient variables
687 for (std::size_t i=0 ; i<_convSet.size() ; i++) {
688 auto cvars = coefVars(i);
689 std::string name = std::string{GetName()} + "_coefVar_" + std::to_string(i);
690 varList.addOwned(std::make_unique<RooConvCoefVar>(name.c_str(),"coefVar",*this,i,&*cvars));
691 }
692
693}
694
695
696////////////////////////////////////////////////////////////////////////////////
697/// Return set of parameters with are used exclusively by the coefficient functions
698
700{
701 std::unique_ptr<RooArgSet> cVars{getParameters(static_cast<RooArgSet*>(nullptr))};
702 std::vector<RooAbsArg*> tmp;
703 for (auto arg : *cVars) {
704 for (auto convSetArg : _convSet) {
705 if (convSetArg->dependsOn(*arg)) {
706 tmp.push_back(arg);
707 }
708 }
709 }
710
711 cVars->remove(tmp.begin(), tmp.end(), true, true);
712
713 return RooFit::makeOwningPtr(std::move(cVars));
714}
715
716
717
718
719////////////////////////////////////////////////////////////////////////////////
720/// Print info about this object to the specified stream. In addition to the info
721/// from RooAbsPdf::printStream() we add:
722///
723/// Verbose : detailed information on convolution integrals
724
725void RooAbsAnaConvPdf::printMultiline(ostream& os, Int_t contents, bool verbose, TString indent) const
726{
727 RooAbsPdf::printMultiline(os,contents,verbose,indent);
728
729 os << indent << "--- RooAbsAnaConvPdf ---" << std::endl;
730 for (RooAbsArg * conv : _convSet) {
731 conv->printMultiline(os,contents,verbose,indent) ;
732 }
733}
734
735
736///////////////////////////////////////////////////////////////////////////////
737/// Label OK'ed components with cache-and-track
739{
740 for (auto const* carg : static_range_cast<RooAbsArg*>(_convSet)) {
741 if (carg->canNodeBeCached()==Always) {
742 trackNodes.add(*carg) ;
743 //cout << "tracking node RooAddPdf component " << carg->ClassName() << "::" << carg->GetName() << std::endl ;
744 }
745 }
746}
747
748std::unique_ptr<RooAbsArg>
750{
751 // If there is only one component in the linear sum of convolutions, we can
752 // just return that one, normalized.
753 if(_convSet.size() == 1) {
754 if (normSet.empty()) {
755 return _convSet[0].compileForNormSet(normSet, ctx);
756 }
757 std::unique_ptr<RooAbsPdf> pdfClone(static_cast<RooAbsPdf *>(_convSet[0].Clone()));
759
760 auto newArg = std::make_unique<RooFit::Detail::RooNormalizedPdf>(*pdfClone, normSet);
761
762 // The direct servers are this pdf and the normalization integral, which
763 // don't need to be compiled further.
764 for (RooAbsArg *server : newArg->servers()) {
765 server->setAttribute("_COMPILED");
766 }
767 newArg->setAttribute("_COMPILED");
768 newArg->addOwnedComponents(std::move(pdfClone));
769 return newArg;
770 }
771
772 // Here, we can't use directly the function from the RooAbsPdf base class,
773 // because the convolution argument servers need to be evaluated
774 // unnormalized, even if they are pdfs.
775
776 if (normSet.empty()) {
778 }
779 std::unique_ptr<RooAbsAnaConvPdf> pdfClone(static_cast<RooAbsAnaConvPdf *>(this->Clone()));
780
781 // The other servers will be compiled with the original normSet, but the
782 // _convSet has to be evaluated unnormalized.
784 for (RooAbsArg *convArg : _convSet) {
785 if (auto convArgClone = ctx.compile(*convArg, *pdfClone, {})) {
787 }
788 }
789 pdfClone->redirectServers(convArgClones, false, true);
790
791 // Compile remaining servers that are evaluated normalized
793
794 // Finally, this RooAbsAnaConvPdf needs to be normalized
795 auto newArg = std::make_unique<RooFit::Detail::RooNormalizedPdf>(*pdfClone, normSet);
796
797 // The direct servers are this pdf and the normalization integral, which
798 // don't need to be compiled further.
799 for (RooAbsArg *server : newArg->servers()) {
800 server->setAttribute("_COMPILED");
801 }
802 newArg->setAttribute("_COMPILED");
803 newArg->addOwnedComponents(std::move(pdfClone));
804 return newArg;
805}
#define c(i)
Definition RSha256.hxx:101
#define coutI(a)
#define cxcoutD(a)
#define coutE(a)
static void indent(ostringstream &buf, int indent_level)
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
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:148
const_iterator begin() const
const_iterator end() const
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...
Base class for PDFs that represent a physics model that can be analytically convolved with a resoluti...
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.
RooResolutionModel * _model
Original resolution model (not a server)
virtual double coefficient(Int_t basisIndex) const =0
void ioStreamerPass2() override
Second-pass schema evolution, called by the RooWorkspace after reading.
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.
bool _ownModel
Flag indicating ownership of _model.
~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.
void printMultiline(std::ostream &stream, Int_t contents, bool verbose=false, TString indent="") const override
Print info about this object to the specified stream.
RooRealProxy _convVar
Convolution variable.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursive) override
Forward server redirection to the original resolution model.
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.
virtual RooFit::OwningPtr< RooArgSet > coefVars(Int_t coefIdx) const
Return set of parameters with are used exclusively by the coefficient functions.
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
Common abstract base class for objects that represent a value and a "shape" in RooFit.
Definition RooAbsArg.h:76
RooFit::OwningPtr< 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...
void removeServer(RooAbsArg &server, bool force=false)
Unregister another RooAbsArg as a server to us, ie, declare that we no longer depend on its value and...
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
bool redirectServers(const RooAbsCollection &newServerList, bool mustReplaceAll=false, bool nameChange=false, bool isRecursionStep=false)
Replace all direct servers of this object with the new servers in newServerList.
virtual void ioStreamerPass2()
Method called by workspace container to finalize schema evolution issues that cannot be handled in a ...
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:88
RooAbsArg * findServer(const char *name) const
Return server of this with name name. Returns nullptr if not found.
Definition RooAbsArg.h:147
OperMode operMode() const
Query the operation mode of this node.
Definition RooAbsArg.h:419
Abstract container object that can hold multiple RooAbsArg objects.
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
Storage_t::size_type size() const
virtual bool addOwned(RooAbsArg &var, bool silent=false)
Add an argument and transfer the ownership to the collection.
const_iterator begin() const
Abstract base class for generator contexts of RooAbsPdf objects.
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:32
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
virtual bool isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print multi line detailed information of this RooAbsPdf.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursiveStep) override
Hook function intercepting redirectServer calls.
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...
Abstract base class for objects that represent a real value that may appear on the left hand side of ...
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:63
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:107
bool _forceNumInt
Force numerical integration if flag set.
Definition RooAbsReal.h:542
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Create an object that represents the integral of the function over one or more observables listed in ...
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:46
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
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.
Container class to hold unbinned data.
Definition RooDataSet.h:32
void compileServers(RooAbsArg &arg, RooArgSet const &normSet)
T * compile(T &arg, RooAbsArg &owner, RooArgSet const &normSet)
A RooFormulaVar is a generic implementation of a real-valued object, which takes a RooArgList of serv...
Implements a universal generator context for all RooAbsPdf classes that do not have or need a special...
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition RooNameReg.h:39
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Variable that can be changed from the outside.
Definition RooRealVar.h:37
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.
bool isBasisSupported(const char *name) const
virtual RooAbsGenContext * modelGenContext(const RooAbsAnaConvPdf &, const RooArgSet &, const RooDataSet *, const RooArgSet *, bool) const
RooAbsRealLValue & convVar() const
Return the convolution variable of the resolution model.
const T & arg() const
Return reference to object held in proxy.
Implements a RooResolution model that corresponds to a delta function.
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:49
Basic string class.
Definition TString.h:138
int CompareTo(const char *cs, ECaseCompare cmp=kExact) const
Compare a string to char *cs2.
Definition TString.cxx:464
T * OwningPtr
An alias for raw pointers for indicating that the return type of a RooFit function is an owning point...
Definition Config.h:35
OwningPtr< T > makeOwningPtr(std::unique_ptr< T > &&ptr)
Internal helper to turn a std::unique_ptr<T> into an OwningPtr.
Definition Config.h:40