Rf 3 1 1_Rangeplot

Multidimensional models: projecting p.d.f and data ranges in continuous observables

Author: Clemens Lange, Wouter Verkerke (C++ version)
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Sunday, July 05, 2020 at 08:24 AM.

In [ ]:
import ROOT

Create 3D pdf and data

Create observables

In [ ]:
x = ROOT.RooRealVar("x", "x", -5, 5)
y = ROOT.RooRealVar("y", "y", -5, 5)
z = ROOT.RooRealVar("z", "z", -5, 5)

Create signal pdf gauss(x)gauss(y)gauss(z)

In [ ]:
gx = ROOT.RooGaussian(
    "gx", "gx", x, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
gy = ROOT.RooGaussian(
    "gy", "gy", y, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
gz = ROOT.RooGaussian(
    "gz", "gz", z, ROOT.RooFit.RooConst(0), ROOT.RooFit.RooConst(1))
sig = ROOT.RooProdPdf("sig", "sig", ROOT.RooArgList(gx, gy, gz))

Create background pdf poly(x)poly(y)poly(z)

In [ ]:
px = ROOT.RooPolynomial("px", "px", x, ROOT.RooArgList(
    ROOT.RooFit.RooConst(-0.1), ROOT.RooFit.RooConst(0.004)))
py = ROOT.RooPolynomial("py", "py", y, ROOT.RooArgList(
    ROOT.RooFit.RooConst(0.1), ROOT.RooFit.RooConst(-0.004)))
pz = ROOT.RooPolynomial("pz", "pz", z)
bkg = ROOT.RooProdPdf("bkg", "bkg", ROOT.RooArgList(px, py, pz))

Create composite pdf sig+bkg

In [ ]:
fsig = ROOT.RooRealVar("fsig", "signal fraction", 0.1, 0., 1.)
model = ROOT.RooAddPdf(
    "model", "model", ROOT.RooArgList(
        sig, bkg), ROOT.RooArgList(fsig))

data = model.generate(ROOT.RooArgSet(x, y, z), 20000)

Project pdf and data on x

Make plain projection of data and pdf on x observable

In [ ]:
frame = x.frame(ROOT.RooFit.Title(
    "Projection of 3D data and pdf on X"), ROOT.RooFit.Bins(40))
data.plotOn(frame)
model.plotOn(frame)

Project pdf and data on x in signal range

Define signal region in y and z observables

In [ ]:
y.setRange("sigRegion", -1, 1)
z.setRange("sigRegion", -1, 1)

Make plot frame

In [ ]:
frame2 = x.frame(ROOT.RooFit.Title(
    "Same projection on X in signal range of (Y,Z)"), ROOT.RooFit.Bins(40))

Plot subset of data in which all observables are inside "sigRegion" For observables that do not have an explicit "sigRegion" range defined (e.g. observable) an implicit definition is used that is identical to the full range (i.e. [-5,5] for x)

In [ ]:
data.plotOn(frame2, ROOT.RooFit.CutRange("sigRegion"))

Project model on x, projected observables (y,z) only in "sigRegion"

In [ ]:
model.plotOn(frame2, ROOT.RooFit.ProjectionRange("sigRegion"))

c = ROOT.TCanvas("rf311_rangeplot", "rf310_rangeplot", 800, 400)
c.Divide(2)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
frame.GetYaxis().SetTitleOffset(1.4)
frame.Draw()
c.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
frame2.GetYaxis().SetTitleOffset(1.4)
frame2.Draw()

c.SaveAs("rf311_rangeplot.png")

Draw all canvases

In [ ]:
from ROOT import gROOT 
gROOT.GetListOfCanvases().Draw()