Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf703_effpdfprod.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Special pdf's: using a product of an (acceptance) efficiency and a pdf as pdf
5##
6## \macro_image
7## \macro_code
8## \macro_output
9##
10## \date February 2018
11## \authors Clemens Lange, Wouter Verkerke (C++ version)
12
13import ROOT
14
15
16# Define observables and decay pdf
17# ---------------------------------------------------------------
18
19# Declare observables
20t = ROOT.RooRealVar("t", "t", 0, 5)
21
22# Make pdf
23tau = ROOT.RooRealVar("tau", "tau", -1.54, -4, -0.1)
24model = ROOT.RooExponential("model", "model", t, tau)
25
26# Define efficiency function
27# ---------------------------------------------------
28
29# Use error function to simulate turn-on slope
30eff = ROOT.RooFormulaVar("eff", "0.5*(TMath::Erf((t-1)/0.5)+1)", [t])
31
32# Define decay pdf with efficiency
33# ---------------------------------------------------------------
34
35# Multiply pdf(t) with efficiency in t
36modelEff = ROOT.RooEffProd("modelEff", "model with efficiency", model, eff)
37
38# Plot efficiency, pdf
39# ----------------------------------------
40
41frame1 = t.frame(Title="Efficiency")
42eff.plotOn(frame1, LineColor="r")
43
44frame2 = t.frame(Title="Pdf with and without efficiency")
45
46model.plotOn(frame2, LineStyle="--")
47modelEff.plotOn(frame2)
48
49# Generate toy data, fit model eff to data
50# ------------------------------------------------------------------------------
51
52# Generate events. If the input pdf has an internal generator, internal generator
53# is used and an accept/reject sampling on the efficiency is applied.
54data = modelEff.generate({t}, 10000)
55
56# Fit pdf. The normalization integral is calculated numerically.
57modelEff.fitTo(data, PrintLevel=-1)
58
59# Plot generated data and overlay fitted pdf
60frame3 = t.frame(Title="Fitted pdf with efficiency")
61data.plotOn(frame3)
62modelEff.plotOn(frame3)
63
64c = ROOT.TCanvas("rf703_effpdfprod", "rf703_effpdfprod", 1200, 400)
65c.Divide(3)
66c.cd(1)
67ROOT.gPad.SetLeftMargin(0.15)
68frame1.GetYaxis().SetTitleOffset(1.4)
69frame1.Draw()
70c.cd(2)
71ROOT.gPad.SetLeftMargin(0.15)
72frame2.GetYaxis().SetTitleOffset(1.6)
73frame2.Draw()
74c.cd(3)
75ROOT.gPad.SetLeftMargin(0.15)
76frame3.GetYaxis().SetTitleOffset(1.6)
77frame3.Draw()
78
79c.SaveAs("rf703_effpdfprod.png")