Implement splitting chopper pulses into neutron pulses

eos/event_analysis.py

@@ -9,7 +9,7 @@ from typing import Tuple
 
 from . import const
 from .event_data_types import EventDataAction, EventDatasetProtocol, append_fields, EVENT_BITMASKS
-from .helpers import filter_project_x, merge_frames, extract_walltime, add_log_to_pulses
+from .helpers import filter_project_x, merge_frames, extract_walltime, add_log_to_pulses, merge_frames_w_index
 from .instrument import Detector
 from .options import IncidentAngle
 from .header import Header
@@ -25,8 +25,9 @@ class ExtractWalltime(EventDataAction):
         dataset.data.events = new_events
 
 class MergeFrames(EventDataAction):
-    def __init__(self, lamdaCut=None):
+    def __init__(self, lamdaCut=None, preciseTime=False):
         self.lamdaCut=lamdaCut
+        self.preciseTime=preciseTime
 
     def perform_action(self, dataset: EventDatasetProtocol)->None:
         if self.lamdaCut is None:
@@ -36,7 +37,33 @@ class MergeFrames(EventDataAction):
         tofCut = lamdaCut*dataset.geometry.chopperDetectorDistance/const.hdm*1e-13
         total_offset = (tofCut +
                         dataset.timing.tau * (dataset.timing.ch1TriggerPhase + dataset.timing.chopperPhase/2)/180)
-        dataset.data.events.tof = merge_frames(dataset.data.events.tof, tofCut, dataset.timing.tau, total_offset)
+        if self.preciseTime and 'wallTime' in dataset.data.events.dtype.names:
+            d = dataset.data
+            # put events into precise sub-frame
+            d.events.tof, subframes = merge_frames_w_index(d.events.tof, tofCut, dataset.timing.tau, total_offset)
+            # add a sub-pulse time 1-tau before and after each existing time
+            utimes, uidxs = np.unique(d.events.wallTime, return_inverse=True)
+            inter_times = np.empty(2*utimes.shape[0]+1, dtype=d.events.wallTime.dtype)
+            tauns = dataset.timing.tau*1e9
+            inter_times[0] = utimes[0]-tauns
+            inter_times[1::2] = utimes
+            inter_times[2::2] = utimes+tauns
+            # use subframe indices to sort existing events into new times
+            d.events.wallTime = inter_times[2*uidxs+subframes+1]
+            # expand pulses array with additional sub-frames
+            new_pulses = np.recarray(2*d.pulses.shape[0]+1, dtype=d.pulses.dtype)
+            new_pulses[0] = d.pulses[0]
+            new_pulses[1::2] = d.pulses
+            new_pulses[2::2] = d.pulses
+            new_pulses.time[0] = d.pulses.time[0]-tauns
+            new_pulses.time[2::2] = d.pulses.time+tauns
+            new_pulses.monitor /= 2.0 # ~preserve total monitor counts
+            d.pulses = new_pulses
+        else:
+            if self.preciseTime:
+                logging.error("    Trying to use precise time frame separation while wallTime is not extracted, yet!")
+            dataset.data.events.tof = merge_frames(dataset.data.events.tof, tofCut, dataset.timing.tau, total_offset)
+
 
 
 class AnalyzePixelIDs(EventDataAction):

eos/helpers.py

@@ -8,7 +8,10 @@ from .event_data_types import EventDatasetProtocol, append_fields
 try:
     from .helpers_numba import merge_frames, extract_walltime, filter_project_x, calculate_derived_properties_focussing
 except ImportError:
+    import logging
+    logging.warning('Cannot import numba enhanced functions, is it installed?')
     from .helpers_fallback import merge_frames, extract_walltime, filter_project_x, calculate_derived_properties_focussing
+from .helpers_fallback import merge_frames_w_index
 
 def add_log_to_pulses(key, dataset: EventDatasetProtocol):
     """

eos/helpers_fallback.py

@@ -8,6 +8,17 @@ def merge_frames(tof_e, tofCut, tau, total_offset):
     # tof shifted to 1 frame
     return np.remainder(tof_e-(tofCut-tau), tau)+total_offset
 
+def merge_frames_w_index(tof_e, tofCut, tau, total_offset):
+    """
+    Version of merge frames that also returns a frame index for each pulse:
+    0  - belongs to the frame it was measured in
+    -1 - arrived in this frame but belongs to the previous neutron pulse
+    1  - belongs to the second neutron pulse of the original frame
+    """
+    new_tof = merge_frames(tof_e, tofCut, tau, total_offset)
+    frame_idx = np.floor_divide(tof_e-tofCut, tau).astype(int)
+    return new_tof, frame_idx
+
 def extract_walltime(tof_e, dataPacket_p, dataPacketTime_p):
     output = np.empty(np.shape(tof_e)[0], dtype=np.int64)
     for i in range(len(dataPacket_p)-1):

eos/options.py

@@ -480,6 +480,16 @@ class ReflectivityReductionConfig(ArgParsable):
                 },
             )
 
+    preciseTiming: bool = field(
+            default=False,
+            metadata={
+                'short': 'pt',
+                'group': 'data manicure',
+                'help':  're-assign events to actual neutron pulses => '
+                         '2xbetter filter resolution but extra computation (experimental)',
+                },
+            )
+
 
 class OutputFomatOption(StrEnum):
     Rqz_ort = "Rqz.ort"

eos/reduction_reflectivity.py

@@ -65,7 +65,7 @@ class ReflectivityReduction:
             # the filtering only makes sense if using actual monitor data, not time
             self.dataevent_actions |= eh.FilterMonitorThreshold(self.config.experiment.lowCurrentThreshold)
         self.dataevent_actions |= eh.FilterStrangeTimes()
-        self.dataevent_actions |= ea.MergeFrames()
+        self.dataevent_actions |= ea.MergeFrames(preciseTime=self.config.reduction.preciseTiming)
         self.dataevent_actions |= ea.AnalyzePixelIDs(self.config.experiment.yRange)
         self.dataevent_actions |= eh.TofTimeCorrection(self.config.experiment.incidentAngle==IncidentAngle.alphaF)
         self.dataevent_actions |= ea.CalculateWavelength(self.config.experiment.lambdaRange)