Logo ROOT  
Reference Guide
rf802_mcstudy_addons.py
Go to the documentation of this file.
1## \ingroup tutorial_roofit
2## \notebook
3##
4## 'VALIDATION AND MC STUDIES' RooFit tutorial macro #802
5##
6## RooMCStudy: using separate fit and generator models, the chi^2 calculator model
7##
8## \macro_code
9##
10## \date February 2018
11## \author Clemens Lange
12
13
14import ROOT
15
16
17# Create model
18# -----------------------
19
20# Observables, parameters
21x = ROOT.RooRealVar("x", "x", -10, 10)
22x.setBins(10)
23mean = ROOT.RooRealVar("mean", "mean of gaussian", 0, -2.0, 1.8)
24sigma = ROOT.RooRealVar("sigma", "width of gaussian", 5, 1, 10)
25
26# Create Gaussian pdf
27gauss = ROOT.RooGaussian("gauss", "gaussian PDF", x, mean, sigma)
28
29# Create manager with chi^2 add-on module
30# ----------------------------------------------------------------------------
31
32# Create study manager for binned likelihood fits of a Gaussian pdf in 10
33# bins
34mcs = ROOT.RooMCStudy(gauss, {x}, ROOT.RooFit.Silence(), ROOT.RooFit.Binned())
35
36# Add chi^2 calculator module to mcs
37chi2mod = ROOT.RooChi2MCSModule()
38mcs.addModule(chi2mod)
39
40# Generate 1000 samples of 1000 events
41mcs.generateAndFit(2000, 1000)
42
43# Fill histograms with distributions chi2 and prob(chi2,ndf) that
44# are calculated by ROOT.RooChiMCSModule
45hist_chi2 = ROOT.RooAbsData.createHistogram(mcs.fitParDataSet(), "chi2")
46hist_prob = ROOT.RooAbsData.createHistogram(mcs.fitParDataSet(), "prob")
47
48# Create manager with separate fit model
49# ----------------------------------------------------------------------------
50
51# Create alternate pdf with shifted mean
52mean2 = ROOT.RooRealVar("mean2", "mean of gaussian 2", 2.0)
53gauss2 = ROOT.RooGaussian("gauss2", "gaussian PDF2", x, mean2, sigma)
54
55# Create study manager with separate generation and fit model. ROOT.This configuration
56# is set up to generate bad fits as the fit and generator model have different means
57# and the mean parameter is not floating in the fit
58mcs2 = ROOT.RooMCStudy(gauss2, {x}, ROOT.RooFit.FitModel(gauss), ROOT.RooFit.Silence(), ROOT.RooFit.Binned())
59
60# Add chi^2 calculator module to mcs
61chi2mod2 = ROOT.RooChi2MCSModule()
62mcs2.addModule(chi2mod2)
63
64# Generate 1000 samples of 1000 events
65mcs2.generateAndFit(2000, 1000)
66
67# Request a the pull plot of mean. The pulls will be one-sided because
68# `mean` is limited to 1.8.
69# Note that RooFit will have trouble to compute the pulls because the parameters
70# are called `mean` in the fit, but `mean2` in the generator model. It is not obvious
71# that these are related. RooFit will nevertheless compute pulls, but complain that
72# this is risky.
73pullMeanFrame = mcs2.plotPull(mean)
74
75# Fill histograms with distributions chi2 and prob(chi2,ndf) that
76# are calculated by ROOT.RooChiMCSModule
77hist2_chi2 = ROOT.RooAbsData.createHistogram(mcs2.fitParDataSet(), "chi2")
78hist2_prob = ROOT.RooAbsData.createHistogram(mcs2.fitParDataSet(), "prob")
79hist2_chi2.SetLineColor(ROOT.kRed)
80hist2_prob.SetLineColor(ROOT.kRed)
81
82c = ROOT.TCanvas("rf802_mcstudy_addons", "rf802_mcstudy_addons", 800, 400)
83c.Divide(3)
84c.cd(1)
85ROOT.gPad.SetLeftMargin(0.15)
86hist_chi2.GetYaxis().SetTitleOffset(1.4)
87hist_chi2.Draw()
88hist2_chi2.Draw("esame")
89c.cd(2)
90ROOT.gPad.SetLeftMargin(0.15)
91hist_prob.GetYaxis().SetTitleOffset(1.4)
92hist_prob.Draw()
93hist2_prob.Draw("esame")
94c.cd(3)
95pullMeanFrame.Draw()
96
97c.SaveAs("rf802_mcstudy_addons.png")
98
99# Make RooMCStudy object available on command line after
100# macro finishes
101ROOT.gDirectory.Add(mcs)