Logo ROOT  
Reference Guide
rf306_condpereventerrors.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Multidimensional models: complete example with use of conditional pdf with per-event errors
5##
6## \macro_code
7##
8## \date February 2018
9## \authors Clemens Lange, Wouter Verkerke (C++ version)
10
11import ROOT
12
13# B-physics pdf with per-event Gaussian resolution
14# ----------------------------------------------------------------------------------------------
15
16# Observables
17dt = ROOT.RooRealVar("dt", "dt", -10, 10)
18dterr = ROOT.RooRealVar("dterr", "per-event error on dt", 0.01, 10)
19
20# Build a gaussian resolution model scaled by the per-error =
21# gauss(dt,bias,sigma*dterr)
22bias = ROOT.RooRealVar("bias", "bias", 0, -10, 10)
23sigma = ROOT.RooRealVar(
24 "sigma", "per-event error scale factor", 1, 0.1, 10)
25gm = ROOT.RooGaussModel(
26 "gm1", "gauss model scaled bt per-event error", dt, bias, sigma, dterr)
27
28# Construct decay(dt) (x) gauss1(dt|dterr)
29tau = ROOT.RooRealVar("tau", "tau", 1.548)
30decay_gm = ROOT.RooDecay("decay_gm", "decay", dt,
31 tau, gm, ROOT.RooDecay.DoubleSided)
32
33# Construct fake 'external' data with per-event error
34# ------------------------------------------------------------------------------------------------------
35
36# Use landau pdf to get somewhat realistic distribution with long tail
37pdfDtErr = ROOT.RooLandau("pdfDtErr", "pdfDtErr", dterr, ROOT.RooFit.RooConst(
38 1), ROOT.RooFit.RooConst(0.25))
39expDataDterr = pdfDtErr.generate(ROOT.RooArgSet(dterr), 10000)
40
41# Sample data from conditional decay_gm(dt|dterr)
42# ---------------------------------------------------------------------------------------------
43
44# Specify external dataset with dterr values to use decay_dm as
45# conditional pdf
46data = decay_gm.generate(ROOT.RooArgSet(
47 dt), ROOT.RooFit.ProtoData(expDataDterr))
48
49# Fit conditional decay_dm(dt|dterr)
50# ---------------------------------------------------------------------
51
52# Specify dterr as conditional observable
53decay_gm.fitTo(data, ROOT.RooFit.ConditionalObservables(
54 ROOT.RooArgSet(dterr)))
55
56# Plot conditional decay_dm(dt|dterr)
57# ---------------------------------------------------------------------
58
59# Make two-dimensional plot of conditional pdf in (dt,dterr)
60hh_decay = decay_gm.createHistogram("hh_decay", dt, ROOT.RooFit.Binning(
61 50), ROOT.RooFit.YVar(dterr, ROOT.RooFit.Binning(50)))
62hh_decay.SetLineColor(ROOT.kBlue)
63
64# Plot decay_gm(dt|dterr) at various values of dterr
65frame = dt.frame(ROOT.RooFit.Title(
66 "Slices of decay(dt|dterr) at various dterr"))
67for ibin in range(0, 100, 20):
68 dterr.setBin(ibin)
69 decay_gm.plotOn(frame, ROOT.RooFit.Normalization(5.))
70
71# Make projection of data an dt
72frame2 = dt.frame(ROOT.RooFit.Title("Projection of decay(dt|dterr) on dt"))
73data.plotOn(frame2)
74
75# Make projection of decay(dt|dterr) on dt.
76#
77# Instead of integrating out dterr, a weighted average of curves
78# at values dterr_i as given in the external dataset.
79# (The kTRUE argument bins the data before projection to speed up the process)
80decay_gm.plotOn(frame2, ROOT.RooFit.ProjWData(expDataDterr, ROOT.kTRUE))
81
82# Draw all frames on canvas
83c = ROOT.TCanvas("rf306_condpereventerrors",
84 "rf306_condperventerrors", 1200, 400)
85c.Divide(3)
86c.cd(1)
87ROOT.gPad.SetLeftMargin(0.20)
88hh_decay.GetZaxis().SetTitleOffset(2.5)
89hh_decay.Draw("surf")
90c.cd(2)
91ROOT.gPad.SetLeftMargin(0.15)
92frame.GetYaxis().SetTitleOffset(1.6)
93frame.Draw()
94c.cd(3)
95ROOT.gPad.SetLeftMargin(0.15)
96frame2.GetYaxis().SetTitleOffset(1.6)
97frame2.Draw()
98
99c.SaveAs("rf306_condpereventerrors.png")