Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf302_utilfuncs.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Multidimensional models: utility functions classes available for use in tailoring of
5## composite (multidimensional) pdfs
6##
7## \macro_image
8## \macro_code
9## \macro_output
10##
11## \date February 2018
12## \authors Clemens Lange, Wouter Verkerke (C++ version)
13
14import ROOT
15
16# Create observables, parameters
17# -----------------------------------------------------------
18
19# Create observables
20x = ROOT.RooRealVar("x", "x", -5, 5)
21y = ROOT.RooRealVar("y", "y", -5, 5)
22
23# Create parameters
24a0 = ROOT.RooRealVar("a0", "a0", -1.5, -5, 5)
25a1 = ROOT.RooRealVar("a1", "a1", -0.5, -1, 1)
26sigma = ROOT.RooRealVar("sigma", "width of gaussian", 0.5)
27
28# Using RooFormulaVar to tailor pdf
29# -----------------------------------------------------------------------
30
31# Create interpreted function f(y) = a0 - a1*sqrt(10*abs(y))
32fy_1 = ROOT.RooFormulaVar("fy_1", "a0-a1*sqrt(10*abs(y))", [y, a0, a1])
33
34# Create gauss(x,f(y),s)
35model_1 = ROOT.RooGaussian("model_1", "Gaussian with shifting mean", x, fy_1, sigma)
36
37# Using RooPolyVar to tailor pdf
38# -----------------------------------------------------------------------
39
40# Create polynomial function f(y) = a0 + a1*y
41fy_2 = ROOT.RooPolyVar("fy_2", "fy_2", y, [a0, a1])
42
43# Create gauss(x,f(y),s)
44model_2 = ROOT.RooGaussian("model_2", "Gaussian with shifting mean", x, fy_2, sigma)
45
46# Using RooAddition to tailor pdf
47# -----------------------------------------------------------------------
48
49# Create sum function f(y) = a0 + y
50fy_3 = ROOT.RooAddition("fy_3", "a0+y", [a0, y])
51
52# Create gauss(x,f(y),s)
53model_3 = ROOT.RooGaussian("model_3", "Gaussian with shifting mean", x, fy_3, sigma)
54
55# Using RooProduct to tailor pdf
56# -----------------------------------------------------------------------
57
58# Create product function f(y) = a1*y
59fy_4 = ROOT.RooProduct("fy_4", "a1*y", [a1, y])
60
61# Create gauss(x,f(y),s)
62model_4 = ROOT.RooGaussian("model_4", "Gaussian with shifting mean", x, fy_4, sigma)
63
64# Plot all pdfs
65# ----------------------------
66
67# Make two-dimensional plots in x vs y
68hh_model_1 = model_1.createHistogram("hh_model_1", x, Binning=50, YVar=dict(var=y, Binning=50))
69hh_model_2 = model_2.createHistogram("hh_model_2", x, Binning=50, YVar=dict(var=y, Binning=50))
70hh_model_3 = model_3.createHistogram("hh_model_3", x, Binning=50, YVar=dict(var=y, Binning=50))
71hh_model_4 = model_4.createHistogram("hh_model_4", x, Binning=50, YVar=dict(var=y, Binning=50))
72hh_model_1.SetLineColor(ROOT.kBlue)
73hh_model_2.SetLineColor(ROOT.kBlue)
74hh_model_3.SetLineColor(ROOT.kBlue)
75hh_model_4.SetLineColor(ROOT.kBlue)
76
77# Make canvas and draw ROOT.RooPlots
78c = ROOT.TCanvas("rf302_utilfuncs", "rf302_utilfuncs", 800, 800)
79c.Divide(2, 2)
80c.cd(1)
81ROOT.gPad.SetLeftMargin(0.20)
82hh_model_1.GetZaxis().SetTitleOffset(2.5)
83hh_model_1.Draw("surf")
84c.cd(2)
85ROOT.gPad.SetLeftMargin(0.20)
86hh_model_2.GetZaxis().SetTitleOffset(2.5)
87hh_model_2.Draw("surf")
88c.cd(3)
89ROOT.gPad.SetLeftMargin(0.20)
90hh_model_3.GetZaxis().SetTitleOffset(2.5)
91hh_model_3.Draw("surf")
92c.cd(4)
93ROOT.gPad.SetLeftMargin(0.20)
94hh_model_4.GetZaxis().SetTitleOffset(2.5)
95hh_model_4.Draw("surf")
96
97c.SaveAs("rf302_utilfuncs.png")