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