add error events to deadtime corr

reduction/lr_reduction/DeadTimeCorrection.py

@@ -1,22 +1,20 @@
 """
     Dead time correction algorithm for single-readout detectors.
 """
-import time
-import math
-import os
 from mantid.api import *
 from mantid.simpleapi import *
 from mantid.kernel import *
 import numpy as np
 import scipy
 
-def call(InputWorkspace, DeadTime=4.2, TOFStep=100, Paralyzable=False, TOFRange=[0, 0], OutputWorkspace='correction'):
+def call(InputWorkspace, InputErrorEventsWorkspace=None, DeadTime=4.2, TOFStep=100., Paralyzable=False, TOFRange=[0, 0], OutputWorkspace='correction'):
     """
         Function to make the algorithm call similar to a normal Mantid call
     """
     algo = SingleReadoutDeadTimeCorrection()
     algo.PyInit()
     algo.setProperty("InputWorkspace", InputWorkspace)
+    algo.setProperty("InputErrorEventsWorkspace", InputErrorEventsWorkspace)
     algo.setProperty("DeadTime", DeadTime)
     algo.setProperty("TOFStep", TOFStep)
     algo.setProperty("Paralyzable", Paralyzable)
@@ -43,8 +41,10 @@ def summary(self):
     def PyInit(self):
         self.declareProperty(IEventWorkspaceProperty("InputWorkspace", "", Direction.Input),
                              "Input workspace use to compute dead time correction")
+        self.declareProperty(IEventWorkspaceProperty("InputErrorEventsWorkspace", "", Direction.Input, PropertyMode.Optional),
+                             "Input workspace with error events use to compute dead time correction")
         self.declareProperty("DeadTime", 4.2, doc="Dead time in microseconds")
-        self.declareProperty("TOFStep", 100,
+        self.declareProperty("TOFStep", 100.0,
                              doc="TOF bins to compute deadtime correction for, in microseconds")
         self.declareProperty("Paralyzable", False,
                              doc="If true, paralyzable correction will be applied, non-paralyzing otherwise")
@@ -56,6 +56,7 @@ def PyInit(self):
     def PyExec(self):
         # Event data must include error events (all triggers on the detector)
         ws_event_data = self.getProperty("InputWorkspace").value
+        ws_error_events = self.getProperty("InputErrorEventsWorkspace").value
         dead_time = self.getProperty("DeadTime").value
         tof_step = self.getProperty("TOFStep").value
         paralyzing = self.getProperty("Paralyzable").value
@@ -70,13 +71,18 @@ def PyExec(self):
 
         # Get the total number of counts on the detector for each TOF bin per pulse
         counts_ws = SumSpectra(_ws_sc)
+
+        # If we have error events, add them since those are also detector triggers
+        if ws_error_events is not None:
+            _errors = Rebin(InputWorkspace=ws_error_events, Params="%s,%s,%s" % (tof_min, tof_step, tof_max), PreserveEvents=False)
+            counts_ws += _errors
+
         t_series = np.asarray(_ws_sc.getRun()['proton_charge'].value)
         non_zero = t_series > 0
         n_pulses = np.count_nonzero(non_zero)
         rate = counts_ws.readY(0) / n_pulses
-        tof_bins = counts_ws.readX(0)
 
-         # Compute the dead time correction for each TOF bin
+        # Compute the dead time correction for each TOF bin
         if paralyzing:
             true_rate = -scipy.special.lambertw(-rate * dead_time / tof_step).real / dead_time
             corr = true_rate / (rate / tof_step)

reduction/lr_reduction/event_reduction.py

@@ -237,12 +237,17 @@ def to_dict(self):
     def get_dead_time_correction(self):
         """
             Compute dead time correction to be applied to the reflectivity curve.
+            The method will also try to load the error events from each of the
+            data files to ensure that we properly estimate the dead time correction.
         """
         # Scattering workspace
         tof_min = self._ws_sc.getTofMin()
         tof_max = self._ws_sc.getTofMax()
 
+        run_number = self._ws_sc.getRun().getProperty("run_number").value
+        error_ws = api.LoadErrorEventsNexus(run_number)
         corr_ws = DeadTimeCorrection.call(InputWorkspace=self._ws_sc,
+                                          InputErrorEventsWorkspace=error_ws,
                                           DeadTime=self.DEAD_TIME,
                                           TOFStep=self.DEAD_TIME_TOF_STEP,
                                           Paralyzable=self.paralyzable,
@@ -252,7 +257,10 @@ def get_dead_time_correction(self):
         wl_bins = corr_ws.readX(0) / self.constant
 
         # Direct beam workspace
+        run_number = self._ws_db.getRun().getProperty("run_number").value
+        error_ws = api.LoadErrorEventsNexus(run_number)
         corr_ws = DeadTimeCorrection.call(InputWorkspace=self._ws_db,
+                                          InputErrorEventsWorkspace=error_ws,
                                           DeadTime=self.DEAD_TIME,
                                           TOFStep=self.DEAD_TIME_TOF_STEP,
                                           Paralyzable=self.paralyzable,
@@ -270,11 +278,15 @@ def get_dead_time_correction(self):
 
         # Interpolate to estimate the dead time correction at the Q values we measured
         q_middle = (self.q_bins[1:] + self.q_bins[:-1]) / 2
-        
+
         dead_time_corr = np.interp(q_middle, q_values, dead_time_per_tof)
-
+
+        # Cleanup
+        api.DeleteWorkspace(corr_ws)
+        api.DeleteWorkspace(error_ws)
+
         return dead_time_corr
-        
+
     def specular(self, q_summing=False, tof_weighted=False, bck_in_q=False,
                  clean=False, normalize=True):
         """