- Changed hits_time dtype to float64. float32 is not enough for ORCA times!! (7 decimal places!)

- Add possibility to throw away a certain percentage of events in a file

file_to_hits.py

@@ -79,9 +79,9 @@ def get_event_data(event_blob, geo, do_mc_hits, use_calibrated_file, data_cuts,
         hits = hits.__array__[hits.triggered.astype(bool)]
         #hits = hits.triggered_hits # alternative, though it only works for the triggered condition!
 
-    pos_x, pos_y, pos_z = hits.pos_x.astype('float32'), hits.pos_y.astype('float32'), hits.pos_z.astype('float32')
-    hits_time = hits.time.astype('float32') # enough for the hit times in KM3NeT
-    triggered = hits.triggered.astype('float32')
+    pos_x, pos_y, pos_z = hits.pos_x, hits.pos_y, hits.pos_z
+    hits_time = hits.time
+    triggered = hits.triggered
 
     time_residual_vertex = get_time_residual_nu_interaction_mean_triggered_hits(time_interaction, hits_time, triggered)
 
@@ -90,10 +90,10 @@ def get_event_data(event_blob, geo, do_mc_hits, use_calibrated_file, data_cuts,
                             run_id, vertex_pos_x, vertex_pos_y, vertex_pos_z, time_residual_vertex], dtype=np.float64)
 
     ax = np.newaxis
-    event_hits = np.concatenate([pos_x[:, ax], pos_y[:, ax], pos_z[:, ax], hits_time[:, ax], triggered[:, ax]], axis=1)
+    event_hits = np.concatenate([pos_x[:, ax], pos_y[:, ax], pos_z[:, ax], hits_time[:, ax], triggered[:, ax]], axis=1) # dtype: np.float64
 
     if do4d[0] is True and do4d[1] == 'channel_id' or do4d[1] == 'xzt-c':
-        channel_id = hits.channel_id.astype('float32')
+        channel_id = hits.channel_id
         event_hits = np.concatenate([event_hits, channel_id[:, ax]], axis=1)
 
     return event_hits, event_track
\ No newline at end of file

h5_data_to_h5_input.py

