Logo ROOT  
Reference Guide
RooParamHistFunc.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * *
4 * This code was autogenerated by RooClassFactory *
5 *****************************************************************************/
6
7/** \class RooParamHistFunc
8 * \ingroup Roofit
9 * A histogram function that assigns scale parameters to every bin. Instead of the bare bin contents,
10 * it therefore yields:
11 * \f[
12 * \gamma_{i} * \mathrm{bin}_i
13 * \f]
14 *
15 * The \f$ \gamma_i \f$ can therefore be used to parametrise statistical uncertainties of the histogram
16 * template. In conjunction with a constraint term, this can be used to implement the Barlow-Beeston method.
17 * The constraint can be implemented using RooHistConstraint.
18 *
19 * See also the tutorial rf709_BarlowBeeston.C
20 */
21
22#include "Riostream.h"
23#include "RooParamHistFunc.h"
24#include "RooAbsReal.h"
25#include "RooAbsCategory.h"
26#include "RooRealVar.h"
27#include "RooHelpers.h"
28#include <math.h>
29#include "TMath.h"
30
31
32using namespace std ;
33
35
36////////////////////////////////////////////////////////////////////////////////
37/// Populate x with observables
38
39RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, RooDataHist& dh, bool paramRelative) :
40 RooAbsReal(name,title),
41 _x("x","x",this),
42 _p("p","p",this),
43 _dh(dh),
44 _relParam(paramRelative)
45{
46 _x.add(*_dh.get()) ;
47
48 // Now populate p with parameters
49 RooArgSet allVars ;
50 for (Int_t i=0 ; i<_dh.numEntries() ; i++) {
51 _dh.get(i) ;
52
53 const char* vname = Form("%s_gamma_bin_%i",GetName(),i) ;
54 RooRealVar* var = new RooRealVar(vname,vname,0,1000) ;
55 var->setVal(_relParam ? 1 : _dh.weight()) ;
56 var->setError(_relParam ? 1 / sqrt(_dh.weight()) : sqrt(_dh.weight()));
57 var->setConstant(true) ;
58 allVars.add(*var) ;
59 _p.add(*var) ;
60 }
61 addOwnedComponents(allVars) ;
62}
63
64////////////////////////////////////////////////////////////////////////////////
65/// Populate x with observables
66
67RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, const RooAbsArg& /*x*/, RooDataHist& dh, bool paramRelative) :
68 RooAbsReal(name,title),
69 _x("x","x",this),
70 _p("p","p",this),
71 _dh(dh),
72 _relParam(paramRelative)
73{
74 _x.add(*_dh.get()) ;
75
76 // Now populate p with parameters
77 RooArgSet allVars ;
78 for (Int_t i=0 ; i<_dh.numEntries() ; i++) {
79 _dh.get(i) ;
80 const char* vname = Form("%s_gamma_bin_%i",GetName(),i) ;
81 RooRealVar* var = new RooRealVar(vname,vname,0,1000) ;
82 var->setVal(_relParam ? 1 : _dh.weight()) ;
83 var->setError(_relParam ? 1 / sqrt(_dh.weight()) : sqrt(_dh.weight()));
84 var->setConstant(true) ;
85 allVars.add(*var) ;
86 _p.add(*var) ;
87 }
88 addOwnedComponents(allVars) ;
89}
90
91////////////////////////////////////////////////////////////////////////////////
92
93RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, RooDataHist& dh, const RooParamHistFunc& paramSource, bool paramRelative) :
94 RooAbsReal(name,title),
95 _x("x","x",this),
96 _p("p","p",this),
97 _dh(dh),
98 _relParam(paramRelative)
99{
100 // Populate x with observables
101 _x.add(*_dh.get()) ;
102
103 // Now populate p with existing parameters
104 _p.add(paramSource._p) ;
105}
106
107////////////////////////////////////////////////////////////////////////////////
108
110 RooAbsReal(other,name),
111 _x("x",this,other._x),
112 _p("p",this,other._p),
113 _dh(other._dh),
114 _relParam(other._relParam)
115{
116}
117
118////////////////////////////////////////////////////////////////////////////////
119
121{
122 Int_t idx = ((RooDataHist&)_dh).getIndex(_x,true) ;
123 double ret = ((RooAbsReal*)_p.at(idx))->getVal() ;
124 if (_relParam) {
125 double nom = getNominal(idx) ;
126 ret *= nom ;
127 }
128 return ret ;
129}
130
131////////////////////////////////////////////////////////////////////////////////
132
134{
135 return ((RooAbsReal&)_p[ibin]).getVal() ;
136}
137
138////////////////////////////////////////////////////////////////////////////////
139
140void RooParamHistFunc::setActual(Int_t ibin, double newVal)
141{
142 ((RooRealVar&)_p[ibin]).setVal(newVal) ;
143}
144
145////////////////////////////////////////////////////////////////////////////////
146
148{
149 _dh.get(ibin) ;
150 return _dh.weight() ;
151}
152
153////////////////////////////////////////////////////////////////////////////////
154
156{
157 _dh.get(ibin) ;
158 return _dh.weightError() ;
159}
160
161////////////////////////////////////////////////////////////////////////////////
162/// Return sampling hint for making curves of (projections) of this function
163/// as the recursive division strategy of RooCurve cannot deal efficiently
164/// with the vertical lines that occur in a non-interpolated histogram
165
166list<double>* RooParamHistFunc::plotSamplingHint(RooAbsRealLValue& obs, double xlo, double xhi) const
167{
168 // Check that observable is in dataset, if not no hint is generated
169 RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dh.get()->find(obs.GetName())) ;
170 if (!lvarg) {
171 return 0 ;
172 }
173
174 // Retrieve position of all bin boundaries
175 const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
176 double* boundaries = binning->array() ;
177
178 list<double>* hint = new list<double> ;
179
180 // Widen range slightly
181 xlo = xlo - 0.01*(xhi-xlo) ;
182 xhi = xhi + 0.01*(xhi-xlo) ;
183
184 double delta = (xhi-xlo)*1e-8 ;
185
186 // Construct array with pairs of points positioned epsilon to the left and
187 // right of the bin boundaries
188 for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
189 if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
190 hint->push_back(boundaries[i]-delta) ;
191 hint->push_back(boundaries[i]+delta) ;
192 }
193 }
194
195 return hint ;
196}
197
198////////////////////////////////////////////////////////////////////////////////
199/// Return sampling hint for making curves of (projections) of this function
200/// as the recursive division strategy of RooCurve cannot deal efficiently
201/// with the vertical lines that occur in a non-interpolated histogram
202
203std::list<double>* RooParamHistFunc::binBoundaries(RooAbsRealLValue& obs, double xlo, double xhi) const
204{
205 // Check that observable is in dataset, if not no hint is generated
206 RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dh.get()->find(obs.GetName())) ;
207 if (!lvarg) {
208 return 0 ;
209 }
210
211 // Retrieve position of all bin boundaries
212 const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
213 double* boundaries = binning->array() ;
214
215 list<double>* hint = new list<double> ;
216
217 // Construct array with pairs of points positioned epsilon to the left and
218 // right of the bin boundaries
219 for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
220 if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
221 hint->push_back(boundaries[i]) ;
222 }
223 }
224
225 return hint ;
226}
227
228////////////////////////////////////////////////////////////////////////////////
229/// Advertise that all integrals can be handled internally.
230
232 const RooArgSet* /*normSet*/, const char* /*rangeName*/) const
233{
234 // Simplest scenario, integrate over all dependents
235 RooAbsCollection *allVarsCommon = allVars.selectCommon(_x) ;
236 bool intAllObs = (allVarsCommon->getSize()==_x.getSize()) ;
237 delete allVarsCommon ;
238 if (intAllObs && matchArgs(allVars,analVars,_x)) {
239 return 1 ;
240 }
241
242 return 0 ;
243}
244
245////////////////////////////////////////////////////////////////////////////////
246/// Implement analytical integrations by doing appropriate weighting from component integrals
247/// functions to integrators of components
248
249double RooParamHistFunc::analyticalIntegralWN(Int_t code, const RooArgSet* /*normSet2*/,const char* rangeName) const
250{
251 // Supports only the scenario of integration over all dependents
252 R__ASSERT(code==1) ;
253
254 // The logic for summing over the histogram is borrowed from RooHistPdf with some differences:
255 //
256 // - a lambda function is used to inject the parameters for bin scaling into the RooDataHist::sum method
257 //
258 // - for simplicity, there is no check for the possibility of full-range integration with another overload of
259 // RooDataHist::sum
260 std::map<const RooAbsArg*, std::pair<double, double> > ranges;
261 for (const auto obs : _x) {
262 ranges[obs] = RooHelpers::getRangeOrBinningInterval(obs, rangeName);
263 }
264
265 auto getBinScale = [&](int iBin){ return static_cast<const RooAbsReal&>(_p[iBin]).getVal(); };
266
267 RooArgSet sliceSet{};
268 return const_cast<RooDataHist&>(_dh).sum(_x, sliceSet, true, false, ranges, getBinScale);
269}
#define e(i)
Definition: RSha256.hxx:103
#define ClassImp(name)
Definition: Rtypes.h:375
#define R__ASSERT(e)
Definition: TError.h:118
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:72
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2185
RooAbsBinning is the abstract base class for RooRealVar binning definitions.
Definition: RooAbsBinning.h:26
virtual Int_t numBoundaries() const =0
virtual double * array() const =0
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
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.
Abstract base class for objects that are lvalues, i.e.
Definition: RooAbsLValue.h:26
virtual const RooAbsBinning * getBinningPtr(const char *rangeName) const =0
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
void setConstant(bool value=true)
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
bool matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:57
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...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:45
double weight(std::size_t i) const
Return weight of i-th bin.
Definition: RooDataHist.h:110
void weightError(double &lo, double &hi, ErrorType etype=Poisson) const override
Return the asymmetric errors on the current weight.
Int_t numEntries() const override
Return the number of bins.
const RooArgSet * get() const override
Get bin centre of current bin.
Definition: RooDataHist.h:82
A histogram function that assigns scale parameters to every bin.
void setActual(Int_t ibin, double newVal)
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
double getNominalError(Int_t ibin) const
double getActual(Int_t ibin)
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const override
Implement analytical integrations by doing appropriate weighting from component integrals functions t...
double getNominal(Int_t ibin) const
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const override
Advertise that all integrals can be handled internally.
std::list< double > * binBoundaries(RooAbsRealLValue &, double, double) const override
Return sampling hint for making curves of (projections) of this function as the recursive division st...
std::list< double > * plotSamplingHint(RooAbsRealLValue &obs, double xlo, double xhi) const override
Return sampling hint for making curves of (projections) of this function as the recursive division st...
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:40
void setVal(double value) override
Set value of variable to 'value'.
Definition: RooRealVar.cxx:281
void setError(double value)
Definition: RooRealVar.h:65
const char * GetName() const override
Returns name of object.
Definition: TNamed.h:47
VecExpr< UnaryOp< Sqrt< T >, VecExpr< A, T, D >, T >, T, D > sqrt(const VecExpr< A, T, D > &rhs)
std::pair< double, double > getRangeOrBinningInterval(RooAbsArg const *arg, const char *rangeName)
Get the lower and upper bound of parameter range if arg can be casted to RooAbsRealLValue.
Definition: RooHelpers.cxx:168
static uint64_t sum(uint64_t i)
Definition: Factory.cxx:2345