From 7e5aecf48bbf1c812e23aa75f3dc9b14a55bc315 Mon Sep 17 00:00:00 2001
From: Wouter Deconinck <wdconinc@gmail.com>
Date: Sun, 5 Jun 2022 23:36:55 +0000
Subject: [PATCH] Reconstruction benchmark fixes for ecce:main

---
 benchmarks/clustering/full_cal_clusters.sh    |   2 +-
 .../clustering/options/full_cal_reco.py       | 148 +++++++++++++-----
 2 files changed, 113 insertions(+), 37 deletions(-)

diff --git a/benchmarks/clustering/full_cal_clusters.sh b/benchmarks/clustering/full_cal_clusters.sh
index 9d7f8864..6818e9d1 100644
--- a/benchmarks/clustering/full_cal_clusters.sh
+++ b/benchmarks/clustering/full_cal_clusters.sh
@@ -130,7 +130,7 @@ fi
 # Run analysis scripts
 FULL_CAL_SCRIPT_DIR=benchmarks/clustering/scripts
 python ${FULL_CAL_SCRIPT_DIR}/cluster_plots.py ${JUGGLER_SIM_FILE} ${JUGGLER_REC_FILE} -o results \
-    --collections "EcalEndcapNClusters, EcalEndcapPClusters, EcalBarrelImagingClusters,
+    --collections "EcalEndcapNClusters, EcalEndcapPClusters, EcalBarrelClusters, EcalBarrelImagingClusters,
                    HcalEndcapNClusters, HcalEndcapPClusters, HcalBarrelClusters"
 
 root_filesize=$(stat --format=%s "${JUGGLER_REC_FILE}")
diff --git a/benchmarks/clustering/options/full_cal_reco.py b/benchmarks/clustering/options/full_cal_reco.py
index 9dfd7a15..5565c40b 100644
--- a/benchmarks/clustering/options/full_cal_reco.py
+++ b/benchmarks/clustering/options/full_cal_reco.py
@@ -10,6 +10,8 @@ from Configurables import ApplicationMgr, EICDataSvc, PodioInput, PodioOutput, G
 from GaudiKernel.SystemOfUnits import MeV, GeV, mm, cm, mrad
 
 detector_name = str(os.environ.get("JUGGLER_DETECTOR", "athena"))
+detector_version = str(os.environ.get("JUGGLER_DETECTOR", "master"))
+
 detector_path = str(os.environ.get("DETECTOR_PATH", "."))
 compact_path = os.path.join(detector_path, detector_name)
 
