Logo ROOT  
Reference Guide
RooProjectedPdf.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2 * Project: RooFit *
3 * *
4 * Copyright (c) 2000-2005, Regents of the University of California *
5 * and Stanford University. All rights reserved. *
6 * *
7 * Redistribution and use in source and binary forms, *
8 * with or without modification, are permitted according to the terms *
9 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
10 *****************************************************************************/
11
12/**
13\file RooProjectedPdf.cxx
14\class RooProjectedPdf
15\ingroup Roofitcore
16
17Class RooProjectedPdf is a RooAbsPdf implementation that represent a projection
18of a given input p.d.f and the object returned by RooAbsPdf::createProjection.
19The actual projection integral for it value and normalization are
20calculated on the fly in getVal() once the normalization observables are known.
21Class RooProjectedPdf can cache projected p.d.f.s for multiple normalization
22observables simultaneously.
23The createProjection() method of RooProjectedPdf is overloaded and will
24return a new RooProjectedPdf that performs the projection of itself
25and the requested additional projections in one integration step
26The performance of <pre>f->createProjection(x)->createProjection(y)</pre>
27is therefore identical to that of <pre>f->createProjection(RooArgSet(x,y))</pre>
28**/
29
30#include "Riostream.h"
31
32#include "RooProjectedPdf.h"
33#include "RooMsgService.h"
34#include "RooAbsReal.h"
35#include "RooRealVar.h"
36#include "RooNameReg.h"
37
38using namespace std;
39
41 ;
42
43
44////////////////////////////////////////////////////////////////////////////////
45/// Default constructor
46
48{
49}
50
51
52
53////////////////////////////////////////////////////////////////////////////////
54/// Construct projection of input pdf '_intpdf' over observables 'intObs'
55
56 RooProjectedPdf::RooProjectedPdf(const char *name, const char *title, RooAbsReal& _intpdf, const RooArgSet& intObs) :
57 RooAbsPdf(name,title),
58 intpdf("!IntegratedPdf","intpdf",this,_intpdf,false,false),
59 intobs("!IntegrationObservables","intobs",this,false,false),
60 deps("!Dependents","deps",this,true,true),
61 _cacheMgr(this,10)
62 {
63 intobs.add(intObs) ;
64
65 // Add all other dependens of projected p.d.f. directly
66 RooArgSet* tmpdeps = _intpdf.getParameters(intObs) ;
67 deps.add(*tmpdeps) ;
68 delete tmpdeps ;
69 }
70
71
72
73////////////////////////////////////////////////////////////////////////////////
74/// Copy constructor
75
77 RooAbsPdf(other,name),
78 intpdf("!IntegratedPdf",this,other.intpdf),
79 intobs("!IntegrationObservable",this,other.intobs),
80 deps("!Dependents",this,other.deps),
81 _cacheMgr(other._cacheMgr,this)
82{
83 }
84
85
86
87////////////////////////////////////////////////////////////////////////////////
88/// Evaluate projected p.d.f
89
91{
92 // Calculate current unnormalized value of object
93 int code ;
94 const RooAbsReal* proj = getProjection(&intobs, _normSet, 0, code);
95
96 return proj->getVal() ;
97}
98
99
100
101////////////////////////////////////////////////////////////////////////////////
102/// Retrieve object representing projection integral of input p.d.f
103/// over observables iset, while normalizing over observables
104/// nset. The code argument returned by reference is the unique code
105/// defining this particular projection configuration
106
107const RooAbsReal* RooProjectedPdf::getProjection(const RooArgSet* iset, const RooArgSet* nset, const char* rangeName, int& code) const
108{
109
110 // Check if this configuration was created before
111 Int_t sterileIdx(-1) ;
112 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(iset,nset,&sterileIdx,RooNameReg::ptr(rangeName)) ;
113 if (cache) {
114 code = _cacheMgr.lastIndex() ;
115 return static_cast<const RooAbsReal*>(cache->_projection);
116 }
117
118 RooArgSet* nset2 = intpdf.arg().getObservables(*nset) ;
119
120 if (iset) {
121 nset2->add(*iset) ;
122 }
123 RooAbsReal* proj = intpdf.arg().createIntegral(iset?*iset:RooArgSet(),nset2,0,rangeName) ;
124 delete nset2 ;
125
126 cache = new CacheElem ;
127 cache->_projection = proj ;
128
129 code = _cacheMgr.setObj(iset,nset,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
130
131 coutI(Integration) << "RooProjectedPdf::getProjection(" << GetName() << ") creating new projection " << proj->GetName() << " with code " << code << endl ;
132
133 return proj ;
134}
135
136
137
138////////////////////////////////////////////////////////////////////////////////
139/// Special version of RooAbsReal::createProjection that deals with
140/// projections of projections. Instead of integrating twice, a new
141/// RooProjectedPdf is returned that is configured to perform the
142/// complete integration in one step
143
145{
146 RooArgSet combiset(iset) ;
147 combiset.add(intobs) ;
148 return static_cast<RooAbsPdf&>( const_cast<RooAbsReal&>(intpdf.arg()) ).createProjection(combiset) ;
149}
150
151
152
153////////////////////////////////////////////////////////////////////////////////
154/// Force RooRealIntegral to relegate integration of all observables to internal logic
155
157{
158 return true ;
159}
160
161
162
163////////////////////////////////////////////////////////////////////////////////
164/// Mark all requested variables as internally integrated
165
166Int_t RooProjectedPdf::getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName) const
167{
168 analVars.add(allVars) ;
169
170 // Create the appropriate integral
171 int code ;
172 RooArgSet allVars2(allVars) ;
173 allVars2.add(intobs) ;
174 getProjection(&allVars2,normSet,rangeName,code) ;
175
176 return code+1 ;
177}
178
179
180
181////////////////////////////////////////////////////////////////////////////////
182/// Return analytical integral represent by appropriate element of projection cache
183
184double RooProjectedPdf::analyticalIntegralWN(Int_t code, const RooArgSet* /*normSet*/, const char* rangeName) const
185{
186 CacheElem *cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1) ;
187
188 if (cache) {
189 return cache->_projection->getVal() ;
190 } else {
191
192 std::unique_ptr<RooArgSet> vars{getParameters(RooArgSet())} ;
193 vars->add(intobs) ;
194 RooArgSet iset = _cacheMgr.selectFromSet1(*vars, code-1) ;
195 RooArgSet nset = _cacheMgr.selectFromSet2(*vars, code-1) ;
196
197 int code2 = -1 ;
198
199 return getProjection(&iset,&nset,rangeName,code2)->getVal() ;
200 }
201
202}
203
204
205
206////////////////////////////////////////////////////////////////////////////////
207/// No internal generator is implemented
208
209Int_t RooProjectedPdf::getGenerator(const RooArgSet& /*directVars*/, RooArgSet& /*generateVars*/, bool /*staticInitOK*/) const
210 {
211 return 0 ;
212 }
213
214
215
216////////////////////////////////////////////////////////////////////////////////
217/// No internal generator is implemented
218
220 {
221 return;
222 }
223
224
225
226////////////////////////////////////////////////////////////////////////////////
227/// Intercept a server redirection all and update list of dependents if necessary
228/// Specifically update the set proxy 'deps' which introduces the dependency
229/// on server value dirty flags of ourselves
230
231bool RooProjectedPdf::redirectServersHook(const RooAbsCollection& newServerList, bool /*mustReplaceAll*/,
232 bool /*nameChange*/, bool /*isRecursive*/)
233{
234 // Redetermine explicit list of dependents if intPdf is being replaced
235 RooAbsArg* newPdf = newServerList.find(intpdf.arg().GetName()) ;
236 if (newPdf) {
237
238 // Determine if set of dependens of new p.d.f is different from old p.d.f.
239 RooArgSet olddeps(deps) ;
240 RooArgSet* newdeps = newPdf->getParameters(intobs) ;
241 RooArgSet* common = (RooArgSet*) newdeps->selectCommon(deps) ;
242 newdeps->remove(*common,true,true) ;
243 olddeps.remove(*common,true,true) ;
244
245 // If so, adjust composition of deps Listproxy
246 if (newdeps->getSize()>0) {
247 deps.add(*newdeps) ;
248 }
249 if (olddeps.getSize()>0) {
250 deps.remove(olddeps,true,true) ;
251 }
252
253 delete common ;
254 delete newdeps ;
255 }
256
257 return false ;
258}
259
260
261
262////////////////////////////////////////////////////////////////////////////////
263/// Return RooAbsArg elements contained in projection cache element.
264
266{
267 RooArgList ret(*_projection) ;
268 return ret ;
269}
270
271
272
273////////////////////////////////////////////////////////////////////////////////
274/// Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the
275/// integration operation
276
277void RooProjectedPdf::printMetaArgs(ostream& os) const
278{
279 os << "Int " << intpdf.arg().GetName() ;
280
281 os << " d" ;
282 os << intobs ;
283 os << " " ;
284
285}
286
287
288
289
290////////////////////////////////////////////////////////////////////////////////
291/// Print contents of cache when printing self as part of object tree
292
293void RooProjectedPdf::CacheElem::printCompactTreeHook(ostream& os, const char* indent, Int_t curElem, Int_t maxElem)
294{
295 if (curElem==0) {
296 os << indent << "RooProjectedPdf begin projection cache" << endl ;
297 }
298
299 TString indent2(indent) ;
300 indent2 += Form("[%d] ",curElem) ;
301
302 _projection->printCompactTree(os,indent2) ;
303
304 if(curElem==maxElem) {
305 os << indent << "RooProjectedPdf end projection cache" << endl ;
306 }
307}
308
309
#define coutI(a)
Definition: RooMsgService.h:34
#define ClassImp(name)
Definition: Rtypes.h:375
static void indent(ostringstream &buf, int indent_level)
char name[80]
Definition: TGX11.cxx:110
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition: TString.cxx:2447
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:77
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:317
friend class RooArgSet
Definition: RooAbsArg.h:648
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:569
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
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.
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
RooAbsArg * find(const char *name) const
Find object with given name in list.
friend class CacheElem
The cache manager.
Definition: RooAbsPdf.h:366
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:354
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Return a p.d.f that represent a projection of this p.d.f integrated over given observables.
Definition: RooAbsPdf.cxx:3281
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
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:553
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Getter function without integration set.
RooArgSet selectFromSet1(RooArgSet const &argSet, int index) const
Create RooArgSet contatining the objects that are both in the cached set 1.
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
RooArgSet selectFromSet2(RooArgSet const &argSet, int index) const
Create RooArgSet contatining the objects that are both in the cached set 2.
Int_t lastIndex() const
Return index of slot used in last get or set operation.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Setter function without integration set.
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.
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:92
RooArgList containedArgs(Action) override
Return RooAbsArg elements contained in projection cache element.
void printCompactTreeHook(std::ostream &, const char *, Int_t, Int_t) override
Print contents of cache when printing self as part of object tree.
Class RooProjectedPdf is a RooAbsPdf implementation that represent a projection of a given input p....
RooRealProxy intpdf
p.d.f that is integrated
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const override
Return analytical integral represent by appropriate element of projection cache.
RooObjCacheManager _cacheMgr
! The cache manager
double evaluate() const override
Evaluate projected p.d.f.
bool forceAnalyticalInt(const RooAbsArg &dep) const override
Force RooRealIntegral to relegate integration of all observables to internal logic.
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const override
Mark all requested variables as internally integrated.
bool redirectServersHook(const RooAbsCollection &newServerList, bool, bool, bool) override
Intercept a server redirection all and update list of dependents if necessary Specifically update the...
void generateEvent(Int_t code) override
No internal generator is implemented.
RooAbsPdf * createProjection(const RooArgSet &iset) override
Special version of RooAbsReal::createProjection that deals with projections of projections.
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool staticInitOK=true) const override
No internal generator is implemented.
void printMetaArgs(std::ostream &os) const override
Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the...
RooProjectedPdf()
Default constructor.
const RooAbsReal * getProjection(const RooArgSet *iset, const RooArgSet *nset, const char *rangeName, int &code) const
Retrieve object representing projection integral of input p.d.f over observables iset,...
RooSetProxy intobs
observables that p.d.f is integrated over
RooSetProxy deps
dependents of this p.d.f
const T & arg() const
Return reference to object held in proxy.
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
Basic string class.
Definition: TString.h:136
@ Integration
Definition: RooGlobalFunc.h:63