Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rf708_bphysics.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_roofit
3## \notebook
4## Special pdf's: special decay pdf for B physics with mixing and/or CP violation
5##
6## \macro_code
7##
8## \date February 2018
9## \authors Clemens Lange, Wouter Verkerke (C++ version)
10
11import ROOT
12
13
14# B-decay with mixing
15# -------------------------
16
17# Construct pdf
18# -------------------------
19
20# Observable
21dt = ROOT.RooRealVar("dt", "dt", -10, 10)
22dt.setBins(40)
23
24# Parameters
25dm = ROOT.RooRealVar("dm", "delta m(B0)", 0.472)
26tau = ROOT.RooRealVar("tau", "tau (B0)", 1.547)
27w = ROOT.RooRealVar("w", "flavour mistag rate", 0.1)
28dw = ROOT.RooRealVar("dw", "delta mistag rate for B0/B0bar", 0.1)
29
30mixState = ROOT.RooCategory("mixState", "B0/B0bar mixing state")
31mixState.defineType("mixed", -1)
32mixState.defineType("unmixed", 1)
33
34tagFlav = ROOT.RooCategory("tagFlav", "Flavour of the tagged B0")
35tagFlav.defineType("B0", 1)
36tagFlav.defineType("B0bar", -1)
37
38# Use delta function resolution model
39tm = ROOT.RooTruthModel("tm", "truth model", dt)
40
41# Construct Bdecay with mixing
42bmix = ROOT.RooBMixDecay(
43 "bmix",
44 "decay",
45 dt,
46 mixState,
47 tagFlav,
48 tau,
49 dm,
50 w,
51 dw,
52 tm,
53 ROOT.RooBMixDecay.DoubleSided)
54
55# Plot pdf in various slices
56# ---------------------------------------------------
57
58# Generate some data
59data = bmix.generate(ROOT.RooArgSet(dt, mixState, tagFlav), 10000)
60
61# Plot B0 and B0bar tagged data separately
62# For all plots below B0 and B0 tagged data will look somewhat differently
63# if the flavor tagging mistag rate for B0 and B0 is different (i.e. dw!=0)
64frame1 = dt.frame(ROOT.RooFit.Title(
65 "B decay distribution with mixing (B0/B0bar)"))
66
67data.plotOn(frame1, ROOT.RooFit.Cut("tagFlav==tagFlav::B0"))
68bmix.plotOn(frame1, ROOT.RooFit.Slice(tagFlav, "B0"))
69
70data.plotOn(frame1, ROOT.RooFit.Cut("tagFlav==tagFlav::B0bar"),
71 ROOT.RooFit.MarkerColor(ROOT.kCyan))
72bmix.plotOn(frame1, ROOT.RooFit.Slice(tagFlav, "B0bar"),
73 ROOT.RooFit.LineColor(ROOT.kCyan))
74
75# Plot mixed slice for B0 and B0bar tagged data separately
76frame2 = dt.frame(ROOT.RooFit.Title(
77 "B decay distribution of mixed events (B0/B0bar)"))
78
79data.plotOn(frame2, ROOT.RooFit.Cut(
80 "mixState==mixState::mixed&&tagFlav==tagFlav::B0"))
81bmix.plotOn(frame2, ROOT.RooFit.Slice(tagFlav, "B0"),
82 ROOT.RooFit.Slice(mixState, "mixed"))
83
84data.plotOn(
85 frame2,
86 ROOT.RooFit.Cut("mixState==mixState::mixed&&tagFlav==tagFlav::B0bar"),
87 ROOT.RooFit.MarkerColor(
88 ROOT.kCyan))
89bmix.plotOn(frame2, ROOT.RooFit.Slice(tagFlav, "B0bar"), ROOT.RooFit.Slice(
90 mixState, "mixed"), ROOT.RooFit.LineColor(ROOT.kCyan))
91
92# Plot unmixed slice for B0 and B0bar tagged data separately
93frame3 = dt.frame(ROOT.RooFit.Title(
94 "B decay distribution of unmixed events (B0/B0bar)"))
95
96data.plotOn(frame3, ROOT.RooFit.Cut(
97 "mixState==mixState::unmixed&&tagFlav==tagFlav::B0"))
98bmix.plotOn(frame3, ROOT.RooFit.Slice(tagFlav, "B0"),
99 ROOT.RooFit.Slice(mixState, "unmixed"))
100
101data.plotOn(
102 frame3,
103 ROOT.RooFit.Cut("mixState==mixState::unmixed&&tagFlav==tagFlav::B0bar"),
104 ROOT.RooFit.MarkerColor(
105 ROOT.kCyan))
106bmix.plotOn(frame3, ROOT.RooFit.Slice(tagFlav, "B0bar"), ROOT.RooFit.Slice(
107 mixState, "unmixed"), ROOT.RooFit.LineColor(ROOT.kCyan))
108
109# B-decay with CP violation
110# -------------------------
111
112# Construct pdf
113# -------------------------
114
115# Additional parameters needed for B decay with CPV
116CPeigen = ROOT.RooRealVar("CPeigen", "CP eigen value", -1)
117absLambda = ROOT.RooRealVar("absLambda", "|lambda|", 1, 0, 2)
118argLambda = ROOT.RooRealVar("absLambda", "|lambda|", 0.7, -1, 1)
119effR = ROOT.RooRealVar("effR", "B0/B0bar reco efficiency ratio", 1)
120
121# Construct Bdecay with CP violation
122bcp = ROOT.RooBCPEffDecay(
123 "bcp",
124 "bcp",
125 dt,
126 tagFlav,
127 tau,
128 dm,
129 w,
130 CPeigen,
131 absLambda,
132 argLambda,
133 effR,
134 dw,
135 tm,
136 ROOT.RooBCPEffDecay.DoubleSided)
137
138# Plot scenario 1 - sin(2b)=0.7, |l|=1
139# ---------------------------------------------------------------------------
140
141# Generate some data
142data2 = bcp.generate(ROOT.RooArgSet(dt, tagFlav), 10000)
143
144# Plot B0 and B0bar tagged data separately
145frame4 = dt.frame(ROOT.RooFit.Title(
146 "B decay distribution with CPV(|l|=1,Im(l)=0.7) (B0/B0bar)"))
147
148data2.plotOn(frame4, ROOT.RooFit.Cut("tagFlav==tagFlav::B0"))
149bcp.plotOn(frame4, ROOT.RooFit.Slice(tagFlav, "B0"))
150
151data2.plotOn(frame4, ROOT.RooFit.Cut("tagFlav==tagFlav::B0bar"),
152 ROOT.RooFit.MarkerColor(ROOT.kCyan))
153bcp.plotOn(frame4, ROOT.RooFit.Slice(tagFlav, "B0bar"),
154 ROOT.RooFit.LineColor(ROOT.kCyan))
155
156# # Plot scenario 2 - sin(2b)=0.7, |l|=0.7
157# -------------------------------------------------------------------------------
158
159absLambda.setVal(0.7)
160
161# Generate some data
162data3 = bcp.generate(ROOT.RooArgSet(dt, tagFlav), 10000)
163
164# Plot B0 and B0bar tagged data separately (sin2b = 0.7 plus direct CPV
165# |l|=0.5)
166frame5 = dt.frame(ROOT.RooFit.Title(
167 "B decay distribution with CPV(|l|=0.7,Im(l)=0.7) (B0/B0bar)"))
168
169data3.plotOn(frame5, ROOT.RooFit.Cut("tagFlav==tagFlav::B0"))
170bcp.plotOn(frame5, ROOT.RooFit.Slice(tagFlav, "B0"))
171
172data3.plotOn(frame5, ROOT.RooFit.Cut("tagFlav==tagFlav::B0bar"),
173 ROOT.RooFit.MarkerColor(ROOT.kCyan))
174bcp.plotOn(frame5, ROOT.RooFit.Slice(tagFlav, "B0bar"),
175 ROOT.RooFit.LineColor(ROOT.kCyan))
176
177
178# Generic B-decay with user coefficients
179# -------------------------
180
181# Construct pdf
182# -------------------------
183
184# Model parameters
185DGbG = ROOT.RooRealVar("DGbG", "DGamma/GammaAvg", 0.5, -1, 1)
186Adir = ROOT.RooRealVar("Adir", "-[1-abs(l)**2]/[1+abs(l)**2]", 0)
187Amix = ROOT.RooRealVar("Amix", "2Im(l)/[1+abs(l)**2]", 0.7)
188Adel = ROOT.RooRealVar("Adel", "2Re(l)/[1+abs(l)**2]", 0.7)
189
190# Derived input parameters for pdf
191DG = ROOT.RooFormulaVar("DG", "Delta Gamma", "@1/@0",
192 ROOT.RooArgList(tau, DGbG))
193
194# Construct coefficient functions for sin,cos, modulations of decay
195# distribution
196fsin = ROOT.RooFormulaVar(
197 "fsin", "fsin", "@0*@1*(1-2*@2)", ROOT.RooArgList(Amix, tagFlav, w))
198fcos = ROOT.RooFormulaVar(
199 "fcos", "fcos", "@0*@1*(1-2*@2)", ROOT.RooArgList(Adir, tagFlav, w))
200fsinh = ROOT.RooFormulaVar("fsinh", "fsinh", "@0", ROOT.RooArgList(Adel))
201
202# Construct generic B decay pdf using above user coefficients
203bcpg = ROOT.RooBDecay("bcpg", "bcpg", dt, tau, DG, ROOT.RooFit.RooConst(
204 1), fsinh, fcos, fsin, dm, tm, ROOT.RooBDecay.DoubleSided)
205
206# Plot - Im(l)=0.7, e(l)=0.7 |l|=1, G/G=0.5
207# -------------------------------------------------------------------------------------
208
209# Generate some data
210data4 = bcpg.generate(ROOT.RooArgSet(dt, tagFlav), 10000)
211
212# Plot B0 and B0bar tagged data separately
213frame6 = dt.frame(ROOT.RooFit.Title(
214 "B decay distribution with CPV(Im(l)=0.7,Re(l)=0.7,|l|=1,dG/G=0.5) (B0/B0bar)"))
215
216data4.plotOn(frame6, ROOT.RooFit.Cut("tagFlav==tagFlav::B0"))
217bcpg.plotOn(frame6, ROOT.RooFit.Slice(tagFlav, "B0"))
218
219data4.plotOn(frame6, ROOT.RooFit.Cut("tagFlav==tagFlav::B0bar"),
220 ROOT.RooFit.MarkerColor(ROOT.kCyan))
221bcpg.plotOn(frame6, ROOT.RooFit.Slice(tagFlav, "B0bar"),
222 ROOT.RooFit.LineColor(ROOT.kCyan))
223
224c = ROOT.TCanvas("rf708_bphysics", "rf708_bphysics", 1200, 800)
225c.Divide(3, 2)
226c.cd(1)
227ROOT.gPad.SetLeftMargin(0.15)
228frame1.GetYaxis().SetTitleOffset(1.6)
229frame1.Draw()
230c.cd(2)
231ROOT.gPad.SetLeftMargin(0.15)
232frame2.GetYaxis().SetTitleOffset(1.6)
233frame2.Draw()
234c.cd(3)
235ROOT.gPad.SetLeftMargin(0.15)
236frame3.GetYaxis().SetTitleOffset(1.6)
237frame3.Draw()
238c.cd(4)
239ROOT.gPad.SetLeftMargin(0.15)
240frame4.GetYaxis().SetTitleOffset(1.6)
241frame4.Draw()
242c.cd(5)
243ROOT.gPad.SetLeftMargin(0.15)
244frame5.GetYaxis().SetTitleOffset(1.6)
245frame5.Draw()
246c.cd(6)
247ROOT.gPad.SetLeftMargin(0.15)
248frame6.GetYaxis().SetTitleOffset(1.6)
249frame6.Draw()
250
251c.SaveAs("rf708_bphysics.png")