Making some cleaner paper plots with timing information

python/LIController.py

@@ -1,6 +1,7 @@
 import h5py
 import numpy as np
 import awkward as ak
+import time
 
 from . import utilities as _utilities
 from . import detector as _detector
@@ -27,6 +28,8 @@ def __init__(self, events_to_inject, experiment, seed=0):
         :param int seed: Optional random number generator seed
         """
 
+        self.global_start = time.time()
+
         self.resources_dir = _util.resource_package_dir()
 
         # Initialize a random number generator
@@ -382,10 +385,16 @@ def GenerateEvents(self, N=None, fill_tables_at_exit=True):
         if N is None:
             N = self.events_to_inject
         count = 0
+        self.gen_times,self.global_times = [],[]
+        prev_time = time.time()
         while (self.injector.InjectedEvents() < self.events_to_inject) and (count < N):
             print("Injecting Event %d/%d  " % (count, N), end="\r")
             tree = self.injector.GenerateEvent()
             self.events.append(tree)
+            t = time.time()
+            self.gen_times.append(t-prev_time)
+            self.global_times.append(t-self.global_start)
+            prev_time = t
             count += 1
         if hasattr(self, "DN_processes"):
             self.DN_processes.SaveCrossSectionTables(fill_tables_at_exit=fill_tables_at_exit)
@@ -396,6 +405,8 @@ def SaveEvents(self, filename, fill_tables_at_exit=True, hdf5=True, parquet=True
         # A dictionary containing each dataset we'd like to save
         datasets = {
             "event_weight":[], # weight of entire event
+            "event_gen_time":[], # generation time of each event
+            "event_global_time":[], # global time of each event
             "num_interactions":[], # number of interactions per event
             "vertex":[], # vertex of each interaction in an event
             "in_fiducial":[], # whether or not each vertex is in the fiducial volume
@@ -410,6 +421,8 @@ def SaveEvents(self, filename, fill_tables_at_exit=True, hdf5=True, parquet=True
         for ie, event in enumerate(self.events):
             print("Saving Event %d/%d  " % (ie, len(self.events)), end="\r")
             datasets["event_weight"].append(self.weighter.EventWeight(event))
+            datasets["event_gen_time"].append(self.gen_times[ie])
+            datasets["event_global_time"].append(self.global_times[ie])
             # add empty lists for each per interaction dataset
             for k in ["vertex",
                       "in_fiducial",

python/LIDarkNews.py

@@ -204,7 +204,7 @@ def __init__(
         table_dir=None,  # table to store
         tolerance=1e-6,  # supposed to represent machine epsilon
         interp_tolerance=5e-2,  # relative interpolation tolerance
-        always_interpolate=False, # bool whether to always interpolate the total/differential cross section
+        always_interpolate=True, # bool whether to always interpolate the total/differential cross section
     ):
         DarkNewsCrossSection.__init__(self)  # C++ constructor
 
@@ -380,7 +380,7 @@ def _query_interpolation_table(self, inputs, mode):
             # check if energy is within table range
 
             if interpolator is None or inputs[0] > interp_table[-1,0]:
-                print("Requested interpolation at %2.2f GeV. Either this above the table boundary or the interpolator doesn't yet exist. Filling %s table"%(inputs[0],mode))
+                print("Requested interpolation at %2.2f GeV. Either this is above the table boundary or the interpolator doesn't yet exist. Filling %s table"%(inputs[0],mode))
                 n = self.FillInterpolationTables(total=(mode=="total"),
                                                  diff=(mode=="differential"),
                                                  Emax = (1+self.interp_tolerance)*inputs[0])
@@ -390,9 +390,11 @@ def _query_interpolation_table(self, inputs, mode):
             elif inputs[0] < interp_table[0,0]:
                 print("Requested interpolation at %2.2f GeV below table boundary. Requring calculation"%inputs[0])
                 return 0
-            val = interpolator(inputs)
+            val = max(0,interpolator(inputs))
             if val<0:
-                print("WARNING: negative interpolated value for %s cross section at,"%mode,inputs) 
+                print("WARNING: negative interpolated value for %s-%s %s cross section at,"%(self.ups_case.nuclear_target.name,
+                                                                                             self.ups_case.scattering_regime,
+                                                                                             mode),inputs) 
             return val
 
         UseSinglePoint, Interpolate, closest_idx = self._interpolation_flags(
@@ -407,7 +409,7 @@ def _query_interpolation_table(self, inputs, mode):
         elif Interpolate:
             return interpolator(inputs)
         else:
-            return 0
+            return -1
 
     def FillTableAtEnergy(self, E, total=True, diff=True, factor=0.8):
         num_added_points = 0
@@ -581,7 +583,7 @@ def DifferentialCrossSection(self, arg1, target=None, energy=None, Q2=None):
         if self.always_interpolate:
             # Check if we can interpolate
             val = self._query_interpolation_table([energy, z], mode="differential")
-            if val > 0:
+            if val >= 0:
                 # we have recovered the differential cross section from the interpolation table
                 return val
 
@@ -632,7 +634,7 @@ def TotalCrossSection(self, arg1, energy=None, target=None):
 
         # Check if we can interpolate
         val = self._query_interpolation_table([energy], mode="total")
-        if val > 0:
+        if val >= 0:
             # we have recovered the cross section from the interpolation table
             return val
 

resources/Examples/Example2/DipolePortal_CCM.py

@@ -23,7 +23,7 @@
 }
 
 # Number of events to inject
-events_to_inject = 10000
+events_to_inject = 100000
 
 # Expeirment to run
 experiment = "CCM"
@@ -40,7 +40,8 @@
     xs_path,
     "Dipole_M%2.2e_mu%2.2e" % (model_kwargs["m4"], model_kwargs["mu_tr_mu4"]),
 )
-controller.InputDarkNewsModel(primary_type, table_dir, **model_kwargs, **xs_kwargs)
+controller.InputDarkNewsModel(primary_type, table_dir,
+                              **model_kwargs, **xs_kwargs)
 
 # Primary distributions
 primary_injection_distributions = {}

resources/Examples/Example2/DipolePortal_MINERvA.py

@@ -18,7 +18,7 @@
 
 # cross section class arguments
 xs_kwargs = {
-    "always_interpolate": True
+    "always_interpolate": True,
 }
 
 # Number of events to inject

resources/Examples/Example2/DipolePortal_MiniBooNE.py

@@ -39,8 +39,7 @@
     xs_path,
     "Dipole_M%2.2e_mu%2.2e" % (model_kwargs["m4"], model_kwargs["mu_tr_mu4"]),
 )
-controller.InputDarkNewsModel(primary_type, table_dir,
-                              **model_kwargs, **xs_kwargs)
+controller.InputDarkNewsModel(primary_type, table_dir, **model_kwargs, **xs_kwargs)
 
 # Primary distributions
 primary_injection_distributions = {}