Logo ROOT   6.08/07
Reference Guide
rf110_normintegration.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook -js
4 /// 'BASIC FUNCTIONALITY' RooFit tutorial macro #110
5 ///
6 /// Examples on normalization of p.d.f.s,
7 /// integration of p.d.fs, construction
8 /// of cumulative distribution functions from p.d.f.s
9 /// in one dimension
10 ///
11 /// \macro_image
12 /// \macro_output
13 /// \macro_code
14 /// \author 07/2008 - Wouter Verkerke
15 
16 #include "RooRealVar.h"
17 #include "RooGaussian.h"
18 #include "RooConstVar.h"
19 #include "RooAbsReal.h"
20 #include "RooPlot.h"
21 #include "TCanvas.h"
22 #include "TAxis.h"
23 using namespace RooFit ;
24 
25 
26 void rf110_normintegration()
27 {
28  // S e t u p m o d e l
29  // ---------------------
30 
31  // Create observables x,y
32  RooRealVar x("x","x",-10,10) ;
33 
34  // Create p.d.f. gaussx(x,-2,3)
35  RooGaussian gx("gx","gx",x,RooConst(-2),RooConst(3)) ;
36 
37 
38  // R e t r i e v e r a w & n o r m a l i z e d v a l u e s o f R o o F i t p . d . f . s
39  // --------------------------------------------------------------------------------------------------
40 
41  // Return 'raw' unnormalized value of gx
42  cout << "gx = " << gx.getVal() << endl ;
43 
44  // Return value of gx normalized over x in range [-10,10]
45  RooArgSet nset(x) ;
46  cout << "gx_Norm[x] = " << gx.getVal(&nset) << endl ;
47 
48  // Create object representing integral over gx
49  // which is used to calculate gx_Norm[x] == gx / gx_Int[x]
50  RooAbsReal* igx = gx.createIntegral(x) ;
51  cout << "gx_Int[x] = " << igx->getVal() << endl ;
52 
53 
54  // I n t e g r a t e n o r m a l i z e d p d f o v e r s u b r a n g e
55  // ----------------------------------------------------------------------------
56 
57  // Define a range named "signal" in x from -5,5
58  x.setRange("signal",-5,5) ;
59 
60  // Create an integral of gx_Norm[x] over x in range "signal"
61  // This is the fraction of of p.d.f. gx_Norm[x] which is in the
62  // range named "signal"
63  RooAbsReal* igx_sig = gx.createIntegral(x,NormSet(x),Range("signal")) ;
64  cout << "gx_Int[x|signal]_Norm[x] = " << igx_sig->getVal() << endl ;
65 
66 
67 
68  // C o n s t r u c t c u m u l a t i v e d i s t r i b u t i o n f u n c t i o n f r o m p d f
69  // -----------------------------------------------------------------------------------------------------
70 
71  // Create the cumulative distribution function of gx
72  // i.e. calculate Int[-10,x] gx(x') dx'
73  RooAbsReal* gx_cdf = gx.createCdf(x) ;
74 
75  // Plot cdf of gx versus x
76  RooPlot* frame = x.frame(Title("c.d.f of Gaussian p.d.f")) ;
77  gx_cdf->plotOn(frame) ;
78 
79 
80 
81  // Draw plot on canvas
82  new TCanvas("rf110_normintegration","rf110_normintegration",600,600) ;
83  gPad->SetLeftMargin(0.15) ; frame->GetYaxis()->SetTitleOffset(1.6) ;
84  frame->Draw() ;
85 
86 
87 }
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title Offset is a correction factor with respect to the "s...
Definition: TAttAxis.cxx:262
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg(), const RooCmdArg &arg10=RooCmdArg()) const
Plot (project) PDF on specified frame.
TAxis * GetYaxis() const
Definition: RooPlot.cxx:1118
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
RooCmdArg NormSet(const RooArgSet &nset)
RooCmdArg Title(const char *name)
RooCmdArg Range(const char *rangeName, Bool_t adjustNorm=kTRUE)
Double_t x[n]
Definition: legend1.C:17
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
RooAbsReal * createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Create an object that represents the integral of the function over one or more observables listed in ...
Definition: RooAbsReal.cxx:503
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
Definition: RooPlot.h:41
The Canvas class.
Definition: TCanvas.h:41
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:53
RooConstVar & RooConst(Double_t val)
#define gPad
Definition: TVirtualPad.h:289
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
Definition: RooPlot.cxx:559