You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

535 lines
19 KiB

10 months ago
10 months ago
10 months ago
10 months ago
10 months ago
10 months ago
10 months ago
10 months ago
  1. # flake8: noqa
  2. """
  3. Takes data from Recent_get_resolution_and_eff_data.py and calculates efficiencies
  4. python scripts/MyPrCheckerEfficiency.py --filename data/resolutions_and_effs_B_normal_weights.root data/resolutions_and_effs_B_electron_weights.root
  5. --outfile data_results/PrCheckerNormalElectron.root --label normal electron --trackers Match BestLong
  6. python scripts/MyPrCheckerEfficiency.py --filename data/resolutions_and_effs_D_electron_weights.root --outfile data_results/PrCheckerDElectron.root
  7. """
  8. import argparse
  9. from ROOT import TMultiGraph, TLatex, TCanvas, TFile, TGaxis
  10. from ROOT import (
  11. kGreen,
  12. kBlue,
  13. kBlack,
  14. kAzure,
  15. kOrange,
  16. kMagenta,
  17. kCyan,
  18. kGray,
  19. kViolet,
  20. kTeal,
  21. )
  22. from ROOT import gROOT, gStyle, gPad
  23. from ROOT import TEfficiency
  24. from array import array
  25. gROOT.SetBatch(True)
  26. from utils.components import unique_name_ext_re, findRootObjByName
  27. def getEfficiencyHistoNames():
  28. return ["p", "pt", "phi", "eta", "nPV"]
  29. def getTrackers(trackers):
  30. return trackers
  31. # data/resolutions_and_effs_Bd2KstEE_MDmaster.root:Track/...
  32. def getOriginFolders():
  33. basedict = {
  34. "Velo": {},
  35. "Upstream": {},
  36. "Forward": {},
  37. "Match": {},
  38. "BestLong": {},
  39. "Seed": {},
  40. "MergedMatch": {},
  41. "DefaultMatch": {},
  42. }
  43. # evtl anpassen wenn die folders anders heissen
  44. basedict["Velo"]["folder"] = "VeloTrackChecker/"
  45. basedict["Upstream"]["folder"] = "UpstreamTrackChecker/"
  46. basedict["Forward"]["folder"] = "ForwardTrackChecker" + unique_name_ext_re() + "/"
  47. basedict["Match"]["folder"] = "MatchTrackChecker" + unique_name_ext_re() + "/"
  48. basedict["BestLong"]["folder"] = "BestLongTrackChecker" + unique_name_ext_re() + "/"
  49. basedict["Seed"]["folder"] = "SeedTrackChecker" + unique_name_ext_re() + "/"
  50. basedict["MergedMatch"]["folder"] = (
  51. "MergedMatchTrackChecker" + unique_name_ext_re() + "/"
  52. )
  53. basedict["DefaultMatch"]["folder"] = (
  54. "DefaultMatchTrackChecker" + unique_name_ext_re() + "/"
  55. )
  56. # basedict["Forward"]["folder"] = "ForwardTrackChecker_7a0dbfa7/"
  57. # basedict["Match"]["folder"] = "MatchTrackChecker_29e3152a/"
  58. # basedict["BestLong"]["folder"] = "BestLongTrackChecker_4ddacce1/"
  59. # basedict["Seed"]["folder"] = "SeedTrackChecker_1b1d5575/"
  60. return basedict
  61. def getTrackNames():
  62. basedict = {
  63. "Velo": {},
  64. "Upstream": {},
  65. "Forward": {},
  66. "Match": {},
  67. "BestLong": {},
  68. "Seed": {},
  69. "MergedMatch": {},
  70. "DefaultMatch": {},
  71. }
  72. basedict["Velo"] = "Velo"
  73. basedict["Upstream"] = "VeloUT"
  74. basedict["Forward"] = "Forward"
  75. basedict["Match"] = "Match"
  76. basedict["BestLong"] = "BestLong"
  77. basedict["Seed"] = "Seed"
  78. basedict["MergedMatch"] = "MergedMatch"
  79. basedict["DefaultMatch"] = "DefaultMatch"
  80. return basedict
  81. def get_colors():
  82. # [kBlack, kGreen + 3, kAzure, kMagenta + 2, kOrange, kCyan + 2]
  83. return [kBlack, kAzure, kGreen + 3, kMagenta + 2, kOrange, kCyan + 2]
  84. def get_elec_colors():
  85. # [kBlack, kGreen + 3, kAzure, kMagenta + 2, kOrange, kCyan + 2]
  86. return [kGray + 2, kBlue - 7, kGreen + 1, kViolet, kOrange - 3, kTeal - 1]
  87. def get_markers():
  88. # [20, 24, 21, 22, 23, 25]
  89. return [20, 21, 24, 25, 22, 23, 26, 32]
  90. def get_fillstyles():
  91. return [1003, 3001, 3002, 3325, 3144, 3244, 3444]
  92. def getGhostHistoNames():
  93. basedict = {
  94. "Velo": {},
  95. "Upstream": {},
  96. "Forward": {},
  97. "Match": {},
  98. "MergedMatch": {},
  99. "DefaultMatch": {},
  100. "BestLong": {},
  101. "Seed": {},
  102. }
  103. basedict["Velo"] = ["eta", "nPV"]
  104. basedict["Upstream"] = ["eta", "p", "pt", "nPV"]
  105. basedict["Forward"] = ["eta", "p", "pt", "nPV"]
  106. basedict["Match"] = ["eta", "p", "pt", "nPV"]
  107. basedict["MergedMatch"] = basedict["Match"]
  108. basedict["DefaultMatch"] = basedict["Match"]
  109. basedict["BestLong"] = ["eta", "p", "pt", "nPV"]
  110. basedict["Seed"] = ["eta", "p", "pt", "nPV"]
  111. return basedict
  112. def argument_parser():
  113. parser = argparse.ArgumentParser(description="location of the tuple file")
  114. parser.add_argument(
  115. "--filename",
  116. type=str,
  117. default=["data/resolutions_and_effs_B.root"],
  118. nargs="+",
  119. help="input files, including path",
  120. )
  121. parser.add_argument(
  122. "--outfile",
  123. type=str,
  124. default="data_results/efficiency_plots.root",
  125. help="output file",
  126. )
  127. parser.add_argument(
  128. "--trackers",
  129. type=str,
  130. nargs="+",
  131. default=["Forward", "Match", "BestLong", "Seed"], # DefaultMatch
  132. help="Trackers to plot.",
  133. )
  134. parser.add_argument(
  135. "--label",
  136. nargs="+",
  137. default=["EffChecker"],
  138. help="label for files",
  139. )
  140. parser.add_argument(
  141. "--savepdf",
  142. action="store_true",
  143. help="save plots in pdf format",
  144. )
  145. parser.add_argument(
  146. "--plot-electrons",
  147. action="store_true",
  148. default=True,
  149. help="plot electrons",
  150. )
  151. parser.add_argument(
  152. "--plot-electrons-only",
  153. action="store_true",
  154. help="plot only electrons",
  155. )
  156. return parser
  157. def get_files(tf, filename, label):
  158. for i, f in enumerate(filename):
  159. tf[label[i]] = TFile(f, "read")
  160. return tf
  161. def get_nicer_var_string(var: str):
  162. nice_vars = dict(pt="p_{T}", eta="#eta", phi="#phi")
  163. try:
  164. return nice_vars[var]
  165. except KeyError:
  166. return var
  167. def get_eff(eff, hist, tf, histoName, label, var):
  168. eff = {}
  169. hist = {}
  170. var = get_nicer_var_string(var)
  171. for i, lab in enumerate(label):
  172. numeratorName = histoName + "_reconstructed"
  173. numerator = findRootObjByName(tf[lab], numeratorName)
  174. denominatorName = histoName + "_reconstructible"
  175. denominator = findRootObjByName(tf[lab], denominatorName)
  176. if numerator.GetEntries() == 0 or denominator.GetEntries() == 0:
  177. continue
  178. teff = TEfficiency(numerator, denominator)
  179. teff.SetStatisticOption(7)
  180. eff[lab] = teff.CreateGraph()
  181. eff[lab].SetName(lab)
  182. eff[lab].SetTitle(lab + " not e^{-}")
  183. if histoName.find("strange") != -1:
  184. eff[lab].SetTitle(lab + " from stranges")
  185. if histoName.find("electron") != -1:
  186. eff[lab].SetTitle(lab + " e^{-}")
  187. hist[lab] = denominator.Clone()
  188. hist[lab].SetName("h_numerator_notElectrons")
  189. hist[lab].SetTitle(var + " distribution, not e^{-}")
  190. if histoName.find("strange") != -1:
  191. hist[lab].SetTitle(var + " distribution, stranges")
  192. if histoName.find("electron") != -1:
  193. hist[lab].SetTitle(var + " distribution, e^{-}")
  194. return eff, hist
  195. def get_ghost(eff, hist, tf, histoName, label):
  196. ghost = {}
  197. for i, lab in enumerate(label):
  198. numeratorName = histoName + "_Ghosts"
  199. denominatorName = histoName + "_Total"
  200. numerator = findRootObjByName(tf[lab], numeratorName)
  201. denominator = findRootObjByName(tf[lab], denominatorName)
  202. print("Numerator = " + numeratorName.replace(unique_name_ext_re(), ""))
  203. print("Denominator = " + denominatorName.replace(unique_name_ext_re(), ""))
  204. teff = TEfficiency(numerator, denominator)
  205. teff.SetStatisticOption(7)
  206. ghost[lab] = teff.CreateGraph()
  207. print(lab)
  208. ghost[lab].SetName(lab)
  209. return ghost
  210. def PrCheckerEfficiency(
  211. filename,
  212. outfile,
  213. label,
  214. trackers,
  215. savepdf,
  216. plot_electrons,
  217. plot_electrons_only,
  218. ):
  219. from utils.LHCbStyle import setLHCbStyle, set_style
  220. from utils.ConfigHistos import (
  221. efficiencyHistoDict,
  222. ghostHistoDict,
  223. categoriesDict,
  224. getCuts,
  225. )
  226. from utils.Legend import place_legend
  227. # from utils.components import unique_name_ext_re, findRootObjByName
  228. setLHCbStyle()
  229. markers = get_markers()
  230. colors = get_colors()
  231. elec_colors = get_elec_colors()
  232. styles = get_fillstyles()
  233. tf = {}
  234. tf = get_files(tf, filename, label)
  235. outputfile = TFile(outfile, "recreate")
  236. latex = TLatex()
  237. latex.SetNDC()
  238. latex.SetTextSize(0.05)
  239. efficiencyHistoDict = efficiencyHistoDict()
  240. efficiencyHistos = getEfficiencyHistoNames()
  241. ghostHistos = getGhostHistoNames()
  242. ghostHistoDict = ghostHistoDict()
  243. categories = categoriesDict()
  244. cuts = getCuts()
  245. trackers = getTrackers(trackers)
  246. folders = getOriginFolders()
  247. # names = getTrackNames()
  248. for tracker in trackers:
  249. outputfile.cd()
  250. trackerDir = outputfile.mkdir(tracker)
  251. trackerDir.cd()
  252. for cut in cuts[tracker]:
  253. cutDir = trackerDir.mkdir(cut)
  254. cutDir.cd()
  255. folder = folders[tracker]["folder"]
  256. print("folder: " + folder.replace(unique_name_ext_re(), ""))
  257. histoBaseName = "Track/" + folder + tracker + "/" + cut + "_"
  258. # calculate efficiency
  259. for histo in efficiencyHistos:
  260. canvastitle = (
  261. "efficiency_" + histo + ", " + categories[tracker][cut]["title"]
  262. )
  263. # get efficiency for not electrons category
  264. histoName = histoBaseName + "" + efficiencyHistoDict[histo]["variable"]
  265. print("not electrons: " + histoName.replace(unique_name_ext_re(), ""))
  266. eff = {}
  267. hist_den = {}
  268. eff, hist_den = get_eff(eff, hist_den, tf, histoName, label, histo)
  269. if categories[tracker][cut]["plotElectrons"] and plot_electrons:
  270. histoNameElec = (
  271. "Track/"
  272. + folder
  273. + tracker
  274. + "/"
  275. + categories[tracker][cut]["Electrons"]
  276. )
  277. histoName_e = (
  278. histoNameElec + "_" + efficiencyHistoDict[histo]["variable"]
  279. )
  280. print("electrons: " + histoName_e.replace(unique_name_ext_re(), ""))
  281. eff_elec = {}
  282. hist_elec = {}
  283. eff_elec, hist_elec = get_eff(
  284. eff_elec,
  285. hist_elec,
  286. tf,
  287. histoName_e,
  288. label,
  289. histo,
  290. )
  291. name = "efficiency_" + histo
  292. canvas = TCanvas(name, canvastitle)
  293. canvas.SetRightMargin(0.1)
  294. mg = TMultiGraph()
  295. for i, lab in enumerate(label):
  296. if not plot_electrons_only:
  297. mg.Add(eff[lab])
  298. set_style(eff[lab], colors[i], markers[i], styles[i])
  299. if categories[tracker][cut]["plotElectrons"] and plot_electrons:
  300. mg.Add(eff_elec[lab])
  301. set_style(eff_elec[lab], elec_colors[i], markers[i], styles[i])
  302. mg.Draw("AP")
  303. mg.GetYaxis().SetRangeUser(0, 1.05)
  304. xtitle = efficiencyHistoDict[histo]["xTitle"]
  305. unit_l = xtitle.split("[")
  306. if "]" in unit_l[-1]:
  307. unit = unit_l[-1].replace("]", "")
  308. else:
  309. unit = "a.u."
  310. print(unit)
  311. mg.GetXaxis().SetTitle(xtitle)
  312. mg.GetXaxis().SetTitleSize(0.06)
  313. mg.GetYaxis().SetTitle(
  314. "Efficiency of Long Tracks",
  315. ) # (" + str(round(hist_den[label[0]].GetBinWidth(1), 2)) + f"{unit})"+"^{-1}")
  316. mg.GetYaxis().SetTitleSize(0.06)
  317. mg.GetYaxis().SetTitleOffset(1.1)
  318. mg.GetXaxis().SetRangeUser(*efficiencyHistoDict[histo]["range"])
  319. mg.GetXaxis().SetNdivisions(10, 5, 0)
  320. mygray = 18
  321. myblue = kBlue - 9
  322. for i, lab in enumerate(label):
  323. rightmax = 1.05 * hist_den[lab].GetMaximum()
  324. scale = gPad.GetUymax() / rightmax
  325. hist_den[lab].Scale(scale)
  326. if categories[tracker][cut]["plotElectrons"] and plot_electrons:
  327. rightmax = 1.05 * hist_elec[lab].GetMaximum()
  328. scale = gPad.GetUymax() / rightmax
  329. hist_elec[lab].Scale(scale)
  330. if i == 0:
  331. if not plot_electrons_only:
  332. set_style(hist_den[lab], mygray, markers[i], styles[i])
  333. gStyle.SetPalette(2, array("i", [mygray - 1, myblue + 1]))
  334. hist_den[lab].Draw("HIST PLC SAME")
  335. if categories[tracker][cut]["plotElectrons"] and plot_electrons:
  336. set_style(hist_elec[lab], myblue, markers[i], styles[i])
  337. hist_elec[lab].SetFillColorAlpha(myblue, 0.35)
  338. hist_elec[lab].Draw("HIST PLC SAME")
  339. # else:
  340. # print(
  341. # "No distribution plotted for other labels.",
  342. # "Can be added by uncommenting the code below this print statement.",
  343. # )
  344. # set_style(hist_den[lab], mygray, markers[i], styles[i])
  345. # gStyle.SetPalette(2, array("i", [mygray - 1, myblue + 1]))
  346. # hist_den[lab].Draw("HIST PLC SAME")
  347. if histo == "p":
  348. pos = [0.5, 0.3, 1.0, 0.5] # [0.53, 0.4, 1.01, 0.71]
  349. elif histo == "pt":
  350. pos = [0.5, 0.3, 0.99, 0.5] # [0.5, 0.4, 0.98, 0.71]
  351. elif histo == "phi":
  352. pos = [0.4, 0.3, 0.9, 0.5]
  353. elif histo == "eta":
  354. pos = [0.5, 0.25, 1.0, 0.45]
  355. else:
  356. pos = [0.35, 0.25, 0.85, 0.45]
  357. legend = place_legend(
  358. canvas, *pos, header="LHCb Simulation", option="LPE"
  359. )
  360. for le in legend.GetListOfPrimitives():
  361. if "distribution" in le.GetLabel():
  362. le.SetOption("LF")
  363. legend.SetTextFont(132)
  364. legend.SetTextSize(0.04)
  365. legend.Draw()
  366. for lab in label:
  367. if not plot_electrons_only:
  368. eff[lab].Draw("P SAME")
  369. if categories[tracker][cut]["plotElectrons"] and plot_electrons:
  370. eff_elec[lab].Draw("P SAME")
  371. cutName = categories[tracker][cut]["title"]
  372. latex.DrawLatex(legend.GetX1() + 0.01, legend.GetY1() - 0.05, cutName)
  373. low = 0
  374. high = 1.05
  375. gPad.Update()
  376. axis = TGaxis(
  377. gPad.GetUxmax(),
  378. gPad.GetUymin(),
  379. gPad.GetUxmax(),
  380. gPad.GetUymax(),
  381. low,
  382. high,
  383. 510,
  384. "+U",
  385. )
  386. axis.SetTitleFont(132)
  387. axis.SetTitleSize(0.06)
  388. axis.SetTitleOffset(0.55)
  389. axis.SetTitle(
  390. "# Tracks " + get_nicer_var_string(histo) + " distribution [a.u.]",
  391. )
  392. axis.SetLabelSize(0)
  393. axis.Draw()
  394. canvas.RedrawAxis()
  395. if savepdf:
  396. filestypes = ["pdf"] # , "png", "eps", "C", "ps", "tex"]
  397. for ftype in filestypes:
  398. if not plot_electrons_only:
  399. canvasName = tracker + "_" + cut + "_" + histo + "." + ftype
  400. else:
  401. canvasName = (
  402. tracker + "Electrons_" + cut + "_" + histo + "." + ftype
  403. )
  404. canvas.SaveAs("checks/" + canvasName)
  405. # canvas.SetRightMargin(0.05)
  406. canvas.Write()
  407. # calculate ghost rate
  408. print("\ncalculate ghost rate: ")
  409. histoBaseName = "Track/" + folder + tracker + "/"
  410. for histo in ghostHistos[tracker]:
  411. trackerDir.cd()
  412. title = "ghost_rate_vs_" + histo
  413. gPad.SetTicks()
  414. histoName = histoBaseName + ghostHistoDict[histo]["variable"]
  415. ghost = {}
  416. hist_den = {}
  417. ghost = get_ghost(ghost, hist_den, tf, histoName, label)
  418. canvas = TCanvas(title, title)
  419. mg = TMultiGraph()
  420. for i, lab in enumerate(label):
  421. mg.Add(ghost[lab])
  422. set_style(ghost[lab], colors[i], markers[2 * i], styles[i])
  423. xtitle = ghostHistoDict[histo]["xTitle"]
  424. mg.GetXaxis().SetTitle(xtitle)
  425. mg.GetYaxis().SetTitle("Fraction of fake tracks")
  426. mg.Draw("ap")
  427. mg.GetXaxis().SetTitleSize(0.06)
  428. mg.GetYaxis().SetTitleSize(0.06)
  429. mg.GetYaxis().SetTitleOffset(1.1)
  430. mg.GetXaxis().SetRangeUser(*efficiencyHistoDict[histo]["range"])
  431. mg.GetXaxis().SetNdivisions(10, 5, 0)
  432. # for lab in label:
  433. # ghost[lab].Draw("P SAME")
  434. if histo == "p":
  435. pos = [0.53, 0.4, 1.00, 0.71]
  436. elif histo == "pt":
  437. pos = [0.5, 0.4, 0.98, 0.71]
  438. elif histo == "eta":
  439. pos = [0.35, 0.6, 0.85, 0.9]
  440. elif histo == "phi":
  441. pos = [0.3, 0.3, 0.9, 0.6]
  442. else:
  443. pos = [0.4, 0.37, 0.80, 0.68]
  444. legend = place_legend(canvas, *pos, header="LHCb Simulation", option="LPE")
  445. legend.SetTextFont(132)
  446. legend.SetTextSize(0.04)
  447. legend.Draw()
  448. # if histo != "nPV":
  449. # latex.DrawLatex(0.7, 0.85, "LHCb simulation")
  450. # else:
  451. # latex.DrawLatex(0.2, 0.85, "LHCb simulation")
  452. # mg.GetYaxis().SetRangeUser(0, 0.4)
  453. if histo == "eta":
  454. mg.GetYaxis().SetRangeUser(0, 0.4)
  455. # track_name = names[tracker] + " tracks"
  456. # latex.DrawLatex(0.7, 0.75, track_name)
  457. # canvas.PlaceLegend()
  458. if savepdf:
  459. filestypes = ["pdf"] # , "png", "eps", "C", "ps", "tex"]
  460. for ftype in filestypes:
  461. canvas.SaveAs(
  462. "checks/" + tracker + "ghost_rate_" + histo + "." + ftype,
  463. )
  464. canvas.Write()
  465. outputfile.Write()
  466. outputfile.Close()
  467. if __name__ == "__main__":
  468. parser = argument_parser()
  469. args = parser.parse_args()
  470. PrCheckerEfficiency(**vars(args))