Skip to content
Snippets Groups Projects
Commit d0ce03e6 authored by Johannes Schumann's avatar Johannes Schumann
Browse files

Restructue and fixed filename fields

parent cc8b2059
No related branches found
No related tags found
1 merge request!1Merge python environment
Pipeline #13992 failed
Stage: test
Stage: coverage
Stage: doc
Hide whitespace changes
Inline Side-by-side
Showing
with 56 additions and 53 deletions
km3buu/output.py
+ 56
− 53
View file @ d0ce03e6
......@@ -32,30 +32,38 @@ LHE_PERT_FILENAME = "EventOutput.Pert.[0-9]{8}.lhe"
LHE_REAL_FILENAME = "EventOutput.Real.[0-9]{8}.lhe"
ROOT_PERT_FILENAME = "EventOutput.Pert.[0-9]{8}.root"
ROOT_REAL_FILENAME = "EventOutput.Real.[0-9]{8}.root"
FLUXDESCR_FILENAME = "neutrino_initialized_energyFlux.dat"
XSECTION_FILENAMES = {"all": "neutrino_absorption_cross_section_ALL.dat"}
LHE_NU_INFO_DTYPE = np.dtype([('type', np.int), ('weight', np.float64),
('mom_lepton_in_E', np.float64),
('mom_lepton_in_x', np.float64),
('mom_lepton_in_y', np.float64),
('mom_lepton_in_z', np.float64),
('mom_lepton_out_E', np.float64),
('mom_lepton_out_x', np.float64),
('mom_lepton_out_y', np.float64),
('mom_lepton_out_z', np.float64),
('mom_nucleon_in_E', np.float64),
('mom_nucleon_in_x', np.float64),
('mom_nucleon_in_y', np.float64),
('mom_nucleon_in_z', np.float64)])
LHE_NU_INFO_DTYPE = np.dtype([
("type", np.int),
("weight", np.float64),
("mom_lepton_in_E", np.float64),
("mom_lepton_in_x", np.float64),
("mom_lepton_in_y", np.float64),
("mom_lepton_in_z", np.float64),
("mom_lepton_out_E", np.float64),
("mom_lepton_out_x", np.float64),
("mom_lepton_out_y", np.float64),
("mom_lepton_out_z", np.float64),
("mom_nucleon_in_E", np.float64),
("mom_nucleon_in_x", np.float64),
("mom_nucleon_in_y", np.float64),
("mom_nucleon_in_z", np.float64),
])
FLUX_INFORMATION_DTYPE = np.dtype([("energy", np.float64),
("flux", np.float64),
("events", np.float64)])
EVENT_TYPE = {
1: 'QE',
32: 'pi neutron-background',
33: 'pi proton-background',
34: 'DIS',
35: '2p2h QE',
36: '2p2h Delta',
37: '2pi background'
1: "QE",
32: "pi neutron-background",
33: "pi proton-background",
34: "DIS",
35: "2p2h QE",
36: "2p2h Delta",
37: "2pi background",
}
......@@ -105,10 +113,21 @@ class GiBUUOutput:
f for f in listdir(self._data_path)
if isfile(join(self._data_path, f))
]
if EVENT_FILENAME in self.output_files:
setattr(self, "final_events",
FinalEvents(join(self._data_path, EVENT_FILENAME)))
self._read_xsection_file()
self._read_root_output()
self._read_flux_file()
self._read_jobcard()
def _read_root_output(self):
root_pert_regex = re.compile(ROOT_PERT_FILENAME)
self.root_pert_files = list(
filter(root_pert_regex.match, self.output_files))
root_real_regex = re.compile(ROOT_REAL_FILENAME)
self.root_real_files = list(
filter(root_real_regex.match, self.output_files))
def _read_xsection_file(self):
if XSECTION_FILENAMES["all"] in self.output_files:
setattr(
self,
......@@ -116,36 +135,19 @@ class GiBUUOutput:
read_nu_abs_xsection(
join(self._data_path, XSECTION_FILENAMES["all"])),
)
root_pert_regex = re.compile(ROOT_PERT_FILENAME)
self.lhe_pert_files = list(
filter(lhe_pert_regex.match, self.output_files))
root_real_regex = re.compile(ROOT_REAL_FILENAME)
self.lhe_real_files = list(
filter(lhe_real_regex.match, self.output_files))
jobcard_regex = re.compile('.*.job')
def _read_jobcard(self):
jobcard_regex = re.compile(".*.job")
jobcard_files = list(filter(jobcard_regex.match, self.output_files))
print(self.output_files)
if len(jobcard_files) == 1:
self._jobcard_fname = jobcard_files[0]
self.jobcard = read_jobcard(self._jobcard_fname)
else:
self.jobcard = None
def associate_jobcard(self, jobcard):
"""
Append a jobcard to the wrapped output directory
Parameters
----------
jobcard: Jobcard
Jobcard object instance
"""
self.jobcard = jobcard
# def read_pert_gibuu_output(filenames):
def _read_flux_file(self):
self.flux_data = np.loadtxt(join(self._data_path, FLUXDESCR_FILENAME),
dtype=FLUX_INFORMATION_DTYPE)
def write_detector_file(gibuu_output,
......@@ -153,6 +155,7 @@ def write_detector_file(gibuu_output,
can=(0, 476.5, 403.4),
livetime=3.156e7):
import aa, ROOT
aafile = ROOT.EventFile()
aafile.set_output(ofile)
mc_event_id = 0
......@@ -162,13 +165,13 @@ def write_detector_file(gibuu_output,
if gibuu_output.jobcard is None:
raise EnvironmentError("No jobcard provided within the GiBUU output!")
nu_type = PDGID_LOOKUP[gibuu_output.jobcard['neutrino_induced']
['flavor_id']]
nu_type = PDGID_LOOKUP[gibuu_output.jobcard["neutrino_induced"]
["flavor_id"]]
sec_lep_type = nu_type
if abs(gibuu_output.jobcard['neutrino_induced']['process_id']) == 2:
if abs(gibuu_output.jobcard["neutrino_induced"]["process_id"]) == 2:
is_cc = True
sec_lep_type -= 1
if gibuu_output.jobcard['neutrino_induced']['process_id'] < 0:
if gibuu_output.jobcard["neutrino_induced"]["process_id"] < 0:
nu_type *= -1
sec_lep_type *= -1
......@@ -209,9 +212,9 @@ def write_detector_file(gibuu_output,
nu_in_trk.type = nu_type
nu_in_trk.pos.set(*vtx_pos)
nu_in_trk.dir.set(*direction)
nu_in_trk.E = event_info['mom_lepton_in_E'].item()
nu_in_trk.E = event_info["mom_lepton_in_E"].item()
nu_in_trk.t = timestamp
nu_in_trk.setusr('cc', is_cc)
nu_in_trk.setusr("cc", is_cc)
aafile.evt.mc_trks.push_back(nu_in_trk)
lep_out_trk = ROOT.Trk()
......@@ -220,11 +223,11 @@ def write_detector_file(gibuu_output,
lep_out_trk.type = sec_lep_type
lep_out_trk.pos.set(*vtx_pos)
mom = np.array(event_info[[
'mom_lepton_out_x', 'mom_lepton_out_y', 'mom_lepton_out_z'
"mom_lepton_out_x", "mom_lepton_out_y", "mom_lepton_out_z"
]].item())
p_dir = R.apply(mom / np.linalg.norm(mom))
lep_out_trk.dir.set(*p_dir)
lep_out_trk.E = event_info['mom_lepton_out_E'].item()
lep_out_trk.E = event_info["mom_lepton_out_E"].item()
lep_out_trk.t = timestamp
aafile.evt.mc_trks.push_back(lep_out_trk)
......
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