From d0ce03e60219506c651a77d486ebd3b6c8db7815 Mon Sep 17 00:00:00 2001 From: Johannes Schumann <johannes.schumann@fau.de> Date: Thu, 24 Sep 2020 02:32:15 +0200 Subject: [PATCH] Restructue and fixed filename fields --- km3buu/output.py | 109 ++++++++++++++++++++++++----------------------- 1 file changed, 56 insertions(+), 53 deletions(-) diff --git a/km3buu/output.py b/km3buu/output.py index eda2658..013ce49 100644 --- a/km3buu/output.py +++ b/km3buu/output.py @@ -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) -- GitLab