ROOT  6.06/09
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 // Your description goes here...
8 
9 #include "Riostream.h"
10 #include "RooParamHistFunc.h"
11 #include "RooAbsReal.h"
12 #include "RooAbsCategory.h"
13 #include "RooRealVar.h"
14 #include <math.h>
15 #include "TMath.h"
16 
17 
18 using namespace std ;
19 
20 
22  ;
23 
24 
25 ////////////////////////////////////////////////////////////////////////////////
26 /// Populate x with observables
27 
28 RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, RooDataHist& dh, Bool_t paramRelative) :
29  RooAbsReal(name,title),
30  _x("x","x",this),
31  _p("p","p",this),
32  _dh(dh),
33  _relParam(paramRelative)
34 {
35  _x.add(*_dh.get()) ;
36 
37  // Now populate p with parameters
38  RooArgSet allVars ;
39  for (Int_t i=0 ; i<_dh.numEntries() ; i++) {
40  _dh.get(i) ;
41  const char* vname = Form("%s_gamma_bin_%i",GetName(),i) ;
42  RooRealVar* var = new RooRealVar(vname,vname,0,1000) ;
43  var->setVal(_relParam ? 1 : _dh.weight()) ;
44  var->setConstant(kTRUE) ;
45  allVars.add(*var) ;
46  _p.add(*var) ;
47  }
48  addOwnedComponents(allVars) ;
49 }
50 
51 
52 
53 ////////////////////////////////////////////////////////////////////////////////
54 /// Populate x with observables
55 
56 RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, const RooAbsArg& /*x*/, RooDataHist& dh, Bool_t paramRelative) :
57  RooAbsReal(name,title),
58  _x("x","x",this),
59  _p("p","p",this),
60  _dh(dh),
61  _relParam(paramRelative)
62 {
63  _x.add(*_dh.get()) ;
64 
65  // Now populate p with parameters
66  RooArgSet allVars ;
67  for (Int_t i=0 ; i<_dh.numEntries() ; i++) {
68  _dh.get(i) ;
69  const char* vname = Form("%s_gamma_bin_%i",GetName(),i) ;
70  RooRealVar* var = new RooRealVar(vname,vname,0,1000) ;
71  var->setVal(_relParam ? 1 : _dh.weight()) ;
72  var->setConstant(kTRUE) ;
73  allVars.add(*var) ;
74  _p.add(*var) ;
75  }
76  addOwnedComponents(allVars) ;
77 }
78 
79 
80 
81 ////////////////////////////////////////////////////////////////////////////////
82 
83 RooParamHistFunc::RooParamHistFunc(const char *name, const char *title, RooDataHist& dh, const RooParamHistFunc& paramSource, Bool_t paramRelative) :
84  RooAbsReal(name,title),
85  _x("x","x",this),
86  _p("p","p",this),
87  _dh(dh),
88  _relParam(paramRelative)
89 {
90  // Populate x with observables
91  _x.add(*_dh.get()) ;
92 
93  // Now populate p with existing parameters
94  _p.add(paramSource._p) ;
95 }
96 
97 
98 
99 ////////////////////////////////////////////////////////////////////////////////
100 
102  RooAbsReal(other,name),
103  _x("x",this,other._x),
104  _p("p",this,other._p),
105  _dh(other._dh),
106  _relParam(other._relParam)
107 {
108 }
109 
110 
111 
112 ////////////////////////////////////////////////////////////////////////////////
113 
115 {
116  Int_t idx = ((RooDataHist&)_dh).getIndex(_x,kTRUE) ;
117  Double_t ret = ((RooAbsReal*)_p.at(idx))->getVal() ;
118  if (_relParam) {
119  Double_t nom = getNominal(idx) ;
120  ret *= nom ;
121  }
122  return ret ;
123 }
124 
125 
126 
127 ////////////////////////////////////////////////////////////////////////////////
128 
130 {
131  return ((RooAbsReal&)_p[ibin]).getVal() ;
132 }
133 
134 
135 
136 ////////////////////////////////////////////////////////////////////////////////
137 
139 {
140  ((RooRealVar&)_p[ibin]).setVal(newVal) ;
141 }
142 
143 
144 
145 ////////////////////////////////////////////////////////////////////////////////
146 
148 {
149  _dh.get(ibin) ;
150  return _dh.weight() ;
151 }
152 
153 
154 ////////////////////////////////////////////////////////////////////////////////
155 
157 {
158  _dh.get(ibin) ;
159  return _dh.weightError() ;
160 }
161 
162 
163 
164 ////////////////////////////////////////////////////////////////////////////////
165 /// Return sampling hint for making curves of (projections) of this function
166 /// as the recursive division strategy of RooCurve cannot deal efficiently
167 /// with the vertical lines that occur in a non-interpolated histogram
168 
170 {
171  // Check that observable is in dataset, if not no hint is generated
172  RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dh.get()->find(obs.GetName())) ;
173  if (!lvarg) {
174  return 0 ;
175  }
176 
177  // Retrieve position of all bin boundaries
178  const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
179  Double_t* boundaries = binning->array() ;
180 
181  list<Double_t>* hint = new list<Double_t> ;
182 
183  // Widen range slighty
184  xlo = xlo - 0.01*(xhi-xlo) ;
185  xhi = xhi + 0.01*(xhi-xlo) ;
186 
187  Double_t delta = (xhi-xlo)*1e-8 ;
188 
189  // Construct array with pairs of points positioned epsilon to the left and
190  // right of the bin boundaries
191  for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
192  if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
193  hint->push_back(boundaries[i]-delta) ;
194  hint->push_back(boundaries[i]+delta) ;
195  }
196  }
197 
198  return hint ;
199 }
200 
201 
202 ////////////////////////////////////////////////////////////////////////////////
203 /// Return sampling hint for making curves of (projections) of this function
204 /// as the recursive division strategy of RooCurve cannot deal efficiently
205 /// with the vertical lines that occur in a non-interpolated histogram
206 
207 std::list<Double_t>* RooParamHistFunc::binBoundaries(RooAbsRealLValue& obs, Double_t xlo, Double_t xhi) const
208 {
209  // Check that observable is in dataset, if not no hint is generated
210  RooAbsLValue* lvarg = dynamic_cast<RooAbsLValue*>(_dh.get()->find(obs.GetName())) ;
211  if (!lvarg) {
212  return 0 ;
213  }
214 
215  // Retrieve position of all bin boundaries
216  const RooAbsBinning* binning = lvarg->getBinningPtr(0) ;
217  Double_t* boundaries = binning->array() ;
218 
219  list<Double_t>* hint = new list<Double_t> ;
220 
221  // Construct array with pairs of points positioned epsilon to the left and
222  // right of the bin boundaries
223  for (Int_t i=0 ; i<binning->numBoundaries() ; i++) {
224  if (boundaries[i]>=xlo && boundaries[i]<=xhi) {
225  hint->push_back(boundaries[i]) ;
226  }
227  }
228 
229  return hint ;
230 }
231 
232 
233 
234 ////////////////////////////////////////////////////////////////////////////////
235 /// Advertise that all integrals can be handled internally.
236 
238  const RooArgSet* /*normSet*/, const char* /*rangeName*/) const
239 {
240  // Simplest scenario, integrate over all dependents
241  RooAbsCollection *allVarsCommon = allVars.selectCommon(_x) ;
242  Bool_t intAllObs = (allVarsCommon->getSize()==_x.getSize()) ;
243  delete allVarsCommon ;
244  if (intAllObs && matchArgs(allVars,analVars,_x)) {
245  return 1 ;
246  }
247 
248  return 0 ;
249 }
250 
251 
252 
253 
254 ////////////////////////////////////////////////////////////////////////////////
255 /// Implement analytical integrations by doing appropriate weighting from component integrals
256 /// functions to integrators of components
257 
258 Double_t RooParamHistFunc::analyticalIntegralWN(Int_t code, const RooArgSet* /*normSet2*/,const char* /*rangeName*/) const
259 {
260  R__ASSERT(code==1) ;
261 
263  RooAbsReal* p ;
264  Double_t ret(0) ;
265  Int_t i(0) ;
266  while((p=(RooAbsReal*)iter.next())) {
267  Double_t bin = p->getVal() ;
268  if (_relParam) bin *= getNominal(i++) ;
269  ret += bin ;
270  }
271 
272  // WVE fix this!!! Assume uniform binning for now!
273  RooFIter xiter = _x.fwdIterator() ;
274  RooAbsArg* obs ;
275  Double_t binV(1) ;
276  while((obs=xiter.next())) {
277  RooRealVar* xx = (RooRealVar*) obs ;
278  binV *= (xx->getMax()-xx->getMin())/xx->numBins() ;
279  }
280 
281  return ret*binV ;
282 }
283 
284 
virtual std::list< Double_t > * plotSamplingHint(RooAbsRealLValue &obs, Double_t xlo, Double_t xhi) const
Return sampling hint for making curves of (projections) of this function as the recursive division st...
virtual Int_t numBins(const char *rangeName=0) const
RooAbsCollection * selectCommon(const RooAbsCollection &refColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
RooFIter fwdIterator() const
Double_t getActual(Int_t ibin)
#define R__ASSERT(e)
Definition: TError.h:98
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Implement analytical integrations by doing appropriate weighting from component integrals functions t...
virtual Double_t getMin(const char *name=0) const
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
Double_t evaluate() const
Bool_t addOwnedComponents(const RooArgSet &comps)
Take ownership of the contents of 'comps'.
Definition: RooAbsArg.cxx:2281
STL namespace.
Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &analVars, const RooArgSet *normSet, const char *rangeName=0) const
Advertise that all integrals can be handled internally.
virtual Int_t numBoundaries() const =0
virtual const RooArgSet * get() const
Definition: RooDataHist.h:77
ClassImp(RooParamHistFunc)
virtual const RooAbsBinning * getBinningPtr(const char *rangeName) const =0
std::map< std::string, std::string >::const_iterator iter
Definition: TAlienJob.cxx:54
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
virtual void weightError(Double_t &lo, Double_t &hi, ErrorType etype=Poisson) const
Return the error on current weight.
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:202
RooAbsArg * find(const char *name) const
Find object with given name in list.
virtual Double_t weight() const
Definition: RooDataHist.h:96
void setConstant(Bool_t value=kTRUE)
char * Form(const char *fmt,...)
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
RooAbsArg * next()
RooAbsArg * at(Int_t idx) const
Definition: RooArgList.h:84
double Double_t
Definition: RtypesCore.h:55
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
virtual std::list< Double_t > * binBoundaries(RooAbsRealLValue &, Double_t, Double_t) const
Return sampling hint for making curves of (projections) of this function as the recursive division st...
virtual Int_t numEntries() const
Return the number of bins.
Double_t getNominalError(Int_t ibin) const
#define name(a, b)
Definition: linkTestLib0.cpp:5
virtual Double_t getMax(const char *name=0) const
void setActual(Int_t ibin, Double_t newVal)
Double_t getNominal(Int_t ibin) const
virtual Double_t * array() const =0
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Int_t getSize() const
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
const Bool_t kTRUE
Definition: Rtypes.h:91
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add element to non-owning set.
Definition: RooArgSet.cxx:448