Skip to content
Snippets Groups Projects
Commit 90175e49 authored by Chao Peng's avatar Chao Peng
Browse files

modify clustering to use the latest components

parent a9d56427
No related branches found
No related tags found
1 merge request!107modify clustering to use the latest components
Hide whitespace changes
Inline Side-by-side
benchmarks/clustering/config.yml
+ 1
− 1
View file @ 90175e49
......@@ -3,7 +3,7 @@ clustering:process :
stage: process
timeout: 8 hour
script:
- bash benchmarks/clustering/barrel_clusters.sh
- bash benchmarks/clustering/full_cal_clusters.sh
clustering:results:
extends: .rec_benchmark
......
benchmarks/clustering/barrel_clusters.sh → benchmarks/clustering/full_cal_clusters.sh
+ 5
− 3
View file @ 90175e49
......@@ -14,8 +14,8 @@ print_env.sh
# deprecated
export JUGGLER_DETECTOR_PATH=${DETECTOR_PATH}
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=1
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=100
fi
if [[ ! -n "${CB_EMCAL_SAMP_FRAC}" ]] ; then
......@@ -34,8 +34,10 @@ echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
root -b -q "benchmarks/clustering/scripts/gen_central_pions.cxx(${JUGGLER_N_EVENTS}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
### run geant4 simulations
# part.minimalKineticEnergy sets Ek limit for saving particles in mcparticles branch
# 1*TeV means only saving particles from the event generator
npsim --runType batch \
--part.minimalKineticEnergy 1000*GeV \
--part.minimalKineticEnergy "1*TeV" \
-v WARNING \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${DETECTOR_PATH}/${JUGGLER_DETECTOR}.xml \
......
benchmarks/clustering/options/full_cal_clustering.py 0 → 100644
+ 172
− 0
View file @ 90175e49
from Gaudi.Configuration import *
import os
import ROOT
from GaudiKernel.DataObjectHandleBase import DataObjectHandleBase
from Configurables import ApplicationMgr, EICDataSvc, PodioInput, PodioOutput, GeoSvc
from GaudiKernel.SystemOfUnits import MeV, GeV, cm
detector_name = str(os.environ.get("JUGGLER_DETECTOR", "athena"))
detector_path = str(os.environ.get("DETECTOR_PATH", "."))
compact_path = os.path.join(detector_path, detector_name)
# get sampling fractions from system environment variable, 1.0 by default
cb_ecal_sf = float(os.environ.get("CB_EMCAL_SAMP_FRAC", 1.0))
cb_hcal_sf = float(os.environ.get("CB_HCAL_SAMP_FRAC", 1.0))
# input and output
input_sims = [f.strip() for f in str.split(os.environ["JUGGLER_SIM_FILE"], ",") if f.strip()]
output_rec = str(os.environ["JUGGLER_REC_FILE"])
n_events = int(os.environ["JUGGLER_N_EVENTS"])
# geometry service
geo_service = GeoSvc("GeoSvc", detectors=["{}.xml".format(detector_name)])
# data service
podioevent = EICDataSvc("EventDataSvc", inputs=[input_sims])
# juggler components
from Configurables import Jug__Base__InputCopier_dd4pod__Geant4ParticleCollection_dd4pod__Geant4ParticleCollection_ as MCCopier
from Configurables import Jug__Digi__CalorimeterHitDigi as CalHitDigi
from Configurables import Jug__Reco__CalorimeterHitReco as CalHitReco
from Configurables import Jug__Reco__CalorimeterHitsMerger as CalHitsMerger
from Configurables import Jug__Reco__CalorimeterIslandCluster as IslandCluster
from Configurables import Jug__Reco__ImagingPixelReco as ImCalPixelReco
from Configurables import Jug__Reco__ImagingTopoCluster as ImagingCluster
from Configurables import Jug__Reco__ClusterRecoCoG as RecoCoG
from Configurables import Jug__Reco__ImagingClusterReco as ImagingClusterReco
# branches needed from simulation root file
sim_coll = [
"mcparticles",
"CrystalEcalHits",
"EcalBarrelHits",
"HcalBarrelHits",
]
# input and output
podin = PodioInput("PodioReader", collections=sim_coll)
podout = PodioOutput("out", filename=output_rec_file)
# copier needed to get around input --> output copy bug. So truth (mcparticles) can be saved in output file
copier = MCCopier("MCCopier", inputCollection="mcparticles", outputCollection="mcparticles2")
# Crystal Endcap Ecal
ce_ecal_digi = CalHitDigi("ce_ecal_digi",
inputHitCollection="CrystalEcalHits",
outputHitCollection="CrystalEcalHitsDigi",
energyResolutions=[0., 0.02, 0.],
dynamicRangeADC=5.*GeV, # digi settings must match with reco
capacityADC=32768,
pedestalMean=400,
pedestalSigma=3)
ce_ecal_reco = CalHitReco("ce_ecal_reco",
inputHitCollection="CrystalEcalHitsDigi",
outputHitCollection="CrystalEcalHitsReco",
dynamicRangeADC=5.*GeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=3,
thresholdFactor=4, # 4 sigma
minimumHitEdep=1.0*MeV) # discard low energy hits
ce_ecal_cl = IslandCluster("ce_ecal_cl",
# OutputLevel=DEBUG,
inputHitCollection="CrystalEcalHitsReco",
outputClusterCollection="CrystalEcalClusters",
splitCluster=False,
minClusterCenterEdep=30*units.MeV,
groupRanges=[2.2*cm, 2.2*cm])
ce_ecal_clreco = RecoCoG("ce_ecal_clreco",
clusterCollection="CrystalEcalClusters",
logWeightBase=4.6)
# central barrel Ecal
cb_ecal_digi = CalHitDigi("cb_ecal_digi",
inputHitCollection="EcalBarrelHits",
outputHitCollection="EcalBarrelHitsDigi",
energyResolutions=[0., 0.02, 0.], # 2% flat resolution
dynamicRangeADC=3*MeV,
capacityADC=8192,
pedestalMean=400,
pedestalSigma=20) # about 6 keV
cb_ecal_reco = ImaCalPixelReco("cb_ecal_reco",
inputHitCollection="EcalBarrelHitsDigi",
outputHitCollection="EcalBarrelHitsReco",
dynamicRangeADC=3.*MeV,
capacityADC=8192,
pedestalMean=400,
pedestalSigma=20,
thresholdFactor=3, # about 20 keV
readoutClass="EcalBarrelHits", # readout class
layerField="layer", # field to get layer id
sectorField="module") # field to get sector id
cb_ecal_cl = ImagingCluster("cb_ecal_cl",
inputHitCollection="EcalBarrelHitsReco",
outputClusterCollection="EcalBarrelClusters",
localRanges=[2.*units.mm, 2*units.mm], # same layer
adjLayerRanges=[10*units.mrad, 10*units.mrad], # adjacent layer
adjLayerDiff=2, # id diff for adjacent layer
adjSectorDist=3.*units.cm) # different sector
cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
inputClusterCollection="EcalBarrelClusters",
outputLayerCollection="EcalBarrelLayers")
# central barrel hcal
cb_hcal_digi = CalHitDigi("cb_hcal_digi",
inputHitCollection="HcalBarrelHits",
outputHitCollection="HcalBarrelHitsDigi",
dynamicRangeADC=50.*MeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=10)
cb_hcal_reco = CalHitReco("cb_hcal_reco",
inputHitCollection="HcalBarrelHitsDigi",
outputHitCollection="HcalBarrelHitsReco",
dynamicRangeADC=50.*MeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=10,
thresholdFactor=5.0)
# merge hits in different layer (projection to local x-y plane)
cb_hcal_merger = CalHitsMerger("cb_hcal_merger",
fields=["layer", "slice"],
fieldRefNumbers=[1, 0],
inputHitCollection="HcalBarrelHitsReco",
outputHitCollection="HcalBarrelHitsRecoXY")
cb_hcal_cl = IslandCluster("cb_hcal_cl",
inputHitCollection="HcalBarrelHitsRecoXY",
outputClusterCollection="HcalBarrelClusters",
splitCluster=False,
minClusterCenterEdep=30.*MeV,
groupRanges=[15.*cm, 15.*cm],
OutputLevel=DEBUG)
cb_hcal_clreco = RecoCoG("cb_hcal_clreco",
clusterCollection="HcalBarrelClusters",
logWeightBase=6.2,
samplingFraction=cb_hcal_sf,
OutputLevel=DEBUG)
out.outputCommands = ["keep *"]
ApplicationMgr(
TopAlg = [podioinput, copier,
ce_ecal_digi, ce_ecal_reco, ce_ecal_cl, ce_ecal_clreco,
cb_ecal_digi, cb_ecal_reco, cb_ecal_cl, cb_ecal_clreco,
cb_hcal_digi, cb_hcal_reco, cb_hcal_cl, cb_hcal_clreco,
out],
EvtSel = 'NONE',
EvtMax = n_events,
ExtSvc = [podioevent],
OutputLevel=DEBUG
)
benchmarks/clustering/options/fullcalo_clustering.py deleted 100644 → 0
+ 0
− 135
View file @ a9d56427
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
)
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment