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