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