@@ -176,6 +176,12 @@ def main(n_bins, geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4
             # Cutting events with energy < threshold (default=0) and with energy > threshold (default=200)
             continue
 
+        if data_cuts['throw_away_prob'] > 0:
+            throw_away_prob = data_cuts['throw_away_prob']
+            throw_away = np.random.choice([False, True], p=[1-throw_away_prob, throw_away_prob])
+            if throw_away is True:
+                continue
+
         # event_track: [event_id, particle_type, energy, isCC, bjorkeny, dir_x/y/z, time]
         mc_infos.append(event_track)
 
@@ -220,13 +226,19 @@ def main(n_bins, geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4
 
 
 if __name__ == '__main__':
-    # main(n_bins=(11,13,18,60), geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4d=(True, 'time'),
+    # 3-100GeV
+    main(n_bins=(11,13,18,300), geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4d=(True, 'time'),
+         timecut = ('trigger_cluster', 'all'), do_mc_hits=False, use_calibrated_file=True,
+         data_cuts = {'triggered': False, 'energy_lower_limit': 0, 'energy_upper_limit': 200, 'throw_away_prob': 0})
+
+    # 1-5GeV
+    # main(n_bins=(11,13,18,300), geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4d=(True, 'time'),
     #      timecut = ('trigger_cluster', 'all'), do_mc_hits=False, use_calibrated_file=True,
-    #      data_cuts = {'triggered': False, 'energy_lower_limit': 0, 'energy_upper_limit': 200})
+    #      data_cuts = {'triggered': False, 'energy_lower_limit': 0, 'energy_upper_limit': 3, 'throw_away_prob': 0.25})
 
-    main(n_bins=(11,13,18,60), geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4d=(True, 'channel_id'),
-         timecut = ('trigger_cluster', 'all'), do_mc_hits=False, use_calibrated_file=True,
-         data_cuts = {'triggered': False, 'energy_lower_limit': 0, 'energy_upper_limit': 200})
+    # main(n_bins=(11,13,18,60), geo_fix=True, do2d=False, do2d_pdf=(False, 10), do3d=False, do4d=(True, 'channel_id'),
+    #      timecut = ('trigger_cluster', 'all'), do_mc_hits=False, use_calibrated_file=True,
+    #      data_cuts = {'triggered': False, 'energy_lower_limit': 0, 'energy_upper_limit': 200, 'throw_away_prob': 0})
 
 
 

submit_h5_data_to_h5_input.sh

 #!/bin/bash
 #
-#PBS -l nodes=1:ppn=4:sl32g,walltime=1:30:00 # TODO change
+#PBS -l nodes=1:ppn=4:sl32g,walltime=4:00:00
 #PBS -o /home/woody/capn/mppi033h/logs/submit_h5_to_histo_${PBS_JOBID}_${PBS_ARRAYID}.out -e /home/woody/capn/mppi033h/logs/submit_h5_to_histo_${PBS_JOBID}_${PBS_ARRAYID}.err
 # first non-empty non-comment line ends PBS options
 
@@ -17,18 +17,22 @@ n=${PBS_ARRAYID}
 CodeFolder=/home/woody/capn/mppi033h/Code/OrcaSong
 cd ${CodeFolder}
 
+ParticleType=muon-CC
+#ParticleType=elec-CC
+#ParticleType=elec-NC
+#ParticleType=tau-CC
 # ----- 3-100GeV------
 #FileName=JTE.KM3Sim.gseagen.muon-CC.3-100GeV-9.1E7-1bin-3.0gspec.ORCA115_9m_2016 #muon-CC
-#FileName=JTE.KM3Sim.gseagen.elec-NC.3-100GeV-3.4E6-1bin-3.0gspec.ORCA115_9m_2016 #elec-NC
 #FileName=JTE.KM3Sim.gseagen.elec-CC.3-100GeV-1.1E6-1bin-3.0gspec.ORCA115_9m_2016 #elec-CC
+#FileName=JTE.KM3Sim.gseagen.elec-NC.3-100GeV-3.4E6-1bin-3.0gspec.ORCA115_9m_2016 #elec-NC
 #FileName=JTE.KM3Sim.gseagen.tau-CC.3.4-100GeV-2.0E8-1bin-3.0gspec.ORCA115_9m_2016 #tau-CC
-#HDFFOLDER=/home/woody/capn/mppi033h/Data/ORCA_JTE_NEMOWATER/raw_data/calibrated/without_mc_time_fix/3-100GeV/${ParticleType}
+#HDFFOLDER=/home/woody/capn/mppi033h/Data/ORCA_JTE_NEMOWATER/raw_data/calibrated/with_jte_times/3-100GeV/${ParticleType}
 # ----- 3-100GeV------
 # ----- 1-5GeV------
-#FileName=JTE.KM3Sim.gseagen.muon-CC.1-5GeV-9.2E5-1bin-1.0gspec.ORCA115_9m_2016 #muon-CC
+FileName=JTE.KM3Sim.gseagen.muon-CC.1-5GeV-9.2E5-1bin-1.0gspec.ORCA115_9m_2016 #muon-CC
 #FileName=JTE.KM3Sim.gseagen.elec-NC.1-5GeV-2.2E6-1bin-1.0gspec.ORCA115_9m_2016 #elec-NC
 #FileName=JTE.KM3Sim.gseagen.elec-CC.1-5GeV-2.7E5-1bin-1.0gspec.ORCA115_9m_2016 #elec-CC
-#HDFFOLDER=/home/woody/capn/mppi033h/Data/ORCA_JTE_NEMOWATER/raw_data/calibrated/without_mc_time_fix/1-5GeV/${ParticleType}
+HDFFOLDER=/home/woody/capn/mppi033h/Data/ORCA_JTE_NEMOWATER/raw_data/calibrated/with_jte_times/1-5GeV/${ParticleType}
 # ----- 1-5GeV------
 
 # -- denser detector study --
@@ -38,17 +42,14 @@ cd ${CodeFolder}
 # 9m: 400 m-CC files, 170 e-CC files. time/file: m-CC 40s, e-CC - --> files_per_job_m-CC = 80, files_per_job_e-CC = 36, 5 jobs
 # 6m: 400 m-CC files, 170 e-CC files. time/file: m-CC 53s, e-CC - --> files_per_job_m-CC = 80, files_per_job_e-CC = 36, 5 jobs
 # 4.5m: 400 m-CC files, 170 e-CC files. time/file: m-CC 1m3s, e-CC - --> files_per_job_m-CC = 80, files_per_job_e-CC = 36, 5 jobs
-
-vert_space=4p5m # 15m fpj m: 80, e: 36 ; 12m fpj m: , e: ;
+#vert_space=4p5m # 15m fpj m: 80, e: 36 ; 12m fpj m: , e: ;
 #FileName=JTE.KM3Sim.selectedEventsInCylinder_gseagen.muon-CC.1-20GeV-2.0E8-1bin-3.0gspec.ORCA115_h23v${vert_space}_2016 #muon-CC
-FileName=JTE.KM3Sim.selectedEventsInCylinder_gseagen.elec-CC.1-20GeV-1.0E6-1bin-3.0gspec.ORCA115_h23v${vert_space}_2016 #elec-CC
-
+#FileName=JTE.KM3Sim.selectedEventsInCylinder_gseagen.elec-CC.1-20GeV-1.0E6-1bin-3.0gspec.ORCA115_h23v${vert_space}_2016 #elec-CC
 #HDFFOLDER=/home/woody/capn/mppi033h/Data/ORCA_JTE_NEMOWATER/raw_data/calibrated/without_mc_time_fix/denser_detector_study/${vert_space}
-HDFFOLDER=/home/vault/capn/mppi033h/4p5m # 4p5m
-
+#HDFFOLDER=/home/vault/capn/mppi033h/4p5m # 4p5m
 # -- denser detector study --
 
-files_per_job=36 # total number of files per job, e.g. 10 jobs for 600: 600/10 = 60
+files_per_job=60 # total number of files per job, e.g. 10 jobs for 600: 600/10 = 60
 
 # run
 no_of_loops=$((${files_per_job}/4)) # divide by 4 cores -> e.g, 15 4-core loops needed for files_per_job=60