Logo ROOT   6.10/09
Reference Guide
RooParametricStepFunction.cxx
Go to the documentation of this file.
1 /*****************************************************************************
2  * Project: RooFit *
3  * Package: RooFitBabar *
4  * @(#)root/roofit:$Id$
5  * Authors: *
6  * Aaron Roodman, Stanford Linear Accelerator Center, Stanford University *
7  * *
8  * Copyright (c) 2004, Stanford University. All rights reserved. *
9  *
10  * Redistribution and use in source and binary forms, *
11  * with or without modification, are permitted according to the terms *
12  * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
13  *****************************************************************************/
14 
15 /** \class RooParametricStepFunction
16  \ingroup Roofit
17 
18 The Parametric Step Function PDF is a binned distribution whose parameters
19 are the heights of each bin. This PDF was first used in BaBar's B0->pi0pi0
20 paper BABAR Collaboration (B. Aubert et al.) Phys.Rev.Lett.91:241801,2003.
21 
22 This PDF may be used to describe oddly shaped distributions. It differs
23 from a RooKeysPdf or a RooHistPdf in that a RooParametricStepFunction
24 has free parameters. In particular, any statistical uncertainty in
25 sample used to model this PDF may be understood with these free parameters;
26 this is not possible with non-parametric PDFs.
27 
28 The RooParametricStepFunction has Nbins-1 free parameters. Note that
29 the limits of the dependent variable must match the low and hi bin limits.
30 
31 An example of usage is:
32 
33 ~~~ {.cpp}
34 Int_t nbins(10);
35 TArrayD limits(nbins+1);
36 limits[0] = 0.0; //etc...
37 RooArgList* list = new RooArgList("list");
38 RooRealVar* binHeight0 = new RooRealVar("binHeight0","bin 0 Value",0.1,0.0,1.0);
39 list->add(binHeight0); // up to binHeight8, ie. 9 parameters
40 
41 RooParametricStepFunction aPdf = ("aPdf","PSF",*x,
42  *list,limits,nbins);
43 ~~~
44 */
45 
46 #include "RooFit.h"
47 
48 #include "Riostream.h"
49 #include "TArrayD.h"
50 #include <math.h>
51 
53 #include "RooAbsReal.h"
54 #include "RooRealVar.h"
55 #include "RooArgList.h"
56 
57 #include "TError.h"
58 
59 using namespace std;
60 
62 
63 ////////////////////////////////////////////////////////////////////////////////
64 /// Constructor
65 
67  RooAbsReal& x, const RooArgList& coefList, TArrayD& limits, Int_t nBins) :
68  RooAbsPdf(name, title),
69  _x("x", "Dependent", this, x),
70  _coefList("coefList","List of coefficients",this),
71  _nBins(nBins)
72 {
73  _coefIter = _coefList.createIterator() ;
74 
75  // Check lowest order
76  if (_nBins<0) {
77  cout << "RooParametricStepFunction::ctor(" << GetName()
78  << ") WARNING: nBins must be >=0, setting value to 0" << endl ;
79  _nBins=0 ;
80  }
81 
82  TIterator* coefIter = coefList.createIterator() ;
83  RooAbsArg* coef ;
84  while((coef = (RooAbsArg*)coefIter->Next())) {
85  if (!dynamic_cast<RooAbsReal*>(coef)) {
86  cout << "RooParametricStepFunction::ctor(" << GetName() << ") ERROR: coefficient " << coef->GetName()
87  << " is not of type RooAbsReal" << endl ;
88  R__ASSERT(0) ;
89  }
90  _coefList.add(*coef) ;
91  }
92  delete coefIter ;
93 
94  // Bin limits
95  limits.Copy(_limits);
96 
97 }
98 
99 ////////////////////////////////////////////////////////////////////////////////
100 /// Copy constructor
101 
103  RooAbsPdf(other, name),
104  _x("x", this, other._x),
105  _coefList("coefList",this,other._coefList),
106  _nBins(other._nBins)
107 {
109  (other._limits).Copy(_limits);
110 }
111 
112 ////////////////////////////////////////////////////////////////////////////////
113 /// Destructor
114 
116 {
117  delete _coefIter ;
118 }
119 
120 
121 
122 ////////////////////////////////////////////////////////////////////////////////
123 
124 Int_t RooParametricStepFunction::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
125 {
126  if (matchArgs(allVars, analVars, _x)) return 1;
127  return 0;
128 }
129 
130 ////////////////////////////////////////////////////////////////////////////////
131 
133 {
134  R__ASSERT(code==1) ;
135 
136  // Case without range is trivial: p.d.f is by construction normalized
137  if (!rangeName) {
138  return 1.0 ;
139  }
140 
141  // Case with ranges, calculate integral explicitly
142  Double_t xmin = _x.min(rangeName) ;
143  Double_t xmax = _x.max(rangeName) ;
144 
145  Double_t sum=0 ;
146  Int_t i ;
147  for (i=1 ; i<=_nBins ; i++) {
148  Double_t binVal = (i<_nBins) ? (static_cast<RooRealVar*>(_coefList.at(i-1))->getVal()) : lastBinValue() ;
149  if (_limits[i-1]>=xmin && _limits[i]<=xmax) {
150  // Bin fully in the integration domain
151  sum += (_limits[i]-_limits[i-1])*binVal ;
152  } else if (_limits[i-1]<xmin && _limits[i]>xmax) {
153  // Domain is fully contained in this bin
154  sum += (xmax-xmin)*binVal ;
155  // Exit here, this is the last bin to be processed by construction
156  return sum ;
157  } else if (_limits[i-1]<xmin && _limits[i]<=xmax && _limits[i]>xmin) {
158  // Lower domain boundary is in bin
159  sum += (_limits[i]-xmin)*binVal ;
160  } else if (_limits[i-1]>=xmin && _limits[i]>xmax && _limits[i-1]<xmax) {
161  sum += (xmax-_limits[i-1])*binVal ;
162  // Upper domain boundary is in bin
163  // Exit here, this is the last bin to be processed by construction
164  return sum ;
165  }
166  }
167 
168  return sum;
169 }
170 
171 ////////////////////////////////////////////////////////////////////////////////
172 
174 {
175  Double_t sum(0.);
176  Double_t binSize(0.);
177  for (Int_t j=1;j<_nBins;j++){
178  RooRealVar* tmp = (RooRealVar*) _coefList.at(j-1);
179  binSize = _limits[j] - _limits[j-1];
180  sum = sum + tmp->getVal()*binSize;
181  }
182  binSize = _limits[_nBins] - _limits[_nBins-1];
183  return (1.0 - sum)/binSize;
184 }
185 
186 ////////////////////////////////////////////////////////////////////////////////
187 
189 {
190  Double_t value(0.);
191  if (_x >= _limits[0] && _x < _limits[_nBins]){
192 
193  for (Int_t i=1;i<=_nBins;i++){
194  if (_x < _limits[i]){
195  // in Bin i-1 (starting with Bin 0)
196  if (i<_nBins) {
197  // not in last Bin
198  RooRealVar* tmp = (RooRealVar*) _coefList.at(i-1);
199  value = tmp->getVal();
200  break;
201  } else {
202  // in last Bin
203  Double_t sum(0.);
204  Double_t binSize(0.);
205  for (Int_t j=1;j<_nBins;j++){
206  RooRealVar* tmp = (RooRealVar*) _coefList.at(j-1);
207  binSize = _limits[j] - _limits[j-1];
208  sum = sum + tmp->getVal()*binSize;
209  }
210  binSize = _limits[_nBins] - _limits[_nBins-1];
211  value = (1.0 - sum)/binSize;
212  if (value<=0.0){
213  value = 0.000001;
214  // cout << "RooParametricStepFunction: sum of values gt 1.0 -- beware!!" <<endl;
215  }
216  break;
217  }
218  }
219  }
220 
221  }
222  return value;
223 
224 }
225 
226 ////////////////////////////////////////////////////////////////////////////////
227 
229  return _nBins;
230 }
231 
232 ////////////////////////////////////////////////////////////////////////////////
233 
235  Double_t* limoutput = _limits.GetArray();
236  return limoutput;
237 }
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
std::string GetName(const std::string &scope_name)
Definition: Cppyy.cxx:145
TIterator * createIterator(Bool_t dir=kIterForward) const
static long int sum(long int i)
Definition: Factory.cxx:2162
float xmin
Definition: THbookFile.cxx:93
The Parametric Step Function PDF is a binned distribution whose parameters are the heights of each bi...
Double_t evaluate() const
do not persist
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
const Double_t * GetArray() const
Definition: TArrayD.h:43
#define R__ASSERT(e)
Definition: TError.h:96
int Int_t
Definition: RtypesCore.h:41
STL namespace.
Iterator abstract base class.
Definition: TIterator.h:30
Double_t x[n]
Definition: legend1.C:17
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:36
RooAbsArg * at(Int_t idx) const
Definition: RooArgList.h:84
float xmax
Definition: THbookFile.cxx:93
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
#define ClassImp(name)
Definition: Rtypes.h:336
Double_t min(const char *rname=0) const
Definition: RooRealProxy.h:56
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
Array of doubles (64 bits per element).
Definition: TArrayD.h:27
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57
virtual void Copy(TObject &named) const
Copy this to obj.
Definition: TNamed.cxx:85
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
virtual TObject * Next()=0
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
virtual ~RooParametricStepFunction()
Destructor.
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...