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