Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf309_ndimplot.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Multidimensional models: making 2/3 dimensional plots of pdfs and datasets
5##
6## \macro_image
7## \macro_code
8## \macro_output
9##
10## \date February 2018
11## \authors Clemens Lange, Wouter Verkerke (C++ version)
12
13import ROOT
14
15# Create 2D model and dataset
16# -----------------------------------------------------
17
18# Create observables
19x = ROOT.RooRealVar("x", "x", -5, 5)
20y = ROOT.RooRealVar("y", "y", -5, 5)
21
22# Create parameters
23a0 = ROOT.RooRealVar("a0", "a0", -3.5, -5, 5)
24a1 = ROOT.RooRealVar("a1", "a1", -1.5, -1, 1)
25sigma = ROOT.RooRealVar("sigma", "width of gaussian", 1.5)
26
27# Create interpreted function f(y) = a0 - a1*sqrt(10*abs(y))
28fy = ROOT.RooFormulaVar("fy", "a0-a1*sqrt(10*abs(y))", [y, a0, a1])
29
30# Create gauss(x,f(y),s)
31model = ROOT.RooGaussian("model", "Gaussian with shifting mean", x, fy, sigma)
32
33# Sample dataset from gauss(x,y)
34data = model.generate({x, y}, 10000)
35
36# Make 2D plots of data and model
37# -------------------------------------------------------------
38
39# Create and fill ROOT 2D histogram (20x20 bins) with contents of dataset
40hh_data = data.createHistogram("x,y", x, Binning=20, YVar=dict(var=y, Binning=20))
41
42# Create and fill ROOT 2D histogram (50x50 bins) with sampling of pdf
43# hh_pdf = model.createHistogram("hh_model", x, Binning=50, YVar=dict(var=y, Binning=50))
44hh_pdf = model.createHistogram("x,y", 50, 50)
45hh_pdf.SetLineColor(ROOT.kBlue)
46
47# Create 3D model and dataset
48# -----------------------------------------------------
49
50# Create observables
51z = ROOT.RooRealVar("z", "z", -5, 5)
52
53gz = ROOT.RooGaussian("gz", "gz", z, 0.0, 2.0)
54model3 = ROOT.RooProdPdf("model3", "model3", [model, gz])
55
56data3 = model3.generate({x, y, z}, 10000)
57
58# Make 3D plots of data and model
59# -------------------------------------------------------------
60
61# Create and fill ROOT 2D histogram (8x8x8 bins) with contents of dataset
62# hh_data3 = data3.createHistogram("hh_data3", x, Binning=8, YVar=(var=y, Binning=8), ZVar=(var=z, Binning=8))
63hh_data3 = data3.createHistogram(
64 "hh_data3",
65 x,
66 Binning=8,
67 YVar=dict(var=y, Binning=8),
68 ZVar=dict(var=z, Binning=8),
69)
70
71# Create and fill ROOT 2D histogram (20x20x20 bins) with sampling of pdf
72hh_pdf3 = model3.createHistogram(
73 "hh_model3",
74 x,
75 Binning=20,
76 YVar=dict(var=y, Binning=20),
77 ZVar=dict(var=z, Binning=20),
78)
79hh_pdf3.SetFillColor(ROOT.kBlue)
80
81c1 = ROOT.TCanvas("rf309_2dimplot", "rf309_2dimplot", 800, 800)
82c1.Divide(2, 2)
83c1.cd(1)
84ROOT.gPad.SetLeftMargin(0.15)
85hh_data.GetZaxis().SetTitleOffset(1.4)
86hh_data.Draw("lego")
87c1.cd(2)
88ROOT.gPad.SetLeftMargin(0.20)
89hh_pdf.GetZaxis().SetTitleOffset(2.5)
90hh_pdf.Draw("surf")
91c1.cd(3)
92ROOT.gPad.SetLeftMargin(0.15)
93hh_data.GetZaxis().SetTitleOffset(1.4)
94hh_data.Draw("box")
95c1.cd(4)
96ROOT.gPad.SetLeftMargin(0.15)
97hh_pdf.GetZaxis().SetTitleOffset(2.5)
98hh_pdf.Draw("cont3")
99c1.SaveAs("rf309_2dimplot.png")
100
101c2 = ROOT.TCanvas("rf309_3dimplot", "rf309_3dimplot", 800, 400)
102c2.Divide(2)
103c2.cd(1)
104ROOT.gPad.SetLeftMargin(0.15)
105hh_data3.GetZaxis().SetTitleOffset(1.4)
106hh_data3.Draw("lego")
107c2.cd(2)
108ROOT.gPad.SetLeftMargin(0.15)
109hh_pdf3.GetZaxis().SetTitleOffset(1.4)
110hh_pdf3.Draw("iso")
111c2.SaveAs("rf309_3dimplot.png")