12 Validation of the Invariant Collision Energy calibration
16 from prompt
import ValidationSettings
25 import matplotlib.pyplot
as plt
26 from matplotlib
import rcParams
32 from datetime
import datetime
34 from collections
import OrderedDict
35 from shutil
import copyfile, rmtree
39 settings = ValidationSettings(name=
'Ecms Calibrations',
41 download_files=[
'stdout'],
45 def setPlotRange(df, tag):
47 Function which adjusts the y-axis range of the plot according to tag
51 df4S = df[df[
'Ecms'] > 10560][
'Ecms']
54 plt.ylim(bottom=yMin-3)
56 dfHigh = df[df[
'Ecms'] < 10600][
'Ecms']
62 dfOff = df[df[
'Ecms'] < 10560][
'Ecms']
70 Converts time from UTC to the JST
72 return np.vectorize(
lambda t: datetime.utcfromtimestamp((t + 9) * 3600))(times)
75 def plotSplitLines(dfC):
77 Plot vertical lines in places where run energy type changes
79 indx = np.where(np.diff(dfC[
'pull'] == 0))[0]
81 tt1 = toJST(dfC[
't2'].loc[i]).item()
82 tt2 = toJST(dfC[
't1'].loc[i+1]).item()
83 tt = tt1 + (tt2 - tt1) / 2
84 plt.axvline(x=tt, color=
'g', linestyle=
'--')
87 def get_Ecms_values(path):
89 Load the values of the Ecms properties from the payloads into the list
93 with open(path +
'/database.txt')
as fDB:
99 fN = ll[0].replace(
'/',
'_') +
'_rev_' + ll[1] +
'.root'
102 runDict[fN] = ((int(r[0]), int(r[1])), (int(r[2]), int(r[3])))
106 fList = glob(path +
'/*.root')
108 bName = os.path.basename(fName)
109 runExp = runDict[bName]
111 f = ROOT.TFile.Open(fName)
112 ecmsObj = f.Get(
"CollisionInvariantMass")
113 assert(ecmsObj.ClassName() ==
"Belle2::CollisionInvariantMass")
115 eCMS = ecmsObj.getMass()
116 eCMSe = ecmsObj.getMassError()
117 eCMSs = ecmsObj.getMassSpread()
119 arr.append((runExp, (eCMS, eCMSe, eCMSs)))
123 arr = sorted(arr, key=
lambda x: x[0])
131 Data are loaded from text files to pandas
135 self.
dfBdfB = pd.read_csv(f
'{location}/BonlyEcmsCalib.txt', delim_whitespace=
True)
138 self.
dfCdfC = pd.read_csv(f
'{location}/finalEcmsCalib.txt', delim_whitespace=
True)
141 self.
dfMdfM = pd.read_csv(f
'{location}/mumuEcalib.txt', delim_whitespace=
True)
149 for i
in range(len(dfM)):
150 if (dfM[
'id'][i]) == 0:
153 dfM.at[i,
'type'] = state
154 dfC.at[i,
'type'] = state
159 Plot single line with error band
163 nan = np.full(len(df), np.nan)
164 avg = (df[
't1']+df[
't2']).to_numpy()/2
165 times = np.c_[df[[
't1',
't2']].to_numpy(), avg].ravel()
166 eCMS = np.c_[df[[var, var]].to_numpy(), nan].ravel()
167 eCMSu = np.c_[df[[varUnc, varUnc]].to_numpy(), nan].ravel()
169 timesg = np.c_[df[
't1'].to_numpy(), avg, df[
't2'].to_numpy()].ravel()
170 eCMSg = np.c_[df[var].to_numpy(), nan, df[var].to_numpy()].ravel()
173 timesg = toJST(timesg)
175 plt.fill_between(times, eCMS-eCMSu, eCMS+eCMSu, alpha=0.2, color=color)
176 plt.plot(times, eCMS, linewidth=2, color=color, label=label)
177 plt.plot(timesg, eCMSg, linewidth=2, color=color, alpha=0.35)
181 Plot Ecms obtained from combined, hadB and mumu method
184 Plotter.plotLine(self.
dfBdfB,
'Ecms',
'green', label=
'B decay method')
185 Plotter.plotLine(self.
dfCdfC,
'Ecms',
'red', label=
'Combined method')
187 d = np.nanmedian(self.
dfCdfC[
'Ecms']-self.
dfMdfM[
'Ecms'])
189 dfMc = self.
dfMdfM.copy()
191 Plotter.plotLine(dfMc,
'Ecms',
'blue', label=f
'mumu method (+{round(d,1)} MeV)')
193 plotSplitLines(self.
dfCdfC)
194 setPlotRange(self.
dfCdfC, tag)
197 plt.ylabel(
'Ecms [MeV]')
201 loc =
'plots/allData'
202 if limits
is not None:
203 plt.xlim(datetime.strptime(limits[0],
'%Y-%m-%d'), datetime.strptime(limits[1],
'%Y-%m-%d'))
204 loc =
'plots/' + limits[0] +
'to' + limits[1]
206 os.makedirs(loc, exist_ok=
True)
207 plt.savefig(f
'{loc}/EcmsComparison{tag}.png')
212 Plot Ecms spread estimated from the combined method and B-decay method
215 Plotter.plotLine(self.
dfBdfB,
'spread',
'green', label=
'B decay method')
216 Plotter.plotLine(self.
dfCdfC,
'spread',
'red', label=
'Combined method')
218 plotSplitLines(self.
dfCdfC)
223 plt.ylabel(
'spread [MeV]')
225 loc =
'plots/allData'
226 if limits
is not None:
227 plt.xlim(datetime.strptime(limits[0],
'%Y-%m-%d'), datetime.strptime(limits[1],
'%Y-%m-%d'))
228 loc =
'plots/' + limits[0] +
'to' + limits[1]
230 os.makedirs(loc, exist_ok=
True)
231 plt.savefig(f
'{loc}/SpreadComparison.png')
235 def plotCurve(df, var, label, withBand=True, withCurve=True):
237 Plot curve with possible error band where large intervals are distinguished by color
240 def plotBands(df, var, color, withBand=True, withCurve=True):
242 nan = np.full(len(df), np.nan)
243 avg = (df[
't1']+df[
't2']).to_numpy()/2
244 times = np.c_[df[[
't1',
't2']].to_numpy(), avg].ravel()
245 eCMS = np.c_[df[[var, var]].to_numpy(), nan].ravel()
249 eCMSu = np.c_[df[[varUnc, varUnc]].to_numpy(), nan].ravel()
250 plt.fill_between(times, eCMS-eCMSu, eCMS+eCMSu, alpha=0.15, color=color)
252 plt.plot(times, eCMS, linewidth=2, color=color)
254 nan = np.full(len(df), np.nan)
255 avg = (df[
't1']+df[
't2']).to_numpy()/2
256 timesg = np.c_[df[
't1'].to_numpy(), avg, df[
't2'].to_numpy()].ravel()
257 eCMSg = np.c_[df[var].to_numpy(), nan, df[var].to_numpy()].ravel()
258 timesg = toJST(timesg)
261 plt.plot(timesg, eCMSg, linewidth=2, color=
'gray', alpha=0.35)
263 plotBands(df[df[
'type'] == 1], var,
'red', withBand, withCurve)
264 plotBands(df[df[
'type'] == -1], var,
'blue', withBand, withCurve)
268 Plot shift between Bhad and mumu method
271 Plotter.plotCurve(self.
dfCdfC,
'shift', label=
'Combined method')
273 plotSplitLines(self.
dfCdfC)
276 plt.ylabel(
'shift [MeV]')
278 loc =
'plots/allData'
279 if limits
is not None:
280 plt.xlim(datetime.strptime(limits[0],
'%Y-%m-%d'), datetime.strptime(limits[1],
'%Y-%m-%d'))
281 loc =
'plots/' + limits[0] +
'to' + limits[1]
283 os.makedirs(loc, exist_ok=
True)
284 plt.savefig(f
'{loc}/Shift.png')
290 Plot Ecms of the combined fit, 'tag' can be '4S' means that y-axis is zoomed to 4S mass region, 'Off' to off-resonance
293 Plotter.plotCurve(self.
dfCdfC,
'Ecms', label=
'Combined method')
295 plotSplitLines(self.
dfCdfC)
298 plt.ylabel(
'Ecms [MeV]')
300 setPlotRange(self.
dfCdfC, tag)
302 loc =
'plots/allData'
303 if limits
is not None:
304 plt.xlim(datetime.strptime(limits[0],
'%Y-%m-%d'), datetime.strptime(limits[1],
'%Y-%m-%d'))
305 loc =
'plots/' + limits[0] +
'to' + limits[1]
307 os.makedirs(loc, exist_ok=
True)
308 plt.savefig(f
'{loc}/Ecms{tag}.png')
313 Plot pulls of the combined fit.
316 dfC = self.
dfCdfC.copy()
319 Plotter.plotCurve(dfC,
'pull', label=
'Combined method', withBand=
False)
320 Plotter.plotCurve(dfC,
'ref', label=
'Combined method', withCurve=
False)
322 plotSplitLines(self.
dfCdfC)
328 loc =
'plots/allData'
329 if limits
is not None:
330 plt.xlim(datetime.strptime(limits[0],
'%Y-%m-%d'), datetime.strptime(limits[1],
'%Y-%m-%d'))
331 loc =
'plots/' + limits[0] +
'to' + limits[1]
333 os.makedirs(loc, exist_ok=
True)
334 plt.savefig(f
'{loc}/Pull.png')
339 def read_Combined_data(outputDir):
341 It reads the calibration table from the text file produced by the CAF calibration.
342 This text file includes the results from the final calibration.
343 (combined or mumu-only)
347 with open(outputDir +
'/finalEcmsCalib.txt',
"r")
as text_file:
348 for i, ll
in enumerate(text_file):
351 ll = ll.strip().split()
355 ll[1])), (int(ll[2]), int(ll[3]), int(ll[4]), int(ll[5])), int(
368 def read_Bonly_data(outputDir):
370 It reads the calibration table from the text file produced by the CAF calibration.
371 This text file includes the results from the B-only calibration.
375 with open(outputDir +
'/BonlyEcmsCalib.txt',
"r")
as text_file:
376 for i, ll
in enumerate(text_file):
379 ll = ll.strip().split()
383 ll[1])), (int(ll[2]), int(ll[3]), int(ll[4]), int(ll[5])),
384 (float(ll[2+4]), float(ll[3+4])), (float(ll[4+4]), float(ll[5+4]))))
388 def read_mumu_data(outputDir):
390 It reads the calibration table from the text file produced by the CAF calibration.
391 This text file includes the results from the mumu-based calibration.
395 with open(outputDir +
'/mumuEcalib.txt',
"r")
as text_file:
396 for i, ll
in enumerate(text_file):
399 ll = ll.strip().split()
403 ll[2+1])), (int(ll[2+2]), int(ll[2+3]), int(ll[2+4]), int(ll[2+5])),
404 (float(ll[2+6]), float(ll[2+7])), (int(ll[0]), int(ll[1]))))
409 def create_hadB_fit_plots(outputDir, pdflatex):
411 Create multi-page pdf file with the fit plots for hadronic B decays (with mumu constraint).
412 The file is created using pdflatex
415 arr = read_Combined_data(outputDir)
419 files = glob(outputDir+
'/'+dName +
'/*.pdf')
420 files = list(map(os.path.basename, files))
422 items = OrderedDict()
425 res = re.search(
'B[p0]_([0-9]*)_([0-9]*)\\.pdf', f)
426 t = int(res.group(1))
427 i = int(res.group(2))
431 items[t] = max(items[t], i + 1)
433 times = sorted(items)
435 header =
"""\\documentclass[aspectratio=169]{beamer}
436 \\usepackage{graphicx}
445 for i0Temp
in range(len(arr)):
446 if int(round(arr[i0Temp][0][0], 0)) == t:
449 assert(i0
is not None)
451 frac = 1. / (items[t] + 0.2)
453 frac = 1. / (items[t] + 0.3)
455 body +=
'\\begin{frame}[t]\n'
456 shift, shiftE = str(round(arr[i0][5][0], 1)), str(round(arr[i0][5][1], 1))
457 spread, spreadE = str(round(arr[i0][6][0], 1)), str(round(arr[i0][6][1], 1))
459 body +=
'$E_\\mathrm{shift} = (' + shift +
'\\pm' + shiftE +
')$~MeV \\hspace{2cm} \n'
460 body +=
'$E_\\mathrm{spread} = (' + spread +
'\\pm' + spreadE +
')$~MeV \\\\ \\vspace{0.5cm}\n'
461 for iShort
in range(items[t]):
463 tStart, tEnd = arr[i][0][0], arr[i][0][1]
465 exp1, run1, exp2, run2 = map(str, arr[i][1])
467 eCMS, eCMSe = str(round(arr[i][3][0], 1)), str(round(arr[i][3][1], 1))
468 pull = str(round(arr[i][4], 2))
470 t1 = datetime.utcfromtimestamp((tStart + 9) * 3600).strftime(
'%y-%m-%d %H:%M')
471 t2 = datetime.utcfromtimestamp((tEnd + 9) * 3600).strftime(
'%y-%m-%d %H:%M')
473 body +=
'\\begin{minipage}{' + str(frac) +
'\\textwidth}\n'
474 body +=
'\\begin{center}\n'
475 body +=
'\\scriptsize ' + exp1 +
' '+run1 +
' \\\\\n'
476 body +=
'\\scriptsize ' + exp2 +
' '+run2 +
' \\\\\n'
477 body +=
'\\scriptsize ' + t1 +
' \\\\\n'
478 body +=
'\\scriptsize ' + t2 +
' \\\\\n'
479 body +=
'\\scriptsize ' + str(round(tEnd - tStart, 1)) +
' hours \\\\ \\vspace{0.3cm}\n'
480 body +=
'$(' + eCMS +
'\\pm' + eCMSe +
')$~MeV \\\\\n'
481 body +=
'$\\mathrm{pull} = ' + pull +
'$ \\\\\n'
482 body +=
'\\includegraphics[width=1.0\\textwidth]{' + outputDir + \
483 '/' + dName +
'/B0_' + str(t) +
'_' + str(iShort) +
'.pdf}\\\\\n'
484 body +=
'\\includegraphics[width=1.0\\textwidth]{' + outputDir + \
485 '/' + dName +
'/Bp_' + str(t) +
'_' + str(iShort) +
'.pdf}\n'
486 body +=
'\\end{center}\n'
487 body +=
'\\end{minipage}\n'
489 body +=
'\\end{frame}\n\n'
491 tail =
'\n\\end{document}'
493 whole = header + body + tail
495 os.makedirs(
'tmp', exist_ok=
True)
497 with open(
"tmp/hadBfits.tex",
"w")
as text_file:
498 text_file.write(whole)
500 subprocess.call(f
'{pdflatex} tmp/hadBfits.tex', shell=
True)
502 os.makedirs(
'plots', exist_ok=
True)
503 copyfile(
'hadBfits.pdf',
'plots/hadBfits.pdf')
505 ext = [
'aux',
'log',
'nav',
'out',
'pdf',
'snm',
'toc']
507 if os.path.exists(f
'hadBfits.{e}'):
508 os.remove(f
'hadBfits.{e}')
514 def create_hadBonly_fit_plots(outputDir, pdflatex):
516 Create multi-page pdf file with the fit plots for hadronic B decays (without using mumu constraint).
517 The file is created using pdflatex
520 arr = read_Bonly_data(outputDir)
522 dName =
'plotsHadBonly'
524 files = glob(outputDir+
'/'+dName +
'/*.pdf')
525 files = list(map(os.path.basename, files))
530 res = re.search(
'B[p0]Single_([0-9]*)\\.pdf', f)
531 t = int(res.group(1))
534 items = sorted(items)
536 header =
"""\\documentclass[aspectratio=169]{beamer}
537 \\usepackage{graphicx}
543 for i, t
in enumerate(items):
545 body +=
'\\begin{frame}[t]\n'
547 exp1, run1, exp2, run2 = map(str, arr[i][1])
548 eCMS, eCMSe = str(round(arr[i][2][0], 1)), str(round(arr[i][2][1], 1))
549 spread, spreadE = str(round(arr[i][3][0], 1)), str(round(arr[i][3][1], 1))
551 body +=
'$E_\\mathrm{cms} = (' + eCMS +
'\\pm' + eCMSe +
')$~MeV \\\\\n'
552 body +=
'$E_\\mathrm{spread} = (' + spread +
'\\pm' + spreadE +
')$~MeV \\\\ \\vspace{0.5cm}\n'
554 tStart, tEnd = arr[i][0][0], arr[i][0][1]
556 t1 = datetime.utcfromtimestamp((tStart + 9) * 3600).strftime(
'%y-%m-%d %H:%M')
557 t2 = datetime.utcfromtimestamp((tEnd + 9) * 3600).strftime(
'%y-%m-%d %H:%M')
559 body +=
'\\begin{minipage}{' + str(0.99) +
'\\textwidth}\n'
560 body +=
'\\begin{center}\n'
561 body +=
'\\scriptsize ' + exp1 +
' ' + run1 +
'\\hspace{1cm} ' + t1 +
' \\\\\n'
562 body +=
'\\scriptsize ' + exp2 +
' ' + run2 +
'\\hspace{1cm} ' + t2 +
' \\\\\n'
563 body +=
'\\scriptsize ' + str(round(tEnd - tStart, 1)) +
' hours \\\\ \\vspace{0.3cm}\n'
564 body +=
'\\includegraphics[width=0.48\\textwidth]{'+outputDir +
'/' + dName +
'/B0Single_' + str(t) +
'.pdf}\n'
565 body +=
'\\includegraphics[width=0.48\\textwidth]{'+outputDir +
'/' + dName +
'/BpSingle_' + str(t) +
'.pdf}\n'
566 body +=
'\\end{center}\n'
567 body +=
'\\end{minipage}\n'
569 body +=
'\\end{frame}\n\n'
571 tail =
'\n\\end{document}'
573 whole = header + body + tail
575 os.makedirs(
'tmp', exist_ok=
True)
577 with open(
"tmp/hadBonlyFits.tex",
"w")
as text_file:
578 text_file.write(whole)
580 subprocess.call(f
'{pdflatex} tmp/hadBonlyFits.tex', shell=
True)
582 os.makedirs(
'plots', exist_ok=
True)
583 copyfile(
'hadBonlyFits.pdf',
'plots/hadBonlyFits.pdf')
585 ext = [
'aux',
'log',
'nav',
'out',
'pdf',
'snm',
'toc']
587 if os.path.exists(f
'hadBonlyFits.{e}'):
588 os.remove(f
'hadBonlyFits.{e}')
590 if os.path.exists(
"tmp/hadBonlyFits.tex"):
591 os.remove(
"tmp/hadBonlyFits.tex")
595 def create_mumu_fit_plots(outputDir, pdflatex):
597 Create multi-page pdf file with the fit plots for mumu method.
598 The file is created using pdflatex
601 arr = read_mumu_data(outputDir)
604 for i, a
in enumerate(arr):
605 if a[3][1] == a[3][0] - 1:
606 limits.append((i - a[3][1], a[3][0]))
610 files = glob(outputDir+
'/'+dName +
'/*.pdf')
611 files = list(map(os.path.basename, files))
616 res = re.search(
'mumu_([0-9]*)\\.pdf', f)
617 t = int(res.group(1))
620 items = sorted(items)
622 header =
"""\\documentclass[aspectratio=169]{beamer}
623 \\usepackage{graphicx}
647 body +=
'\\begin{frame}[t]\n'
649 for i
in range(k, k+n):
651 body +=
'\\begin{minipage}{' + str(frac) +
'\\textwidth}\n'
653 exp1, run1, exp2, run2 = map(str, arr[i][1])
654 eCMS, eCMSe = str(round(arr[i][2][0], 1)), str(round(arr[i][2][1], 1))
656 tStart, tEnd = arr[i][0][0], arr[i][0][1]
658 t1 = datetime.utcfromtimestamp((tStart + 9) * 3600).strftime(
'%y-%m-%d %H:%M')
659 t2 = datetime.utcfromtimestamp((tEnd + 9) * 3600).strftime(
'%y-%m-%d %H:%M')
661 body +=
'\\begin{center}\n'
662 body +=
'\\tiny $E_\\mathrm{cms} = (' + eCMS +
'\\pm' + eCMSe +
')$~MeV \\\\\n'
663 body +=
'\\tiny ' + exp1 +
' ' + run1 +
'\\hspace{0.3cm} ' + t1 +
' \\\\\n'
664 body +=
'\\tiny ' + exp2 +
' ' + run2 +
'\\hspace{0.3cm} ' + t2 +
' \\\\\n'
665 body +=
'\\tiny ' + str(round(tEnd - tStart, 1)) +
' hours \\vspace{0.3cm}\n'
666 body +=
'\\includegraphics[trim=0.3cm 0.0cm 1.3cm 0.7cm,clip=true,width=0.99\\textwidth]{' + \
667 outputDir +
'/' + dName +
'/mumu_' + str(items[i]) +
'.pdf}\n'
668 body +=
'\\end{center}\n'
669 body +=
'\\end{minipage}\n'
671 body +=
'\\end{frame}\n\n'
673 tail =
'\n\\end{document}'
675 whole = header + body + tail
677 os.makedirs(
'tmp', exist_ok=
True)
679 with open(
"tmp/mumuFits.tex",
"w")
as text_file:
680 text_file.write(whole)
682 subprocess.call(f
'{pdflatex} tmp/mumuFits.tex', shell=
True)
684 os.makedirs(
'plots', exist_ok=
True)
685 copyfile(
'mumuFits.pdf',
'plots/mumuFits.pdf')
687 ext = [
'aux',
'log',
'nav',
'out',
'pdf',
'snm',
'toc']
689 if os.path.exists(f
'mumuFits.{e}'):
690 os.remove(f
'mumuFits.{e}')
692 if os.path.exists(
"tmp/mumuFits.tex"):
693 os.remove(
"tmp/mumuFits.tex")
697 def run_validation(job_path, input_data_path, requested_iov, expert_config):
699 Create validation plots related to the Ecms calibration
703 result = subprocess.run(
'locate pdflatex | grep "/pdflatex$"', stdout=subprocess.PIPE, shell=
True)
704 pdflatex = result.stdout.strip().decode(
'utf-8')
707 if expert_config !=
'':
708 expert_config = json.loads(expert_config)
711 if expert_config
is not None and 'plotsRanges' in expert_config:
712 allLimits += expert_config[
'plotsRanges']
714 plt.figure(figsize=(18, 9))
715 rcParams[
'axes.formatter.useoffset'] =
False
717 location = f
'{job_path}/eCMS/0/algorithm_output'
721 for limits
in allLimits:
722 plotter.plotEcmsComparison(limits, tag=
'4S')
723 plotter.plotEcmsComparison(limits, tag=
'Off')
724 plotter.plotEcmsComparison(limits)
725 plotter.plotSpreadComparison(limits)
727 plotter.plotEcms(limits, tag=
'4S')
728 plotter.plotEcms(limits, tag=
'Off')
729 plotter.plotEcms(limits, tag=
'')
731 plotter.plotShift(limits)
732 plotter.plotPulls(limits)
735 create_hadB_fit_plots(location, pdflatex)
736 create_hadBonly_fit_plots(location, pdflatex)
737 create_mumu_fit_plots(location, pdflatex)
740 copyfile(f
'{location}/BonlyEcmsCalib.txt',
'plots/BonlyEcmsCalib.txt')
741 copyfile(f
'{location}/finalEcmsCalib.txt',
'plots/finalEcmsCalib.txt')
742 copyfile(f
'{location}/mumuEcalib.txt',
'plots/mumuEcalib.txt')
745 if __name__ ==
"__main__":
def plotEcms(self, limits=None, tag='')
def __init__(self, location)
def plotSpreadComparison(self, limits=None)
def plotLine(df, var, color, label)
def plotEcmsComparison(self, limits=None, tag='')
def plotPulls(self, limits=None)
def plotShift(self, limits=None)
def plotCurve(df, var, label, withBand=True, withCurve=True)