Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
rf209_anaconv.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ## Addition and convolution: decay function p.d.fs with optional B physics effects (mixing and CP violation) that can be analytically convolved with e.g. Gaussian resolution functions
5 ##
6 ## pdf1 = decay(t,tau) (x) delta(t)
7 ## pdf2 = decay(t,tau) (x) gauss(t,m,s)
8 ## pdf3 = decay(t,tau) (x) (f*gauss1(t,m1,s1) + (1-f)*gauss2(t,m1,s1))
9 ##
10 ## \macro_code
11 ##
12 ## \date February 2018
13 ## \author Clemens Lange, Wouter Verkerke (C++ version)
14 
15 import ROOT
16 
17 # B-physics pdf with truth resolution
18 # ---------------------------------------------------------------------
19 
20 # Variables of decay p.d.f.
21 dt = ROOT.RooRealVar("dt", "dt", -10, 10)
22 tau = ROOT.RooRealVar("tau", "tau", 1.548)
23 
24 # Build a truth resolution model (delta function)
25 tm = ROOT.RooTruthModel("tm", "truth model", dt)
26 
27 # Construct decay(t) (x) delta(t)
28 decay_tm = ROOT.RooDecay("decay_tm", "decay", dt,
29  tau, tm, ROOT.RooDecay.DoubleSided)
30 
31 # Plot p.d.f. (dashed)
32 frame = dt.frame(ROOT.RooFit.Title("Bdecay (x) resolution"))
33 decay_tm.plotOn(frame, ROOT.RooFit.LineStyle(ROOT.kDashed))
34 
35 # B-physics pdf with Gaussian resolution
36 # ----------------------------------------------------------------------------
37 
38 # Build a gaussian resolution model
39 bias1 = ROOT.RooRealVar("bias1", "bias1", 0)
40 sigma1 = ROOT.RooRealVar("sigma1", "sigma1", 1)
41 gm1 = ROOT.RooGaussModel("gm1", "gauss model 1", dt, bias1, sigma1)
42 
43 # Construct decay(t) (x) gauss1(t)
44 decay_gm1 = ROOT.RooDecay("decay_gm1", "decay",
45  dt, tau, gm1, ROOT.RooDecay.DoubleSided)
46 
47 # Plot p.d.f.
48 decay_gm1.plotOn(frame)
49 
50 # B-physics pdf with double Gaussian resolution
51 # ------------------------------------------------------------------------------------------
52 
53 # Build another gaussian resolution model
54 bias2 = ROOT.RooRealVar("bias2", "bias2", 0)
55 sigma2 = ROOT.RooRealVar("sigma2", "sigma2", 5)
56 gm2 = ROOT.RooGaussModel("gm2", "gauss model 2", dt, bias2, sigma2)
57 
58 # Build a composite resolution model f*gm1+(1-f)*gm2
59 gm1frac = ROOT.RooRealVar("gm1frac", "fraction of gm1", 0.5)
60 gmsum = ROOT.RooAddModel(
61  "gmsum",
62  "sum of gm1 and gm2",
63  ROOT.RooArgList(
64  gm1,
65  gm2),
66  ROOT.RooArgList(gm1frac))
67 
68 # Construct decay(t) (x) (f*gm1 + (1-f)*gm2)
69 decay_gmsum = ROOT.RooDecay(
70  "decay_gmsum", "decay", dt, tau, gmsum, ROOT.RooDecay.DoubleSided)
71 
72 # Plot p.d.f. (red)
73 decay_gmsum.plotOn(frame, ROOT.RooFit.LineColor(ROOT.kRed))
74 
75 # Draw all frames on canvas
76 c = ROOT.TCanvas("rf209_anaconv", "rf209_anaconv", 600, 600)
77 ROOT.gPad.SetLeftMargin(0.15)
78 frame.GetYaxis().SetTitleOffset(1.6)
79 frame.Draw()
80 
81 c.SaveAs("rf209_anaconv.png")