Special pdf's: using a product of an (acceptance) efficiency and a pdf as pdf
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 Tuesday, March 19, 2024 at 07:17 PM.
import ROOT
Declare observables
t = ROOT.RooRealVar("t", "t", 0, 5)
Make pdf
tau = ROOT.RooRealVar("tau", "tau", -1.54, -4, -0.1)
model = ROOT.RooExponential("model", "model", t, tau)
Use error function to simulate turn-on slope
eff = ROOT.RooFormulaVar("eff", "0.5*(TMath::Erf((t-1)/0.5)+1)", [t])
Multiply pdf(t) with efficiency in t
modelEff = ROOT.RooEffProd("modelEff", "model with efficiency", model, eff)
frame1 = t.frame(Title="Efficiency")
eff.plotOn(frame1, LineColor="r")
frame2 = t.frame(Title="Pdf with and without efficiency")
model.plotOn(frame2, LineStyle="--")
modelEff.plotOn(frame2)
<cppyy.gbl.RooPlot object at 0x94f29d0>
[#1] INFO:NumericIntegration -- RooRealIntegral::init(modelEff_Int[t]) using numeric integrator RooIntegrator1D to calculate Int(t)
Generate events. If the input pdf has an internal generator, internal generator is used and an accept/reject sampling on the efficiency is applied.
data = modelEff.generate({t}, 10000)
Fit pdf. The normalization integral is calculated numerically.
modelEff.fitTo(data, PrintLevel=-1)
<cppyy.gbl.RooFitResult object at 0x(nil)>
[#1] INFO:Fitting -- RooAbsPdf::fitTo(modelEff_over_modelEff_Int[t]) fixing normalization set for coefficient determination to observables in data [#1] INFO:Fitting -- using CPU computation library compiled with -mavx2 [#1] INFO:Fitting -- RooAddition::defaultErrorLevel(nll_modelEff_over_modelEff_Int[t]_modelEffData) Summation contains a RooNLLVar, using its error level [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: activating const optimization [#1] INFO:NumericIntegration -- RooRealIntegral::init(modelEff_Int[t]) using numeric integrator RooIntegrator1D to calculate Int(t) [#1] INFO:Minimization -- RooAbsMinimizerFcn::setOptimizeConst: deactivating const optimization
Plot generated data and overlay fitted pdf
frame3 = t.frame(Title="Fitted pdf with efficiency")
data.plotOn(frame3)
modelEff.plotOn(frame3)
c = ROOT.TCanvas("rf703_effpdfprod", "rf703_effpdfprod", 1200, 400)
c.Divide(3)
c.cd(1)
ROOT.gPad.SetLeftMargin(0.15)
frame1.GetYaxis().SetTitleOffset(1.4)
frame1.Draw()
c.cd(2)
ROOT.gPad.SetLeftMargin(0.15)
frame2.GetYaxis().SetTitleOffset(1.6)
frame2.Draw()
c.cd(3)
ROOT.gPad.SetLeftMargin(0.15)
frame3.GetYaxis().SetTitleOffset(1.6)
frame3.Draw()
c.SaveAs("rf703_effpdfprod.png")
[#1] INFO:NumericIntegration -- RooRealIntegral::init(modelEff_Int[t]) using numeric integrator RooIntegrator1D to calculate Int(t)
Info in <TCanvas::Print>: png file rf703_effpdfprod.png has been created
Draw all canvases
from ROOT import gROOT
gROOT.GetListOfCanvases().Draw()