ROOT   Reference Guide
Searching...
No Matches
rf308_normintegration2d.py
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 February 2018
12## \authors Clemens Lange, Wouter Verkerke (C++ version)
13
14from __future__ import print_function
15import ROOT
16
17# Set up model
18# ---------------------
19
20# Create observables x,y
21x = ROOT.RooRealVar("x", "x", -10, 10)
22y = ROOT.RooRealVar("y", "y", -10, 10)
23
24# Create pdf gaussx(x,-2,3), gaussy(y,2,2)
25gx = ROOT.RooGaussian("gx", "gx", x, -2.0, 3.0)
26gy = ROOT.RooGaussian("gy", "gy", y, +2.0, 2.0)
27
28# gxy = gx(x)*gy(y)
29gxy = ROOT.RooProdPdf("gxy", "gxy", [gx, gy])
30
31# Retrieve raw & normalized values of RooFit pdfs
32# --------------------------------------------------------------------------------------------------
33
34# Return 'raw' unnormalized value of gx
35print("gxy = ", gxy.getVal())
36
37# Return value of gxy normalized over x _and_ y in range [-10,10]
38nset_xy = {x, y}
39print("gx_Norm[x,y] = ", gxy.getVal(nset_xy))
40
41# Create object representing integral over gx
42# which is used to calculate gx_Norm[x,y] == gx / gx_Int[x,y]
43x_and_y = {x, y}
44igxy = gxy.createIntegral(x_and_y)
45print("gx_Int[x,y] = ", igxy.getVal())
46
47# NB: it is also possible to do the following
48
49# Return value of gxy normalized over x in range [-10,10] (i.e. treating y
50# as parameter)
51nset_x = {x}
52print("gx_Norm[x] = ", gxy.getVal(nset_x))
53
54# Return value of gxy normalized over y in range [-10,10] (i.e. treating x
55# as parameter)
56nset_y = {y}
57print("gx_Norm[y] = ", gxy.getVal(nset_y))
58
59# Integrate normalized pdf over subrange
60# ----------------------------------------------------------------------------
61
62# Define a range named "signal" in x from -5,5
63x.setRange("signal", -5, 5)
64y.setRange("signal", -3, 3)
65
66# Create an integral of gxy_Norm[x,y] over x and y in range "signal"
67# ROOT.This is the fraction of of pdf gxy_Norm[x,y] which is in the
68# range named "signal"
69
70igxy_sig = gxy.createIntegral(x_and_y, NormSet=x_and_y, Range="signal")
71print("gx_Int[x,y|signal]_Norm[x,y] = ", igxy_sig.getVal())
72
73# Construct cumulative distribution function from pdf
74# -----------------------------------------------------------------------------------------------------
75
76# Create the cumulative distribution function of gx
77# i.e. calculate Int[-10,x] gx(x') dx'
78gxy_cdf = gxy.createCdf({x, y})
79
80# Plot cdf of gx versus x
81hh_cdf = gxy_cdf.createHistogram("hh_cdf", x, Binning=40, YVar=dict(var=y, Binning=40))
82hh_cdf.SetLineColor(ROOT.kBlue)
83
84c = ROOT.TCanvas("rf308_normintegration2d", "rf308_normintegration2d", 600, 600)