Logo ROOT  
Reference Guide
rf110_normintegration.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Basic functionality: examples on normalization and integration of p.d.fs, construction of cumulative distribution functions from monodimensional p.d.f.s
5##
6## \macro_code
7##
8## \date February 2018
9## \author Clemens Lange, Wouter Verkerke (C++ version)
10
11from __future__ import print_function
12import ROOT
13
14# Set up model
15# ---------------------
16
17# Create observables x,y
18x = ROOT.RooRealVar("x", "x", -10, 10)
19
20# Create p.d.f. gaussx(x,-2,3)
21gx = ROOT.RooGaussian(
22 "gx", "gx", x, ROOT.RooFit.RooConst(-2), ROOT.RooFit.RooConst(3))
23
24# Retrieve raw & normalized values of RooFit p.d.f.s
25# --------------------------------------------------------------------------------------------------
26
27# Return 'raw' unnormalized value of gx
28print("gx = ", gx.getVal())
29
30# Return value of gx normalized over x in range [-10,10]
31nset = ROOT.RooArgSet(x)
32print("gx_Norm[x] = ", gx.getVal(nset))
33
34# Create object representing integral over gx
35# which is used to calculate gx_Norm[x] == gx / gx_Int[x]
36igx = gx.createIntegral(ROOT.RooArgSet(x))
37print("gx_Int[x] = ", igx.getVal())
38
39# Integrate normalized pdf over subrange
40# ----------------------------------------------------------------------------
41
42# Define a range named "signal" in x from -5,5
43x.setRange("signal", -5, 5)
44
45# Create an integral of gx_Norm[x] over x in range "signal"
46# ROOT.This is the fraction of of p.d.f. gx_Norm[x] which is in the
47# range named "signal"
48xset = ROOT.RooArgSet(x)
49igx_sig = gx.createIntegral(xset, ROOT.RooFit.NormSet(xset), ROOT.RooFit.Range("signal"))
50print("gx_Int[x|signal]_Norm[x] = ", igx_sig.getVal())
51
52# Construct cumulative distribution function from pdf
53# -----------------------------------------------------------------------------------------------------
54
55# Create the cumulative distribution function of gx
56# i.e. calculate Int[-10,x] gx(x') dx'
57gx_cdf = gx.createCdf(ROOT.RooArgSet(x))
58
59# Plot cdf of gx versus x
60frame = x.frame(ROOT.RooFit.Title("c.d.f of Gaussian p.d.f"))
61gx_cdf.plotOn(frame)
62
63# Draw plot on canvas
64c = ROOT.TCanvas("rf110_normintegration",
65 "rf110_normintegration", 600, 600)
66ROOT.gPad.SetLeftMargin(0.15)
67frame.GetYaxis().SetTitleOffset(1.6)
68frame.Draw()
69
70c.SaveAs("rf110_normintegration.png")