From c8d5494085a8d1eb4755301a2a3dda246e5b07ca Mon Sep 17 00:00:00 2001
From: Miguel Arratia <miguela@ucr.edu>
Date: Sat, 5 Jun 2021 20:32:17 +0000
Subject: [PATCH] adding fullcalo_clustering, which is based on
 calo_clustering.py but adds HCal hits and clusters

---
 benchmarks/clustering/barrel_clusters.sh      |  22 +--
 .../clustering/options/fullcalo_clustering.py | 135 ++++++++++++++++++
 2 files changed, 146 insertions(+), 11 deletions(-)
 create mode 100644 benchmarks/clustering/options/fullcalo_clustering.py

diff --git a/benchmarks/clustering/barrel_clusters.sh b/benchmarks/clustering/barrel_clusters.sh
index 4223d3b6..6e701cc3 100644
--- a/benchmarks/clustering/barrel_clusters.sh
+++ b/benchmarks/clustering/barrel_clusters.sh
@@ -59,21 +59,21 @@ if [[ "$?" -ne "0" ]] ; then
 fi
 
 # Need to figure out how to pass file name to juggler from the commandline
-xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv gaudirun.py benchmarks/clustering/options/calorimeter_clustering.py
+xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv gaudirun.py benchmarks/clustering/options/fullcaloclustering.py
 if [[ "$?" -ne "0" ]] ; then
   echo "ERROR running juggler"
   exit 1
 fi
 
-pwd
-mkdir -p results/clustering
+#pwd
+#mkdir -p results/clustering
 
-root -b -q "benchmarks/clustering/scripts/barrel_clusters.cxx(\"${JUGGLER_REC_FILE}\")"
+#root -b -q "benchmarks/clustering/scripts/barrel_clusters.cxx(\"${JUGGLER_REC_FILE}\")"
 
