Logo ROOT  
Reference Guide
rf308_normintegration2d.C
Go to the documentation of this file.
1 /// \file
2 /// \ingroup tutorial_roofit
3 /// \notebook
4 /// Multidimensional models: normalization and integration of pdfs, construction of
5 /// cumulative distribution functions from pdfs in two dimensions
6 ///
7 /// \macro_image
8 /// \macro_output
9 /// \macro_code
10 ///
11 /// \date July 2008
12 /// \author Wouter Verkerke
13 
14 #include "RooRealVar.h"
15 #include "RooGaussian.h"
16 #include "RooConstVar.h"
17 #include "RooProdPdf.h"
18 #include "RooAbsReal.h"
19 #include "RooPlot.h"
20 #include "TCanvas.h"
21 #include "TAxis.h"
22 #include "TH1.h"
23 using namespace RooFit;
24 
26 {
27  // S e t u p m o d e l
28  // ---------------------
29 
30  // Create observables x,y
31  RooRealVar x("x", "x", -10, 10);
32  RooRealVar y("y", "y", -10, 10);
33 
34  // Create pdf gaussx(x,-2,3), gaussy(y,2,2)
35  RooGaussian gx("gx", "gx", x, RooConst(-2), RooConst(3));
36  RooGaussian gy("gy", "gy", y, RooConst(+2), RooConst(2));
37 
38  // Create gxy = gx(x)*gy(y)
39  RooProdPdf gxy("gxy", "gxy", RooArgSet(gx, gy));
40 
41  // 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
42  // --------------------------------------------------------------------------------------------------
43 
44  // Return 'raw' unnormalized value of gx
45  cout << "gxy = " << gxy.getVal() << endl;
46 
47  // Return value of gxy normalized over x _and_ y in range [-10,10]
48  RooArgSet nset_xy(x, y);
49  cout << "gx_Norm[x,y] = " << gxy.getVal(&nset_xy) << endl;
50 
51  // Create object representing integral over gx
52  // which is used to calculate gx_Norm[x,y] == gx / gx_Int[x,y]
53  RooAbsReal *igxy = gxy.createIntegral(RooArgSet(x, y));
54  cout << "gx_Int[x,y] = " << igxy->getVal() << endl;
55 
56  // NB: it is also possible to do the following
57 
58  // Return value of gxy normalized over x in range [-10,10] (i.e. treating y as parameter)
59  RooArgSet nset_x(x);
60  cout << "gx_Norm[x] = " << gxy.getVal(&nset_x) << endl;
61 
62  // Return value of gxy normalized over y in range [-10,10] (i.e. treating x as parameter)
63  RooArgSet nset_y(y);
64  cout << "gx_Norm[y] = " << gxy.getVal(&nset_y) << endl;
65 
66  // 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
67  // ----------------------------------------------------------------------------
68 
69  // Define a range named "signal" in x from -5,5
70  x.setRange("signal", -5, 5);
71  y.setRange("signal", -3, 3);
72 
73  // Create an integral of gxy_Norm[x,y] over x and y in range "signal"
74  // This is the fraction of of pdf gxy_Norm[x,y] which is in the
75  // range named "signal"
76  RooAbsReal *igxy_sig = gxy.createIntegral(RooArgSet(x, y), NormSet(RooArgSet(x, y)), Range("signal"));
77  cout << "gx_Int[x,y|signal]_Norm[x,y] = " << igxy_sig->getVal() << endl;
78 
79  // 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
80  // -----------------------------------------------------------------------------------------------------
81 
82  // Create the cumulative distribution function of gx
83  // i.e. calculate Int[-10,x] gx(x') dx'
84  RooAbsReal *gxy_cdf = gxy.createCdf(RooArgSet(x, y));
85 
86  // Plot cdf of gx versus x
87  TH1 *hh_cdf = gxy_cdf->createHistogram("hh_cdf", x, Binning(40), YVar(y, Binning(40)));
88  hh_cdf->SetLineColor(kBlue);
89 
90  new TCanvas("rf308_normintegration2d", "rf308_normintegration2d", 600, 600);
91  gPad->SetLeftMargin(0.15);
92  hh_cdf->GetZaxis()->SetTitleOffset(1.8);
93  hh_cdf->Draw("surf");
94 }
RooAbsReal.h
RooAbsReal::createIntegral
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:596
RooGaussian.h
RooAbsReal::getVal
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:92
x
Double_t x[n]
Definition: legend1.C:17
RooGaussian
Plain Gaussian p.d.f.
Definition: RooGaussian.h:25
RooAbsReal
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:62
TAttLine::SetLineColor
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition: TAttLine.h:40
TCanvas.h
RooFit::YVar
RooCmdArg YVar(const RooAbsRealLValue &var, const RooCmdArg &arg=RooCmdArg::none())
Definition: RooGlobalFunc.cxx:246
RooFit::Binning
RooCmdArg Binning(const RooAbsBinning &binning)
Definition: RooGlobalFunc.cxx:83
TH1::GetZaxis
TAxis * GetZaxis()
Definition: TH1.h:322
RooProdPdf.h
RooFit
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition: RooCFunction1Binding.h:29
RooPlot.h
y
Double_t y[n]
Definition: legend1.C:17
RooRealVar.h
RooConstVar.h
RooFit::NormSet
RooCmdArg NormSet(const RooArgSet &nset)
Definition: RooGlobalFunc.cxx:262
RooAbsReal::createHistogram
TH1 * createHistogram(const char *varNameList, Int_t xbins=0, Int_t ybins=0, Int_t zbins=0) const
Create and fill a ROOT histogram TH1, TH2 or TH3 with the values of this function for the variables w...
Definition: RooAbsReal.cxx:1326
TCanvas
The Canvas class.
Definition: TCanvas.h:23
TAxis.h
TH1
TH1 is the base class of all histogramm classes in ROOT.
Definition: TH1.h:58
RooFit::Range
RooCmdArg Range(const char *rangeName, Bool_t adjustNorm=kTRUE)
Definition: RooGlobalFunc.cxx:53
kBlue
@ kBlue
Definition: Rtypes.h:66
gPad
#define gPad
Definition: TVirtualPad.h:287
TAttAxis::SetTitleOffset
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition: TAttAxis.cxx:293
RooRealVar
RooRealVar represents a variable that can be changed from the outside.
Definition: RooRealVar.h:37
RooProdPdf
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:37
TH1.h
rf308_normintegration2d
Definition: rf308_normintegration2d.py:1
RooArgSet
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:29
TH1::Draw
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:3050
RooFit::RooConst
RooConstVar & RooConst(Double_t val)
Definition: RooGlobalFunc.cxx:347