Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
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_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# B-physics pdf with per-event Gaussian resolution
16# ----------------------------------------------------------------------------------------------
17
18# Observables
19dt = ROOT.RooRealVar("dt", "dt", -10, 10)
20dterr = ROOT.RooRealVar("dterr", "per-event error on dt", 0.01, 10)
21
22# Build a gaussian resolution model scaled by the per-error =
23# gauss(dt,bias,sigma*dterr)
24bias = ROOT.RooRealVar("bias", "bias", 0, -10, 10)
25sigma = ROOT.RooRealVar("sigma", "per-event error scale factor", 1, 0.1, 10)
26gm = ROOT.RooGaussModel("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, tau, gm, type="DoubleSided")
31
32# Construct fake 'external' data with per-event error
33# ------------------------------------------------------------------------------------------------------
34
35# Use landau pdf to get somewhat realistic distribution with long tail
36pdfDtErr = ROOT.RooLandau("pdfDtErr", "pdfDtErr", dterr, 1.0, 0.25)
37expDataDterr = pdfDtErr.generate({dterr}, 10000)
38
39# Sample data from conditional decay_gm(dt|dterr)
40# ---------------------------------------------------------------------------------------------
41
42# Specify external dataset with dterr values to use decay_dm as
43# conditional pdf
44data = decay_gm.generate({dt}, ProtoData=expDataDterr)
45
46# Fit conditional decay_dm(dt|dterr)
47# ---------------------------------------------------------------------
48
49# Specify dterr as conditional observable
50decay_gm.fitTo(data, ConditionalObservables={dterr}, PrintLevel=-1)
51
52# Plot conditional decay_dm(dt|dterr)
53# ---------------------------------------------------------------------
54
55# Make two-dimensional plot of conditional pdf in (dt,dterr)
56hh_decay = decay_gm.createHistogram("hh_decay", dt, Binning=50, YVar=dict(var=dterr, Binning=50))
57hh_decay.SetLineColor(ROOT.kBlue)
58
59# Plot decay_gm(dt|dterr) at various values of dterr
60frame = dt.frame(Title="Slices of decay(dt|dterr) at various dterr")
61for ibin in range(0, 100, 20):
62 dterr.setBin(ibin)
63 decay_gm.plotOn(frame, Normalization=5.0)
64
65# Make projection of data an dt
66frame2 = dt.frame(Title="Projection of decay(dt|dterr) on dt")
67data.plotOn(frame2)
68
69# Make projection of decay(dt|dterr) on dt.
70#
71# Instead of integrating out dterr, a weighted average of curves
72# at values dterr_i as given in the external dataset.
73# (The kTRUE argument bins the data before projection to speed up the process)
74decay_gm.plotOn(frame2, ProjWData=(expDataDterr, True))
75
76# Draw all frames on canvas
77c = ROOT.TCanvas("rf306_condpereventerrors", "rf306_condperventerrors", 1200, 400)
78c.Divide(3)
79c.cd(1)
80ROOT.gPad.SetLeftMargin(0.20)
81hh_decay.GetZaxis().SetTitleOffset(2.5)
82hh_decay.Draw("surf")
83c.cd(2)
84ROOT.gPad.SetLeftMargin(0.15)
85frame.GetYaxis().SetTitleOffset(1.6)
86frame.Draw()
87c.cd(3)
88ROOT.gPad.SetLeftMargin(0.15)
89frame2.GetYaxis().SetTitleOffset(1.6)
90frame2.Draw()
91
92c.SaveAs("rf306_condpereventerrors.png")