Logo ROOT  
Reference Guide
rf310_sliceplot.py
Go to the documentation of this file.
1 ## \file
2 ## \ingroup tutorial_roofit
3 ## \notebook
4 ## Multidimensional models: projecting pdf and data slices in discrete observables
5 ##
6 ## \macro_code
7 ##
8 ## \date February 2018
9 ## \authors Clemens Lange, Wouter Verkerke (C++ version)
10 
11 import ROOT
12 
13 # Create B decay pdf with mixing
14 # ----------------------------------------------------------
15 
16 # Decay time observables
17 dt = ROOT.RooRealVar("dt", "dt", -20, 20)
18 
19 # Discrete observables mixState (B0tag==B0reco?) and tagFlav
20 # (B0tag==B0(bar)?)
21 mixState = ROOT.RooCategory("mixState", "B0/B0bar mixing state")
22 tagFlav = ROOT.RooCategory("tagFlav", "Flavour of the tagged B0")
23 
24 # Define state labels of discrete observables
25 mixState.defineType("mixed", -1)
26 mixState.defineType("unmixed", 1)
27 tagFlav.defineType("B0", 1)
28 tagFlav.defineType("B0bar", -1)
29 
30 # Model parameters
31 dm = ROOT.RooRealVar("dm", "delta m(B)", 0.472, 0., 1.0)
32 tau = ROOT.RooRealVar("tau", "B0 decay time", 1.547, 1.0, 2.0)
33 w = ROOT.RooRealVar("w", "Flavor Mistag rate", 0.03, 0.0, 1.0)
34 dw = ROOT.RooRealVar(
35  "dw", "Flavor Mistag rate difference between B0 and B0bar", 0.01)
36 
37 # Build a gaussian resolution model
38 bias1 = ROOT.RooRealVar("bias1", "bias1", 0)
39 sigma1 = ROOT.RooRealVar("sigma1", "sigma1", 0.01)
40 gm1 = ROOT.RooGaussModel("gm1", "gauss model 1", dt, bias1, sigma1)
41 
42 # Construct a decay pdf, with single gaussian resolution model
43 bmix_gm1 = ROOT.RooBMixDecay(
44  "bmix",
45  "decay",
46  dt,
47  mixState,
48  tagFlav,
49  tau,
50  dm,
51  w,
52  dw,
53  gm1,
54  ROOT.RooBMixDecay.DoubleSided)
55 
56 # Generate BMixing data with above set of event errors
57 data = bmix_gm1.generate(ROOT.RooArgSet(dt, tagFlav, mixState), 20000)
58 
59 # Plot full decay distribution
60 # ----------------------------------------------------------
61 
62 # Create frame, data and pdf projection (integrated over tagFlav and
63 # mixState)
64 frame = dt.frame(ROOT.RooFit.Title("Inclusive decay distribution"))
65 data.plotOn(frame)
66 bmix_gm1.plotOn(frame)
67 
68 # Plot decay distribution for mixed and unmixed slice of mixState
69 # -------------------------------------------------------------------------------------------
70 
71 # Create frame, data (mixed only)
72 frame2 = dt.frame(ROOT.RooFit.Title("Decay distribution of mixed events"))
73 data.plotOn(frame2, ROOT.RooFit.Cut("mixState==mixState::mixed"))
74 
75 # Position slice in mixState at "mixed" and plot slice of pdf in mixstate
76 # over data (integrated over tagFlav)
77 bmix_gm1.plotOn(frame2, ROOT.RooFit.Slice(mixState, "mixed"))
78 
79 # Create frame, data (unmixed only)
80 frame3 = dt.frame(ROOT.RooFit.Title(
81  "Decay distribution of unmixed events"))
82 data.plotOn(frame3, ROOT.RooFit.Cut("mixState==mixState::unmixed"))
83 
84 # Position slice in mixState at "unmixed" and plot slice of pdf in
85 # mixstate over data (integrated over tagFlav)
86 bmix_gm1.plotOn(frame3, ROOT.RooFit.Slice(mixState, "unmixed"))
87 
88 c = ROOT.TCanvas("rf310_sliceplot", "rf310_sliceplot", 1200, 400)
89 c.Divide(3)
90 c.cd(1)
91 ROOT.gPad.SetLeftMargin(0.15)
92 frame.GetYaxis().SetTitleOffset(1.4)
93 ROOT.gPad.SetLogy()
94 frame.Draw()
95 c.cd(2)
96 ROOT.gPad.SetLeftMargin(0.15)
97 frame2.GetYaxis().SetTitleOffset(1.4)
98 ROOT.gPad.SetLogy()
99 frame2.Draw()
100 c.cd(3)
101 ROOT.gPad.SetLeftMargin(0.15)
102 frame3.GetYaxis().SetTitleOffset(1.4)
103 ROOT.gPad.SetLogy()
104 frame3.Draw()
105 
106 c.SaveAs("rf310_sliceplot.png")