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