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

Resolve "Use the latest algorithms for calorimeter clustering"

parent 623c07b6
No related branches found
No related tags found
1 merge request!117Resolve "Use the latest algorithms for calorimeter clustering"
Hide whitespace changes
Inline Side-by-side
Showing
with 554 additions and 585 deletions
benchmarks/clustering/options/full_cal_reco.py
+ 96
− 88
View file @ 5fae2c48
''' '''
An example script to digitize/reconstruct calorimeter hits An example option file to digitize/reconstruct/clustering calorimeter hits
''' '''
from Gaudi.Configuration import * from Gaudi.Configuration import *
import os import os
...@@ -14,7 +14,7 @@ detector_path = str(os.environ.get("JUGGLER_DETECTOR_PATH", ".")) ...@@ -14,7 +14,7 @@ detector_path = str(os.environ.get("JUGGLER_DETECTOR_PATH", "."))
compact_path = os.path.join(detector_path, detector_name) compact_path = os.path.join(detector_path, detector_name)
# get sampling fractions from system environment variable, 1.0 by default # get sampling fractions from system environment variable, 1.0 by default
ce_ecal_sf = float(os.environ.get("CE_ECAL_SAMP_FRAC", 0.253)) ci_ecal_sf = float(os.environ.get("CI_ECAL_SAMP_FRAC", 0.253))
cb_ecal_sf = float(os.environ.get("CB_ECAL_SAMP_FRAC", 0.01324)) cb_ecal_sf = float(os.environ.get("CB_ECAL_SAMP_FRAC", 0.01324))
cb_hcal_sf = float(os.environ.get("CB_HCAL_SAMP_FRAC", 0.038)) cb_hcal_sf = float(os.environ.get("CB_HCAL_SAMP_FRAC", 0.038))
ci_hcal_sf = float(os.environ.get("CI_HCAL_SAMP_FRAC", 0.025)) ci_hcal_sf = float(os.environ.get("CI_HCAL_SAMP_FRAC", 0.025))
...@@ -63,60 +63,62 @@ copier = MCCopier("MCCopier", inputCollection="mcparticles", outputCollection="m ...@@ -63,60 +63,62 @@ copier = MCCopier("MCCopier", inputCollection="mcparticles", outputCollection="m
# Crystal Endcap Ecal # Crystal Endcap Ecal
ce_ecal_daq = dict(
dynamicRangeADC=5.*GeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=3)
ce_ecal_digi = CalHitDigi("ce_ecal_digi", ce_ecal_digi = CalHitDigi("ce_ecal_digi",
inputHitCollection="EcalEndcapNHits", inputHitCollection="EcalEndcapNHits",
outputHitCollection="EcalEndcapNHitsDigi", outputHitCollection="EcalEndcapNHitsDigi",
energyResolutions=[0., 0.02, 0.], energyResolutions=[0., 0.02, 0.],
dynamicRangeADC=5.*GeV, # digi settings must match with reco **ce_ecal_daq)
capacityADC=32768,
pedestalMean=400,
pedestalSigma=3)
ce_ecal_reco = CalHitReco("ce_ecal_reco", ce_ecal_reco = CalHitReco("ce_ecal_reco",
inputHitCollection="EcalEndcapNHitsDigi", inputHitCollection="EcalEndcapNHitsDigi",
outputHitCollection="EcalEndcapNHitsReco", outputHitCollection="EcalEndcapNHitsReco",
dynamicRangeADC=5.*GeV, thresholdFactor=4, # 4 sigma cut on pedestal sigma
capacityADC=32768,
pedestalMean=400,
pedestalSigma=3,
thresholdFactor=4, # 4 sigma
minimumHitEdep=1.0*MeV, # discard low energy hits
readoutClass="EcalEndcapNHits", readoutClass="EcalEndcapNHits",
sectorField="sector") sectorField="sector",
**ce_ecal_daq)
ce_ecal_cl = IslandCluster("ce_ecal_cl", ce_ecal_cl = IslandCluster("ce_ecal_cl",
# OutputLevel=DEBUG, # OutputLevel=DEBUG,
inputHitCollection="EcalEndcapNHitsReco", inputHitCollection="EcalEndcapNHitsReco",
outputHitCollection="EcalEndcapNClusterHits", outputHitCollection="EcalEndcapNClusterHits",
splitCluster=False, splitCluster=False,
minClusterHitEdep=1.0*MeV, # discard low energy hits
minClusterCenterEdep=30*MeV, minClusterCenterEdep=30*MeV,
groupRanges=[2.2*cm, 2.2*cm]) dimScaledDist=1.8) # dimension scaled dist is good for hybrid sectors with different module size
ce_ecal_clreco = RecoCoG("ce_ecal_clreco", ce_ecal_clreco = RecoCoG("ce_ecal_clreco",
inputHitCollection="EcalEndcapNClusterHits", inputHitCollection="EcalEndcapNClusterHits",
outputClusterCollection="EcalEndcapNClusters", outputClusterCollection="EcalEndcapNClusters",
samplingFraction=0.998, # this accounts for a small fraction of leakage samplingFraction=0.998, # this accounts for a small fraction of leakage
logWeightBase=4.6) logWeightBase=4.6)
# Endcap Sampling Ecal # Endcap Sampling Ecal
ci_ecal_daq = dict(
dynamicRangeADC=50.*MeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=10)
ci_ecal_digi = CalHitDigi("ci_ecal_digi", ci_ecal_digi = CalHitDigi("ci_ecal_digi",
inputHitCollection="EcalEndcapPHits", inputHitCollection="EcalEndcapPHits",
outputHitCollection="EcalEndcapPHitsDigi", outputHitCollection="EcalEndcapPHitsDigi",
dynamicRangeADC=50.*MeV, **ci_ecal_daq)
capacityADC=32768,
pedestalMean=400,
pedestalSigma=10)
ci_ecal_reco = CalHitReco("ci_ecal_reco", ci_ecal_reco = CalHitReco("ci_ecal_reco",
inputHitCollection="EcalEndcapPHitsDigi", inputHitCollection="EcalEndcapPHitsDigi",
outputHitCollection="EcalEndcapPHitsReco", outputHitCollection="EcalEndcapPHitsReco",
dynamicRangeADC=50.*MeV, thresholdFactor=5.0,
capacityADC=32768, **ci_ecal_daq)
pedestalMean=400,
pedestalSigma=10,
thresholdFactor=5.0)
# merge hits in different layer (projection to local x-y plane) # merge hits in different layer (projection to local x-y plane)
ci_ecal_merger = CalHitsMerger("ci_ecal_merger", ci_ecal_merger = CalHitsMerger("ci_ecal_merger",
# OutputLevel=DEBUG,
inputHitCollection="EcalEndcapPHitsReco", inputHitCollection="EcalEndcapPHitsReco",
outputHitCollection="EcalEndcapPHitsRecoXY", outputHitCollection="EcalEndcapPHitsRecoXY",
fields=["layer", "slice"], fields=["layer", "slice"],
...@@ -124,45 +126,49 @@ ci_ecal_merger = CalHitsMerger("ci_ecal_merger", ...@@ -124,45 +126,49 @@ ci_ecal_merger = CalHitsMerger("ci_ecal_merger",
readoutClass="EcalEndcapPHits") readoutClass="EcalEndcapPHits")
ci_ecal_cl = IslandCluster("ci_ecal_cl", ci_ecal_cl = IslandCluster("ci_ecal_cl",
# OutputLevel=DEBUG,
inputHitCollection="EcalEndcapPHitsRecoXY", inputHitCollection="EcalEndcapPHitsRecoXY",
outputHitCollection="EcalEndcapPClusterHits", outputHitCollection="EcalEndcapPClusterHits",
splitCluster=False, splitCluster=False,
minClusterCenterEdep=30.*MeV, minClusterCenterEdep=10.*MeV,
groupRanges=[5*mm, 5*mm]) localDistXY=[10*mm, 10*mm])
ci_ecal_clreco = RecoCoG("ci_ecal_clreco", ci_ecal_clreco = RecoCoG("ci_ecal_clreco",
inputHitCollection="EcalEndcapPClusterHits", inputHitCollection="EcalEndcapPClusterHits",
outputClusterCollection="EcalEndcapPClusters", outputClusterCollection="EcalEndcapPClusters",
logWeightBase=6.2, logWeightBase=6.2,
samplingFraction=ce_ecal_sf) samplingFraction=ci_ecal_sf)
# Central Barrel Ecal (Imaging Cal.) # Central Barrel Ecal (Imaging Cal.)
cb_ecal_digi = CalHitDigi("cb_ecal_digi", cb_ecal_daq = dict(
inputHitCollection="EcalBarrelHits",
outputHitCollection="EcalBarrelHitsDigi",
energyResolutions=[0., 0.02, 0.], # 2% flat resolution
dynamicRangeADC=3*MeV, dynamicRangeADC=3*MeV,
capacityADC=8192, capacityADC=8192,
pedestalMean=400, pedestalMean=400,
pedestalSigma=20) # about 6 keV pedestalSigma=20) # about 6 keV
cb_ecal_digi = CalHitDigi("cb_ecal_digi",
inputHitCollection="EcalBarrelHits",
outputHitCollection="EcalBarrelHitsDigi",
energyResolutions=[0., 0.02, 0.], # 2% flat resolution
**cb_ecal_daq)
cb_ecal_reco = ImCalPixelReco("cb_ecal_reco", cb_ecal_reco = ImCalPixelReco("cb_ecal_reco",
inputHitCollection="EcalBarrelHitsDigi", inputHitCollection="EcalBarrelHitsDigi",
outputHitCollection="EcalBarrelHitsReco", outputHitCollection="EcalBarrelHitsReco",
dynamicRangeADC=3.*MeV,
capacityADC=8192,
pedestalMean=400,
pedestalSigma=20,
thresholdFactor=3, # about 20 keV thresholdFactor=3, # about 20 keV
readoutClass="EcalBarrelHits", # readout class readoutClass="EcalBarrelHits", # readout class
layerField="layer", # field to get layer id layerField="layer", # field to get layer id
sectorField="module") # field to get sector id sectorField="module", # field to get sector id
**cb_ecal_daq)
cb_ecal_cl = ImagingCluster("cb_ecal_cl", cb_ecal_cl = ImagingCluster("cb_ecal_cl",
inputHitCollection="EcalBarrelHitsReco", inputHitCollection="EcalBarrelHitsReco",
outputHitCollection="EcalBarrelClusterHits", outputHitCollection="EcalBarrelClusterHits",
localRanges=[2.*mm, 2*mm], # same layer localDistXY=[2.*mm, 2*mm], # same layer
adjLayerRanges=[10*mrad, 10*mrad], # adjacent layer layerDistEtaPhi=[10*mrad, 10*mrad], # adjacent layer
adjLayerDiff=2, # id diff for adjacent layer neighbourLayersRange=2, # id diff for adjacent layer
adjSectorDist=3.*cm) # different sector sectorDist=3.*cm) # different sector
cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco", cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
samplingFraction=cb_ecal_sf, samplingFraction=cb_ecal_sf,
inputHitCollection="EcalBarrelClusterHits", inputHitCollection="EcalBarrelClusterHits",
...@@ -170,25 +176,26 @@ cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco", ...@@ -170,25 +176,26 @@ cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
outputLayerCollection="EcalBarrelLayers") outputLayerCollection="EcalBarrelLayers")
# Central Barrel Hcal # Central Barrel Hcal
cb_hcal_digi = CalHitDigi("cb_hcal_digi", cb_hcal_daq = dict(
inputHitCollection="HcalBarrelHits",
outputHitCollection="HcalBarrelHitsDigi",
dynamicRangeADC=50.*MeV, dynamicRangeADC=50.*MeV,
capacityADC=32768, capacityADC=32768,
pedestalMean=400, pedestalMean=400,
pedestalSigma=10) pedestalSigma=10)
cb_hcal_digi = CalHitDigi("cb_hcal_digi",
inputHitCollection="HcalBarrelHits",
outputHitCollection="HcalBarrelHitsDigi",
**cb_hcal_daq)
cb_hcal_reco = CalHitReco("cb_hcal_reco", cb_hcal_reco = CalHitReco("cb_hcal_reco",
inputHitCollection="HcalBarrelHitsDigi", inputHitCollection="HcalBarrelHitsDigi",
outputHitCollection="HcalBarrelHitsReco", outputHitCollection="HcalBarrelHitsReco",
dynamicRangeADC=50.*MeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=10,
thresholdFactor=5.0, thresholdFactor=5.0,
readoutClass="HcalBarrelHits", readoutClass="HcalBarrelHits",
layerField="layer", layerField="layer",
sectorField="module") sectorField="module",
**cb_hcal_daq)
cb_hcal_merger = CalHitsMerger("cb_hcal_merger", cb_hcal_merger = CalHitsMerger("cb_hcal_merger",
inputHitCollection="HcalBarrelHitsReco", inputHitCollection="HcalBarrelHitsReco",
outputHitCollection="HcalBarrelHitsRecoXY", outputHitCollection="HcalBarrelHitsRecoXY",
...@@ -201,32 +208,32 @@ cb_hcal_cl = IslandCluster("cb_hcal_cl", ...@@ -201,32 +208,32 @@ cb_hcal_cl = IslandCluster("cb_hcal_cl",
outputHitCollection="HcalBarrelClusterHits", outputHitCollection="HcalBarrelClusterHits",
splitCluster=False, splitCluster=False,
minClusterCenterEdep=30.*MeV, minClusterCenterEdep=30.*MeV,
groupRanges=[15.*cm, 15.*cm]) localDistXY=[15.*cm, 15.*cm])
cb_hcal_clreco = RecoCoG("cb_hcal_clreco", cb_hcal_clreco = RecoCoG("cb_hcal_clreco",
inputHitCollection="HcalBarrelClusterHits", inputHitCollection="HcalBarrelClusterHits",
outputClusterCollection="HcalBarrelClusters", outputClusterCollection="HcalBarrelClusters",
logWeightBase=6.2, logWeightBase=6.2,
samplingFraction=cb_hcal_sf) samplingFraction=cb_hcal_sf)
# Hcal Hadron Endcap # Hcal Hadron Endcap
ci_hcal_digi = CalHitDigi("ci_hcal_digi", ci_hcal_daq = dict(
inputHitCollection="HcalHadronEndcapHits",
outputHitCollection="HcalHadronEndcapHitsDigi",
dynamicRangeADC=50.*MeV, dynamicRangeADC=50.*MeV,
capacityADC=32768, capacityADC=32768,
pedestalMean=400, pedestalMean=400,
pedestalSigma=10) pedestalSigma=10)
ci_hcal_digi = CalHitDigi("ci_hcal_digi",
inputHitCollection="HcalHadronEndcapHits",
outputHitCollection="HcalHadronEndcapHitsDigi",
**ci_hcal_daq)
ci_hcal_reco = CalHitReco("ci_hcal_reco", ci_hcal_reco = CalHitReco("ci_hcal_reco",
inputHitCollection="HcalHadronEndcapHitsDigi", inputHitCollection="HcalHadronEndcapHitsDigi",
outputHitCollection="HcalHadronEndcapHitsReco", outputHitCollection="HcalHadronEndcapHitsReco",
dynamicRangeADC=50.*MeV, thresholdFactor=5.0,
capacityADC=32768, **ci_hcal_daq)
pedestalMean=400,
pedestalSigma=10,
thresholdFactor=5.0)
ci_hcal_merger = CalHitsMerger("ci_hcal_merger", ci_hcal_merger = CalHitsMerger("ci_hcal_merger",
inputHitCollection="HcalHadronEndcapHitsReco", inputHitCollection="HcalHadronEndcapHitsReco",
outputHitCollection="HcalHadronEndcapHitsRecoXY", outputHitCollection="HcalHadronEndcapHitsRecoXY",
...@@ -239,31 +246,32 @@ ci_hcal_cl = IslandCluster("ci_hcal_cl", ...@@ -239,31 +246,32 @@ ci_hcal_cl = IslandCluster("ci_hcal_cl",
outputHitCollection="HcalHadronEndcapClusterHits", outputHitCollection="HcalHadronEndcapClusterHits",
splitCluster=False, splitCluster=False,
minClusterCenterEdep=30.*MeV, minClusterCenterEdep=30.*MeV,
groupRanges=[15.*cm, 15.*cm]) localDistXY=[15.*cm, 15.*cm])
ci_hcal_clreco = RecoCoG("ci_hcal_clreco", ci_hcal_clreco = RecoCoG("ci_hcal_clreco",
inputHitCollection="HcalHadronEndcapClusterHits", inputHitCollection="HcalHadronEndcapClusterHits",
outputClusterCollection="HcalHadronEndcapClusters", outputClusterCollection="HcalHadronEndcapClusters",
logWeightBase=6.2, logWeightBase=6.2,
samplingFraction=ci_hcal_sf) samplingFraction=ci_hcal_sf)
# Hcal Electron Endcap # Hcal Electron Endcap
ce_hcal_digi = CalHitDigi("ce_hcal_digi", ce_hcal_daq = dict(
inputHitCollection="HcalElectronEndcapHits",
outputHitCollection="HcalElectronEndcapHitsDigi",
dynamicRangeADC=50.*MeV, dynamicRangeADC=50.*MeV,
capacityADC=32768, capacityADC=32768,
pedestalMean=400, pedestalMean=400,
pedestalSigma=10) pedestalSigma=10)
ce_hcal_digi = CalHitDigi("ce_hcal_digi",
inputHitCollection="HcalElectronEndcapHits",
outputHitCollection="HcalElectronEndcapHitsDigi",
**ce_hcal_daq)
ce_hcal_reco = CalHitReco("ce_hcal_reco", ce_hcal_reco = CalHitReco("ce_hcal_reco",
inputHitCollection="HcalElectronEndcapHitsDigi", inputHitCollection="HcalElectronEndcapHitsDigi",
outputHitCollection="HcalElectronEndcapHitsReco", outputHitCollection="HcalElectronEndcapHitsReco",
dynamicRangeADC=50.*MeV, thresholdFactor=5.0,
capacityADC=32768, **ce_hcal_daq)
pedestalMean=400,
pedestalSigma=10,
thresholdFactor=5.0)
ce_hcal_merger = CalHitsMerger("ce_hcal_merger", ce_hcal_merger = CalHitsMerger("ce_hcal_merger",
inputHitCollection="HcalElectronEndcapHitsReco", inputHitCollection="HcalElectronEndcapHitsReco",
outputHitCollection="HcalElectronEndcapHitsRecoXY", outputHitCollection="HcalElectronEndcapHitsRecoXY",
...@@ -276,34 +284,34 @@ ce_hcal_cl = IslandCluster("ce_hcal_cl", ...@@ -276,34 +284,34 @@ ce_hcal_cl = IslandCluster("ce_hcal_cl",
outputHitCollection="HcalElectronEndcapClusterHits", outputHitCollection="HcalElectronEndcapClusterHits",
splitCluster=False, splitCluster=False,
minClusterCenterEdep=30.*MeV, minClusterCenterEdep=30.*MeV,
groupRanges=[15.*cm, 15.*cm]) localDistXY=[15.*cm, 15.*cm])
ce_hcal_clreco = RecoCoG("ce_hcal_clreco", ce_hcal_clreco = RecoCoG("ce_hcal_clreco",
inputHitCollection="HcalElectronEndcapClusterHits", inputHitCollection="HcalElectronEndcapClusterHits",
outputClusterCollection="HcalElectronEndcapClusters", outputClusterCollection="HcalElectronEndcapClusters",
logWeightBase=6.2, logWeightBase=6.2,
samplingFraction=ce_hcal_sf) samplingFraction=ce_hcal_sf)
podout.outputCommands = ['drop *', 'keep mcparticles*', 'keep *Reco', 'keep *Digi', 'keep *Clusters', 'keep *Layers', 'drop *ProtoClusters'] podout.outputCommands = ['drop *',
'keep mcparticles2',
'keep *Digi',
'keep *Reco*',
'keep *ClusterHits',
'keep *Clusters',
'keep *Layers']
ApplicationMgr( ApplicationMgr(
TopAlg = [podin, copier, TopAlg = [podin, copier,
ce_ecal_digi, ce_ecal_reco, ce_ecal_digi, ce_ecal_reco, ce_ecal_cl, ce_ecal_clreco,
ci_ecal_digi, ci_ecal_reco, ci_ecal_digi, ci_ecal_reco, ci_ecal_merger, ci_ecal_cl, ci_ecal_clreco,
cb_ecal_digi, cb_ecal_reco, cb_ecal_digi, cb_ecal_reco, cb_ecal_cl, cb_ecal_clreco,
cb_hcal_digi, cb_hcal_reco, cb_hcal_digi, cb_hcal_reco, cb_hcal_merger, cb_hcal_cl, cb_hcal_clreco,
ce_hcal_digi, ce_hcal_reco, ce_hcal_digi, ce_hcal_reco, ce_hcal_merger, ce_hcal_cl, ce_hcal_clreco,
ci_hcal_digi, ci_hcal_reco, ci_hcal_digi, ci_hcal_reco, ci_hcal_merger, ci_hcal_cl, ci_hcal_clreco,
ce_ecal_cl, ce_ecal_clreco,
ci_ecal_merger, ci_ecal_cl, ci_ecal_clreco,
cb_ecal_cl, cb_ecal_clreco,
cb_hcal_merger, cb_hcal_cl, cb_hcal_clreco,
ce_hcal_merger, ce_hcal_cl, ce_hcal_clreco,
ci_hcal_merger, ci_hcal_cl, ci_hcal_clreco,
podout], podout],
EvtSel = 'NONE', EvtSel = 'NONE',
EvtMax = n_events, EvtMax = n_events,
ExtSvc = [podioevent], ExtSvc = [podioevent],
OutputLevel=DEBUG OutputLevel=DEBUG
) )
benchmarks/clustering/scripts/cluster_plots.py
+ 1
− 1
View file @ 5fae2c48
...@@ -97,7 +97,7 @@ def general_clusters_figure(df, collection, save, min_nhits=3): ...@@ -97,7 +97,7 @@ def general_clusters_figure(df, collection, save, min_nhits=3):
fig, axs = plt.subplots(2, 2, figsize=(16, 12), dpi=120) fig, axs = plt.subplots(2, 2, figsize=(16, 12), dpi=120)
labels = [ labels = [
('Number of Hits', 'Counts'), ('Number of Hits', 'Counts'),
(r'$E$ (MeV)', 'Counts'), (r'$E$ (GeV)', 'Counts'),
(r'$\theta$ (rad)', 'Counts'), (r'$\theta$ (rad)', 'Counts'),
(r'$\phi$ (rad)', 'Counts'), (r'$\phi$ (rad)', 'Counts'),
] ]
......
benchmarks/ecal/config.yml
+ 33
− 19
View file @ 5fae2c48
emcal_crystal_electrons: ecal_endcap_e_electrons:
extends: .rec_benchmark extends: .rec_benchmark
stage: process stage: process
timeout: 8 hours timeout: 8 hours
script: script:
- bash benchmarks/ecal/emcal_electrons.sh - python benchmarks/ecal/run_emcal_benchmarks.py endcap_e --particle "electron" -n 100 -t "endcap_e_electron"
--pmin 5.0 --pmax 5.0 --amin 135 --amax 175
emcal_crystal_pi0s: ecal_endcap_i_electrons:
extends: .rec_benchmark extends: .rec_benchmark
stage: run stage: process
timeout: 12 hours timeout: 8 hours
script: script:
- bash benchmarks/ecal/emcal_pi0s.sh - python benchmarks/ecal/run_emcal_benchmarks.py endcap_i --particle "electron" -n 100 -t "endcap_i_electron"
--pmin 5.0 --pmax 5.0 --amin 5 --amax 45
emcal_barrel_electrons: ecal_barrel_electrons:
extends: .rec_benchmark extends: .rec_benchmark
timeout: 48 hours stage: process
stage: run timeout: 8 hours
script: script:
- bash benchmarks/ecal/run_emcal_barrel_electrons.sh - python benchmarks/ecal/run_emcal_benchmarks.py barrel --particle "electron" -n 100 -t "barrel_electron"
--pmin 5.0 --pmax 5.0 --amin 45 --amax 135
emcal_barrel_pions: ecal_endcap_e_pions:
extends: .rec_benchmark extends: .rec_benchmark
timeout: 48 hours stage: process
stage: run timeout: 8 hours
script: script:
- bash benchmarks/ecal/run_emcal_barrel_pions.sh - python benchmarks/ecal/run_emcal_benchmarks.py endcap_e --particle "pion-" -n 100 -t "endcap_e_pion"
--pmin 5.0 --pmax 5.0 --amin 135 --amax 175
full_emcal_barrel_electrons: ecal_endcap_i_pions:
extends: .rec_benchmark extends: .rec_benchmark
timeout: 48 hours stage: process
stage: run timeout: 8 hours
script:
- python benchmarks/ecal/run_emcal_benchmarks.py endcap_i --particle "pion-" -n 100 -t "endcap_i_pion"
--pmin 5.0 --pmax 5.0 --amin 5 --amax 45
ecal_barrel_pions:
extends: .rec_benchmark
stage: process
timeout: 8 hours
script: script:
- bash benchmarks/ecal/full_emcal_electrons.sh - python benchmarks/ecal/run_emcal_benchmarks.py barrel --particle "pion-" -n 100 -t "barrel_pion"
--pmin 5.0 --pmax 5.0 --amin 45 --amax 135
ecal_collect: ecal_collect:
stage: collect stage: collect
needs: needs:
- ["emcal_crystal_electrons", "emcal_crystal_pi0s", "emcal_barrel_electrons", "emcal_barrel_pions", "full_emcal_barrel_electrons"] - ["ecal_endcap_e_electrons", "ecal_endcap_e_pions", "ecal_endcap_i_electrons", "ecal_endcap_i_pions",
"ecal_barrel_electrons", "ecal_barrel_pions"]
script: script:
- echo "Done collecting artifacts." - echo "Done collecting artifacts."
benchmarks/ecal/emcal_electrons.sh deleted 100644 → 0
+ 0
− 116
View file @ 623c07b6
#!/bin/bash
print_env.sh
## To run the reconstruction, we need the following global variables:
## - JUGGLER_INSTALL_PREFIX: Install prefix for Juggler (simu/recon)
## - JUGGLER_DETECTOR: the detector package we want to use for this benchmark
## - JUGGLER_DETECTOR_VERSION: the detector package we want to use for this benchmark
## - DETECTOR_PATH: full path to the detector definitions
##
## You can ready options/env.sh for more in-depth explanations of the variables
## and how they can be controlled.
source options/env.sh
if [[ ! -n "${JUGGLER_DETECTOR_PATH}" ]] ; then
export JUGGLER_DETECTOR_PATH=${DETECTOR_PATH}
fi
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=100
fi
if [[ ! -n "${E_start}" ]] ; then
export E_start=0.0
fi
if [[ ! -n "${E_end}" ]] ; then
export E_end=30.0
fi
export JUGGLER_FILE_NAME_TAG="emcal_uniform_electrons"
export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.root"
export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
# generate the input events
# note datasets is now only used to develop datasets.
#git clone https://eicweb.phy.anl.gov/EIC/datasets.git datasets
root -b -q "benchmarks/ecal/scripts/emcal_electrons.cxx(${JUGGLER_N_EVENTS}, ${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script"
exit 1
fi
root -b -q "benchmarks/ecal/scripts/emcal_electrons_reader.cxx(${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script"
exit 1
fi
# run geant4 simulations
npsim --runType batch \
--part.minimalKineticEnergy 1000*GeV \
-v WARNING \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${JUGGLER_DETECTOR_PATH}/${JUGGLER_DETECTOR}.xml \
--inputFiles ${JUGGLER_FILE_NAME_TAG}.hepmc \
--outputFile ${JUGGLER_SIM_FILE}
# Need to figure out how to pass file name to juggler from the commandline
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running npdet"
exit 1
fi
xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv \
gaudirun.py benchmarks/ecal/options/crystal_calorimeter_reco.py
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running juggler"
exit 1
fi
mkdir -p results
#rootls ${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}
root -b -q "benchmarks/ecal/scripts/rec_emcal_electrons_reader.C(${E_start}, ${E_end}, \"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
#root -b -q "ecal/scripts/makeplot.C(${E_start}, ${E_end}, \"${JUGGLER_DETECTOR}/${JUGGLER_REC_FILE}\", \"results/rec_${JUGGLER_FILE_NAME_TAG}.txt\")"
#root -b -q "ecal/scripts/makeplot_input.C(\"${JUGGLER_DETECTOR}/${JUGGLER_SIM_FILE}\", \"results/sim_${JUGGLER_FILE_NAME_TAG}.txt\")"
root -b -q "benchmarks/ecal/scripts/crystal_cal_electrons.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
root -b -q "benchmarks/ecal/scripts/emcal_electrons_analysis.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
# Script to generate energy resolution plots
root -b -q "benchmarks/ecal/scripts/rec_emcal_resolution_analysis.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
#paste results/sim_${JUGGLER_FILE_NAME_TAG}.txt results/rec_${JUGGLER_FILE_NAME_TAG}.txt > results/eng_${JUGGLER_FILE_NAME_TAG}.txt
#root -b -q "ecal/scripts/read_eng.C(\"results/eng_${JUGGLER_FILE_NAME_TAG}.root\", \"results/eng_${JUGGLER_FILE_NAME_TAG}.txt\")"
#root -b -q "ecal/scripts/cal_eng_res.C(\"results/eng_${JUGGLER_FILE_NAME_TAG}.root\")"
root_filesize=$(stat --format=%s "${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/.
fi
fi
benchmarks/ecal/emcal_pi0s.sh deleted 100644 → 0
+ 0
− 86
View file @ 623c07b6
#!/bin/bash
# Based on emcal_electrons.sh script
print_env.sh
## To run the reconstruction, we need the following global variables:
## - JUGGLER_INSTALL_PREFIX: Install prefix for Juggler (simu/recon)
## - JUGGLER_DETECTOR: the detector package we want to use for this benchmark
## - JUGGLER_DETECTOR_VERSION: the detector package we want to use for this benchmark
## - DETECTOR_PATH: full path to the detector definitions
##
## You can read .local/bin/env.sh for more in-depth explanations of the variables
## and how they can be controlled.
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=100
fi
if [[ ! -n "${JUGGLER_DETECTOR_PATH}" ]] ; then
export JUGGLER_DETECTOR_PATH=${DETECTOR_PATH}
fi
# File names
export JUGGLER_FILE_NAME_TAG="emcal_uniform_pi0s"
export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.root"
export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
# Datasets
#git clone https://eicweb.phy.anl.gov/EIC/datasets.git datasets
root -b -q "benchmarks/ecal/scripts/emcal_pi0.cxx(${JUGGLER_N_EVENTS}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
root -b -q "benchmarks/ecal/scripts/emcal_pi0_reader.cxx(\"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
echo "CHECK POINT FOR GEANT4 SIMULATION"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
echo "JUGGLER_FILE_NAME_TAG = ${JUGGLER_FILE_NAME_TAG}"
echo "JUGGLER_SIM_FILE = ${JUGGLER_SIM_FILE}"
# run geant4 simulations
npsim --runType batch \
--part.minimalKineticEnergy 1000*GeV \
-v WARNING \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${JUGGLER_DETECTOR_PATH}/${JUGGLER_DETECTOR}.xml \
--inputFiles ${JUGGLER_FILE_NAME_TAG}.hepmc \
--outputFile ${JUGGLER_SIM_FILE}
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running npdet"
exit 1
fi
# Need to figure out how to pass file name to juggler from the commandline
gaudirun.py benchmarks/ecal/options/crystal_calorimeter_reco.py
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running juggler"
exit 1
fi
mkdir -p results
root -b -q "benchmarks/ecal/scripts/emcal_pion_anlaysis.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
root_filesize=$(stat --format=%s "${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/.
fi
fi
benchmarks/ecal/full_emcal_electrons.sh deleted 100644 → 0
+ 0
− 82
View file @ 623c07b6
#!/bin/bash
print_env.sh
## To run the reconstruction, we need the following global variables:
## - JUGGLER_INSTALL_PREFIX: Install prefix for Juggler (simu/recon)
## - JUGGLER_DETECTOR: the detector package we want to use for this benchmark
## - JUGGLER_DETECTOR_VERSION: the detector package we want to use for this benchmark
## - DETECTOR_PATH: full path to the detector definitions
##
## You can ready options/env.sh for more in-depth explanations of the variables
## and how they can be controlled.
source options/env.sh
if [[ ! -n "${JUGGLER_DETECTOR_PATH}" ]] ; then
export JUGGLER_DETECTOR_PATH=${DETECTOR_PATH}
fi
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=100
fi
if [[ ! -n "${E_start}" ]] ; then
export E_start=0.5
fi
if [[ ! -n "${E_end}" ]] ; then
export E_end=5.0
fi
export JUGGLER_FILE_NAME_TAG="emcal_uniform_electrons"
export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.root"
export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
# generate the input events
# note datasets is now only used to develop datasets.
#git clone https://eicweb.phy.anl.gov/EIC/datasets.git datasets
root -b -q "benchmarks/ecal/scripts/full_emcal_electrons.cxx(${JUGGLER_N_EVENTS}, ${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script"
exit 1
fi
# run geant4 simulations
npsim --runType batch \
-v WARNING \
--part.minimalKineticEnergy 1*GeV \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${JUGGLER_DETECTOR_PATH}/${JUGGLER_DETECTOR}.xml \
--inputFiles ${JUGGLER_FILE_NAME_TAG}.hepmc \
--outputFile ${JUGGLER_SIM_FILE}
# Need to figure out how to pass file name to juggler from the commandline
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running npdet"
exit 1
fi
xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv \
gaudirun.py benchmarks/ecal/options/full_em_calorimeter_reco.py
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running juggler"
exit 1
fi
root -b -q "benchmarks/ecal/scripts/full_emcal_electrons_analysis.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
root_filesize=$(stat --format=%s "${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/.
fi
fi
benchmarks/ecal/options/barrel.py 0 → 100644
+ 99
− 0
View file @ 5fae2c48
'''
An example script to digitize/reconstruct/clustering endcap ecal hits
'''
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, mm, cm, mrad
detector_name = str(os.environ.get("JUGGLER_DETECTOR", "athena"))
detector_path = str(os.environ.get("JUGGLER_DETECTOR_PATH", "."))
compact_path = str(os.environ.get("JUGGLER_COMPACT_PATH", "{}.xml".format(os.path.join(detector_path, detector_name))))
# 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"])
cb_ecal_sf = float(os.environ.get("CB_ECAL_SAMP_FRAC", 0.01324))
# geometry service
geo_service = GeoSvc("GeoSvc", detectors=[compact_path], OutputLevel=INFO)
# 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__ImagingPixelReco as ImCalPixelReco
from Configurables import Jug__Reco__ImagingTopoCluster as ImagingCluster
from Configurables import Jug__Reco__ImagingClusterReco as ImagingClusterReco
# branches needed from simulation root file
sim_coll = [
"mcparticles",
"EcalBarrelHits",
]
# input and output
podin = PodioInput("PodioReader", collections=sim_coll)
podout = PodioOutput("out", filename=output_rec)
# 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")
# Central Barrel Ecal (Imaging Cal.)
cb_ecal_daq = dict(
dynamicRangeADC=3*MeV,
capacityADC=8192,
pedestalMean=400,
pedestalSigma=20) # about 6 keV
cb_ecal_digi = CalHitDigi("cb_ecal_digi",
inputHitCollection="EcalBarrelHits",
outputHitCollection="EcalBarrelHitsDigi",
energyResolutions=[0., 0.02, 0.], # 2% flat resolution
**cb_ecal_daq)
cb_ecal_reco = ImCalPixelReco("cb_ecal_reco",
inputHitCollection="EcalBarrelHitsDigi",
outputHitCollection="EcalBarrelHitsReco",
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_daq)
cb_ecal_cl = ImagingCluster("cb_ecal_cl",
inputHitCollection="EcalBarrelHitsReco",
outputHitCollection="EcalBarrelClusterHits",
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
cb_ecal_clreco = ImagingClusterReco("cb_ecal_clreco",
samplingFraction=cb_ecal_sf,
inputHitCollection="EcalBarrelClusterHits",
outputClusterCollection="EcalBarrelClusters",
outputLayerCollection="EcalBarrelLayers")
podout.outputCommands = ['drop *',
'keep mcparticles2',
'keep *HitsReco',
'keep *HitsDigi',
'keep *ClusterHits',
'keep *Clusters']
ApplicationMgr(
TopAlg = [podin, copier,
cb_ecal_digi, cb_ecal_reco, cb_ecal_cl, cb_ecal_clreco,
podout],
EvtSel = 'NONE',
EvtMax = n_events,
ExtSvc = [podioevent],
OutputLevel=DEBUG
)
benchmarks/ecal/options/endcap_e.py 0 → 100644
+ 98
− 0
View file @ 5fae2c48
'''
An example script to digitize/reconstruct/clustering endcap ecal hits
'''
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, mm, cm, mrad
detector_name = str(os.environ.get("JUGGLER_DETECTOR", "athena"))
detector_path = str(os.environ.get("JUGGLER_DETECTOR_PATH", "."))
compact_path = str(os.environ.get("JUGGLER_COMPACT_PATH", "{}.xml".format(os.path.join(detector_path, detector_name))))
# 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=[compact_path], OutputLevel=INFO)
# 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__CalorimeterIslandCluster as IslandCluster
from Configurables import Jug__Reco__ClusterRecoCoG as RecoCoG
# branches needed from simulation root file
sim_coll = [
"mcparticles",
"EcalEndcapNHits",
]
# input and output
podin = PodioInput("PodioReader", collections=sim_coll)
podout = PodioOutput("out", filename=output_rec)
# 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_daq = dict(
dynamicRangeADC=5.*GeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=3)
ce_ecal_digi = CalHitDigi("ce_ecal_digi",
inputHitCollection="EcalEndcapNHits",
outputHitCollection="EcalEndcapNHitsDigi",
energyResolutions=[0., 0.02, 0.],
**ce_ecal_daq)
ce_ecal_reco = CalHitReco("ce_ecal_reco",
inputHitCollection="EcalEndcapNHitsDigi",
outputHitCollection="EcalEndcapNHitsReco",
thresholdFactor=4, # 4 sigma cut on pedestal sigma
readoutClass="EcalEndcapNHits",
sectorField="sector",
**ce_ecal_daq)
ce_ecal_cl = IslandCluster("ce_ecal_cl",
# OutputLevel=DEBUG,
inputHitCollection="EcalEndcapNHitsReco",
outputHitCollection="EcalEndcapNClusterHits",
splitCluster=False,
minClusterHitEdep=1.0*MeV, # discard low energy hits
minClusterCenterEdep=30*MeV,
dimScaledDist=1.8) # hybrid calorimeter with different dimensions, using a scaled dist
ce_ecal_clreco = RecoCoG("ce_ecal_clreco",
inputHitCollection="EcalEndcapNClusterHits",
outputClusterCollection="EcalEndcapNClusters",
samplingFraction=0.998, # this accounts for a small fraction of leakage
logWeightBase=4.6)
podout.outputCommands = ['drop *',
'keep mcparticles2',
'keep *HitsReco',
'keep *HitsDigi',
'keep *ClusterHits',
'keep *Clusters']
ApplicationMgr(
TopAlg = [podin, copier,
ce_ecal_digi, ce_ecal_reco, ce_ecal_cl, ce_ecal_clreco,
podout],
EvtSel = 'NONE',
EvtMax = n_events,
ExtSvc = [podioevent],
OutputLevel=DEBUG
)
benchmarks/ecal/options/endcap_i.py 0 → 100644
+ 106
− 0
View file @ 5fae2c48
'''
An example script to digitize/reconstruct/clustering endcap ecal hits
'''
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, mm, cm, mrad
detector_name = str(os.environ.get("JUGGLER_DETECTOR", "athena"))
detector_path = str(os.environ.get("JUGGLER_DETECTOR_PATH", "."))
compact_path = str(os.environ.get("JUGGLER_COMPACT_PATH", "{}.xml".format(os.path.join(detector_path, detector_name))))
ci_ecal_sf = float(os.environ.get("CI_ECAL_SAMP_FRAC", 0.253))
# 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=[compact_path], OutputLevel=INFO)
# 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__ClusterRecoCoG as RecoCoG
# branches needed from simulation root file
sim_coll = [
"mcparticles",
"EcalEndcapPHits",
]
# input and output
podin = PodioInput("PodioReader", collections=sim_coll)
podout = PodioOutput("out", filename=output_rec)
# 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")
# Endcap sampling Ecal
ci_ecal_daq = dict(
dynamicRangeADC=50.*MeV,
capacityADC=32768,
pedestalMean=400,
pedestalSigma=10)
ci_ecal_digi = CalHitDigi("ci_ecal_digi",
inputHitCollection="EcalEndcapPHits",
outputHitCollection="EcalEndcapPHitsDigi",
**ci_ecal_daq)
ci_ecal_reco = CalHitReco("ci_ecal_reco",
inputHitCollection="EcalEndcapPHitsDigi",
outputHitCollection="EcalEndcapPHitsReco",
thresholdFactor=5.0,
**ci_ecal_daq)
# merge hits in different layer (projection to local x-y plane)
ci_ecal_merger = CalHitsMerger("ci_ecal_merger",
# OutputLevel=DEBUG,
inputHitCollection="EcalEndcapPHitsReco",
outputHitCollection="EcalEndcapPHitsRecoXY",
fields=["layer", "slice"],
fieldRefNumbers=[1, 0],
readoutClass="EcalEndcapPHits")
ci_ecal_cl = IslandCluster("ci_ecal_cl",
# OutputLevel=DEBUG,
inputHitCollection="EcalEndcapPHitsRecoXY",
outputHitCollection="EcalEndcapPClusterHits",
splitCluster=False,
minClusterCenterEdep=10.*MeV,
localDistXY=[10*mm, 10*mm])
ci_ecal_clreco = RecoCoG("ci_ecal_clreco",
inputHitCollection="EcalEndcapPClusterHits",
outputClusterCollection="EcalEndcapPClusters",
logWeightBase=6.2,
samplingFraction=ci_ecal_sf)
podout.outputCommands = ['drop *',
'keep mcparticles2',
'keep *HitsReco*',
'keep *HitsDigi',
'keep *ClusterHits',
'keep *Clusters']
ApplicationMgr(
TopAlg = [podin, copier,
ci_ecal_digi, ci_ecal_reco, ci_ecal_merger, ci_ecal_cl, ci_ecal_clreco,
podout],
EvtSel = 'NONE',
EvtMax = n_events,
ExtSvc = [podioevent],
OutputLevel=DEBUG
)
benchmarks/ecal/run_emcal_barrel_electrons.sh deleted 100755 → 0
+ 0
− 98
View file @ 623c07b6
#!/bin/bash
print_env.sh
## To run the reconstruction, we need the following global variables:
## - JUGGLER_INSTALL_PREFIX: Install prefix for Juggler (simu/recon)
## - JUGGLER_DETECTOR: the detector package we want to use for this benchmark
## - JUGGLER_DETECTOR_VERSION: the detector package we want to use for this benchmark
## - DETECTOR_PATH: full path to the detector definitions
##
## You can ready options/env.sh for more in-depth explanations of the variables
## and how they can be controlled.
source options/env.sh
if [[ ! -n "${JUGGLER_DETECTOR_PATH}" ]] ; then
export JUGGLER_DETECTOR_PATH=${DETECTOR_PATH}
fi
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=1000
fi
if [[ ! -n "${E_start}" ]] ; then
export E_start=5.0
fi
if [[ ! -n "${E_end}" ]] ; then
export E_end=5.0
fi
if [[ ! -n "${CB_EMCAL_SAMP_FRAC}" ]] ; then
export CB_EMCAL_SAMP_FRAC=0.014
fi
export JUGGLER_FILE_NAME_TAG="emcal_barrel_uniform_electrons"
export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.root"
export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
# Generate the input events
root -b -q "benchmarks/ecal/scripts/emcal_barrel_electrons.cxx(${JUGGLER_N_EVENTS}, ${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script: generating input events"
exit 1
fi
# Plot the input events
root -b -q "benchmarks/ecal/scripts/emcal_barrel_electrons_reader.cxx(${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script: plotting input events"
exit 1
fi
# Run geant4 simulations
npsim --runType batch \
-v WARNING \
--part.minimalKineticEnergy 0.5*GeV \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${JUGGLER_DETECTOR_PATH}/${JUGGLER_DETECTOR}.xml \
--inputFiles ${JUGGLER_FILE_NAME_TAG}.hepmc \
--outputFile ${JUGGLER_SIM_FILE}
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running npdet"
exit 1
fi
rootls -t "${JUGGLER_SIM_FILE}"
# Run Juggler
xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv \
gaudirun.py benchmarks/ecal/options/emcal_barrel_reco.py
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running juggler"
exit 1
fi
# Directory for plots
mkdir -p results
# Run analysis script
root -b -q "benchmarks/ecal/scripts/emcal_barrel_electrons_analysis.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
root_filesize=$(stat --format=%s "${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/.
fi
fi
benchmarks/ecal/run_emcal_barrel_pions.sh deleted 100755 → 0
+ 0
− 95
View file @ 623c07b6
#!/bin/bash
print_env.sh
## To run the reconstruction, we need the following global variables:
## - JUGGLER_INSTALL_PREFIX: Install prefix for Juggler (simu/recon)
## - JUGGLER_DETECTOR: the detector package we want to use for this benchmark
## - JUGGLER_DETECTOR_VERSION: the detector package we want to use for this benchmark
## - DETECTOR_PATH: full path to the detector definitions
##
## You can ready options/env.sh for more in-depth explanations of the variables
## and how they can be controlled.
source options/env.sh
if [[ ! -n "${JUGGLER_DETECTOR_PATH}" ]] ; then
export JUGGLER_DETECTOR_PATH=${DETECTOR_PATH}
fi
if [[ ! -n "${JUGGLER_N_EVENTS}" ]] ; then
export JUGGLER_N_EVENTS=1000
fi
if [[ ! -n "${E_start}" ]] ; then
export E_start=5.0
fi
if [[ ! -n "${E_end}" ]] ; then
export E_end=5.0
fi
if [[ ! -n "${CB_EMCAL_SAMP_FRAC}" ]] ; then
export CB_EMCAL_SAMP_FRAC=0.01
fi
export JUGGLER_FILE_NAME_TAG="emcal_barrel_uniform_pions"
export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.root"
export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
echo "JUGGLER_DETECTOR = ${JUGGLER_DETECTOR}"
# Generate the input events
root -b -q "benchmarks/ecal/scripts/emcal_barrel_pions.cxx(${JUGGLER_N_EVENTS}, ${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script: generating input events"
exit 1
fi
# Plot the input events
root -b -q "benchmarks/ecal/scripts/emcal_barrel_pions_reader.cxx(${E_start}, ${E_end}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running script: plotting input events"
exit 1
fi
# Run geant4 simulations
npsim --runType batch \
-v WARNING \
--part.minimalKineticEnergy 0.5*GeV \
--numberOfEvents ${JUGGLER_N_EVENTS} \
--compactFile ${JUGGLER_DETECTOR_PATH}/${JUGGLER_DETECTOR}.xml \
--inputFiles ${JUGGLER_FILE_NAME_TAG}.hepmc \
--outputFile ${JUGGLER_SIM_FILE}
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running npdet"
exit 1
fi
# Run Juggler
xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv \
gaudirun.py benchmarks/ecal/options/emcal_barrel_reco.py
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running juggler"
exit 1
fi
# Directory for plots
mkdir -p results
# Run analysis script
root -b -q "benchmarks/ecal/scripts/emcal_barrel_pions_analysis.cxx(\"${JUGGLER_REC_FILE}\")"
if [[ "$?" -ne "0" ]] ; then
echo "ERROR running analysis script"
exit 1
fi
root_filesize=$(stat --format=%s "${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/.
fi
fi
benchmarks/ecal/run_emcal_benchmarks.py 0 → 100755
+ 121
− 0
View file @ 5fae2c48
#! /usr/local/bin/python3
'''
A python script to help run emcal benchmarks
simulation -> reconstruction -> analysis
Author: Chao Peng, Jihee Kim (ANL)
Date: 07/04/2021
'''
import os
import sys
import subprocess
import argparse
# type: option files, analysis files
default_type = {
'endcap_e': [['endcap_e.py'], []],
'endcap_i': [['endcap_i.py'], []],
'barrel': [['barrel.py'], []],
# 'all': [['all_ecal.py'], []],
}
default_compact = os.path.join(os.environ.get('JUGGLER_DETECTOR_PATH', os.environ.get('DETECTOR_PATH', '')),
'{}.xml'.format(os.environ.get('JUGGLER_DETECTOR', '')))
parser = argparse.ArgumentParser()
parser.add_argument('run_type', help='Run type, support {}'.format(list(default_type.keys())))
parser.add_argument('-n', '--numberOfEvents', dest='nev', type=int, default=100, help='Number of events to process.')
parser.add_argument('-t', '--nametag', type=str, default='IMCAL_ML', help='Name tag for output files.')
parser.add_argument('--seed', type=int, default=-1, help='Random seed to scripts.')
parser.add_argument('--process', type=str, default='sim, rec, ana', help='Processes to be executed (sim, rec, ana).')
parser.add_argument('--numberOfLayers', dest='nlayer', type=int, default=20, help='number of layers in ML data.')
parser.add_argument('--numberOfHits', dest='nhit', type=int, default=20, help='number of hits in ML data.')
parser.add_argument('--particles', type=str, default='electron', help='Partcile names, separated by \",\".')
parser.add_argument('--amin', type=float, default=5, help='Minimum polar angle of particles (degree).')
parser.add_argument('--amax', type=float, default=175, help='Maximum polar angle of particles (degree).')
parser.add_argument('--pmin', type=float, default=0.5, help='Minimum momentum of particles (GeV).')
parser.add_argument('--pmax', type=float, default=10, help='Maximum momentum of particles (GeV).')
parser.add_argument('--compact', type=str, default=default_compact, help='Path to detector compact file.')
parser.add_argument('--uploadSizeLimit', type=float, default=10, help='Upper limit of file size (Mb) to be uploaded.')
args = parser.parse_args()
kwargs = vars(args)
if args.run_type not in default_type:
print('ERROR unsupported run type {}, must in {}'.format(args.run_type, list(default_type.keys())))
exit(1)
opt_scripts, ana_scripts = default_type[args.run_type]
gen_file = 'gen_{nametag}_{pmin}_{pmax}.hepmc'.format(**kwargs)
sim_file = 'sim_{nametag}_{pmin}_{pmax}.root'.format(**kwargs)
rec_file = 'rec_{nametag}_{pmin}_{pmax}.root'.format(**kwargs)
procs = [p.strip() for p in args.process.split(',')]
sdir = os.path.dirname(os.path.realpath(__file__))
if 'sim' in procs:
# generate particles
gen_cmd = ['python', os.path.join(sdir, 'scripts', 'gen_particles.py'), gen_file,
'-n', '{}'.format(args.nev),
'-s', '{}'.format(args.seed),
'--angmin', '{}'.format(args.amin), '--angmax', '{}'.format(args.amax),
'--pmin', '{}'.format(args.pmin), '--pmax', '{}'.format(args.pmax),
'--particles', args.particles]
subprocess.run(gen_cmd)
# simulation
sim_cmd = ['npsim',
'--part.minimalKineticEnergy', '1*TeV',
'--numberOfEvents', '{}'.format(args.nev),
'--runType', 'batch',
'--inputFiles', gen_file,
'--outputFile', sim_file,
'--compact', args.compact,
'-v', 'WARNING']
if args.seed > 0:
sim_cmd += ['--random.seed', args.seed]
print(sim_cmd)
return_code = subprocess.run(sim_cmd).returncode
if return_code is not None and return_code < 0:
print("ERROR running simulation!")
exit(1)
subprocess.run(['rootls', '-t', sim_file])
if 'rec' in procs:
# export to environment variables (used to pass arguments to the option file)
os.environ['JUGGLER_SIM_FILE'] = sim_file
os.environ['JUGGLER_REC_FILE'] = rec_file
os.environ['JUGGLER_COMPACT_PATH'] = args.compact
os.environ['JUGGLER_N_EVENTS'] = '{}'.format(args.nev)
juggler_xenv = os.path.join(os.environ.get('JUGGLER_INSTALL_PREFIX', '../local'), 'Juggler.xenv')
for opt in opt_scripts:
rec_cmd = ['xenv', '-x', juggler_xenv, 'gaudirun.py', os.path.join(sdir, 'options', opt)]
return_code = subprocess.run(rec_cmd).returncode
if return_code is not None and return_code < 0:
print('ERROR running juggler ({})!'.format(opt))
exit(1)
process = subprocess.run(['rootls', '-t', rec_file])
if 'ana' in procs:
os.makedirs('results', exist_ok=True)
for ana in ana_scripts:
ana_cmd = ['python', os.path.join(sdir, 'scripts', ana), rec_file, '-o', 'results']
return_code = subprocess.run(ana_cmd).returncode
if return_code is not None and return_code < 0:
print('ERROR running analysis ({})!'.format(ana))
exit(1)
# upload generated data files if it was small enough (< 10 Mb)
for upload in [rec_file]:
process = subprocess.run(['stat', '--format=%s', upload], stdout=subprocess.PIPE)
size = int(process.stdout)/1024./1024.
if size > args.uploadSizeLimit:
print('Skip uploading file \"{}\" because its size ({:.2f} Mb) > {:.2f} Mb'\
.format(upload, size, args.uploadSizeLimit))
else:
os.system('cp {} results/.'.format(upload))
print('Upload file \"{}\", size = {:.2f} Mb'.format(upload, size))
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment