Logo ROOT  
Reference Guide
rf307_fullpereventerrors.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ##
5 ## Multidimensional models: usage of full 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 empirical pdf for per-event error
35 # -----------------------------------------------------------------
36 
37 # Use landau p.d.f to get empirical 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 # Construct a histogram pdf to describe the shape of the dtErr distribution
43 expHistDterr = expDataDterr.binnedClone()
44 pdfErr = ROOT.RooHistPdf(
45  "pdfErr", "pdfErr", ROOT.RooArgSet(dterr), expHistDterr)
46 
47 # Construct conditional product decay_dm(dt|dterr)*pdf(dterr)
48 # ----------------------------------------------------------------------------------------------------------------------
49 
50 # Construct production of conditional decay_dm(dt|dterr) with empirical
51 # pdfErr(dterr)
52 model = ROOT.RooProdPdf(
53  "model",
54  "model",
55  ROOT.RooArgSet(pdfErr),
56  ROOT.RooFit.Conditional(
57  ROOT.RooArgSet(decay_gm),
58  ROOT.RooArgSet(dt)))
59 
60 # (Alternatively you could also use the landau shape pdfDtErr)
61 # ROOT.RooProdPdf model("model", "model",pdfDtErr,
62 # ROOT.RooFit.Conditional(decay_gm,dt))
63 
64 # Sample, fit and plot product model
65 # ------------------------------------------------------------------
66 
67 # Specify external dataset with dterr values to use model_dm as
68 # conditional p.d.f.
69 data = model.generate(ROOT.RooArgSet(dt, dterr), 10000)
70 
71 # Fit conditional decay_dm(dt|dterr)
72 # ---------------------------------------------------------------------
73 
74 # Specify dterr as conditional observable
75 model.fitTo(data)
76 
77 # Plot conditional decay_dm(dt|dterr)
78 # ---------------------------------------------------------------------
79 
80 # Make two-dimensional plot of conditional p.d.f in (dt,dterr)
81 hh_model = model.createHistogram("hh_model", dt, ROOT.RooFit.Binning(
82  50), ROOT.RooFit.YVar(dterr, ROOT.RooFit.Binning(50)))
83 hh_model.SetLineColor(ROOT.kBlue)
84 
85 # Make projection of data an dt
86 frame = dt.frame(ROOT.RooFit.Title("Projection of model(dt|dterr) on dt"))
87 data.plotOn(frame)
88 model.plotOn(frame)
89 
90 # Draw all frames on canvas
91 c = ROOT.TCanvas("rf307_fullpereventerrors",
92  "rf307_fullpereventerrors", 800, 400)
93 c.Divide(2)
94 c.cd(1)
95 ROOT.gPad.SetLeftMargin(0.20)
96 hh_model.GetZaxis().SetTitleOffset(2.5)
97 hh_model.Draw("surf")
98 c.cd(2)
99 ROOT.gPad.SetLeftMargin(0.15)
100 frame.GetYaxis().SetTitleOffset(1.6)
101 frame.Draw()
102 
103 c.SaveAs("rf307_fullpereventerrors.png")