Belle II Software  release-06-00-14
deltaTVertexTagVResolution.py
1 #!/usr/bin/env python
2 # -*- coding: utf-8 -*-
3 
4 
11 
12 
26 
27 
28 from basf2 import B2FATAL
29 import ROOT
30 import sys
31 import glob
32 import math
33 from operator import itemgetter
34 
35 if len(sys.argv) != 4:
36  sys.exit(
37  "Must provide 3 arguments: [input_sim_file] or ['input_sim_file*'] wildcards, [treeName] " +
38  "and [VXD_requirement (PXD or SVD)]")
39 
40 PATH = "."
41 
42 workingFile = sys.argv[1]
43 workingFiles = glob.glob(str(workingFile))
44 
45 limDeltaT = 5
46 limDeltaTErr = 3.0
47 limZSig = 0.03
48 limZTag = 0.03
49 
50 treeName = str(sys.argv[2])
51 
52 # No PXD hit equired: PXD0. At least one PXD (SVD) hit for one of the muon tracks: PXD1 (SVD1).
53 # Hit required for both muon tracks: PXD2 (SVD2)"
54 VXDReqs = []
55 if str(sys.argv[3]) == "PXD":
56  VXDReqs = ["PXD0", "PXD2"]
57 elif str(sys.argv[3]) == "SVD":
58  VXDReqs = ["SVD0", "SVD2"]
59 else:
60  B2FATAL('Not available VXD requirement " + str(sys.argv[3]) + ". Available are "PXD" and "SVD".')
61 
62 ROOT.RooMsgService.instance().setGlobalKillBelow(ROOT.RooFit.WARNING)
63 
64 fitResults = []
65 numberOfEntries = []
66 
67 for VXDReq in VXDReqs:
68  iResult = []
69 
70  tdat = ROOT.TChain(treeName)
71 
72  for iFile in workingFiles:
73  tdat.AddFile(iFile)
74 
75  B0_DeltaT = ROOT.RooRealVar("DeltaT", "DeltaT", 0.)
76  deltaTErr = ROOT.RooRealVar("DeltaTErr", "DeltaTErr", 0, limDeltaTErr, "ps")
77  B0_TruthDeltaT = ROOT.RooRealVar("MCDeltaT", "MCDeltaT", 0.)
78 
79  B0_qrMC = ROOT.RooRealVar("qrMC", "qrMC", 0., -100, 100)
80  B0_isSignal = ROOT.RooRealVar("isSignal", "isSignal", 0., 0., 512.)
81 
82  B0_TagVz = ROOT.RooRealVar("TagVz", "TagVz", 0., -100, 100, "cm")
83  B0_TruthTagVz = ROOT.RooRealVar("mcTagVz", "mcTagVz", 0., -100, 100, "cm")
84 
85  B0_Z = ROOT.RooRealVar("z", "z", 0., -100, 100, "cm")
86  B0_TruthZ = ROOT.RooRealVar("mcZ", "mcZ", 0., -100, 100, "cm")
87 
88  B0_Jpsi_mu0_nPXDHits = ROOT.RooRealVar("Jpsi_mu_0_nPXDHits", "Jpsi_mu_0_nPXDHits", 0., -10., 100.)
89  B0_Jpsi_mu1_nPXDHits = ROOT.RooRealVar("Jpsi_mu_1_nPXDHits", "Jpsi_mu_1_nPXDHits", 0., -10., 100.)
90  B0_Jpsi_mu0_nSVDHits = ROOT.RooRealVar("Jpsi_mu_0_nSVDHits", "Jpsi_mu_0_nSVDHits", 0., -10., 100.)
91  B0_Jpsi_mu1_nSVDHits = ROOT.RooRealVar("Jpsi_mu_1_nSVDHits", "Jpsi_mu_1_nSVDHits", 0., -10., 100.)
92 
93  # B0_Jpsi_mu0_firstSVDLayer = ROOT.RooRealVar("B0_Jpsi_mu0__firstSVDLayer", "B0_pi0_e1__firstSVDLayer", 1., -10., 100.)
94  # B0_Jpsi_mu1_firstSVDLayer = ROOT.RooRealVar("B0_Jpsi_mu1__firstSVDLayer", "B0_pi0_e0__firstSVDLayer", 1., -10., 100.)
95 
96  DT = ROOT.RooRealVar("DT", "DT", 0., -limDeltaT, limDeltaT, "ps")
97  DSigZ = ROOT.RooRealVar("DSigZ", "DSigZ", 0., -limZSig, limZSig, "cm")
98  # DSigZ = ROOT.RooFormulaVar("DSigZ", "DSigZ", "@@0-@@1", ROOT.RooArgList(B0_Z, B0_TruthZ))
99  DTagZ = ROOT.RooRealVar("DSigZ", "DSigZ", 0., -limZTag, limZTag, "cm")
100  # DTagZ = ROOT.RooFormulaVar("DTagZ", "DTagZ", "@@0-@@1", ROOT.RooArgList(B0_TagVz, B0_TruthTagVz))
101 
102  histo_DeltaT = ROOT.TH1F('B0_DeltaT_' + VXDReq, 'Residual of DeltaT',
103  100, -limDeltaT, limDeltaT)
104  histo_DeltaTErr = ROOT.TH1F('B0_DeltaTErr_' + VXDReq, 'Residual of DeltaZsig',
105  100, 0, 30)
106  histo_DeltaZSig = ROOT.TH1F('B0_DeltaZsig_' + VXDReq, 'Residual of DeltaZsig',
107  100, -limZSig, limZSig)
108  histo_DeltaZTag = ROOT.TH1F('B0_DeltaZtag_' + VXDReq, 'Residual of DeltaZsig',
109  100, -limZTag, limZTag)
110 
111  cut = "isSignal == 1 " # + "&& abs(B0_DeltaTErr)< " + str(limDeltaTErr) + " "
112 
113  if VXDReq == 'PXD1':
114  cut = cut + "&& (Jpsi_mu_0_nPXDHits> 0 || Jpsi_mu_1_nPXDHits> 0) "
115  if VXDReq == 'PXD2':
116  cut = cut + "&& Jpsi_mu_0_nPXDHits> 0 && Jpsi_mu_1_nPXDHits> 0 "
117 
118  if VXDReq == 'SVD1':
119  cut = cut + "&& (Jpsi_mu_0_SVDHits> 0 || Jpsi_mu_1_nSVDHits> 0) "
120  if VXDReq == 'SVD2':
121  cut = cut + "&& Jpsi_mu_0_nSVDHits> 0 && Jpsi_mu_1_nSVDHits> 0 "
122 
123  argSet = ROOT.RooArgSet(
124  B0_DeltaT,
125  B0_TruthDeltaT,
126  B0_TagVz,
127  B0_TruthTagVz,
128  B0_Z,
129  B0_TruthZ)
130 
131  argSet.add(B0_Jpsi_mu0_nPXDHits)
132  argSet.add(B0_Jpsi_mu1_nPXDHits)
133  argSet.add(B0_Jpsi_mu0_nSVDHits)
134  argSet.add(B0_Jpsi_mu1_nSVDHits)
135 
136  argSet.add(B0_isSignal)
137 
138  data = ROOT.RooDataSet(
139  "data",
140  "data",
141  tdat,
142  argSet,
143  cut)
144 
145  if VXDReq == 'PXD1' or VXDReq == 'PXD2':
146  fitDataDTErr = ROOT.RooDataSet("data", "data", tdat, ROOT.RooArgSet(
147  B0_isSignal, B0_Jpsi_mu0_nPXDHits, B0_Jpsi_mu1_nPXDHits, deltaTErr), cut)
148  elif VXDReq == 'SVD1' or VXDReq == 'SVD2':
149  fitDataDTErr = ROOT.RooDataSet("data", "data", tdat, ROOT.RooArgSet(
150  B0_isSignal, B0_Jpsi_mu0_nSVDHits, B0_Jpsi_mu1_nSVDHits, deltaTErr), cut)
151  else:
152  fitDataDTErr = ROOT.RooDataSet("data", "data", tdat, ROOT.RooArgSet(B0_isSignal, deltaTErr), cut)
153 
154  # fitData.append(data)
155 
156  fitDataDT = ROOT.RooDataSet("fitDataDT", "fitDataDT", ROOT.RooArgSet(DT))
157  fitDataSigZ = ROOT.RooDataSet("fitDataSigZ", "fitDataSigZ", ROOT.RooArgSet(DSigZ))
158  fitDataTagZ = ROOT.RooDataSet("fitDataTagZ", "fitDataTagZ", ROOT.RooArgSet(DTagZ))
159 
160  for i in range(data.numEntries()):
161 
162  row = data.get(i)
163 
164  tDT = row.getRealValue("DeltaT", 0, ROOT.kTRUE) - row.getRealValue("MCDeltaT", 0, ROOT.kTRUE)
165  if abs(tDT) < limDeltaT:
166  DT.setVal(tDT)
167  fitDataDT.add(ROOT.RooArgSet(DT))
168 
169  tDSigZ = row.getRealValue("z", 0, ROOT.kTRUE) - row.getRealValue("mcZ", 0, ROOT.kTRUE)
170 
171  if abs(tDSigZ) < limZSig:
172  DSigZ.setVal(tDSigZ)
173  fitDataSigZ.add(ROOT.RooArgSet(DSigZ))
174 
175  tDTagZ = row.getRealValue("TagVz", 0, ROOT.kTRUE) - row.getRealValue("mcTagVz", 0, ROOT.kTRUE)
176 
177  if abs(tDTagZ) < limZTag:
178  DTagZ.setVal(tDTagZ)
179  fitDataTagZ.add(ROOT.RooArgSet(DTagZ))
180 
181  fitDataDTErr.Print()
182  fitDataDT.Print()
183  fitDataSigZ.Print()
184  fitDataTagZ.Print()
185  numberOfEntries.append(data.numEntries())
186 
187 # Fit and plot of the DeltaT Error and DeltaTRECO - DeltaTMC
188 
189  Mu1 = ROOT.RooRealVar("Mu1", "Mu1", 0., -limDeltaT, limDeltaT)
190  Mu2 = ROOT.RooRealVar("Mu2", "Mu2", 0., -limDeltaT, limDeltaT)
191  Mu3 = ROOT.RooRealVar("Mu3", "Mu3", 0., -limDeltaT, limDeltaT)
192  Sigma1 = ROOT.RooRealVar("Sigma1", "Sigma1", 1.88046e+00, 0., limDeltaT)
193  Sigma2 = ROOT.RooRealVar("Sigma2", "Sigma2", 3.40331e+00, 0., limDeltaT)
194  Sigma3 = ROOT.RooRealVar("Sigma3", "Sigma3", 8.23171e-01, 0., limDeltaT)
195  frac1 = ROOT.RooRealVar("frac1", "frac1", 5.48703e-01, 0.0, 1.)
196  frac2 = ROOT.RooRealVar("frac2", "frac2", 2.60604e-01, 0.0, 1.)
197 
198  g1 = ROOT.RooGaussModel("g1", "g1", DT, Mu1, Sigma1)
199  g2 = ROOT.RooGaussModel("g2", "g2", DT, Mu2, Sigma2)
200  g3 = ROOT.RooGaussModel("g3", "g3", DT, Mu3, Sigma3)
201 
202  argset1 = ROOT.RooArgSet(g1)
203  argset2 = ROOT.RooArgSet(g2)
204  argset3 = ROOT.RooArgSet(g3)
205 
206  model = ROOT.RooAddModel("model", "model", ROOT.RooArgList(g1, g2, g3), ROOT.RooArgList(frac1, frac2))
207 
208  DT.setRange("fitRange", -limDeltaT, limDeltaT)
209 
210  fitRes = model.fitTo(
211  fitDataDT,
212  ROOT.RooFit.Minos(ROOT.kFALSE), ROOT.RooFit.Extended(ROOT.kFALSE),
213  ROOT.RooFit.NumCPU(8), ROOT.RooFit.Save())
214 
215  fitRes.Print("v")
216 
217  deltaTErr.Print("v")
218 
219  resFrame = DT.frame()
220  fitDataDT.plotOn(resFrame)
221 
222  meanCBS = ROOT.RooRealVar("meanCBS", "meanCBS", 5.37602e-01, 0.5, 1, "ps")
223  sigmaCBS = ROOT.RooRealVar("sigmaCBS", "sigmaCBS", 8.16334e-02, 0, 0.1, "ps")
224  alphaCBS = ROOT.RooRealVar("alphaCBS", "alphaCBS", -4.85571e-01, -1, 0, "")
225  nCBS = ROOT.RooRealVar("nCBS", "nCBS", 1.86325e+00, 0.5, 3, "")
226  dtErrCBS = ROOT.RooCBShape("dtErrGen", "dtErrGen", deltaTErr, meanCBS, sigmaCBS, alphaCBS, nCBS)
227 
228  MuErr1 = ROOT.RooRealVar("MuErr1", "MuErr1", 4.12399e-01, 0.3, 0.6, "ps")
229  SigmaErr1 = ROOT.RooRealVar("SigmaErr1", "SigmaErr1", 5.41152e-02, 0. - 1, 0.07, "ps")
230  gErr1 = ROOT.RooGaussModel("gErr1", "gErr1", deltaTErr, MuErr1, SigmaErr1)
231  fracErr1 = ROOT.RooRealVar("fracErr1", "fracErr1", 7.50810e-01, 0.0, 1.)
232 
233  MuErr2 = ROOT.RooRealVar("MuErr2", "MuErr2", 3.26658e-01, 0.2, 0.4, "ps")
234  SigmaErr2 = ROOT.RooRealVar("SigmaErr2", "SigmaErr2", 3.66794e-02, 0.01, 0.08, "ps")
235  gErr2 = ROOT.RooGaussModel("gErr2", "gErr2", deltaTErr, MuErr2, SigmaErr2)
236  fracErr2 = ROOT.RooRealVar("fracErr2", "fracErr2", 1.82254e-01, 0.0, 1.)
237 
238  modelTErr = ROOT.RooAddModel(
239  "modelErr", "modelErr", ROOT.RooArgList(
240  dtErrCBS, gErr1, gErr2), ROOT.RooArgList(
241  fracErr1, fracErr2))
242 
243  if VXDReq == 'PXD0' or VXDReq == 'PXD1' or VXDReq == 'PXD2':
244  CBSFitRes = modelTErr.fitTo(
245  fitDataDTErr,
246  ROOT.RooFit.Minos(ROOT.kFALSE), ROOT.RooFit.Extended(ROOT.kFALSE),
247  ROOT.RooFit.NumCPU(8), ROOT.RooFit.Save())
248 
249  dtErrCBS.Print()
250 
251  argset1TErr = ROOT.RooArgSet(dtErrCBS)
252  argset2TErr = ROOT.RooArgSet(gErr1)
253  argset3TErr = ROOT.RooArgSet(gErr2)
254 
255  resFrameDtErr = deltaTErr.frame()
256  fitDataDTErr.plotOn(resFrameDtErr)
257  modelTErr.plotOn(resFrameDtErr)
258  modelTErr.plotOn(
259  resFrameDtErr,
260  ROOT.RooFit.Components(argset1TErr),
261  ROOT.RooFit.LineColor(ROOT.kRed + 2),
262  ROOT.RooFit.LineWidth(4))
263  modelTErr.plotOn(
264  resFrameDtErr,
265  ROOT.RooFit.Components(argset2TErr),
266  ROOT.RooFit.LineColor(ROOT.kGreen + 3),
267  ROOT.RooFit.LineWidth(4))
268  modelTErr.plotOn(
269  resFrameDtErr,
270  ROOT.RooFit.Components(argset3TErr),
271  ROOT.RooFit.LineColor(ROOT.kMagenta + 2),
272  ROOT.RooFit.LineWidth(4))
273 
274  f1 = frac1.getVal()
275  f2 = frac2.getVal()
276  f3 = 1 - f1 - f2
277 
278  argsetList = []
279  argsetList.append([argset1, f1])
280  argsetList.append([argset2, f2])
281  argsetList.append([argset3, f3])
282 
283  argsetList = sorted(argsetList, key=itemgetter(1))
284 
285  model.plotOn(resFrame, ROOT.RooFit.LineWidth(3))
286  model.plotOn(resFrame, ROOT.RooFit.Components(argsetList[2][0]), ROOT.RooFit.LineColor(ROOT.kRed + 2), ROOT.RooFit.LineWidth(4))
287  model.plotOn(
288  resFrame, ROOT.RooFit.Components(
289  argsetList[1][0]), ROOT.RooFit.LineColor(
290  ROOT.kGreen + 3), ROOT.RooFit.LineWidth(4))
291  model.plotOn(
292  resFrame, ROOT.RooFit.Components(
293  argsetList[0][0]), ROOT.RooFit.LineColor(
294  ROOT.kMagenta + 2), ROOT.RooFit.LineWidth(4))
295 
296  resFrame.SetTitle("")
297  sXtitle = "#Deltat - Gen. #Delta t / ps"
298  resFrame.GetXaxis().SetTitle(sXtitle)
299  resFrame.GetXaxis().SetTitleSize(0.05)
300  resFrame.GetXaxis().SetLabelSize(0.045)
301  resFrame.GetYaxis().SetTitleSize(0.05)
302  resFrame.GetYaxis().SetTitleOffset(1.5)
303  resFrame.GetYaxis().SetLabelSize(0.045)
304 
305  shift = Mu1.getVal() * f1 + Mu2.getVal() * f2 + Mu3.getVal() * f3
306  resolution = Sigma1.getVal() * f1 + Sigma2.getVal() * f2 + Sigma3.getVal() * f3
307 
308  shiftErr = math.sqrt((Mu1.getError() * f1)**2 + (Mu2.getError() * f2)**2 + (Mu3.getError() * f3)**2)
309  resolutionErr = math.sqrt((Sigma1.getError() * f1)**2 + (Sigma2.getError() * f2)**2 + (Sigma3.getError() * f3)**2)
310 
311  if shiftErr < 0.01:
312  shiftErr = 0.01
313 
314  if resolutionErr < 0.01:
315  resolutionErr = 0.01
316 
317  Numbr = '{:d}'.format(int((f1 + f2) * fitDataDT.numEntries()))
318 
319  c1 = ROOT.TCanvas("c1", "c1", 1400, 1100)
320  c1.cd()
321  Pad = ROOT.TPad("p1", "p1", 0, 0, 1, 1, 0, 0, 0)
322  Pad.SetLeftMargin(0.15)
323  Pad.SetBottomMargin(0.15)
324  Pad.Draw()
325  Pad.cd()
326  resFrame.Draw()
327  legend = ROOT.TLegend(0.59, 0.6, 0.9, 0.9)
328  # legend.AddEntry(0, 'Entries' + '{:>11}'.format(Numbr))
329  legend.AddEntry(
330  0,
331  '#splitline{#mu_{#Delta t} =' +
332  '{: 4.2f}'.format(shift) +
333  '}{ #pm ' +
334  '{:4.2f}'.format(shiftErr) +
335  ' ps}')
336  legend.AddEntry(0, '#splitline{#sigma_{#Delta t} =' + '{: 4.2f}'.format(resolution) + '}{ #pm ' +
337  '{:4.2f}'.format(resolutionErr) + ' ps}')
338  legend.SetTextSize(0.054)
339  legend.SetFillColorAlpha(ROOT.kWhite, 0)
340  legend.Draw()
341  Pad.Update()
342  nPlot = PATH + "/test6_CPVResDeltaT" + VXDReq + ".pdf"
343  c1.SaveAs(nPlot)
344  c1.Clear()
345 
346  iResult.append(['mu = ' + '{: 4.3f}'.format(shift) + ' +- ' + '{:4.3f}'.format(shiftErr) + ' ps',
347  'sigma =' + '{: 4.3f}'.format(resolution) + ' +- ' + '{:4.3f}'.format(resolutionErr) + ' ps'])
348 
349  resFrameDtErr.SetTitle("")
350  sXtitleLandau = "#sigma_{#Deltat} / ps"
351  resFrameDtErr.GetXaxis().SetTitle(sXtitleLandau)
352  resFrameDtErr.GetXaxis().SetTitleSize(0.05)
353  resFrameDtErr.GetXaxis().SetLabelSize(0.045)
354  resFrameDtErr.GetYaxis().SetTitleSize(0.05)
355  resFrameDtErr.GetYaxis().SetTitleOffset(1.5)
356  resFrameDtErr.GetYaxis().SetLabelSize(0.045)
357 
358  c1 = ROOT.TCanvas("c1", "c1", 1400, 1100)
359  c1.cd()
360  Pad = ROOT.TPad("p1", "p1", 0, 0, 1, 1, 0, 0, 0)
361  Pad.SetLeftMargin(0.15)
362  Pad.SetBottomMargin(0.15)
363  Pad.Draw()
364  Pad.cd()
365  resFrameDtErr.Draw()
366  legend = ROOT.TLegend(0.59, 0.6, 0.9, 0.9)
367  # legend.AddEntry(0, 'Entries' + '{:>11}'.format(Numbr))
368  legend.AddEntry(0, '#splitline{#mu_{#Delta t} =' + '{: 4.2f}'.format(meanCBS.getVal()) + '}{ #pm ' +
369  '{:4.2f}'.format(meanCBS.getError()) + ' ps}') # '{:>6}'.format(Shift)
370  legend.AddEntry(0, '#splitline{#sigma_{#Delta t} =' + '{: 4.2f}'.format(sigmaCBS.getVal()) +
371  '}{ #pm ' + '{:4.2f}'.format(sigmaCBS.getError()) + ' ps}') # '{:>4}'.format(Resol)
372  legend.SetTextSize(0.054)
373  legend.SetFillColorAlpha(ROOT.kWhite, 0)
374  # legend.Draw()
375  Pad.Update()
376  nPlot = PATH + "/test6_CPVResDeltaTError" + VXDReq + ".pdf"
377  c1.SaveAs(nPlot)
378  c1.Clear()
379 
380  fitRes.Clear()
381  model.Clear()
382 
383  # Fit of Delta Z for B0_sig
384 
385  Mu1SigZ = ROOT.RooRealVar("Mu1SigZ", "Mu1SigZ", -6.03806e-06, -limZSig, limZSig)
386  Mu2SigZ = ROOT.RooRealVar("Mu2SigZ", "Mu2SigZ", 1.45755e-05, -limZSig, limZSig)
387  Mu3SigZ = ROOT.RooRealVar("Mu3SigZ", "Mu3SigZ", -1.84464e-04, -limZSig, limZSig)
388  Sigma1SigZ = ROOT.RooRealVar("Sigma1SigZ", "Sigma1SigZ", 4.03530e-03, 0., limZSig)
389  Sigma2SigZ = ROOT.RooRealVar("Sigma2SigZ", "Sigma2SigZ", 1.73995e-03, 0., limZSig)
390  Sigma3SigZ = ROOT.RooRealVar("Sigma3SigZ", "Sigma3SigZ", 2.18176e-02, 0., limZSig)
391  frac1SigZ = ROOT.RooRealVar("frac1SigZ", "frac1SigZ", 2.08032e-01, 0.0, 1.)
392  frac2SigZ = ROOT.RooRealVar("frac2SigZ", "frac2SigZ", 7.80053e-01, 0.0, 1.)
393 
394  g1SigZ = ROOT.RooGaussian("g1", "g1", DSigZ, Mu1SigZ, Sigma1SigZ)
395  g2SigZ = ROOT.RooGaussian("g2", "g2", DSigZ, Mu2SigZ, Sigma2SigZ)
396  g3SigZ = ROOT.RooGaussian("g3", "g3", DSigZ, Mu3SigZ, Sigma3SigZ)
397 
398  argset1SigZ = ROOT.RooArgSet(g1SigZ)
399  argset2SigZ = ROOT.RooArgSet(g2SigZ)
400  argset3SigZ = ROOT.RooArgSet(g3SigZ)
401 
402  modelSigZ = ROOT.RooAddPdf(
403  "modelSigZ", "modelSigZ", ROOT.RooArgList(g1SigZ, g2SigZ, g3SigZ),
404  ROOT.RooArgList(frac1SigZ, frac2SigZ))
405 
406  DSigZ.setRange("fitRange", -limZSig, limZSig)
407 
408  fitResSigZ = modelSigZ.fitTo(
409  fitDataSigZ, ROOT.RooFit.Minos(ROOT.kFALSE),
410  ROOT.RooFit.Extended(ROOT.kFALSE), ROOT.RooFit.NumCPU(1), ROOT.RooFit.Save())
411 
412  resFrameSigZ = DSigZ.frame()
413 
414  f1SigZ = frac1SigZ.getVal()
415  f2SigZ = frac2SigZ.getVal()
416  f3SigZ = 1 - f1SigZ - f2SigZ
417 
418  argsetList = []
419  argsetList.append([argset1SigZ, f1SigZ])
420  argsetList.append([argset2SigZ, f2SigZ])
421  argsetList.append([argset3SigZ, f3SigZ])
422 
423  argsetList = sorted(argsetList, key=itemgetter(1))
424 
425  fitDataSigZ.plotOn(resFrameSigZ)
426  modelSigZ.plotOn(resFrameSigZ, ROOT.RooFit.LineWidth(4))
427  modelSigZ.plotOn(resFrameSigZ, ROOT.RooFit.Components(argsetList[2][0]), ROOT.RooFit.LineColor(ROOT.kRed + 2),
428  ROOT.RooFit.LineWidth(4))
429  modelSigZ.plotOn(resFrameSigZ, ROOT.RooFit.Components(argsetList[1][0]), ROOT.RooFit.LineColor(ROOT.kGreen + 3),
430  ROOT.RooFit.LineWidth(4))
431  modelSigZ.plotOn(
432  resFrameSigZ,
433  ROOT.RooFit.Components(argsetList[0][0]),
434  ROOT.RooFit.LineColor(
435  ROOT.kMagenta + 2),
436  ROOT.RooFit.LineWidth(4))
437  resFrameSigZ.SetTitle("")
438 
439  sXtitleSigZ = "z - mcZ / cm"
440 
441  resFrameSigZ.GetXaxis().SetTitle(sXtitleSigZ)
442  resFrameSigZ.GetXaxis().SetLimits(-limZSig, limZSig)
443  resFrameSigZ.GetXaxis().SetTitleSize(0.05)
444  resFrameSigZ.GetXaxis().SetLabelSize(0.045)
445  resFrameSigZ.GetYaxis().SetTitleSize(0.05)
446  resFrameSigZ.GetYaxis().SetTitleOffset(1.5)
447  resFrameSigZ.GetYaxis().SetLabelSize(0.045)
448 
449  shiftSigZ = Mu1SigZ.getVal() * f1SigZ + Mu2SigZ.getVal() * f2SigZ + Mu3SigZ.getVal() * f3SigZ
450  resolutionSigZ = Sigma1SigZ.getVal() * f1SigZ + Sigma2SigZ.getVal() * f2SigZ + Sigma3SigZ.getVal() * f3SigZ
451 
452  shiftSigZ = shiftSigZ * 10000
453  resolutionSigZ = resolutionSigZ * 10000
454 
455  shiftErrSigZ = math.sqrt((Mu1SigZ.getError() * f1SigZ)**2 + (Mu2SigZ.getError() * f2SigZ) ** 2 +
456  (Mu3SigZ.getError() * f3SigZ)**2) * 10000
457  resolutionErrSigZ = math.sqrt((Sigma1SigZ.getError() * f1SigZ)**2 + (Sigma2SigZ.getError() * f2SigZ)**2 +
458  (Sigma3SigZ.getError() * f3SigZ)**2) * 10000
459 
460  cSig = ROOT.TCanvas("c1", "c1", 1400, 1100)
461  cSig.cd()
462  Pad = ROOT.TPad("p1", "p1", 0, 0, 1, 1, 0, 0, 0)
463  Pad.SetLeftMargin(0.15)
464  Pad.SetBottomMargin(0.15)
465  Pad.Draw()
466  Pad.cd()
467  resFrameSigZ.Draw()
468  legend = ROOT.TLegend(0.63, 0.65, 0.9, 0.9)
469  # NumbrSigZ = '{:d}'.format(int((f1+f2)*fitDataSigZ.numEntries()))
470  # legend.AddEntry(0, 'Entries' + '{:>11}'.format(NumbrSigZ))
471 
472  legend.AddEntry(
473  0,
474  '#splitline{#mu_{#Delta z} =' +
475  '{: 1.1f}'.format(shiftSigZ) +
476  '}{ #pm ' +
477  '{:1.1f}'.format(shiftErrSigZ) +
478  ' #mum}')
479  legend.AddEntry(
480  0,
481  '#splitline{#sigma_{#Delta z} =' +
482  '{: 1.1f}'.format(resolutionSigZ) +
483  '}{ #pm ' +
484  '{:1.1f}'.format(resolutionErrSigZ) +
485  ' #mum}')
486 
487  legend.SetTextSize(0.05)
488  legend.SetFillColorAlpha(ROOT.kWhite, 0)
489  legend.Draw()
490  Pad.Update()
491  nPlot = PATH + "/test6_CPVResDeltaZsig" + VXDReq + ".pdf"
492  cSig.SaveAs(nPlot)
493  cSig.Clear()
494 
495  fitResSigZ.Clear()
496  modelSigZ.Clear()
497 
498  iResult.append(['mu = ' + '{:^5.3f}'.format(shiftSigZ) + ' +- ' + '{:^4.3f}'.format(shiftErrSigZ) + ' mum',
499  'sigma = ' + '{:^4.3f}'.format(resolutionSigZ) + ' +- ' + '{:^4.3f}'.format(resolutionErrSigZ) + ' mum'])
500 
501  # Fit of Delta z for B0_Tag
502 
503  Mu1TagZ = ROOT.RooRealVar("Mu1TagZ", "Mu1TagZ", 0., -limZTag, limZTag)
504  Mu2TagZ = ROOT.RooRealVar("Mu2TagZ", "Mu2TagZ", 0., -limZTag, limZTag)
505  Mu3TagZ = ROOT.RooRealVar("Mu3TagZ", "Mu3TagZ", 0., -limZTag, limZTag)
506  Sigma1TagZ = ROOT.RooRealVar("Sigma1TagZ", "Sigma1TagZ", 2.51877e-02, 0., limZTag)
507  Sigma2TagZ = ROOT.RooRealVar("Sigma2TagZ", "Sigma2TagZ", 1.54011e-02, 0., limZTag)
508  Sigma3TagZ = ROOT.RooRealVar("Sigma3TagZ", "Sigma3TagZ", 1.61081e-02, 0., limZTag)
509  frac1TagZ = ROOT.RooRealVar("frac1TagZ", "frac1TagZ", 1.20825e-01, 0.0, 1.)
510  frac2TagZ = ROOT.RooRealVar("frac2TagZ", "frac2TagZ", 1.10840e-01, 0.0, 1.)
511 
512  g1TagZ = ROOT.RooGaussian("g1", "g1", DTagZ, Mu1TagZ, Sigma1TagZ)
513  g2TagZ = ROOT.RooGaussian("g2", "g2", DTagZ, Mu2TagZ, Sigma2TagZ)
514  g3TagZ = ROOT.RooGaussian("g3", "g3", DTagZ, Mu3TagZ, Sigma3TagZ)
515 
516  argset1TagZ = ROOT.RooArgSet(g1TagZ)
517  argset2TagZ = ROOT.RooArgSet(g2TagZ)
518  argset3TagZ = ROOT.RooArgSet(g3TagZ)
519 
520  modelTagZ = ROOT.RooAddPdf("modelTagZ", "modelTagZ", ROOT.RooArgList(
521  g1TagZ, g2TagZ, g3TagZ), ROOT.RooArgList(frac1TagZ, frac2TagZ))
522 
523  DTagZ.setRange("fitRange", -limZTag, limZTag)
524 
525  fitResTagZ = modelTagZ.fitTo(
526  fitDataTagZ, ROOT.RooFit.Minos(
527  ROOT.kFALSE), ROOT.RooFit.Extended(
528  ROOT.kFALSE), ROOT.RooFit.NumCPU(1), ROOT.RooFit.Save())
529 
530  fitResTagZ.Print("v")
531  resFrameTagZ = DTagZ.frame()
532 
533  f1TagZ = frac1TagZ.getVal()
534  f2TagZ = frac2TagZ.getVal()
535  f3TagZ = 1 - f1TagZ - f2TagZ
536 
537  argsetList = []
538  argsetList.append([argset1TagZ, f1TagZ])
539  argsetList.append([argset2TagZ, f2TagZ])
540  argsetList.append([argset3TagZ, f3TagZ])
541 
542  argsetList = sorted(argsetList, key=itemgetter(1))
543 
544  fitDataTagZ.plotOn(resFrameTagZ)
545  modelTagZ.plotOn(resFrameTagZ, ROOT.RooFit.LineWidth(4))
546  modelTagZ.plotOn(
547  resFrameTagZ, ROOT.RooFit.Components(
548  argsetList[2][0]), ROOT.RooFit.LineColor(
549  ROOT.kRed + 2), ROOT.RooFit.LineWidth(4))
550  modelTagZ.plotOn(
551  resFrameTagZ,
552  ROOT.RooFit.Components(argsetList[1][0]),
553  ROOT.RooFit.LineColor(ROOT.kGreen + 3),
554  ROOT.RooFit.LineWidth(4))
555  modelTagZ.plotOn(
556  resFrameTagZ,
557  ROOT.RooFit.Components(argsetList[0][0]),
558  ROOT.RooFit.LineColor(ROOT.kMagenta + 2),
559  ROOT.RooFit.LineWidth(4))
560  resFrameTagZ.SetTitle("")
561 
562  sXtitleTagZ = "TagVz - mcTagVz / cm"
563 
564  resFrameTagZ.GetXaxis().SetTitle(sXtitleTagZ)
565  resFrameTagZ.GetXaxis().SetTitleSize(0.05)
566  resFrameTagZ.GetXaxis().SetLabelSize(0.045)
567  resFrameTagZ.GetYaxis().SetTitleSize(0.05)
568  resFrameTagZ.GetYaxis().SetTitleOffset(1.5)
569  resFrameTagZ.GetYaxis().SetLabelSize(0.045)
570 
571  shiftTagZ = Mu1TagZ.getVal() * f1TagZ + Mu2TagZ.getVal() * f2TagZ + Mu3TagZ.getVal() * f3TagZ
572  resolutionTagZ = Sigma1TagZ.getVal() * f1TagZ + Sigma2TagZ.getVal() * f2TagZ + Sigma3TagZ.getVal() * f3TagZ
573 
574  shiftTagZ = shiftTagZ * 10000
575  resolutionTagZ = resolutionTagZ * 10000
576 
577  shiftErrTagZ = math.sqrt((Mu1TagZ.getError() * f1TagZ)**2 + (Mu2TagZ.getError() * f2TagZ) ** 2 +
578  (Mu3TagZ.getError() * f3TagZ)**2) * 10000
579  resolutionErrTagZ = math.sqrt((Sigma1TagZ.getError() * f1TagZ)**2 + (Sigma2TagZ.getError() * f2TagZ)**2 +
580  (Sigma3TagZ.getError() * f3TagZ)**2) * 10000
581 
582  cTag = ROOT.TCanvas("c1", "c1", 1400, 1100)
583  cTag.cd()
584  Pad = ROOT.TPad("p1", "p1", 0, 0, 1, 1, 0, 0, 0)
585  Pad.SetLeftMargin(0.15)
586  Pad.SetBottomMargin(0.15)
587  Pad.Draw()
588  Pad.cd()
589  resFrameTagZ.Draw()
590  legend = ROOT.TLegend(0.64, 0.65, 0.9, 0.9)
591  # NumbrTagZ = '{:d}'.format(int((f1+f2)*fitDataTagZ.numEntries()))
592  # legend.AddEntry(0, 'Entries' + '{:>11}'.format(NumbrTagZ))
593 
594  legend.AddEntry(0, '#splitline{#mu_{#Delta z} =' + '{: 1.1f}'.format(shiftTagZ) +
595  '}{ #pm ' + '{: 1.1f}'.format(shiftErrTagZ) + ' #mum}')
596  legend.AddEntry(0, '#splitline{#sigma_{#Delta z} =' + '{: 1.1f}'.format(resolutionTagZ) +
597  '}{ #pm ' + '{: 1.1f}'.format(resolutionErrTagZ) + ' #mum}')
598  legend.SetTextSize(0.05)
599  legend.SetFillColorAlpha(ROOT.kWhite, 0)
600  legend.Draw()
601  Pad.Update()
602  nPlot = PATH + "/test6_CPVResDeltaZtag" + VXDReq + ".pdf"
603  cTag.SaveAs(nPlot)
604  cTag.Clear()
605 
606  iResult.append(['mu = ' + '{:^5.3f}'.format(shiftTagZ) + ' +- ' + '{:^4.3f}'.format(shiftErrTagZ) + ' mum',
607  'sigma = ' + '{:^4.3f}'.format(resolutionTagZ) + ' +- ' + '{:^4.3f}'.format(resolutionErrTagZ) + ' mum'])
608 
609  fitResults.append(iResult)
610 
611 
612 print('*********************** FIT RESULTS ***************************')
613 print('* *')
614 print('* WEIGHTED AVERAGES OF THE PARAMETERS *')
615 print('* OF THE FITTED 3 GAUSSIAN *')
616 print('* *')
617 print('**************** WITHOUT VXD HIT REQUIREMENT ******************')
618 print('* *')
619 print('* DeltaT - Gen. DeltaT *')
620 print('* *')
621 print('* ' + fitResults[0][0][0] + ' ' + fitResults[0][0][1] + ' *')
622 print('* *')
623 print('* SigZ - Gen. SigZ *')
624 print('* *')
625 print('* ' + fitResults[0][1][0] + ' ' + fitResults[0][1][1] + ' *')
626 print('* *')
627 print('* TagZ - Gen. TagZ *')
628 print('* *')
629 print('* ' + fitResults[0][2][0] + ' ' + fitResults[0][2][1] + ' *')
630 print('* *')
631 print('********REQUIRING BOTH MUON TRACKS TO HAVE A ' + sys.argv[3] + ' HIT***********')
632 print('* *')
633 print('* Efficiency *')
634 print('* *')
635 print('* N_' + VXDReqs[1] + '/N_' + VXDReqs[0] + ' = ' + str(numberOfEntries[1]) + "/" + str(numberOfEntries[0]) + ' = ' +
636  '{:^3.2f}'.format(float((numberOfEntries[1] / numberOfEntries[0]) * 100)) + '% *')
637 print('* *')
638 print('* DeltaT - Gen. DeltaT *')
639 print('* *')
640 print('* ' + fitResults[1][0][0] + ' ' + fitResults[1][0][1] + ' *')
641 print('* *')
642 print('* SigZ - Gen. SigZ *')
643 print('* *')
644 print('* ' + fitResults[1][1][0] + ' ' + fitResults[1][1][1] + ' *')
645 print('* *')
646 print('* TagZ - Gen. TagZ *')
647 print('* *')
648 print('* ' + fitResults[1][2][0] + ' ' + fitResults[1][2][1] + ' *')
649 print('* *')
650 print('***************************************************************')