@@ -60,8 +62,6 @@ sim_coll = [
     "EcalEndcapPHitsContributions",
     "EcalBarrelHits",
     "EcalBarrelHitsContributions",
-    "EcalBarrelScFiHits",
-    "EcalBarrelScFiHitsContributions",
     "HcalBarrelHits",
     "HcalBarrelHitsContributions",
     "HcalEndcapPHits",
@@ -70,10 +70,18 @@ sim_coll = [
     "HcalEndcapNHitsContributions",
 ]
 
-# input and output
-podin = PodioInput("PodioReader", collections=sim_coll)
-podout = PodioOutput("out", filename=output_rec)
+if 'athena' in detector_name:
+    sim_coll += [
+        "EcalBarrelScFiHits",
+        "EcalBarrelScFiHitsContributions",
+    ]
 
+# list of algorithms
+algs = []
+
+# input
+podin = PodioInput("PodioReader", collections=sim_coll)
+algs.append(podin)
 
 # Crystal Endcap Ecal
 ce_ecal_daq = dict(
@@ -87,6 +95,7 @@ ce_ecal_digi = CalHitDigi("ce_ecal_digi",
         outputHitCollection="EcalEndcapNHitsDigi",
         energyResolutions=[0., 0.02, 0.],
         **ce_ecal_daq)
+algs.append(ce_ecal_digi)
 
 ce_ecal_reco = CalHitReco("ce_ecal_reco",
         inputHitCollection=ce_ecal_digi.outputHitCollection,
@@ -96,6 +105,7 @@ ce_ecal_reco = CalHitReco("ce_ecal_reco",
         sectorField="sector",
         samplingFraction=0.998,      # this accounts for a small fraction of leakage
         **ce_ecal_daq)
+algs.append(ce_ecal_reco)
 
 ce_ecal_cl = IslandCluster("ce_ecal_cl",
         # OutputLevel=DEBUG,
@@ -106,11 +116,13 @@ ce_ecal_cl = IslandCluster("ce_ecal_cl",
         minClusterCenterEdep=30*MeV,
         sectorDist=5.0*cm,
         dimScaledLocalDistXY=[1.8, 1.8])          # dimension scaled dist is good for hybrid sectors with different module size
+algs.append(ce_ecal_cl)
 
 ce_ecal_clreco = RecoCoG("ce_ecal_clreco",
         inputProtoClusterCollection=ce_ecal_cl.outputProtoClusterCollection,
         outputClusterCollection="EcalEndcapNClusters",
         logWeightBase=4.6)
+algs.append(ce_ecal_clreco)
 
 # Endcap Sampling Ecal
 ci_ecal_daq = dict(
@@ -123,6 +135,7 @@ ci_ecal_digi = CalHitDigi("ci_ecal_digi",
         inputHitCollection="EcalEndcapPHits",
         outputHitCollection="EcalEndcapPHitsDigi",
         **ci_ecal_daq)
+algs.append(ci_ecal_digi)
 
 ci_ecal_reco = CalHitReco("ci_ecal_reco",
         inputHitCollection=ci_ecal_digi.outputHitCollection,
@@ -130,6 +143,7 @@ ci_ecal_reco = CalHitReco("ci_ecal_reco",
         thresholdFactor=5.0,
         samplingFraction=ci_ecal_sf,
         **ci_ecal_daq)
+algs.append(ci_ecal_reco)
 
 # merge hits in different layer (projection to local x-y plane)
 ci_ecal_merger = CalHitsMerger("ci_ecal_merger",
@@ -141,6 +155,7 @@ ci_ecal_merger = CalHitsMerger("ci_ecal_merger",
         fields=["fiber_x", "fiber_y"],
         fieldRefNumbers=[1, 1],
         readoutClass="EcalEndcapPHits")
+algs.append(ci_ecal_merger)
 
 ci_ecal_cl = IslandCluster("ci_ecal_cl",
         # OutputLevel=DEBUG,
@@ -149,26 +164,31 @@ ci_ecal_cl = IslandCluster("ci_ecal_cl",
         splitCluster=False,
         minClusterCenterEdep=10.*MeV,
         localDistXY=[10*mm, 10*mm])
+algs.append(ci_ecal_cl)
 
 ci_ecal_clreco = RecoCoG("ci_ecal_clreco",
         inputProtoClusterCollection=ci_ecal_cl.outputProtoClusterCollection,
         outputClusterCollection="EcalEndcapPClusters",
         logWeightBase=6.2)
+algs.append(ci_ecal_clreco)
 
-# Central Barrel Ecal (Imaging Cal.)
-cb_ecal_daq = dict(
+# Central Barrel Ecal
+if 'athena' in detector_name:
+    # Imaging calorimeter
+    cb_ecal_daq = dict(
         dynamicRangeADC=3*MeV,
         capacityADC=8192,
         pedestalMean=400,
         pedestalSigma=20)   # about 6 keV
 
-cb_ecal_digi = CalHitDigi("cb_ecal_digi",
+    cb_ecal_digi = CalHitDigi("cb_ecal_digi",
         inputHitCollection="EcalBarrelHits",
         outputHitCollection="EcalBarrelImagingHitsDigi",
         energyResolutions=[0., 0.02, 0.],   # 2% flat resolution
         **cb_ecal_daq)
+    algs.append(cb_ecal_digi)
 
-cb_ecal_reco = ImCalPixelReco("cb_ecal_reco",
+    cb_ecal_reco = ImCalPixelReco("cb_ecal_reco",
         inputHitCollection=cb_ecal_digi.outputHitCollection,
         outputHitCollection="EcalBarrelImagingHitsReco",
         thresholdFactor=3,  # about 20 keV
@@ -177,35 +197,79 @@ cb_ecal_reco = ImCalPixelReco("cb_ecal_reco",
         sectorField="module",           # field to get sector id
         samplingFraction=cb_ecal_sf,
         **cb_ecal_daq)
+    algs.append(cb_ecal_reco)
 
-cb_ecal_cl = ImagingCluster("cb_ecal_cl",
+    cb_ecal_cl = ImagingCluster("cb_ecal_cl",
         inputHitCollection=cb_ecal_reco.outputHitCollection,
         outputProtoClusterCollection="EcalBarrelImagingProtoClusters",
         localDistXY=[2.*mm, 2*mm],              # same layer
         layerDistEtaPhi=[10*mrad, 10*mrad],     # adjacent layer
         neighbourLayersRange=2,                 # id diff for adjacent layer
         sectorDist=3.*cm)                       # different sector
+    algs.append(cb_ecal_cl)
 
-cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
+    cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
         inputProtoClusters=cb_ecal_cl.outputProtoClusterCollection,
         mcHits="EcalBarrelHits",
         outputClusters="EcalBarrelImagingClusters",
         outputLayers="EcalBarrelImagingLayers")
+    algs.append(cb_ecal_clreco)
+else:
+    # SciGlass calorimeter
+    cb_ecal_daq = dict(
+        dynamicRangeADC=5.*GeV,
+        capacityADC=32768,
+        pedestalMean=400,
+        pedestalSigma=3)
+
+    cb_ecal_digi = CalHitDigi("cb_ecal_digi",
+        inputHitCollection="EcalBarrelHits",
+        outputHitCollection="EcalBarrelHitsDigi",
+        energyResolutions=[0., 0.02, 0.],   # 2% flat resolution
+        **cb_ecal_daq)
+    algs.append(cb_ecal_digi)
+
+    cb_ecal_reco = CalHitReco("cb_ecal_reco",
+        inputHitCollection=cb_ecal_digi.outputHitCollection,
+        outputHitCollection="EcalBarrelHitsReco",
+        thresholdFactor=3,  # about 20 keV
+        readoutClass="EcalBarrelHits",  # readout class
+        sectorField="sector",           # field to get sector id
+        samplingFraction=0.998,         # this accounts for a small fraction of leakage
+        **cb_ecal_daq)
+    algs.append(cb_ecal_reco)
 
-#Central ECAL SciFi
-# use the same daq_setting for digi/reco pair
-scfi_barrel_daq = dict(
+    cb_ecal_cl = IslandCluster("cb_ecal_cl",
+        inputHitCollection=cb_ecal_reco.outputHitCollection,
+        outputProtoClusterCollection="EcalBarrelProtoClusters",
+        splitCluster=False,
+        minClusterHitEdep=1.0*MeV,  # discard low energy hits
+        minClusterCenterEdep=30*MeV,
+        sectorDist=5.0*cm)
+    algs.append(cb_ecal_cl)
+
+    cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
+        inputProtoClusters=cb_ecal_cl.outputProtoClusterCollection,
+        mcHits="EcalBarrelHits",
+        outputClusters="EcalBarrelClusters",
+        outputLayers="EcalBarrelLayers")
+    algs.append(cb_ecal_clreco)
+
+# Central Barrel Ecal SciFi
+if 'athena' in detector_name:
+    scfi_barrel_daq = dict(
         dynamicRangeADC=50.*MeV,
         capacityADC=32768,
         pedestalMean=400,
         pedestalSigma=10)
 
-scfi_barrel_digi = CalHitDigi("scfi_barrel_digi",
+    scfi_barrel_digi = CalHitDigi("scfi_barrel_digi",
         inputHitCollection="EcalBarrelScFiHits",
         outputHitCollection="EcalBarrelScFiHitsDigi",
         **scfi_barrel_daq)
+    algs.append(scfi_barrel_digi)
 
-scfi_barrel_reco = CalHitReco("scfi_barrel_reco",
+    scfi_barrel_reco = CalHitReco("scfi_barrel_reco",
         inputHitCollection=scfi_barrel_digi.outputHitCollection,
         outputHitCollection="EcalBarrelScFiHitsReco",
         thresholdFactor=5.0,
@@ -215,29 +279,32 @@ scfi_barrel_reco = CalHitReco("scfi_barrel_reco",
         localDetFields=["system", "module"], # use local coordinates in each module (stave)
         samplingFraction=scifi_barrel_sf,
         **scfi_barrel_daq)
+    algs.append(scfi_barrel_reco)
 
-# merge hits in different layer (projection to local x-y plane)
-scfi_barrel_merger = CalHitsMerger("scfi_barrel_merger",
+    # merge hits in different layer (projection to local x-y plane)
+    scfi_barrel_merger = CalHitsMerger("scfi_barrel_merger",
         # OutputLevel=DEBUG,
         inputHitCollection=scfi_barrel_reco.outputHitCollection,
         outputHitCollection="EcalBarrelScFiGridReco",
         fields=["fiber"],
         fieldRefNumbers=[1],
         readoutClass="EcalBarrelScFiHits")
+    algs.append(scfi_barrel_merger)
 
-scfi_barrel_cl = IslandCluster("scfi_barrel_cl",
+    scfi_barrel_cl = IslandCluster("scfi_barrel_cl",
         # OutputLevel=DEBUG,
         inputHitCollection=scfi_barrel_merger.outputHitCollection,
         outputProtoClusterCollection="EcalBarrelScFiProtoClusters",
         splitCluster=False,
         minClusterCenterEdep=10.*MeV,
         localDistXZ=[30*mm, 30*mm])
+    algs.append(scfi_barrel_cl)
 
-scfi_barrel_clreco = RecoCoG("scfi_barrel_clreco",
-       inputProtoClusterCollection=scfi_barrel_cl.outputProtoClusterCollection,
-       outputClusterCollection="EcalBarrelScFiClusters",
-       logWeightBase=6.2)
-
+    scfi_barrel_clreco = RecoCoG("scfi_barrel_clreco",
+        inputProtoClusterCollection=scfi_barrel_cl.outputProtoClusterCollection,
+        outputClusterCollection="EcalBarrelScFiClusters",
+        logWeightBase=6.2)
+    algs.append(scfi_barrel_clreco)
 
 # Central Barrel Hcal
 cb_hcal_daq = dict(
@@ -250,6 +317,7 @@ cb_hcal_digi = CalHitDigi("cb_hcal_digi",
          inputHitCollection="HcalBarrelHits",
          outputHitCollection="HcalBarrelHitsDigi",
          **cb_hcal_daq)
+algs.append(cb_hcal_digi)
 
 cb_hcal_reco = CalHitReco("cb_hcal_reco",
         inputHitCollection=cb_hcal_digi.outputHitCollection,
@@ -260,6 +328,7 @@ cb_hcal_reco = CalHitReco("cb_hcal_reco",
         sectorField="module",
         samplingFraction=cb_hcal_sf,
         **cb_hcal_daq)
+algs.append(cb_hcal_reco)
 
 cb_hcal_merger = CalHitsMerger("cb_hcal_merger",
         inputHitCollection=cb_hcal_reco.outputHitCollection,
@@ -267,6 +336,7 @@ cb_hcal_merger = CalHitsMerger("cb_hcal_merger",
         readoutClass="HcalBarrelHits",
         fields=["layer", "slice"],
         fieldRefNumbers=[1, 0])
+algs.append(cb_hcal_merger)
 
 cb_hcal_cl = IslandCluster("cb_hcal_cl",
         inputHitCollection=cb_hcal_merger.outputHitCollection,
@@ -274,12 +344,13 @@ cb_hcal_cl = IslandCluster("cb_hcal_cl",
         splitCluster=False,
         minClusterCenterEdep=30.*MeV,
         localDistXY=[15.*cm, 15.*cm])
+algs.append(cb_hcal_cl)
 
 cb_hcal_clreco = RecoCoG("cb_hcal_clreco",
         inputProtoClusterCollection=cb_hcal_cl.outputProtoClusterCollection,
         outputClusterCollection="HcalBarrelClusters",
         logWeightBase=6.2)
-
+algs.append(cb_hcal_clreco)
 
 # Hcal Hadron Endcap
 ci_hcal_daq = dict(
@@ -287,10 +358,12 @@ ci_hcal_daq = dict(
          capacityADC=32768,
          pedestalMean=400,
          pedestalSigma=10)
+
 ci_hcal_digi = CalHitDigi("ci_hcal_digi",
          inputHitCollection="HcalEndcapPHits",
          outputHitCollection="HcalEndcapPHitsDigi",
          **ci_hcal_daq)
+algs.append(ci_hcal_digi)
 
 ci_hcal_reco = CalHitReco("ci_hcal_reco",
         inputHitCollection=ci_hcal_digi.outputHitCollection,
@@ -298,6 +371,7 @@ ci_hcal_reco = CalHitReco("ci_hcal_reco",
         thresholdFactor=5.0,
         samplingFraction=ci_hcal_sf,
         **ci_hcal_daq)
+algs.append(ci_hcal_reco)
 
 ci_hcal_merger = CalHitsMerger("ci_hcal_merger",
         inputHitCollection=ci_hcal_reco.outputHitCollection,
@@ -305,6 +379,7 @@ ci_hcal_merger = CalHitsMerger("ci_hcal_merger",
         readoutClass="HcalEndcapPHits",
         fields=["layer", "slice"],
         fieldRefNumbers=[1, 0])
+algs.append(ci_hcal_merger)
 
 ci_hcal_cl = IslandCluster("ci_hcal_cl",
         inputHitCollection=ci_hcal_merger.outputHitCollection,
@@ -312,11 +387,13 @@ ci_hcal_cl = IslandCluster("ci_hcal_cl",
         splitCluster=False,
         minClusterCenterEdep=30.*MeV,
         localDistXY=[15.*cm, 15.*cm])
+algs.append(ci_hcal_cl)
 
 ci_hcal_clreco = RecoCoG("ci_hcal_clreco",
         inputProtoClusterCollection=ci_hcal_cl.outputProtoClusterCollection,
         outputClusterCollection="HcalEndcapPClusters",
         logWeightBase=6.2)
+algs.append(ci_hcal_clreco)
 
 # Hcal Electron Endcap
 ce_hcal_daq = dict(
@@ -329,6 +406,7 @@ ce_hcal_digi = CalHitDigi("ce_hcal_digi",
         inputHitCollection="HcalEndcapNHits",
         outputHitCollection="HcalEndcapNHitsDigi",
         **ce_hcal_daq)
+algs.append(ce_hcal_digi)
 
 ce_hcal_reco = CalHitReco("ce_hcal_reco",
         inputHitCollection=ce_hcal_digi.outputHitCollection,
@@ -336,6 +414,7 @@ ce_hcal_reco = CalHitReco("ce_hcal_reco",
         thresholdFactor=5.0,
         samplingFraction=ce_hcal_sf,
         **ce_hcal_daq)
+algs.append(ce_hcal_reco)
 
 ce_hcal_merger = CalHitsMerger("ce_hcal_merger",
         inputHitCollection=ce_hcal_reco.outputHitCollection,
@@ -343,6 +422,7 @@ ce_hcal_merger = CalHitsMerger("ce_hcal_merger",
         readoutClass="HcalEndcapNHits",
         fields=["layer", "slice"],
         fieldRefNumbers=[1, 0])
+algs.append(ce_hcal_merger)
 
 ce_hcal_cl = IslandCluster("ce_hcal_cl",
         inputHitCollection=ce_hcal_merger.outputHitCollection,
@@ -350,32 +430,28 @@ ce_hcal_cl = IslandCluster("ce_hcal_cl",
         splitCluster=False,
         minClusterCenterEdep=30.*MeV,
         localDistXY=[15.*cm, 15.*cm])
+algs.append(ce_hcal_cl)
 
 ce_hcal_clreco = RecoCoG("ce_hcal_clreco",
         inputProtoClusterCollection=ce_hcal_cl.outputProtoClusterCollection,
         outputClusterCollection="HcalEndcapNClusters",
         logWeightBase=6.2)
+algs.append(ce_hcal_clreco)
 
-
+# output
+podout = PodioOutput("out", filename=output_rec)
 podout.outputCommands = ['drop *',
         'keep MCParticles',
         'keep *Digi',
         'keep *Reco*',
         'keep *Cluster*',
         'keep *Layers']
+algs.append(podout)
 
 ApplicationMgr(
-    TopAlg = [podin,
-              ce_ecal_digi, ce_ecal_reco, ce_ecal_cl, ce_ecal_clreco,
-              ci_ecal_digi, ci_ecal_reco, ci_ecal_merger, ci_ecal_cl, ci_ecal_clreco,
-              cb_ecal_digi, cb_ecal_reco, cb_ecal_cl, cb_ecal_clreco,
-              scfi_barrel_digi, scfi_barrel_reco, scfi_barrel_merger, scfi_barrel_cl, scfi_barrel_clreco,
-              cb_hcal_digi, cb_hcal_reco, cb_hcal_merger, cb_hcal_cl, cb_hcal_clreco,
-              ce_hcal_digi, ce_hcal_reco, ce_hcal_merger, ce_hcal_cl, ce_hcal_clreco,
-              ci_hcal_digi, ci_hcal_reco, ci_hcal_merger, ci_hcal_cl, ci_hcal_clreco,
-              podout],
+    TopAlg = algs,
     EvtSel = 'NONE',
     EvtMax = n_events,
     ExtSvc = [podioevent],
-    OutputLevel=DEBUG
+    OutputLevel=WARNING
 )
-- 
GitLab