feat: photon and hit multiplicities

benchmarks/rich/draw_benchmark.py

@@ -61,14 +61,18 @@ canv_dict = {
 
 # draw photon spectra
 canv = canv_dict["photon_spectra"]
-canv.Divide(1,2)
-for i in range(2):
+canv.Divide(2,2)
+for i in range(4):
     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("phot/nphot_dist").Draw()
+canv.cd(3)
 ana_file.Get("digi/phot_spectrum_rec").Draw()
+canv.cd(4)
+ana_file.Get("digi/nhits_dist").Draw()
 
 # draw digitization
 canv = canv_dict["digitization"]

benchmarks/rich/include/RawHitAnalysis.h

@@ -42,6 +42,7 @@ namespace benchmarks {
       TH1D *m_tdc_dist;
       TH2D *m_tdc_vs_adc;
       TH1D *m_phot_spectrum;
+      TH1D *m_nhits_dist;
 
       // logging
       std::shared_ptr<spdlog::logger> m_log;

benchmarks/rich/include/SimHitAnalysis.h

@@ -33,6 +33,7 @@ namespace benchmarks {
     private:
 
       // histograms
+      TH1D *m_nphot_dist;
       TH2D *m_nphot_vs_p;
       TH1D *m_nphot_vs_p__transient; // transient (not written)
       TH1D *m_phot_spectrum;

benchmarks/rich/src/RawHitAnalysis.cc

@@ -22,9 +22,12 @@ namespace benchmarks {
         );
     m_phot_spectrum = new TH1D(
         "phot_spectrum_rec",
-        "Photon wavelength at emission point -- only photons from digitized hits;#lambda [nm]",
+        "Photon wavelength for digitized hits;#lambda [nm], at emission point",
         Tools::wl_bins, Tools::wl_min, Tools::wl_max
         );
+    m_nhits_dist = new TH1D("nhits_dist", "Number of digitized hits;N_{hits}",
+        Tools::npe_max, 0, Tools::npe_max
+        );
 
     // format histograms
     auto format1D = [] (auto h) {
@@ -35,6 +38,8 @@ namespace benchmarks {
     format1D(m_tdc_dist);
     m_phot_spectrum->SetLineColor(kViolet+2);
     m_phot_spectrum->SetFillColor(kViolet+2);
+    m_nhits_dist->SetLineColor(kCyan+2);
+    m_nhits_dist->SetFillColor(kCyan+2);
   }
 
 
@@ -46,6 +51,9 @@ namespace benchmarks {
       )
   {
 
+    // fill nhits
+    m_nhits_dist->Fill(raw_hits.size());
+
     // loop over all raw hits (including noise)
     for(const auto& raw_hit : raw_hits) {
       auto adc = raw_hit.getCharge();

benchmarks/rich/src/SimHitAnalysis.cc

@@ -14,7 +14,10 @@ namespace benchmarks {
     m_log = logger;
 
     // initialize histograms
-    m_nphot_vs_p = new TH2D("nphot_vs_p", "N_{photons} vs. Thrown Momentum;p [GeV];N_{photons}",
+    m_nphot_dist = new TH1D("nphot_dist", "Number of incident photons;N_{photons}",
+        Tools::nphot_max, 0, Tools::nphot_max
+        );
+    m_nphot_vs_p = new TH2D("nphot_vs_p", "Number of incident photons vs. Thrown Momentum;p [GeV];N_{photons}",
         Tools::momentum_bins, 0, Tools::momentum_max,
         Tools::nphot_max,     0, Tools::nphot_max
         );
@@ -25,11 +28,13 @@ namespace benchmarks {
         );
     m_phot_spectrum = new TH1D(
         "phot_spectrum_sim",
-        "Photon wavelength at emission point -- all photons incident on sensors;#lambda [nm]",
+        "Incident photon wavelength;#lambda [nm], at emission point",
         Tools::wl_bins, Tools::wl_min, Tools::wl_max
         );
     m_phot_spectrum->SetLineColor(kViolet+2);
     m_phot_spectrum->SetFillColor(kViolet+2);
+    m_nphot_dist->SetLineColor(kCyan+2);
+    m_nphot_dist->SetFillColor(kCyan+2);
   }
 
 
@@ -46,6 +51,9 @@ namespace benchmarks {
     part->SetSingleParticle(mc_parts);
     auto thrown_momentum = part->GetMomentum();
 
+    // fill nphot distribution
+    m_nphot_dist->Fill(sim_hits.size());
+
     // loop over hits
     for(const auto& sim_hit : sim_hits) {