feat: dRICH benchmarks

Add PID reconstruction benchmarks for the dRICH.

• This follows the modularity paradigm by implementing the benchmarks themselves as independent algorithms; they depend only on the input and output collections
• The executable code benchmark.cc opens the input file(s) and runs the benchmark algorithms; it could be generalized to run any such algorithm

installation tree:

CMAKE_INSTALL_PREFIX
├── bin
│   └── benchmark_rich_reconstruction
├── include
│   └── reconstruction_benchmarks
│       └── rich
│           ├── ChargedParticle.h
│           ├── CherenkovPIDAnalysis.h
│           ├── RawHitAnalysis.h
│           ├── ReconstructedParticleAnalysis.h
│           ├── SimHitAnalysis.h
│           └── Tools.h
└── lib
    └── libreconstruction_benchmarks_rich.so

benchmarks/rich/draw_benchmark.py

+#!/usr/bin/env python
+# Copyright 2023, Christopher Dilks
+# Subject to the terms in the LICENSE file found in the top-level directory.
+
+import ROOT as r
+import sys, getopt, pathlib
+
+# suppress graphics
+r.gROOT.SetBatch(True)
+
+# ARGUMENTS
+################################################################
+
+ana_file_name = 'out/ana.edm4hep.root'
+output_dir    = 'out/ana_plots'
+
+helpStr = f'''
+{sys.argv[0]} [OPTIONS]
+
+    -i <input file>: specify an input file, e.g., hepmc
+       default: {ana_file_name}
+
+    -o <output dir>: specify an output directory
+       default: {output_dir}
+    
+    -h: show this usage guide
+
+    '''
+
+try:
+    opts, args = getopt.getopt(sys.argv[1:], 'i:o:h')
+except getopt.GetoptError:
+    print('\n\nERROR: invalid argument\n', helpStr, file=sys.stderr)
+    sys.exit(2)
+for opt, arg in opts:
+    if(opt == '-i'): ana_file_name = arg.lstrip()
+    if(opt == '-o'): output_dir    = arg.lstrip()
+    if(opt == '-h'):
+        print(helpStr)
+        sys.exit(2)
+print(f'''
+ana_file_name = {ana_file_name}
+output_dir    = {output_dir}
+''')
+
+
+# PLOTTING
+################################################################
+
+# make canvases
+ana_file = r.TFile.Open(ana_file_name, "READ")
+def make_canv(name, dimx=1200, dimy=800):
+    return r.TCanvas(name, name, dimx, dimy)
+canv_dict = {
+    "photon_spectra": make_canv("photon_spectra"),
+    "digitization":   make_canv("digitization"),
+    "pidAerogel":     make_canv("pidAerogel"),
+    "pidGas":         make_canv("pidGas"),
+    "pidMerged":      make_canv("pidMerged"),
+}
+
+# draw photon spectra
+canv = canv_dict["photon_spectra"]
+canv.Divide(1,2)
+for i in range(2):
+    canv.GetPad(i+1).SetGrid(1,1)
+    canv.GetPad(i+1).SetLogy()
+canv.cd(1)
+ana_file.Get("phot/phot_spectrum_sim").Draw()
+canv.cd(2)
+ana_file.Get("digi/phot_spectrum_rec").Draw()
+
+# draw digitization
+canv = canv_dict["digitization"]
+canv.Divide(2,2)
+for i in range(4):
+    canv.GetPad(i+1).SetGrid(1,1)
+    if(i<2):
+        canv.GetPad(i+1).SetLogy()
+    else:
+        canv.GetPad(i+1).SetLogz()
+canv.cd(1)
+ana_file.Get("digi/adc_dist").Draw()
+canv.cd(2)
+ana_file.Get("digi/tdc_dist").Draw()
+canv.cd(3)
+ana_file.Get("digi/tdc_vs_adc").Draw("COLZ")
+
+# draw CherenkovPID
+for rad in ["Aerogel", "Gas", "Merged"]:
+    pid_name = f'pid{rad}'
+    canv = canv_dict[pid_name]
+    canv.Divide(2,2)
+    for i in range(4):
+        canv.GetPad(i+1).SetGrid(1,1)
+    canv.cd(1)
+    ana_file.Get(f'{pid_name}/npe_dist_{rad}').Draw()
+    canv.cd(2)
+    ana_file.Get(f'{pid_name}/theta_dist_{rad}').Draw()
+    canv.cd(3)
+    ana_file.Get(f'{pid_name}/thetaResid_dist_{rad}').Draw()
+    canv.cd(4)
+    ana_file.Get(f'{pid_name}/highestWeight_dist_{rad}').Draw()
+
+
+# FINISH
+################################################################
+
+pathlib.Path(output_dir).mkdir(parents=True, exist_ok=True)
+for name, canvas in canv_dict.items():
+    canvas.SaveAs(f'{output_dir}/{name}.png')
+ana_file.Close()