-root_filesize=$(stat --format=%s "${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}")
-if [[ "${JUGGLER_N_EVENTS}" -lt "500" ]] ; then 
-  # file must be less than 10 MB to upload
-  if [[ "${root_filesize}" -lt "10000000" ]] ; then 
-    cp ${JUGGLER_REC_FILE} results/clustering/.
-  fi
-fi
+#root_filesize=$(stat --format=%s "${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}")
+#if [[ "${JUGGLER_N_EVENTS}" -lt "500" ]] ; then 
+#  # file must be less than 10 MB to upload
+#  if [[ "${root_filesize}" -lt "10000000" ]] ; then 
+#    cp ${JUGGLER_REC_FILE} results/clustering/.
+#  fi
+#fi
diff --git a/benchmarks/clustering/options/fullcalo_clustering.py b/benchmarks/clustering/options/fullcalo_clustering.py
new file mode 100644
index 00000000..d2aa4e68
--- /dev/null
+++ b/benchmarks/clustering/options/fullcalo_clustering.py
@@ -0,0 +1,135 @@
+from Gaudi.Configuration import *
+import os
+  
+from GaudiKernel.DataObjectHandleBase import DataObjectHandleBase
+from Configurables import ApplicationMgr, EICDataSvc, PodioOutput, GeoSvc
+from GaudiKernel import SystemOfUnits as units
+
+detector_name = "topside"
+if "JUGGLER_DETECTOR" in os.environ :
+  detector_name = str(os.environ["JUGGLER_DETECTOR"])
+
+if "JUGGLER_DETECTOR_PATH" in os.environ :
+  detector_name = str(os.environ["JUGGLER_DETECTOR_PATH"]) + "/" + detector_name
+
+# get sampling fraction from system environment variable, 1.0 by default
+sf = 1.0
+if "CB_EMCAL_SAMP_FRAC" in os.environ :
+  sf = str(os.environ["CB_EMCAL_SAMP_FRAC"])
+
+# todo add checks
+input_sim_file  = str(os.environ["JUGGLER_SIM_FILE"])
+output_rec_file = str(os.environ["JUGGLER_REC_FILE"])
+n_events = str(os.environ["JUGGLER_N_EVENTS"])
+
+geo_service  = GeoSvc("GeoSvc", detectors=["{}.xml".format(detector_name)])
+podioevent   = EICDataSvc("EventDataSvc", inputs=[input_sim_file])
+
+from Configurables import PodioInput
+from Configurables import Jug__Base__InputCopier_dd4pod__Geant4ParticleCollection_dd4pod__Geant4ParticleCollection_ as MCCopier
+from Configurables import Jug__Base__InputCopier_dd4pod__CalorimeterHitCollection_dd4pod__CalorimeterHitCollection_ as CalCopier
+
+from Configurables import Jug__Digi__HadronicCalDigi as HadCalorimeterDigi
+from Configurables import Jug__Digi__EMCalorimeterDigi as EMCalorimeterDigi
+from Configurables import Jug__Digi__CrystalEndcapsDigi as CrystalEndcapsDigi
+
+from Configurables import Jug__Reco__CrystalEndcapsReco as CrystalEndcapsReco
+from Configurables import Jug__Reco__EMCalReconstruction as EMCalReconstruction
+from Configurables import Jug__Reco__HCalReconstruction as HCalReconstruction
+
+from Configurables import Jug__Reco__SimpleClustering as SimpleClustering
+from Configurables import Jug__Reco__CalorimeterIslandCluster as IslandCluster
+from Configurables import Jug__Reco__ClusterRecoCoG as RecoCoG
+
+podioinput = PodioInput("PodioReader", collections=["mcparticles","CrystalEcalHits","EcalBarrelHits","HcalBarrelHits"], OutputLevel=DEBUG)
+
+## copiers to get around input --> output copy bug. Note the "2" appended to the output collection.
+copier = MCCopier("MCCopier", 
+        inputCollection="mcparticles", 
+        outputCollection="mcparticles2") 
+calcopier = CalCopier("CalCopier", 
+        inputCollection="CrystalEcalHits", 
+        outputCollection="CrystalEcalHits2") 
+
+
+##raw hits
+emcaldigi = CrystalEndcapsDigi("ecal_digi",
+        inputHitCollection="CrystalEcalHits", 
+        outputHitCollection="RawDigiEcalHits")
+ecdigi = EMCalorimeterDigi("ec_barrel_digi", 
+        inputHitCollection="EcalBarrelHits", 
+        outputHitCollection="RawEcalBarrelHits")
+
+hcaldigi = HadCalorimeterDigi("hcal_barrel_digi",
+         inputHitCollection="HcalBarrelHits",
+         outputHitCollection="RawHcalBarrelHits")
+
+
+
+#digitized hits
+crystal_ec_reco = CrystalEndcapsReco("crystal_ec_reco", 
+        inputHitCollection="RawDigiEcalHits", 
+        outputHitCollection="RecoEcalHits",
+        minModuleEdep=1.0*units.MeV)
+
+ecal_reco = EMCalReconstruction("ecal_reco", 
+        inputHitCollection="RawEcalBarrelHits", 
+        outputHitCollection="RecEcalBarrelHits",
+        samplingFraction=0.25,
+        minModuleEdep=0.0*units.MeV)
+
+hcal_reco = HCalReconstruction("hcal_reco",
+        inputHitCollection="RawHcalBarrelHits",
+        outputHitCollection="RecHcalBarrelHits",
+        samplingFraction=0.25,
+        minModuleEdep=0.0*units.MeV)
+#clusters
+ec_barrel_cluster = IslandCluster("ec_barrel_cluster", 
+        inputHitCollection="RecEcalBarrelHits", 
+        outputClusterCollection="EcalBarrelClusters",
+        splitHitCollection="splitEcalBarrelHitCollection",
+        minClusterCenterEdep=1*units.MeV, 
+        groupRange=2.0,
+        OutputLevel=DEBUG)
+
+crystal_ec_cluster = IslandCluster("crystal_ec_cluster", 
+        inputHitCollection="RecoEcalHits", 
+        outputClusterCollection="EcalClusters",
+        minClusterCenterEdep=30*units.MeV, 
+        groupRange=2.0,
+        OutputLevel=DEBUG)
+
+simple_cluster = SimpleClustering("simple_cluster", 
+        inputHitCollection="RecEcalBarrelHits", 
+        outputClusters="SimpleClusters",
+        OutputLevel=DEBUG)
+
+# finalizing clustering (center-of-gravity step)
+ec_barrel_clusterreco = RecoCoG("ec_barrel_clusterreco",
+        clusterCollection="EcalBarrelClusters", 
+        logWeightBase=6.2,
+        samplingFraction=sf) 
+
+clusterreco = RecoCoG("cluster_reco", 
+        clusterCollection="EcalClusters", 
+        logWeightBase=4.2, 
+        moduleDimZName="CrystalBox_z_length",
+        samplingFraction=sf, 
+        OutputLevel=DEBUG)
+
+out = PodioOutput("out", filename=output_rec_file)
+
+out.outputCommands = ["keep *"]
+
+ApplicationMgr(
+    TopAlg = [podioinput, copier, calcopier,
+              ecdigi, emcaldigi,hcaldigi,
+              crystal_ec_reco, ecal_reco, hcal_reco,  
+              #ec_barrel_cluster, crystal_ec_cluster, ec_barrel_clusterreco, clusterreco,
+              #simple_cluster,
+              out],
+    EvtSel = 'NONE',
+    EvtMax   = n_events,
+    ExtSvc = [podioevent],
+    OutputLevel=DEBUG
+ )
-- 
GitLab