From a8fe4408e1b54123906afedba1293eb40be34886 Mon Sep 17 00:00:00 2001
From: Johannes Schumann <johannes.schumann@fau.de>
Date: Fri, 31 Mar 2023 20:13:07 +0000
Subject: [PATCH] Cleanups
---
km3buu/output.py | 130 +++++++++--------------------------------------
1 file changed, 25 insertions(+), 105 deletions(-)
diff --git a/km3buu/output.py b/km3buu/output.py
index 9e34538..357cc10 100644
--- a/km3buu/output.py
+++ b/km3buu/output.py
@@ -26,6 +26,7 @@ import scipy.constants as constants
from scipy.optimize import curve_fit
import mendeleev
from datetime import datetime
+import km3io
from .physics import visible_energy_fraction
from .jobcard import Jobcard, read_jobcard, PDGID_LOOKUP
@@ -178,88 +179,7 @@ EMPTY_KM3NET_HEADER_DICT = {
"primary": "0",
}
-PARTICLE_MC_STATUS = {
- "TRK_MOTHER_UNDEFINED":
- # mother id was not defined for this MC track (all reco tracks have this value)
- -1,
- "TRK_MOTHER_NONE":
- -2, # mother id of a particle if it has no parent
- "TRK_ST_UNDEFINED":
- # status was not defined for this MC track (all reco tracks have this value)
- 0,
- "TRK_ST_FINALSTATE":
- # particle to be tracked by detector-level MC ('track_in' tag in evt files from gseagen, genhen, mupage).
- 1,
- "TRK_ST_PRIMARYNEUTRINO":
- # initial state neutrino ('neutrino' tag in evt files from gseagen and genhen).
- 100,
- "TRK_ST_PRIMARYCOSMIC":
- # initial state cosmic ray ('track_primary' tag in evt files from corant).
- 200,
- "TRK_ST_ININUCLEI":
- 5, # Initial state nuclei (gseagen)
- "TRK_ST_INTERSTATE":
- 2, # Intermidiate state particles produced in hadronic showers (gseagen)
- "TRK_ST_DECSTATE":
- 3, # Short-lived particles that are forced to decay, like pi0 (gseagen)
- "TRK_ST_NUCTGT":
- 11, # Nucleon target (gseagen)
- "TRK_ST_PREHAD":
- 12, # DIS pre-fragmentation hadronic state (gseagen)
- "TRK_ST_PRERES":
- 13, # resonant pre-decayed state (gseagen)
- "TRK_ST_HADNUC":
- 14, # Hadrons inside the nucleus before FSI (gseagen)
- "TRK_ST_NUCLREM":
- 15, # Low energy nuclear fragments (gseagen)
- "TRK_ST_NUCLCLT":
- 16, # For composite nucleons before phase space decay (gseagen)
- "TRK_ST_FAKECORSIKA":
- # fake particle from corant/CORSIKA to add parent information (gseagen)
- 21,
- "TRK_ST_FAKECORSIKA_DEC_MU_START":
- 22, # fake particle from CORSIKA: decaying mu at start (gseagen)
- "TRK_ST_FAKECORSIKA_DEC_MU_END":
- 23, # fake particle from CORSIKA: decaying mu at end (gseagen)
- "TRK_ST_FAKECORSIKA_ETA_2GAMMA":
- 24, # fake particle from CORSIKA: eta -> 2 gamma (gseagen)
- "TRK_ST_FAKECORSIKA_ETA_3PI0":
- 25, # fake particle from CORSIKA: eta -> 3 pi0 (gseagen)
- "TRK_ST_FAKECORSIKA_ETA_PIP_PIM_PI0":
- 26, # fake particle from CORSIKA: eta -> pi+ pi- pi0 (gseagen)
- "TRK_ST_FAKECORSIKA_ETA_2PI_GAMMA":
- 27, # fake particle from CORSIKA: eta -> pi+ pi- gamma (gseagen)
- "TRK_ST_FAKECORSIKA_CHERENKOV_GAMMA":
- # fake particle from CORSIKA: Cherenkov photons on particle output file (gseagen)
- 28,
- "TRK_ST_PROPLEPTON":
- 1001, # lepton propagated that reaches the can (gseagen)
- "TRK_ST_PROPDECLEPTON":
- 2001 # lepton propagated and decayed before got to the can (gseagen)
-}
-
-W2LIST_LOOKUP = {
- "PS": 0,
- "EG": 1,
- "XSEC_MEAN": 2,
- "COLUMN_DEPTH": 3,
- "P_EARTH": 4,
- "WATER_INT_LEN": 5,
- "BX": 7,
- "BY": 8,
- "ICHAN": 9,
- "CC": 10,
- "XSEC": 13,
- "DXSEC": 14,
- "TARGETA": 15,
- "TARGETZ": 16,
- "VERINCAN": 17,
- "LEPINCAN": 18,
- "GIBUU_WEIGHT": 23,
- "GIBUU_SCAT_TYPE": 24
-}
-
-W2LIST_LENGTH = max(W2LIST_LOOKUP.values()) + 1
+W2LIST_LENGTH = max(km3io.definitions.w2list_km3buu.values()) + 1
GIBUU_FIELDNAMES = [
'weight', 'barcode', 'Px', 'Py', 'Pz', 'E', 'evType', 'lepIn_E',
@@ -769,27 +689,27 @@ def write_detector_file(gibuu_output,
evt.w.push_back(-1.0) # w3 (= w2*flux)
# Event Information (w2list)
evt.w2list.resize(W2LIST_LENGTH)
- evt.w2list[W2LIST_LOOKUP["PS"]] = global_generation_weight
- evt.w2list[W2LIST_LOOKUP["EG"]] = gibuu_output.flux_index
- evt.w2list[W2LIST_LOOKUP["XSEC_MEAN"]] = mean_xsec_func(
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_PS"]] = global_generation_weight
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_EG"]] = gibuu_output.flux_index
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_XSEC_MEAN"]] = mean_xsec_func(
event.lepIn_E)
- evt.w2list[W2LIST_LOOKUP["XSEC"]] = event.xsec
- evt.w2list[W2LIST_LOOKUP["TARGETA"]] = gibuu_output.A
- evt.w2list[W2LIST_LOOKUP["TARGETZ"]] = gibuu_output.Z
- evt.w2list[W2LIST_LOOKUP["BX"]] = bjorkenx[mc_event_id]
- evt.w2list[W2LIST_LOOKUP["BY"]] = bjorkeny[mc_event_id]
- evt.w2list[W2LIST_LOOKUP["CC"]] = ichan
- evt.w2list[W2LIST_LOOKUP["ICHAN"]] = SCATTERING_TYPE_TO_GENIE[
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_XSEC"]] = event.xsec
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_TARGETA"]] = gibuu_output.A
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_TARGETZ"]] = gibuu_output.Z
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_BX"]] = bjorkenx[mc_event_id]
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_BY"]] = bjorkeny[mc_event_id]
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_CC"]] = ichan
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_ICHAN"]] = SCATTERING_TYPE_TO_GENIE[
event.evType]
- evt.w2list[W2LIST_LOOKUP["VERINCAN"]] = 1
- evt.w2list[W2LIST_LOOKUP["LEPINCAN"]] = 1
- evt.w2list[W2LIST_LOOKUP["GIBUU_WEIGHT"]] = event.weight
- evt.w2list[W2LIST_LOOKUP["GIBUU_SCAT_TYPE"]] = event.evType
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_VERINCAN"]] = 1
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_LEPINCAN"]] = 1
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_GIBUU_WEIGHT"]] = event.weight
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_GIBUU_SCAT_TYPE"]] = event.evType
# TODO
- evt.w2list[W2LIST_LOOKUP["DXSEC"]] = np.nan
- evt.w2list[W2LIST_LOOKUP["COLUMN_DEPTH"]] = np.nan
- evt.w2list[W2LIST_LOOKUP["P_EARTH"]] = np.nan
- evt.w2list[W2LIST_LOOKUP["WATER_INT_LEN"]] = np.nan
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_DXSEC"]] = np.nan
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_COLUMN_DEPTH"]] = np.nan
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_P_EARTH"]] = np.nan
+ evt.w2list[km3io.definitions.w2list_km3buu["W2LIST_KM3BUU_WATER_INT_LEN"]] = np.nan
# Vertex Position
vtx_pos = np.array(geometry.random_pos())
@@ -816,7 +736,7 @@ def write_detector_file(gibuu_output,
nu_in_trk.dir.set(*direction)
nu_in_trk.E = event.lepIn_E
nu_in_trk.t = timestamp
- nu_in_trk.status = PARTICLE_MC_STATUS["TRK_ST_PRIMARYNEUTRINO"]
+ nu_in_trk.status = km3io.definitions.trkmembers["TRK_ST_PRIMARYNEUTRINO"]
evt.mc_trks.push_back(nu_in_trk)
lep_out_trk = ROOT.Trk()
@@ -832,20 +752,20 @@ def write_detector_file(gibuu_output,
lep_out_trk.t = timestamp
if tau_secondaries is not None:
- lep_out_trk.status = PARTICLE_MC_STATUS["TRK_ST_UNDEFINED"]
+ lep_out_trk.status = km3io.definitions.trkmembers["TRK_ST_UNDEFINED"]
else:
- lep_out_trk.status = PARTICLE_MC_STATUS["TRK_ST_FINALSTATE"]
+ lep_out_trk.status = km3io.definitions.trkmembers["TRK_ST_FINALSTATE"]
evt.mc_trks.push_back(lep_out_trk)
if tau_secondaries is not None:
event_tau_sec = tau_secondaries[mc_event_id]
add_particles(event_tau_sec, vtx_pos, R, mc_trk_id, timestamp,
- PARTICLE_MC_STATUS["TRK_ST_FINALSTATE"])
+ km3io.definitions.trkmembers["TRK_ST_FINALSTATE"])
mc_trk_id += len(event_tau_sec.E)
add_particles(event, vtx_pos, R, mc_trk_id, timestamp,
- PARTICLE_MC_STATUS["TRK_ST_FINALSTATE"])
+ km3io.definitions.trkmembers["TRK_ST_FINALSTATE"])
tree.Fill()
for k, v in geometry.header_entries(mc_event_id + 1).items():
--
GitLab