ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
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 //////////////////////////////////////////////////////////////////////////////
16 //
17 // The Parametric Step Function PDF is a binned distribution whose parameters
18 // are the heights of each bin. This PDF was first used in BaBar's B0->pi0pi0
19 // paper BABAR Collaboration (B. Aubert et al.) Phys.Rev.Lett.91:241801,2003.
20 //
21 // This PDF may be used to describe oddly shaped distributions. It differs
22 // from a RooKeysPdf or a RooHistPdf in that a RooParametricStepFunction
23 // has free parameters. In particular, any statistical uncertainty in
24 // sample used to model this PDF may be understood with these free parameters;
25 // this is not possible with non-parametric PDFs.
26 //
27 // The RooParametricStepFunction has Nbins-1 free parameters. Note that
28 // the limits of the dependent varaible must match the low and hi bin limits.
29 //
30 // An example of usage is:
31 //
32 // Int_t nbins(10);
33 // TArrayD limits(nbins+1);
34 // limits[0] = 0.0; //etc...
35 // RooArgList* list = new RooArgList("list");
36 // RooRealVar* binHeight0 = new RooRealVar("binHeight0","bin 0 Value",0.1,0.0,1.0);
37 // list->add(binHeight0); // up to binHeight8, ie. 9 parameters
38 //
39 // RooParametricStepFunction aPdf = ("aPdf","PSF",*x,
40 // *list,limits,nbins);
41 //
42 
43 
44 #include "RooFit.h"
45 
46 #include "Riostream.h"
47 #include "TArrayD.h"
48 #include <math.h>
49 
51 #include "RooAbsReal.h"
52 #include "RooRealVar.h"
53 #include "RooArgList.h"
54 
55 #include "TError.h"
56 
57 using namespace std;
58 
60 ;
61 
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  }
90  _coefList.add(*coef) ;
91  }
92  delete coefIter ;
93 
94  // Bin limits
95  limits.Copy(_limits);
96 
97 }
98 
99 
100 
101 ////////////////////////////////////////////////////////////////////////////////
102 /// Copy constructor
103 
105  RooAbsPdf(other, name),
106  _x("x", this, other._x),
107  _coefList("coefList",this,other._coefList),
108  _nBins(other._nBins)
109 {
111  (other._limits).Copy(_limits);
112 }
113 
114 
115 
116 ////////////////////////////////////////////////////////////////////////////////
117 /// Destructor
118 
120 {
121  delete _coefIter ;
122 }
123 
124 
125 
126 ////////////////////////////////////////////////////////////////////////////////
127 
128 Int_t RooParametricStepFunction::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
129 {
130  if (matchArgs(allVars, analVars, _x)) return 1;
131  return 0;
132 }
133 
134 
135 
136 ////////////////////////////////////////////////////////////////////////////////
137 
139 {
140  R__ASSERT(code==1) ;
141 
142  // Case without range is trivial: p.d.f is by construction normalized
143  if (!rangeName) {
144  return 1.0 ;
145  }
146 
147  // Case with ranges, calculate integral explicitly
148  Double_t xmin = _x.min(rangeName) ;
149  Double_t xmax = _x.max(rangeName) ;
150 
151  Double_t sum=0 ;
152  Int_t i ;
153  for (i=1 ; i<=_nBins ; i++) {
154  Double_t binVal = (i<_nBins) ? (static_cast<RooRealVar*>(_coefList.at(i-1))->getVal()) : lastBinValue() ;
155  if (_limits[i-1]>=xmin && _limits[i]<=xmax) {
156  // Bin fully in the integration domain
157  sum += (_limits[i]-_limits[i-1])*binVal ;
158  } else if (_limits[i-1]<xmin && _limits[i]>xmax) {
159  // Domain is fully contained in this bin
160  sum += (xmax-xmin)*binVal ;
161  // Exit here, this is the last bin to be processed by construction
162  return sum ;
163  } else if (_limits[i-1]<xmin && _limits[i]<=xmax && _limits[i]>xmin) {
164  // Lower domain boundary is in bin
165  sum += (_limits[i]-xmin)*binVal ;
166  } else if (_limits[i-1]>=xmin && _limits[i]>xmax && _limits[i-1]<xmax) {
167  sum += (xmax-_limits[i-1])*binVal ;
168  // Upper domain boundary is in bin
169  // Exit here, this is the last bin to be processed by construction
170  return sum ;
171  }
172  }
173 
174  return sum;
175 }
176 
177 
178 
179 ////////////////////////////////////////////////////////////////////////////////
180 
182 {
183  Double_t sum(0.);
184  Double_t binSize(0.);
185  for (Int_t j=1;j<_nBins;j++){
186  RooRealVar* tmp = (RooRealVar*) _coefList.at(j-1);
187  binSize = _limits[j] - _limits[j-1];
188  sum = sum + tmp->getVal()*binSize;
189  }
190  binSize = _limits[_nBins] - _limits[_nBins-1];
191  return (1.0 - sum)/binSize;
192 }
193 
194 
195 
196 ////////////////////////////////////////////////////////////////////////////////
197 
199 {
200  Double_t value(0.);
201  if (_x >= _limits[0] && _x < _limits[_nBins]){
202 
203  for (Int_t i=1;i<=_nBins;i++){
204  if (_x < _limits[i]){
205  // in Bin i-1 (starting with Bin 0)
206  if (i<_nBins) {
207  // not in last Bin
208  RooRealVar* tmp = (RooRealVar*) _coefList.at(i-1);
209  value = tmp->getVal();
210  break;
211  } else {
212  // in last Bin
213  Double_t sum(0.);
214  Double_t binSize(0.);
215  for (Int_t j=1;j<_nBins;j++){
216  RooRealVar* tmp = (RooRealVar*) _coefList.at(j-1);
217  binSize = _limits[j] - _limits[j-1];
218  sum = sum + tmp->getVal()*binSize;
219  }
220  binSize = _limits[_nBins] - _limits[_nBins-1];
221  value = (1.0 - sum)/binSize;
222  if (value<=0.0){
223  value = 0.000001;
224  // cout << "RooParametricStepFunction: sum of values gt 1.0 -- beware!!" <<endl;
225  }
226  break;
227  }
228  }
229  }
230 
231  }
232  return value;
233 
234 }
235 
236 
237 ////////////////////////////////////////////////////////////////////////////////
238 
240  return _nBins;
241 }
242 
243 
244 ////////////////////////////////////////////////////////////////////////////////
245 
247  Double_t* limoutput = _limits.GetArray();
248  return limoutput;
249 }
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported...
float xmin
Definition: THbookFile.cxx:93
#define R__ASSERT(e)
Definition: TError.h:98
ClassImp(RooParametricStepFunction)
int Int_t
Definition: RtypesCore.h:41
void Copy(TArrayD &array) const
Definition: TArrayD.h:44
Iterator abstract base class.
Definition: TIterator.h:32
Double_t x[n]
Definition: legend1.C:17
Double_t evaluate() const
do not persist
virtual void Copy(TObject &named) const
Copy this to obj.
Definition: TNamed.cxx:83
TIterator * createIterator(Bool_t dir=kIterForward) const
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral. ...
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Reimplementation of standard RooArgList::add()
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
TPaveLabel title(3, 27.1, 15, 28.7,"ROOT Environment and Tools")
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
float xmax
Definition: THbookFile.cxx:93
const Double_t * GetArray() const
Definition: TArrayD.h:45
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
Array of doubles (64 bits per element).
Definition: TArrayD.h:29
#define name(a, b)
Definition: linkTestLib0.cpp:5
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
Double_t min(const char *rname=0) const
Definition: RooRealProxy.h:56
virtual TObject * Next()=0
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooAbsArg is the common abstract base class for objects that represent a value (of arbitrary type) an...
Definition: RooAbsArg.h:66
virtual ~RooParametricStepFunction()
Destructor.
float value
Definition: math.cpp:443
Double_t max(const char *rname=0) const
Definition: RooRealProxy.h:57