diff --git a/ecal/emcal_electrons.sh b/ecal/emcal_electrons.sh index 5a8b182330abde55d5830d074e8802f2b11516f1..eac567922c74eb37afab830aee48df824274695b 100644 --- a/ecal/emcal_electrons.sh +++ b/ecal/emcal_electrons.sh @@ -70,7 +70,8 @@ if [[ "$?" -ne "0" ]] ; then echo "ERROR running npdet" exit 1 fi -xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv gaudirun.py ../ecal/options/crystal_calorimeter_reco.py +xenv -x ${JUGGLER_INSTALL_PREFIX}/Juggler.xenv \ + gaudirun.py ../ecal/options/crystal_calorimeter_reco.py if [[ "$?" -ne "0" ]] ; then echo "ERROR running juggler" exit 1 diff --git a/ecal/options/crystal_calorimeter_reco.py b/ecal/options/crystal_calorimeter_reco.py index d6ef2b61ad72d25f2d1924359dce66351b870076..fa85b808b7f964c69c4b4d0122cf10fdc648224d 100644 --- a/ecal/options/crystal_calorimeter_reco.py +++ b/ecal/options/crystal_calorimeter_reco.py @@ -10,9 +10,22 @@ if "JUGGLER_DETECTOR" in os.environ : detector_name = str(os.environ["JUGGLER_DETECTOR"]) # 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"]) +input_sim_file = "jug_input.root" +if "JUGGLER_SIM_FILE" in os.environ : + input_sim_file = str(os.environ["JUGGLER_SIM_FILE"]) +else : + print(" ERROR : JUGGLER_SIM_FILE not set" ) + +output_rec_file = "jug_rec.root" +if "JUGGLER_REC_FILE" in os.environ : + output_rec_file = str(os.environ["JUGGLER_REC_FILE"]) +else : + print(" ERROR : JUGGLER_REC_FILE not set" ) + + +n_events = 100 +if "JUGGLER_N_EVENTS" in os.environ : + n_events = str(os.environ["JUGGLER_N_EVENTS"]) geo_service = GeoSvc("GeoSvc", detectors=["{}.xml".format(detector_name)]) podioevent = EICDataSvc("EventDataSvc", inputs=[input_sim_file], OutputLevel=DEBUG) @@ -31,15 +44,28 @@ from Configurables import Jug__Reco__ClusterRecoCoG as RecoCoG podioinput = PodioInput("PodioReader", collections=["mcparticles","CrystalEcalHits"], 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",OutputLevel=DEBUG) +copier = MCCopier("MCCopier", inputCollection="mcparticles", outputCollection="mcparticles2",OutputLevel=DEBUG) calcopier = CalCopier("CalCopier", inputCollection="CrystalEcalHits", outputCollection="CrystalEcalHits2",OutputLevel=DEBUG) -emcaldigi = CrystalEndcapsDigi("ecal_digi", inputHitCollection="CrystalEcalHits", outputHitCollection="RawDigiEcalHits") -emcalreco = CrystalEndcapsReco("ecal_reco", inputHitCollection="RawDigiEcalHits", outputHitCollection="RecoEcalHits", - minModuleEdep=1.0*units.MeV) -emcalcluster = IslandCluster("emcal_cluster", inputHitCollection="RecoEcalHits", outputClusterCollection="EcalClusters", - minClusterCenterEdep=30*units.MeV, groupRange=2.0) -clusterreco = RecoCoG("cluster_reco", clusterCollection="EcalClusters", logWeightBase=4.2, moduleDimZName="CrystalBox_z_length") +emcaldigi = CrystalEndcapsDigi("ecal_digi", + inputHitCollection="CrystalEcalHits", + outputHitCollection="RawDigiEcalHits") + +emcalreco = CrystalEndcapsReco("ecal_reco", + inputHitCollection="RawDigiEcalHits", + outputHitCollection="RecoEcalHits", + minModuleEdep=0.00001*units.MeV) + +emcalcluster = IslandCluster("emcal_cluster", + inputHitCollection="RecoEcalHits", + outputClusterCollection="EcalClusters", + minClusterCenterEdep=0.00001*units.MeV, + groupRange=2.0) + +clusterreco = RecoCoG("cluster_reco", + clusterCollection="EcalClusters", + logWeightBase=4.2, + moduleDimZName="CrystalBox_z_length") out = PodioOutput("out", filename=output_rec_file) diff --git a/ecal/scripts/emcal_electrons.cxx b/ecal/scripts/emcal_electrons.cxx index 9ef64b872ee75ed1ac32915aff44241f972e6d28..7707f24a0207a19dd075d8c8f0f1a5778519461f 100644 --- a/ecal/scripts/emcal_electrons.cxx +++ b/ecal/scripts/emcal_electrons.cxx @@ -22,7 +22,7 @@ using namespace HepMC3; void emcal_electrons(int n_events = 1e2, double e_start = 1.0, double e_end = 1.0, const char* out_fname = "./data/emcal_electron_0GeVto30GeV_100kEvt.hepmc") { - double cos_theta_min = std::cos(M_PI * (120.0 / 180.0)); + double cos_theta_min = std::cos(M_PI * (155.0 / 180.0)); double cos_theta_max = std::cos(M_PI); WriterAscii hepmc_output(out_fname);