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 "RooFit.h"
33#include "RooProjectedPdf.h"
34#include "RooMsgService.h"
35#include "RooAbsReal.h"
36#include "RooRealVar.h"
37#include "RooNameReg.h"
38
39using namespace std;
40
42 ;
43
44
45////////////////////////////////////////////////////////////////////////////////
46/// Default constructor
47
49{
50}
51
52
53
54////////////////////////////////////////////////////////////////////////////////
55/// Construct projection of input pdf '_intpdf' over observables 'intObs'
56
57 RooProjectedPdf::RooProjectedPdf(const char *name, const char *title, RooAbsReal& _intpdf, const RooArgSet& intObs) :
58 RooAbsPdf(name,title),
59 intpdf("!IntegratedPdf","intpdf",this,_intpdf,kFALSE,kFALSE),
60 intobs("!IntegrationObservables","intobs",this,kFALSE,kFALSE),
61 deps("!Dependents","deps",this,kTRUE,kTRUE),
62 _cacheMgr(this,10)
63 {
64 intobs.add(intObs) ;
65
66 // Add all other dependens of projected p.d.f. directly
67 RooArgSet* tmpdeps = _intpdf.getParameters(intObs) ;
68 deps.add(*tmpdeps) ;
69 delete tmpdeps ;
70 }
71
72
73
74////////////////////////////////////////////////////////////////////////////////
75/// Copy constructor
76
78 RooAbsPdf(other,name),
79 intpdf("!IntegratedPdf",this,other.intpdf),
80 intobs("!IntegrationObservable",this,other.intobs),
81 deps("!Dependents",this,other.deps),
82 _cacheMgr(other._cacheMgr,this)
83{
84 }
85
86
87
88////////////////////////////////////////////////////////////////////////////////
89/// Special version of getVal() overrides RooAbsReal::getValF() to save value of current normalization set
90
92{
93 _curNormSet = (RooArgSet*)set ;
94
95 return RooAbsPdf::getValV(set) ;
96}
97
98
99
100////////////////////////////////////////////////////////////////////////////////
101/// Evaluate projected p.d.f
102
104{
105 // Calculate current unnormalized value of object
106 int code ;
107 const RooAbsReal* proj = getProjection(&intobs,_curNormSet,0,code) ;
108
109 return proj->getVal() ;
110}
111
112
113
114////////////////////////////////////////////////////////////////////////////////
115/// Retrieve object representing projection integral of input p.d.f
116/// over observables iset, while normalizing over observables
117/// nset. The code argument returned by reference is the unique code
118/// defining this particular projection configuration
119
120const RooAbsReal* RooProjectedPdf::getProjection(const RooArgSet* iset, const RooArgSet* nset, const char* rangeName, int& code) const
121{
122
123 // Check if this configuration was created before
124 Int_t sterileIdx(-1) ;
125 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(iset,nset,&sterileIdx,RooNameReg::ptr(rangeName)) ;
126 if (cache) {
127 code = _cacheMgr.lastIndex() ;
128 return static_cast<const RooAbsReal*>(cache->_projection);
129 }
130
131 RooArgSet* nset2 = intpdf.arg().getObservables(*nset) ;
132
133 if (iset) {
134 nset2->add(*iset) ;
135 }
136 RooAbsReal* proj = intpdf.arg().createIntegral(iset?*iset:RooArgSet(),nset2,0,rangeName) ;
137 delete nset2 ;
138
139 cache = new CacheElem ;
140 cache->_projection = proj ;
141
142 code = _cacheMgr.setObj(iset,nset,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
143
144 coutI(Integration) << "RooProjectedPdf::getProjection(" << GetName() << ") creating new projection " << proj->GetName() << " with code " << code << endl ;
145
146 return proj ;
147}
148
149
150
151////////////////////////////////////////////////////////////////////////////////
152/// Special version of RooAbsReal::createProjection that deals with
153/// projections of projections. Instead of integrating twice, a new
154/// RooProjectedPdf is returned that is configured to perform the
155/// complete integration in one step
156
158{
159 RooArgSet combiset(iset) ;
160 combiset.add(intobs) ;
161 return static_cast<RooAbsPdf&>( const_cast<RooAbsReal&>(intpdf.arg()) ).createProjection(combiset) ;
162}
163
164
165
166////////////////////////////////////////////////////////////////////////////////
167/// Force RooRealIntegral to relegate integration of all observables to internal logic
168
170{
171 return kTRUE ;
172}
173
174
175
176////////////////////////////////////////////////////////////////////////////////
177/// Mark all requested variables as internally integrated
178
179Int_t RooProjectedPdf::getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName) const
180{
181 analVars.add(allVars) ;
182
183 // Create the appropriate integral
184 int code ;
185 RooArgSet allVars2(allVars) ;
186 allVars2.add(intobs) ;
187 getProjection(&allVars2,normSet,rangeName,code) ;
188
189 return code+1 ;
190}
191
192
193
194////////////////////////////////////////////////////////////////////////////////
195/// Return analytical integral represent by appropriate element of projection cache
196
197Double_t RooProjectedPdf::analyticalIntegralWN(Int_t code, const RooArgSet* /*normSet*/, const char* rangeName) const
198{
199 CacheElem *cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1) ;
200
201 if (cache) {
202 Double_t ret= cache->_projection->getVal() ;
203 return ret ;
204 } else {
205
207 vars->add(intobs) ;
208 RooArgSet* iset = _cacheMgr.nameSet1ByIndex(code-1)->select(*vars) ;
209 RooArgSet* nset = _cacheMgr.nameSet2ByIndex(code-1)->select(*vars) ;
210
211 Int_t code2(-1) ;
212 const RooAbsReal* proj = getProjection(iset,nset,rangeName,code2) ;
213
214 delete vars ;
215 delete nset ;
216 delete iset ;
217
218 Double_t ret = proj->getVal() ;
219 return ret ;
220 }
221
222}
223
224
225
226////////////////////////////////////////////////////////////////////////////////
227/// No internal generator is implemented
228
229Int_t RooProjectedPdf::getGenerator(const RooArgSet& /*directVars*/, RooArgSet& /*generateVars*/, Bool_t /*staticInitOK*/) const
230 {
231 return 0 ;
232 }
233
234
235
236////////////////////////////////////////////////////////////////////////////////
237/// No internal generator is implemented
238
240 {
241 return;
242 }
243
244
245
246////////////////////////////////////////////////////////////////////////////////
247/// Intercept a server redirection all and update list of dependents if necessary
248/// Specifically update the set proxy 'deps' which introduces the dependency
249/// on server value dirty flags of ourselves
250
251Bool_t RooProjectedPdf::redirectServersHook(const RooAbsCollection& newServerList, Bool_t /*mustReplaceAll*/,
252 Bool_t /*nameChange*/, Bool_t /*isRecursive*/)
253{
254 // Redetermine explicit list of dependents if intPdf is being replaced
255 RooAbsArg* newPdf = newServerList.find(intpdf.arg().GetName()) ;
256 if (newPdf) {
257
258 // Determine if set of dependens of new p.d.f is different from old p.d.f.
259 RooArgSet olddeps(deps) ;
260 RooArgSet* newdeps = newPdf->getParameters(intobs) ;
261 RooArgSet* common = (RooArgSet*) newdeps->selectCommon(deps) ;
262 newdeps->remove(*common,kTRUE,kTRUE) ;
263 olddeps.remove(*common,kTRUE,kTRUE) ;
264
265 // If so, adjust composition of deps Listproxy
266 if (newdeps->getSize()>0) {
267 deps.add(*newdeps) ;
268 }
269 if (olddeps.getSize()>0) {
270 deps.remove(olddeps,kTRUE,kTRUE) ;
271 }
272
273 delete common ;
274 delete newdeps ;
275 }
276
277 return kFALSE ;
278}
279
280
281
282////////////////////////////////////////////////////////////////////////////////
283/// Return RooAbsArg elements contained in projection cache element.
284
286{
287 RooArgList ret(*_projection) ;
288 return ret ;
289}
290
291
292
293////////////////////////////////////////////////////////////////////////////////
294/// Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the
295/// integration operation
296
297void RooProjectedPdf::printMetaArgs(ostream& os) const
298{
299 os << "Int " << intpdf.arg().GetName() ;
300
301 os << " d" ;
302 os << intobs ;
303 os << " " ;
304
305}
306
307
308
309
310////////////////////////////////////////////////////////////////////////////////
311/// Print contents of cache when printing self as part of object tree
312
313void RooProjectedPdf::CacheElem::printCompactTreeHook(ostream& os, const char* indent, Int_t curElem, Int_t maxElem)
314{
315 if (curElem==0) {
316 os << indent << "RooProjectedPdf begin projection cache" << endl ;
317 }
318
319 TString indent2(indent) ;
320 indent2 += Form("[%d] ",curElem) ;
321
322 _projection->printCompactTree(os,indent2) ;
323
324 if(curElem==maxElem) {
325 os << indent << "RooProjectedPdf end projection cache" << endl ;
326 }
327}
328
329
#define coutI(a)
Definition: RooMsgService.h:31
int Int_t
Definition: RtypesCore.h:41
const Bool_t kFALSE
Definition: RtypesCore.h:88
bool Bool_t
Definition: RtypesCore.h:59
double Double_t
Definition: RtypesCore.h:55
const Bool_t kTRUE
Definition: RtypesCore.h:87
#define ClassImp(name)
Definition: Rtypes.h:365
static void indent(ostringstream &buf, int indent_level)
char name[80]
Definition: TGX11.cxx:109
char * Form(const char *fmt,...)
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:71
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:240
friend class RooArgSet
Definition: RooAbsArg.h:551
RooArgSet * getParameters(const RooAbsData *data, Bool_t stripDisconnected=kTRUE) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
Definition: RooAbsArg.cxx:548
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.
RooAbsCollection * selectCommon(const RooAbsCollection &refColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
Int_t getSize() const
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
RooAbsArg * find(const char *name) const
Find object with given name in list.
friend class CacheElem
Definition: RooAbsPdf.h:334
virtual Double_t getValV(const RooArgSet *set=0) const
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:281
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:3362
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:59
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:87
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:562
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
virtual Bool_t add(const RooAbsCollection &col, Bool_t silent=kFALSE)
Add a collection of arguments to this collection by calling add() for each element in the source coll...
Definition: RooArgSet.h:88
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
const RooNameSet * nameSet2ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Int_t lastIndex() const
const RooNameSet * nameSet1ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooArgSet * select(const RooArgSet &list) const
Construct a RooArgSet of objects in input 'list' whose names match to those in the internal name list...
Definition: RooNameSet.cxx:187
virtual RooArgList containedArgs(Action)
Return RooAbsArg elements contained in projection cache element.
virtual void printCompactTreeHook(std::ostream &, const char *, Int_t, Int_t)
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....
void generateEvent(Int_t code)
No internal generator is implemented.
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the...
RooRealProxy intpdf
virtual Bool_t forceAnalyticalInt(const RooAbsArg &dep) const
Force RooRealIntegral to relegate integration of all observables to internal logic.
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
No internal generator is implemented.
virtual Double_t getValV(const RooArgSet *set=0) const
Special version of getVal() overrides RooAbsReal::getValF() to save value of current normalization se...
RooObjCacheManager _cacheMgr
RooArgSet * _curNormSet
virtual Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return analytical integral represent by appropriate element of projection cache.
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Mark all requested variables as internally integrated.
Double_t evaluate() const
Evaluate projected p.d.f.
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,...
Bool_t redirectServersHook(const RooAbsCollection &newServerList, Bool_t, Bool_t, Bool_t)
The cache manager.
RooSetProxy intobs
RooSetProxy deps
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Special version of RooAbsReal::createProjection that deals with projections of projections.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Overloaded RooArgSet::add() method inserts 'var' into set and registers 'var' as server to owner with...
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove object 'var' from set and deregister 'var' as server to owner.
const T & arg() const
Return reference to object held in proxy.
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Basic string class.
Definition: TString.h:131
@ Integration
Definition: RooGlobalFunc.h:67