Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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_code
9/// \macro_output
10///
11/// \date July 2008
12/// \author Wouter Verkerke
13
14#include "RooRealVar.h"
15#include "RooGaussian.h"
16#include "RooProdPdf.h"
17#include "RooAbsReal.h"
18#include "RooPlot.h"
19#include "TCanvas.h"
20#include "TAxis.h"
21#include "TH1.h"
22using namespace RooFit;
23
25{
26 // S e t u p m o d e l
27 // ---------------------
28
29 // Create observables x,y
30 RooRealVar x("x", "x", -10, 10);
31 RooRealVar y("y", "y", -10, 10);
32
33 // Create pdf gaussx(x,-2,3), gaussy(y,2,2)
34 RooGaussian gx("gx", "gx", x, -2.0, 3.0);
35 RooGaussian gy("gy", "gy", y, +2.0, 2.0);
36
37 // Create gxy = gx(x)*gy(y)
38 RooProdPdf gxy("gxy", "gxy", RooArgSet(gx, gy));
39
40 // 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
41 // --------------------------------------------------------------------------------------------------
42
43 // Return 'raw' unnormalized value of gx
44 cout << "gxy = " << gxy.getVal() << endl;
45
46 // Return value of gxy normalized over x _and_ y in range [-10,10]
47 RooArgSet nset_xy(x, y);
48 cout << "gx_Norm[x,y] = " << gxy.getVal(&nset_xy) << endl;
49
50 // Create object representing integral over gx
51 // which is used to calculate gx_Norm[x,y] == gx / gx_Int[x,y]
52 std::unique_ptr<RooAbsReal> igxy{gxy.createIntegral(RooArgSet(x, y))};
53 cout << "gx_Int[x,y] = " << igxy->getVal() << endl;
54
55 // NB: it is also possible to do the following
56
57 // Return value of gxy normalized over x in range [-10,10] (i.e. treating y as parameter)
58 RooArgSet nset_x(x);
59 cout << "gx_Norm[x] = " << gxy.getVal(&nset_x) << endl;
60
61 // Return value of gxy normalized over y in range [-10,10] (i.e. treating x as parameter)
62 RooArgSet nset_y(y);
63 cout << "gx_Norm[y] = " << gxy.getVal(&nset_y) << endl;
64
65 // 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
66 // ----------------------------------------------------------------------------
67
68 // Define a range named "signal" in x from -5,5
69 x.setRange("signal", -5, 5);
70 y.setRange("signal", -3, 3);
71
72 // Create an integral of gxy_Norm[x,y] over x and y in range "signal"
73 // This is the fraction of of pdf gxy_Norm[x,y] which is in the
74 // range named "signal"
75 std::unique_ptr<RooAbsReal> igxy_sig{gxy.createIntegral({x, y}, NormSet(RooArgSet(x, y)), Range("signal"))};
76 cout << "gx_Int[x,y|signal]_Norm[x,y] = " << igxy_sig->getVal() << endl;
77
78 // 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
79 // -----------------------------------------------------------------------------------------------------
80
81 // Create the cumulative distribution function of gx
82 // i.e. calculate Int[-10,x] gx(x') dx'
83 std::unique_ptr<RooAbsReal> gxy_cdf{gxy.createCdf(RooArgSet(x, y))};
84
85 // Plot cdf of gx versus x
86 TH1 *hh_cdf = gxy_cdf->createHistogram("hh_cdf", x, Binning(40), YVar(y, Binning(40)));
87 hh_cdf->SetLineColor(kBlue);
88
89 new TCanvas("rf308_normintegration2d", "rf308_normintegration2d", 600, 600);
90 gPad->SetLeftMargin(0.15);
91 hh_cdf->GetZaxis()->SetTitleOffset(1.8);
92 hh_cdf->Draw("surf");
93}
@ kBlue
Definition Rtypes.h:66
#define gPad
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:24
Plain Gaussian p.d.f.
Definition RooGaussian.h:24
Efficient implementation of a product of PDFs of the form.
Definition RooProdPdf.h:39
Variable that can be changed from the outside.
Definition RooRealVar.h:37
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition TAttAxis.cxx:298
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:40
The Canvas class.
Definition TCanvas.h:23
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:59
TAxis * GetZaxis()
Definition TH1.h:327
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3068
RooCmdArg YVar(const RooAbsRealLValue &var, const RooCmdArg &arg={})
RooCmdArg NormSet(Args_t &&... argsOrArgSet)
RooCmdArg Binning(const RooAbsBinning &binning)
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition CodegenImpl.h:64
Ta Range(0, 0, 1, 1)