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