Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooEffProd.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, NIKHEF
7 * GR, Gerhard Raven, NIKHEF/VU *
8 * *
9 * Redistribution and use in source and binary forms, *
10 * with or without modification, are permitted according to the terms *
11 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
12 *****************************************************************************/
13
14
15/////////////////////////////////////////////////////////////////////////////////////
16/// \class RooEffProd
17/// The class RooEffProd implements the product of a PDF with an efficiency function.
18/// The normalization integral of the product is calculated numerically, but the
19/// event generation is handled by a specialized generator context that implements
20/// the event generation in a more efficient for cases where the PDF has an internal
21/// generator that is smarter than accept reject.
22///
23
24#include "RooFit.h"
25#include "RooEffProd.h"
26#include "RooEffGenContext.h"
27#include "RooNameReg.h"
28#include "RooRealVar.h"
29
30using namespace std;
31
33 ;
34
35
36
37////////////////////////////////////////////////////////////////////////////////
38/// Constructor of a a production of p.d.f inPdf with efficiency
39/// function inEff.
40
41RooEffProd::RooEffProd(const char *name, const char *title,
42 RooAbsPdf& inPdf, RooAbsReal& inEff) :
43 RooAbsPdf(name,title),
44 _cacheMgr(this,10),
45 _pdf("pdf","pre-efficiency pdf", this,inPdf),
46 _eff("eff","efficiency function",this,inEff),
47 _nset(0),
48 _fixedNset(0)
49{
50}
51
52
53
54
55////////////////////////////////////////////////////////////////////////////////
56/// Copy constructor
57
58RooEffProd::RooEffProd(const RooEffProd& other, const char* name) :
59 RooAbsPdf(other, name),
60 _cacheMgr(other._cacheMgr,this),
61 _pdf("pdf",this,other._pdf),
62 _eff("acc",this,other._eff),
63 _nset(0),
64 _fixedNset(0)
65{
66}
67
68
69
70
71////////////////////////////////////////////////////////////////////////////////
72/// Destructor
73
75{
76}
77
78
79
80
81////////////////////////////////////////////////////////////////////////////////
82/// Calculate and return 'raw' unnormalized value of p.d.f
83
85{
86 return eff()->getVal() * pdf()->getVal(_fixedNset ? _fixedNset : _normSet);
87}
88
89
90
91////////////////////////////////////////////////////////////////////////////////
92/// Return specialized generator context for RooEffProds that implements generation
93/// in a more efficient way than can be done for generic correlated products
94
96 const RooArgSet* auxProto, Bool_t verbose) const
97{
98 assert(pdf()!=0);
99 assert(eff()!=0);
100 return new RooEffGenContext(*this,*pdf(),*eff(),vars,prototype,auxProto,verbose) ;
101}
102
103
104
105
106
107////////////////////////////////////////////////////////////////////////////////
108/// Return internal integration capabilities of the p.d.f. Given a set 'allVars' for which
109/// integration is requested, returned the largest subset for which internal (analytical)
110/// integration is implemented (in argument analVars). The return value is a unique integer
111/// code that identifies the integration configuration (integrated observables and range name).
112///
113/// This implementation in RooEffProd catches all integrals without normalization and reroutes them
114/// through a custom integration routine that properly accounts for the use of normalized p.d.f.
115/// in the evaluate() expression, which breaks the default RooAbsPdf normalization handling
116
118 const RooArgSet* normSet, const char* rangeName) const
119{
120
121 // No special handling required if a normalization set is given
122 if (normSet && normSet->getSize()>0) {
123 return 0 ;
124 }
125 // No special handling required if running with a fixed normalization set
126 if (_fixedNset) {
127 return 0 ;
128 }
129
130 // Special handling of integrals w/o normalization set. We need to pass _a_ normalization set
131 // to pdf().getVal(nset) in evaluate() because for certain p.d.fs the shape depends on the
132 // chosen normalization set. This functions correctly automatically for plain getVal() calls,
133 // however when the (numeric) normalization is calculated, getVal() is called without any normalization
134 // set causing the normalization to be calculated for a possibly different shape. To fix that
135 // integrals over a RooEffProd without normalization setup are explicitly handled here. These integrals
136 // are calculated using a clone of the RooEffProd that has a fixed normalization set passed to the
137 // underlying p.d.f regardless of the normalization set passed to getVal(). Here the set of observables
138 // over which is integrated is passed.
139
140 // Declare that we can analytically integrate all requested observables
141 analVars.add(allVars) ;
142
143 // Check if cache was previously created
144 Int_t sterileIndex(-1) ;
145 CacheElem* cache = (CacheElem*) _cacheMgr.getObj(&allVars,&allVars,&sterileIndex,RooNameReg::ptr(rangeName)) ;
146 if (cache) {
147 return _cacheMgr.lastIndex()+1;
148 }
149
150 // Construct cache with clone of p.d.f that has fixed normalization set that is passed to input pdf
151 cache = new CacheElem ;
152 cache->_intObs.addClone(allVars) ;
153 cache->_clone = (RooEffProd*) clone(Form("%s_clone",GetName())) ;
154 cache->_clone->_fixedNset = &cache->_intObs ;
155 cache->_int = cache->_clone->createIntegral(cache->_intObs,rangeName) ;
156
157 // Store cache and return index as code
158 Int_t code = _cacheMgr.setObj(&allVars,&allVars,(RooAbsCacheElement*)cache,RooNameReg::ptr(rangeName)) ;
159
160 return code+1 ;
161}
162
163
164
165
166
167////////////////////////////////////////////////////////////////////////////////
168/// Return value of integral identified by code, which should be a return value of getAnalyticalIntegralWN,
169/// Code zero is always handled and signifies no integration (return value is normalized p.d.f. value)
170
171Double_t RooEffProd::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* /*rangeName*/) const
172{
173 // Return analytical integral defined by given scenario code
174
175 // No integration scenario
176 if (code==0) {
177 return getVal(normSet) ;
178 }
179
180 // Partial integration scenarios
181 CacheElem* cache = (CacheElem*) _cacheMgr.getObjByIndex(code-1) ;
182
183 return cache->_int->getVal() ;
184}
185
186
187
188
189////////////////////////////////////////////////////////////////////////////////
190/// Report all RooAbsArg derived objects contained in Cache Element (used in function optimization and
191/// and server redirect management of the cache)
192
194{
195 RooArgList ret(_intObs) ;
196 ret.add(*_int) ;
197 ret.add(*_clone) ;
198 return ret ;
199}
#define ClassImp(name)
Definition Rtypes.h:364
char name[80]
Definition TGX11.cxx:110
char * Form(const char *fmt,...)
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
Int_t getSize() const
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add the specified argument to list.
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition RooAbsPdf.h:322
friend class RooEffGenContext
Definition RooAbsPdf.h:297
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:61
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:91
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
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:29
Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add element to non-owning set.
RooAbsArg * addClone(const RooAbsArg &var, Bool_t silent=kFALSE) override
Add clone of specified element to an owning set.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=0, const TNamed *isetRangeName=0)
T * getObjByIndex(Int_t index) const
Retrieve payload object by slot index.
Int_t lastIndex() const
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=0)
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:33
RooEffProd * _clone
Definition RooEffProd.h:60
RooAbsReal * _int
Definition RooEffProd.h:61
virtual RooArgList containedArgs(Action)
Report all RooAbsArg derived objects contained in Cache Element (used in function optimization and an...
The class RooEffProd implements the product of a PDF with an efficiency function.
Definition RooEffProd.h:20
virtual ~RooEffProd()
Destructor.
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype, const RooArgSet *auxProto, Bool_t verbose) const
Return specialized generator context for RooEffProds that implements generation in a more efficient w...
virtual Double_t evaluate() const
Calculate and return 'raw' unnormalized value of p.d.f.
const RooAbsPdf * pdf() const
Definition RooEffProd.h:42
RooArgSet * _fixedNset
Normalization set to be used in evaluation.
Definition RooEffProd.h:73
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, const char *rangeName=0) const
Return internal integration capabilities of the p.d.f.
virtual TObject * clone(const char *newname) const
Definition RooEffProd.h:28
RooObjCacheManager _cacheMgr
Definition RooEffProd.h:65
const RooAbsReal * eff() const
Definition RooEffProd.h:46
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Return value of integral identified by code, which should be a return value of getAnalyticalIntegralW...
static const TNamed * ptr(const char *stringPtr)
Return a unique TNamed pointer for given C++ string.
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47