ROOT   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 * *
9 *
10 * Redistribution and use in source and binary forms, *
11 * with or without modification, are permitted according to the terms *
13 *****************************************************************************/
14
15/** \class RooParametricStepFunction
16 \ingroup Roofit
17
18The Parametric Step Function PDF is a binned distribution whose parameters
19are the heights of each bin. This PDF was first used in BaBar's B0->pi0pi0
20paper BABAR Collaboration (B. Aubert et al.) Phys.Rev.Lett.91:241801,2003.
21
22This PDF may be used to describe oddly shaped distributions. It differs
23from a RooKeysPdf or a RooHistPdf in that a RooParametricStepFunction
24has free parameters. In particular, any statistical uncertainty in
25sample used to model this PDF may be understood with these free parameters;
26this is not possible with non-parametric PDFs.
27
28The RooParametricStepFunction has Nbins-1 free parameters. Note that
29the limits of the dependent variable must match the low and hi bin limits.
30
31An example of usage is:
32
33~~~ {.cpp}
34Int_t nbins(10);
35TArrayD limits(nbins+1);
36limits[0] = 0.0; //etc...
37RooArgList* list = new RooArgList("list");
38RooRealVar* binHeight0 = new RooRealVar("binHeight0","bin 0 Value",0.1,0.0,1.0);
39list->add(binHeight0); // up to binHeight8, ie. 9 parameters
40
41RooParametricStepFunction 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
59using 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{
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 }
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
124Int_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}
int Int_t
Definition: RtypesCore.h:41
double Double_t
Definition: RtypesCore.h:55
#define ClassImp(name)
Definition: Rtypes.h:365
#define R__ASSERT(e)
Definition: TError.h:96
char name[80]
Definition: TGX11.cxx:109
float xmin
Definition: THbookFile.cxx:93
float xmax
Definition: THbookFile.cxx:93
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:71
TIterator * createIterator(Bool_t dir=kIterForward) const R__SUGGEST_ALTERNATIVE("begin()
TIterator-style iteration over contained elements.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:59
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:87
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgList.h:21
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition: RooArgList.h:74
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:28
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
The Parametric Step Function PDF is a binned distribution whose parameters are the heights of each bi...
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
virtual ~RooParametricStepFunction()
Destructor.
Double_t evaluate() const
do not persist
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:35
Double_t min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
Double_t max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
Array of doubles (64 bits per element).
Definition: TArrayD.h:27
void Copy(TArrayD &array) const
Definition: TArrayD.h:42
const Double_t * GetArray() const
Definition: TArrayD.h:43
Iterator abstract base class.
Definition: TIterator.h:30
virtual TObject * Next()=0
virtual void Copy(TObject &named) const
Copy this to obj.
Definition: TNamed.cxx:94
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
Double_t x[n]
Definition: legend1.C:17
static long int sum(long int i)
Definition: Factory.cxx:2276