Logo ROOT  
Reference Guide
rf316_llratioplot.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Multidimensional models: using the likelihood ratio techique to construct a signal
5## enhanced one-dimensional projection of a multi-dimensional pdf
6##
7## \macro_code
8##
9## \date February 2018
10## \authors Clemens Lange, Wouter Verkerke (C++ version)
11
12import ROOT
13
14
15# Create 3D pdf and data
16# -------------------------------------------
17
18# Create observables
19x = ROOT.RooRealVar("x", "x", -5, 5)
20y = ROOT.RooRealVar("y", "y", -5, 5)
21z = ROOT.RooRealVar("z", "z", -5, 5)
22
23# Create signal pdf gauss(x)*gauss(y)*gauss(z)
24gx = ROOT.RooGaussian(
25 "gx", "gx", x, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
26gy = ROOT.RooGaussian(
27 "gy", "gy", y, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
28gz = ROOT.RooGaussian(
29 "gz", "gz", z, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
30sig = ROOT.RooProdPdf("sig", "sig", ROOT.RooArgList(gx, gy, gz))
31
32# Create background pdf poly(x)*poly(y)*poly(z)
33px = ROOT.RooPolynomial("px", "px", x, ROOT.RooArgList(
34 ROOT.RooFit.RooConst(-0.1), ROOT.RooFit.RooConst(0.004)))
35py = ROOT.RooPolynomial("py", "py", y, ROOT.RooArgList(
36 ROOT.RooFit.RooConst(0.1), ROOT.RooFit.RooConst(-0.004)))
37pz = ROOT.RooPolynomial("pz", "pz", z)
38bkg = ROOT.RooProdPdf("bkg", "bkg", ROOT.RooArgList(px, py, pz))
39
40# Create composite pdf sig+bkg
41fsig = ROOT.RooRealVar("fsig", "signal fraction", 0.1, 0., 1.)
42model = ROOT.RooAddPdf(
43 "model", "model", ROOT.RooArgList(
44 sig, bkg), ROOT.RooArgList(fsig))
45
46data = model.generate(ROOT.RooArgSet(x, y, z), 20000)
47
48# Project pdf and data on x
49# -------------------------------------------------
50
51# Make plain projection of data and pdf on x observable
52frame = x.frame(ROOT.RooFit.Title(
53 "Projection of 3D data and pdf on X"), ROOT.RooFit.Bins(40))
54data.plotOn(frame)
55model.plotOn(frame)
56
57# Define projected signal likelihood ratio
58# ----------------------------------------------------------------------------------
59
60# Calculate projection of signal and total likelihood on (y,z) observables
61# i.e. integrate signal and composite model over x
62sigyz = sig.createProjection(ROOT.RooArgSet(x))
63totyz = model.createProjection(ROOT.RooArgSet(x))
64
65# Construct the log of the signal / signal+background probability
66llratio_func = ROOT.RooFormulaVar(
67 "llratio", "log10(@0)-log10(@1)", ROOT.RooArgList(sigyz, totyz))
68
69# Plot data with a LL ratio cut
70# -------------------------------------------------------
71
72# Calculate the llratio value for each event in the dataset
73data.addColumn(llratio_func)
74
75# Extract the subset of data with large signal likelihood
76dataSel = data.reduce(ROOT.RooFit.Cut("llratio>0.7"))
77
78# Make plot frame
79frame2 = x.frame(ROOT.RooFit.Title(
80 "Same projection on X with LLratio(y,z)>0.7"), ROOT.RooFit.Bins(40))
81
82# Plot select data on frame
83dataSel.plotOn(frame2)
84
85# Make MC projection of pdf with same LL ratio cut
86# ---------------------------------------------------------------------------------------------
87
88# Generate large number of events for MC integration of pdf projection
89mcprojData = model.generate(ROOT.RooArgSet(x, y, z), 10000)
90
91# Calculate LL ratio for each generated event and select MC events with
92# llratio)0.7
93mcprojData.addColumn(llratio_func)
94mcprojDataSel = mcprojData.reduce(ROOT.RooFit.Cut("llratio>0.7"))
95
96# Project model on x, projected observables (y,z) with Monte Carlo technique
97# on set of events with the same llratio cut as was applied to data
98model.plotOn(frame2, ROOT.RooFit.ProjWData(mcprojDataSel))
99
100c = ROOT.TCanvas("rf316_llratioplot", "rf316_llratioplot", 800, 400)
101c.Divide(2)
102c.cd(1)
103ROOT.gPad.SetLeftMargin(0.15)
104frame.GetYaxis().SetTitleOffset(1.4)
105frame.Draw()
106c.cd(2)
107ROOT.gPad.SetLeftMargin(0.15)
108frame2.GetYaxis().SetTitleOffset(1.4)
109frame2.Draw()
110c.SaveAs("rf316_llratioplot.png")