Logo ROOT  
Reference Guide
rf608_fitresultaspdf.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ##
5 ## Likelihood and minimization: representing the parabolic approximation of the fit as a
6 ## multi-variate Gaussian on the parameters of the fitted p.d.f.
7 ##
8 ## \macro_code
9 ##
10 ## \date February 2018
11 ## \authors Clemens Lange, Wouter Verkerke (C++ version)
12 
13 import ROOT
14 
15 
16 # Create model and dataset
17 # -----------------------------------------------
18 
19 # Observable
20 x = ROOT.RooRealVar("x", "x", -20, 20)
21 
22 # Model (intentional strong correlations)
23 mean = ROOT.RooRealVar("mean", "mean of g1 and g2", 0, -1, 1)
24 sigma_g1 = ROOT.RooRealVar("sigma_g1", "width of g1", 2)
25 g1 = ROOT.RooGaussian("g1", "g1", x, mean, sigma_g1)
26 
27 sigma_g2 = ROOT.RooRealVar("sigma_g2", "width of g2", 4, 3.0, 5.0)
28 g2 = ROOT.RooGaussian("g2", "g2", x, mean, sigma_g2)
29 
30 frac = ROOT.RooRealVar("frac", "frac", 0.5, 0.0, 1.0)
31 model = ROOT.RooAddPdf(
32  "model", "model", ROOT.RooArgList(
33  g1, g2), ROOT.RooArgList(frac))
34 
35 # Generate 1000 events
36 data = model.generate(ROOT.RooArgSet(x), 1000)
37 
38 # Fit model to data
39 # ----------------------------------
40 
41 r = model.fitTo(data, ROOT.RooFit.Save())
42 
43 # Create MV Gaussian pdf of fitted parameters
44 # ------------------------------------------------------------------------------------
45 
46 parabPdf = r.createHessePdf(ROOT.RooArgSet(frac, mean, sigma_g2))
47 
48 # Some exercises with the parameter pdf
49 # -----------------------------------------------------------------------------
50 
51 # Generate 100K points in the parameter space, from the MVGaussian p.d.f.
52 d = parabPdf.generate(ROOT.RooArgSet(mean, sigma_g2, frac), 100000)
53 
54 # Sample a 3-D histogram of the p.d.f. to be visualized as an error
55 # ellipsoid using the GLISO draw option
56 hh_3d = parabPdf.createHistogram("mean,sigma_g2,frac", 25, 25, 25)
57 hh_3d.SetFillColor(ROOT.kBlue)
58 
59 # Project 3D parameter p.d.f. down to 3 permutations of two-dimensional p.d.f.s
60 # The integrations corresponding to these projections are performed analytically
61 # by the MV Gaussian p.d.f.
62 pdf_sigmag2_frac = parabPdf.createProjection(ROOT.RooArgSet(mean))
63 pdf_mean_frac = parabPdf.createProjection(ROOT.RooArgSet(sigma_g2))
64 pdf_mean_sigmag2 = parabPdf.createProjection(ROOT.RooArgSet(frac))
65 
66 # Make 2D plots of the 3 two-dimensional p.d.f. projections
67 hh_sigmag2_frac = pdf_sigmag2_frac.createHistogram("sigma_g2,frac", 50, 50)
68 hh_mean_frac = pdf_mean_frac.createHistogram("mean,frac", 50, 50)
69 hh_mean_sigmag2 = pdf_mean_sigmag2.createHistogram("mean,sigma_g2", 50, 50)
70 hh_mean_frac.SetLineColor(ROOT.kBlue)
71 hh_sigmag2_frac.SetLineColor(ROOT.kBlue)
72 hh_mean_sigmag2.SetLineColor(ROOT.kBlue)
73 
74 # Draw the 'sigar'
75 ROOT.gStyle.SetCanvasPreferGL(True)
76 ROOT.gStyle.SetPalette(1)
77 c1 = ROOT.TCanvas("rf608_fitresultaspdf_1", "rf608_fitresultaspdf_1", 600, 600)
78 hh_3d.Draw("gliso")
79 
80 c1.SaveAs("rf608_fitresultaspdf_1.png")
81 
82 # Draw the 2D projections of the 3D p.d.f.
83 c2 = ROOT.TCanvas("rf608_fitresultaspdf_2",
84  "rf608_fitresultaspdf_2", 900, 600)
85 c2.Divide(3, 2)
86 c2.cd(1)
87 ROOT.gPad.SetLeftMargin(0.15)
88 hh_mean_sigmag2.GetZaxis().SetTitleOffset(1.4)
89 hh_mean_sigmag2.Draw("surf3")
90 c2.cd(2)
91 ROOT.gPad.SetLeftMargin(0.15)
92 hh_sigmag2_frac.GetZaxis().SetTitleOffset(1.4)
93 hh_sigmag2_frac.Draw("surf3")
94 c2.cd(3)
95 ROOT.gPad.SetLeftMargin(0.15)
96 hh_mean_frac.GetZaxis().SetTitleOffset(1.4)
97 hh_mean_frac.Draw("surf3")
98 
99 # Draw the distributions of parameter points sampled from the p.d.f.
100 tmp1 = d.createHistogram(mean, sigma_g2, 50, 50)
101 tmp2 = d.createHistogram(sigma_g2, frac, 50, 50)
102 tmp3 = d.createHistogram(mean, frac, 50, 50)
103 
104 c2.cd(4)
105 ROOT.gPad.SetLeftMargin(0.15)
106 tmp1.GetZaxis().SetTitleOffset(1.4)
107 tmp1.Draw("lego3")
108 c2.cd(5)
109 ROOT.gPad.SetLeftMargin(0.15)
110 tmp2.GetZaxis().SetTitleOffset(1.4)
111 tmp2.Draw("lego3")
112 c2.cd(6)
113 ROOT.gPad.SetLeftMargin(0.15)
114 tmp3.GetZaxis().SetTitleOffset(1.4)
115 tmp3.Draw("lego3")
116 
117 c2.SaveAs("rf608_fitresultaspdf_2.png")