Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf702_efficiencyfit_2D.py File Reference

Namespaces

namespace  rf702_efficiencyfit_2D
 

Detailed Description

View in nbviewer Open in SWAN
Special pdf's: unbinned maximum likelihood fit of an efficiency eff(x) function to a dataset D(x,cut), cut is a category encoding a selection whose efficiency as function of x should be described by eff(x)

import ROOT
flat = False
# Construct efficiency function e(x,y)
# -----------------------------------------------------------------------
# Declare variables x,mean, with associated name, title, value and allowed
# range
x = ROOT.RooRealVar("x", "x", -10, 10)
y = ROOT.RooRealVar("y", "y", -10, 10)
# Efficiency function eff(x;a,b)
ax = ROOT.RooRealVar("ax", "ay", 0.6, 0, 1)
bx = ROOT.RooRealVar("bx", "by", 5)
cx = ROOT.RooRealVar("cx", "cy", -1, -10, 10)
ay = ROOT.RooRealVar("ay", "ay", 0.2, 0, 1)
by = ROOT.RooRealVar("by", "by", 5)
cy = ROOT.RooRealVar("cy", "cy", -1, -10, 10)
effFunc = ROOT.RooFormulaVar(
"effFunc", "((1-ax)+ax*cos((x-cx)/bx))*((1-ay)+ay*cos((y-cy)/by))", [ax, bx, cx, x, ay, by, cy, y]
)
# Acceptance state cut (1 or 0)
cut = ROOT.RooCategory("cut", "cutr", {"accept": 1, "reject": 0})
# Construct conditional efficiency pdf E(cut|x,y)
# ---------------------------------------------------------------------------------------------
# Construct efficiency pdf eff(cut|x)
effPdf = ROOT.RooEfficiency("effPdf", "effPdf", effFunc, cut, "accept")
# Generate data(x,y,cut) from a toy model
# -------------------------------------------------------------------------------
# Construct global shape pdf shape(x) and product model(x,cut) = eff(cut|x)*shape(x)
# (These are _only_ needed to generate some toy MC here to be used later)
shapePdfX = ROOT.RooPolynomial("shapePdfX", "shapePdfX", x, [0 if flat else -0.095])
shapePdfY = ROOT.RooPolynomial("shapePdfY", "shapePdfY", y, [0 if flat else +0.095])
shapePdf = ROOT.RooProdPdf("shapePdf", "shapePdf", [shapePdfX, shapePdfY])
model = ROOT.RooProdPdf("model", "model", {shapePdf}, Conditional=({effPdf}, {cut}))
# Generate some toy data from model
data = model.generate({x, y, cut}, 10000)
# Fit conditional efficiency pdf to data
# --------------------------------------------------------------------------
# Fit conditional efficiency pdf to data
effPdf.fitTo(data, ConditionalObservables={x, y}, PrintLevel=-1)
# Plot fitted, data efficiency
# --------------------------------------------------------
# Make 2D histograms of all data, data and efficiency function
hh_data_all = data.createHistogram("hh_data_all", x, Binning=8, YVar=dict(var=y, Binning=8))
hh_data_sel = data.createHistogram("hh_data_sel", x, Binning=8, YVar=dict(var=y, Binning=8), Cut="cut==cut::accept")
hh_eff = effFunc.createHistogram("hh_eff", x, Binning=50, YVar=dict(var=y, Binning=50))
# Some adjustsment for good visualization
hh_data_all.SetMinimum(0)
hh_data_sel.SetMinimum(0)
hh_eff.SetMinimum(0)
hh_eff.SetLineColor(ROOT.kBlue)
# Draw all frames on a canvas
ca = ROOT.TCanvas("rf702_efficiency_2D", "rf702_efficiency_2D", 1200, 400)
ca.Divide(3)
ca.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
hh_data_all.GetZaxis().SetTitleOffset(1.8)
hh_data_all.Draw("lego")
ca.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
hh_data_sel.GetZaxis().SetTitleOffset(1.8)
hh_data_sel.Draw("lego")
ca.cd(3)
ROOT.gPad.SetLeftMargin(0.15)
hh_eff.GetZaxis().SetTitleOffset(1.8)
hh_eff.Draw("surf")
ca.SaveAs("rf702_efficiency_2D.png")
[#0] WARNING:Generation -- RooAcceptReject::ctor(effPdf_Int[]_Norm[cut]) WARNING: performing accept/reject sampling on a p.d.f in 2 dimensions without prior knowledge on maximum value of p.d.f. Determining maximum value by taking 200000 trial samples. If p.d.f contains sharp peaks smaller than average distance between trial sampling points these may be missed and p.d.f. may be sampled incorrectly.
[#0] WARNING:Generation -- RooAcceptReject::ctor(effPdf_Int[]_Norm[cut]): WARNING: 200000 trial samples requested by p.d.f for 2-dimensional accept/reject sampling, this may take some time
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization
[#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Date
February 2018
Authors
Clemens Lange, Wouter Verkerke (C++ version)

Definition in file rf702_efficiencyfit_2D.py.