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 
17 Class RooProjectedPdf is a RooAbsPdf implementation that represent a projection
18 of a given input p.d.f and the object returned by RooAbsPdf::createProjection.
19 The actual projection integral for it value and normalization are
20 calculated on the fly in getVal() once the normalization observables are known.
21 Class RooProjectedPdf can cache projected p.d.f.s for multiple normalization
22 observables simultaneously.
23 The createProjection() method of RooProjectedPdf is overloaded and will
24 return a new RooProjectedPdf that performs the projection of itself
25 and the requested additional projections in one integration step
26 The performance of <pre>f->createProjection(x)->createProjection(y)</pre>
27 is 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 
39 using 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 /// Evaluate projected p.d.f
90 
92 {
93  // Calculate current unnormalized value of object
94  int code ;
95  const RooAbsReal* proj = getProjection(&intobs, _normSet, 0, code);
96 
97  return proj->getVal() ;
98 }
99 
100 
101 
102 ////////////////////////////////////////////////////////////////////////////////
103 /// Retrieve object representing projection integral of input p.d.f
104 /// over observables iset, while normalizing over observables
105 /// nset. The code argument returned by reference is the unique code
106 /// defining this particular projection configuration
107 
108 const RooAbsReal* RooProjectedPdf::getProjection(const RooArgSet* iset, const RooArgSet* nset, const char* rangeName, int& code) const
109 {
110 
111  // Check if this configuration was created before
112  Int_t sterileIdx(-1) ;
113  CacheElem* cache = (CacheElem*) _cacheMgr.getObj(iset,nset,&sterileIdx,RooNameReg::ptr(rangeName)) ;
114  if (cache) {
115  code = _cacheMgr.lastIndex() ;
116  return static_cast<const RooAbsReal*>(cache->_projection);
117  }
118 
119  RooArgSet* nset2 = intpdf.arg().getObservables(*nset) ;
120 
121  if (iset) {
122  nset2->add(*iset) ;
123  }
124  RooAbsReal* proj = intpdf.arg().createIntegral(iset?*iset:RooArgSet(),nset2,0,rangeName) ;
125  delete nset2 ;
126 
127  cache = new CacheElem ;
128  cache->_projection = proj ;
129 
130  code = _cacheMgr.setObj(iset,nset,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
131 
132  coutI(Integration) << "RooProjectedPdf::getProjection(" << GetName() << ") creating new projection " << proj->GetName() << " with code " << code << endl ;
133 
134  return proj ;
135 }
136 
137 
138 
139 ////////////////////////////////////////////////////////////////////////////////
140 /// Special version of RooAbsReal::createProjection that deals with
141 /// projections of projections. Instead of integrating twice, a new
142 /// RooProjectedPdf is returned that is configured to perform the
143 /// complete integration in one step
144 
146 {
147  RooArgSet combiset(iset) ;
148  combiset.add(intobs) ;
149  return static_cast<RooAbsPdf&>( const_cast<RooAbsReal&>(intpdf.arg()) ).createProjection(combiset) ;
150 }
151 
152 
153 
154 ////////////////////////////////////////////////////////////////////////////////
155 /// Force RooRealIntegral to relegate integration of all observables to internal logic
156 
158 {
159  return kTRUE ;
160 }
161 
162 
163 
164 ////////////////////////////////////////////////////////////////////////////////
165 /// Mark all requested variables as internally integrated
166 
167 Int_t RooProjectedPdf::getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName) const
168 {
169  analVars.add(allVars) ;
170 
171  // Create the appropriate integral
172  int code ;
173  RooArgSet allVars2(allVars) ;
174  allVars2.add(intobs) ;
175  getProjection(&allVars2,normSet,rangeName,code) ;
176 
177  return code+1 ;
178 }
179 
180 
181 
182 ////////////////////////////////////////////////////////////////////////////////
183 /// Return analytical integral represent by appropriate element of projection cache
184 
185 Double_t RooProjectedPdf::analyticalIntegralWN(Int_t code, const RooArgSet* /*normSet*/, const char* rangeName) const
186 {
187  CacheElem *cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1) ;
188 
189  if (cache) {
190  Double_t ret= cache->_projection->getVal() ;
191  return ret ;
192  } else {
193 
194  RooArgSet* vars = getParameters(RooArgSet()) ;
195  vars->add(intobs) ;
196  RooArgSet* iset = _cacheMgr.nameSet1ByIndex(code-1)->select(*vars) ;
197  RooArgSet* nset = _cacheMgr.nameSet2ByIndex(code-1)->select(*vars) ;
198 
199  Int_t code2(-1) ;
200  const RooAbsReal* proj = getProjection(iset,nset,rangeName,code2) ;
201 
202  delete vars ;
203  delete nset ;
204  delete iset ;
205 
206  Double_t ret = proj->getVal() ;
207  return ret ;
208  }
209 
210 }
211 
212 
213 
214 ////////////////////////////////////////////////////////////////////////////////
215 /// No internal generator is implemented
216 
217 Int_t RooProjectedPdf::getGenerator(const RooArgSet& /*directVars*/, RooArgSet& /*generateVars*/, Bool_t /*staticInitOK*/) const
218  {
219  return 0 ;
220  }
221 
222 
223 
224 ////////////////////////////////////////////////////////////////////////////////
225 /// No internal generator is implemented
226 
228  {
229  return;
230  }
231 
232 
233 
234 ////////////////////////////////////////////////////////////////////////////////
235 /// Intercept a server redirection all and update list of dependents if necessary
236 /// Specifically update the set proxy 'deps' which introduces the dependency
237 /// on server value dirty flags of ourselves
238 
239 Bool_t RooProjectedPdf::redirectServersHook(const RooAbsCollection& newServerList, Bool_t /*mustReplaceAll*/,
240  Bool_t /*nameChange*/, Bool_t /*isRecursive*/)
241 {
242  // Redetermine explicit list of dependents if intPdf is being replaced
243  RooAbsArg* newPdf = newServerList.find(intpdf.arg().GetName()) ;
244  if (newPdf) {
245 
246  // Determine if set of dependens of new p.d.f is different from old p.d.f.
247  RooArgSet olddeps(deps) ;
248  RooArgSet* newdeps = newPdf->getParameters(intobs) ;
249  RooArgSet* common = (RooArgSet*) newdeps->selectCommon(deps) ;
250  newdeps->remove(*common,kTRUE,kTRUE) ;
251  olddeps.remove(*common,kTRUE,kTRUE) ;
252 
253  // If so, adjust composition of deps Listproxy
254  if (newdeps->getSize()>0) {
255  deps.add(*newdeps) ;
256  }
257  if (olddeps.getSize()>0) {
258  deps.remove(olddeps,kTRUE,kTRUE) ;
259  }
260 
261  delete common ;
262  delete newdeps ;
263  }
264 
265  return kFALSE ;
266 }
267 
268 
269 
270 ////////////////////////////////////////////////////////////////////////////////
271 /// Return RooAbsArg elements contained in projection cache element.
272 
274 {
275  RooArgList ret(*_projection) ;
276  return ret ;
277 }
278 
279 
280 
281 ////////////////////////////////////////////////////////////////////////////////
282 /// Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the
283 /// integration operation
284 
285 void RooProjectedPdf::printMetaArgs(ostream& os) const
286 {
287  os << "Int " << intpdf.arg().GetName() ;
288 
289  os << " d" ;
290  os << intobs ;
291  os << " " ;
292 
293 }
294 
295 
296 
297 
298 ////////////////////////////////////////////////////////////////////////////////
299 /// Print contents of cache when printing self as part of object tree
300 
301 void RooProjectedPdf::CacheElem::printCompactTreeHook(ostream& os, const char* indent, Int_t curElem, Int_t maxElem)
302 {
303  if (curElem==0) {
304  os << indent << "RooProjectedPdf begin projection cache" << endl ;
305  }
306 
307  TString indent2(indent) ;
308  indent2 += Form("[%d] ",curElem) ;
309 
310  _projection->printCompactTree(os,indent2) ;
311 
312  if(curElem==maxElem) {
313  os << indent << "RooProjectedPdf end projection cache" << endl ;
314  }
315 }
316 
317 
RooAbsPdf::_normSet
RooArgSet * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:321
RooCacheManager::setObj
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:49
RooAbsPdf::CacheElem
friend class CacheElem
Definition: RooAbsPdf.h:333
RooSetProxy::add
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...
Definition: RooSetProxy.cxx:165
RooProjectedPdf::getProjection
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,...
Definition: RooProjectedPdf.cxx:108
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:91
RooAbsReal.h
RooMsgService.h
RooNameReg::ptr
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
Definition: RooNameReg.cxx:93
RooFit.h
RooTemplateProxy::arg
const T & arg() const
Return reference to object held in proxy.
Definition: RooTemplateProxy.h:271
RooProjectedPdf::CacheElem::printCompactTreeHook
virtual void printCompactTreeHook(std::ostream &, const char *, Int_t, Int_t)
Print contents of cache when printing self as part of object tree.
Definition: RooProjectedPdf.cxx:301
RooCacheManager::getObj
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
Definition: RooCacheManager.h:44
RooAbsReal::createIntegral
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:561
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
Form
char * Form(const char *fmt,...)
RooProjectedPdf::createProjection
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Special version of RooAbsReal::createProjection that deals with projections of projections.
Definition: RooProjectedPdf.cxx:145
RooProjectedPdf::evaluate
Double_t evaluate() const
Evaluate projected p.d.f.
Definition: RooProjectedPdf.cxx:91
RooProjectedPdf::getGenerator
Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
No internal generator is implemented.
Definition: RooProjectedPdf.cxx:217
RooProjectedPdf::RooProjectedPdf
RooProjectedPdf()
Default constructor.
Definition: RooProjectedPdf.cxx:48
RooArgList
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooProjectedPdf::CacheElem::containedArgs
virtual RooArgList containedArgs(Action)
Return RooAbsArg elements contained in projection cache element.
Definition: RooProjectedPdf.cxx:273
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:91
RooCacheManager::nameSet2ByIndex
const RooNameSet * nameSet2ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
Definition: RooCacheManager.h:334
RooAbsCollection::remove
virtual Bool_t remove(const RooAbsArg &var, Bool_t silent=kFALSE, Bool_t matchByNameOnly=kFALSE)
Remove the specified argument from our list.
Definition: RooAbsCollection.cxx:584
RooAbsCollection::find
RooAbsArg * find(const char *name) const
Find object with given name in list.
Definition: RooAbsCollection.cxx:813
coutI
#define coutI(a)
Definition: RooMsgService.h:30
indent
static void indent(ostringstream &buf, int indent_level)
Definition: TClingCallFunc.cxx:87
RooArgSet::add
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
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
RooNameSet::select
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:185
RooAbsCacheElement
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
Definition: RooAbsCacheElement.h:26
TString
Basic string class.
Definition: TString.h:136
bool
RooProjectedPdf::generateEvent
void generateEvent(Int_t code)
No internal generator is implemented.
Definition: RooProjectedPdf.cxx:227
RooProjectedPdf.h
RooProjectedPdf::_cacheMgr
RooObjCacheManager _cacheMgr
Definition: RooProjectedPdf.h:62
RooAbsArg::getParameters
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:546
RooSetProxy::remove
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.
Definition: RooSetProxy.cxx:193
RooCacheManager::getObjByIndex
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Definition: RooCacheManager.h:308
RooProjectedPdf::getAnalyticalIntegralWN
virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Mark all requested variables as internally integrated.
Definition: RooProjectedPdf.cxx:167
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:92
RooAbsCollection::selectCommon
RooAbsCollection * selectCommon(const RooAbsCollection &refColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
Definition: RooAbsCollection.cxx:700
RooAbsCollection
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
Definition: RooAbsCollection.h:30
RooAbsArg::RooArgSet
friend class RooArgSet
Definition: RooAbsArg.h:581
RooCacheManager::nameSet1ByIndex
const RooNameSet * nameSet1ByIndex(Int_t index) const
Retrieve RooNameSet associated with slot at given index.
Definition: RooCacheManager.h:321
RooProjectedPdf::forceAnalyticalInt
virtual Bool_t forceAnalyticalInt(const RooAbsArg &dep) const
Force RooRealIntegral to relegate integration of all observables to internal logic.
Definition: RooProjectedPdf.cxx:157
RooProjectedPdf::CacheElem
Definition: RooProjectedPdf.h:53
RooRealVar.h
RooProjectedPdf::deps
RooSetProxy deps
Definition: RooProjectedPdf.h:51
RooProjectedPdf::CacheElem::_projection
RooAbsReal * _projection
Definition: RooProjectedPdf.h:55
RooProjectedPdf::printMetaArgs
void printMetaArgs(std::ostream &os) const
Customized printing of arguments of a RooRealIntegral to more intuitively reflect the contents of the...
Definition: RooProjectedPdf.cxx:285
Double_t
double Double_t
Definition: RtypesCore.h:59
RooAbsArg::getObservables
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Return the observables of this pdf given a set of observables.
Definition: RooAbsArg.h:294
RooAbsPdf::createProjection
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:3382
RooProjectedPdf::intobs
RooSetProxy intobs
Definition: RooProjectedPdf.h:50
RooAbsCacheElement::Action
Action
Definition: RooAbsCacheElement.h:39
RooProjectedPdf::redirectServersHook
Bool_t redirectServersHook(const RooAbsCollection &newServerList, Bool_t, Bool_t, Bool_t)
The cache manager.
Definition: RooProjectedPdf.cxx:239
name
char name[80]
Definition: TGX11.cxx:110
RooProjectedPdf::intpdf
RooRealProxy intpdf
Definition: RooProjectedPdf.h:49
RooCacheManager::lastIndex
Int_t lastIndex() const
Definition: RooCacheManager.h:65
RooAbsArg
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:73
RooAbsPdf
Definition: RooAbsPdf.h:40
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
RooFit::Integration
@ Integration
Definition: RooGlobalFunc.h:67
Riostream.h
RooProjectedPdf
Class RooProjectedPdf is a RooAbsPdf implementation that represent a projection of a given input p....
Definition: RooProjectedPdf.h:21
RooNameReg.h
RooAbsCollection::getSize
Int_t getSize() const
Definition: RooAbsCollection.h:171
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
int
RooProjectedPdf::analyticalIntegralWN
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.
Definition: RooProjectedPdf.cxx:185