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_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
16# B-decay with mixing
17# -------------------------
18
19# Construct pdf
20# -------------------------
21
22# Observable
23dt = ROOT.RooRealVar("dt", "dt", -10, 10)
24dt.setBins(40)
25
26# Parameters
27dm = ROOT.RooRealVar("dm", "delta m(B0)", 0.472)
28tau = ROOT.RooRealVar("tau", "tau (B0)", 1.547)
29w = ROOT.RooRealVar("w", "flavour mistag rate", 0.1)
30dw = ROOT.RooRealVar("dw", "delta mistag rate for B0/B0bar", 0.1)
31
32mixState = ROOT.RooCategory("mixState", "B0/B0bar mixing state", {"mixed": -1, "unmixed": 1})
33
34tagFlav = ROOT.RooCategory("tagFlav", "Flavour of the tagged B0", {"B0": 1, "B0bar": -1})
35
36# Use delta function resolution model
37tm = ROOT.RooTruthModel("tm", "truth model", dt)
38
39# Construct Bdecay with mixing
40bmix = ROOT.RooBMixDecay("bmix", "decay", dt, mixState, tagFlav, tau, dm, w, dw, tm, type="DoubleSided")
41
42# Plot pdf in various slices
43# ---------------------------------------------------
44
45# Generate some data
46data = bmix.generate({dt, mixState, tagFlav}, 10000)
47
48# Plot B0 and B0bar tagged data separately
49# For all plots below B0 and B0 tagged data will look somewhat differently
50# if the flavor tagging mistag rate for B0 and B0 is different (i.e. dw!=0)
51frame1 = dt.frame(Title="B decay distribution with mixing (B0/B0bar)")
52
53data.plotOn(frame1, Cut="tagFlav==tagFlav::B0")
54bmix.plotOn(frame1, Slice=(tagFlav, "B0"))
55
56data.plotOn(frame1, Cut="tagFlav==tagFlav::B0bar", MarkerColor="c")
57bmix.plotOn(frame1, Slice=(tagFlav, "B0bar"), LineColor="c")
58
59# Plot mixed slice for B0 and B0bar tagged data separately
60frame2 = dt.frame(Title="B decay distribution of mixed events (B0/B0bar)")
61
62data.plotOn(frame2, Cut="mixState==mixState::mixed&&tagFlav==tagFlav::B0")
63bmix.plotOn(frame2, Slice={tagFlav: "B0", mixState: "mixed"})
64
65data.plotOn(frame2, Cut="mixState==mixState::mixed&&tagFlav==tagFlav::B0bar", MarkerColor="c")
66bmix.plotOn(frame2, Slice={tagFlav: "B0bar", mixState: "mixed"}, LineColor="c")
67
68# Plot unmixed slice for B0 and B0bar tagged data separately
69frame3 = dt.frame(Title="B decay distribution of unmixed events (B0/B0bar)")
70
71data.plotOn(frame3, Cut="mixState==mixState::unmixed&&tagFlav==tagFlav::B0")
72bmix.plotOn(frame3, Slice={tagFlav: "B0", mixState: "unmixed"})
73
74data.plotOn(frame3, Cut="mixState==mixState::unmixed&&tagFlav==tagFlav::B0bar", MarkerColor="c")
75bmix.plotOn(frame3, Slice={tagFlav: "B0bar", mixState: "unmixed"}, LineColor="c")
76
77# B-decay with CP violation
78# -------------------------
79
80# Construct pdf
81# -------------------------
82
83# Additional parameters needed for B decay with CPV
84CPeigen = ROOT.RooRealVar("CPeigen", "CP eigen value", -1)
85absLambda = ROOT.RooRealVar("absLambda", "|lambda|", 1, 0, 2)
86argLambda = ROOT.RooRealVar("absLambda", "|lambda|", 0.7, -1, 1)
87effR = ROOT.RooRealVar("effR", "B0/B0bar reco efficiency ratio", 1)
88
89# Construct Bdecay with CP violation
90bcp = ROOT.RooBCPEffDecay(
91 "bcp", "bcp", dt, tagFlav, tau, dm, w, CPeigen, absLambda, argLambda, effR, dw, tm, type="DoubleSided"
92)
93
94# Plot scenario 1 - sin(2b)=0.7, |l|=1
95# ---------------------------------------------------------------------------
96
97# Generate some data
98data2 = bcp.generate({dt, tagFlav}, 10000)
99
100# Plot B0 and B0bar tagged data separately
101frame4 = dt.frame(Title="B decay distribution with CPV(|l|=1,Im(l)=0.7) (B0/B0bar)")
102
103data2.plotOn(frame4, Cut="tagFlav==tagFlav::B0")
104bcp.plotOn(frame4, Slice=(tagFlav, "B0"))
105
106data2.plotOn(frame4, Cut="tagFlav==tagFlav::B0bar", MarkerColor="c")
107bcp.plotOn(frame4, Slice=(tagFlav, "B0bar"), LineColor="c")
108
109# # Plot scenario 2 - sin(2b)=0.7, |l|=0.7
110# -------------------------------------------------------------------------------
111
112absLambda.setVal(0.7)
113
114# Generate some data
115data3 = bcp.generate({dt, tagFlav}, 10000)
116
117# Plot B0 and B0bar tagged data separately (sin2b = 0.7 plus direct CPV
118# |l|=0.5)
119frame5 = dt.frame(Title="B decay distribution with CPV(|l|=0.7,Im(l)=0.7) (B0/B0bar)")
120
121data3.plotOn(frame5, Cut="tagFlav==tagFlav::B0")
122bcp.plotOn(frame5, Slice=(tagFlav, "B0"))
123
124data3.plotOn(frame5, Cut="tagFlav==tagFlav::B0bar", MarkerColor="c")
125bcp.plotOn(frame5, Slice=(tagFlav, "B0bar"), LineColor="c")
126
127
128# Generic B-decay with user coefficients
129# -------------------------
130
131# Construct pdf
132# -------------------------
133
134# Model parameters
135DGbG = ROOT.RooRealVar("DGbG", "DGamma/GammaAvg", 0.5, -1, 1)
136Adir = ROOT.RooRealVar("Adir", "-[1-abs(l)**2]/[1+abs(l)**2]", 0)
137Amix = ROOT.RooRealVar("Amix", "2Im(l)/[1+abs(l)**2]", 0.7)
138Adel = ROOT.RooRealVar("Adel", "2Re(l)/[1+abs(l)**2]", 0.7)
139
140# Derived input parameters for pdf
141DG = ROOT.RooFormulaVar("DG", "Delta Gamma", "@1/@0", [tau, DGbG])
142
143# Construct coefficient functions for sin,cos, modulations of decay
144# distribution
145fsin = ROOT.RooFormulaVar("fsin", "fsin", "@0*@1*(1-2*@2)", [Amix, tagFlav, w])
146fcos = ROOT.RooFormulaVar("fcos", "fcos", "@0*@1*(1-2*@2)", [Adir, tagFlav, w])
147fsinh = ROOT.RooFormulaVar("fsinh", "fsinh", "@0", [Adel])
148
149# Construct generic B decay pdf using above user coefficients
150bcpg = ROOT.RooBDecay(
151 "bcpg", "bcpg", dt, tau, DG, ROOT.RooFit.RooConst(1), fsinh, fcos, fsin, dm, tm, type="DoubleSided"
152)
153
154# Plot - Im(l)=0.7, e(l)=0.7 |l|=1, G/G=0.5
155# -------------------------------------------------------------------------------------
156
157# Generate some data
158data4 = bcpg.generate({dt, tagFlav}, 10000)
159
160# Plot B0 and B0bar tagged data separately
161frame6 = dt.frame(Title="B decay distribution with CPV(Im(l)=0.7,Re(l)=0.7,|l|=1,dG/G=0.5) (B0/B0bar)")
162
163data4.plotOn(frame6, Cut="tagFlav==tagFlav::B0")
164bcpg.plotOn(frame6, Slice=(tagFlav, "B0"))
165
166data4.plotOn(frame6, Cut="tagFlav==tagFlav::B0bar", MarkerColor="c")
167bcpg.plotOn(frame6, Slice=(tagFlav, "B0bar"), LineColor="c")
168
169c = ROOT.TCanvas("rf708_bphysics", "rf708_bphysics", 1200, 800)
170c.Divide(3, 2)
171c.cd(1)
172ROOT.gPad.SetLeftMargin(0.15)
173frame1.GetYaxis().SetTitleOffset(1.6)
174frame1.Draw()
175c.cd(2)
176ROOT.gPad.SetLeftMargin(0.15)
177frame2.GetYaxis().SetTitleOffset(1.6)
178frame2.Draw()
179c.cd(3)
180ROOT.gPad.SetLeftMargin(0.15)
181frame3.GetYaxis().SetTitleOffset(1.6)
182frame3.Draw()
183c.cd(4)
184ROOT.gPad.SetLeftMargin(0.15)
185frame4.GetYaxis().SetTitleOffset(1.6)
186frame4.Draw()
187c.cd(5)
188ROOT.gPad.SetLeftMargin(0.15)
189frame5.GetYaxis().SetTitleOffset(1.6)
190frame5.Draw()
191c.cd(6)
192ROOT.gPad.SetLeftMargin(0.15)
193frame6.GetYaxis().SetTitleOffset(1.6)
194frame6.Draw()
195
196c.SaveAs("rf708_bphysics.png")