Merge pull request #102 from cta-observatory/delta_t_option

Add options to disable delta t and spike correction

ctapipe_io_lst/__init__.py

@@ -178,7 +178,17 @@ class LSTEventSource(EventSource):
 
     calibrate_flatfields_and_pedestals = Bool(
         default_value=True,
-        help='To be set to True for calibration processing'
+        help='If True, flat field and pedestal events are also calibrated.'
+    ).tag(config=True)
+
+    apply_drs4_corrections = Bool(
+        default_value=True,
+        help=(
+            'Apply DRS4 corrections.'
+            ' If True, this will fill R1 waveforms with the corrections applied'
+            ' Use the options for the LSTR0Corrections to configure which'
+            ' corrections are applied'
+        ),
     ).tag(config=True)
 
     classes = [PointingSource, EventTimeCalculator, LSTR0Corrections]
@@ -349,7 +359,7 @@ def _generator(self):
             self.fill_pointing_info(array_event)
 
             # apply low level corrections
-            if self.r0_r1_calibrator.drs4_pedestal_path.tel[self.tel_id] is not None:
+            if self.apply_drs4_corrections:
                 self.r0_r1_calibrator.apply_drs4_corrections(array_event)
 
                 # flat field tagging is performed on r1 data, so can only

ctapipe_io_lst/calibration.py

@@ -85,7 +85,10 @@ class LSTR0Corrections(TelescopeComponent):
     """
     offset = IntTelescopeParameter(
         default_value=400,
-        help='Define the offset of the baseline'
+        help=(
+            'Define the offset of the baseline'
+            ', set to 0 to disable offset subtraction'
+        )
     ).tag(config=True)
 
     r1_sample_start = IntTelescopeParameter(
@@ -104,7 +107,10 @@ class LSTR0Corrections(TelescopeComponent):
         trait=Path(exists=True, directory_ok=False),
         allow_none=True,
         default_value=None,
-        help='Path to the LST pedestal file',
+        help=(
+            'Path to the LST pedestal file'
+            ', required when `apply_drs4_pedestal_correction=True`'
+        ),
     ).tag(config=True)
 
     calibration_path = Path(
@@ -134,6 +140,24 @@ class LSTR0Corrections(TelescopeComponent):
         help='Set to False to keep both gains.'
     ).tag(config=True)
 
+    apply_drs4_pedestal_correction = Bool(
+        default_value=True,
+        help=(
+            'Set to False to disable drs4 pedestal correction.'
+            ' Providing the drs4_pedestal_path is required to perform this calibration'
+        ),
+    ).tag(config=True)
+
+    apply_timelapse_correction = Bool(
+        default_value=True,
+        help='Set to False to disable drs4 timelapse correction'
+    ).tag(config=True)
+
+    apply_spike_correction = Bool(
+        default_value=True,
+        help='Set to False to disable drs4 spike correction'
+    ).tag(config=True)
+
     add_calibration_timeshift = Bool(
         default_value=True,
         help=(
@@ -198,9 +222,14 @@ def apply_drs4_corrections(self, event: ArrayEventContainer):
             r1.waveform = event.r0.tel[tel_id].waveform.astype(np.float32)
 
             # apply drs4 corrections
-            self.subtract_pedestal(event, tel_id)
-            self.time_lapse_corr(event, tel_id)
-            self.interpolate_spikes(event, tel_id)
+            if self.apply_drs4_pedestal_correction:
+                self.subtract_pedestal(event, tel_id)
+
+            if self.apply_timelapse_correction:
+                self.time_lapse_corr(event, tel_id)
+
+            if self.apply_spike_correction:
+                self.interpolate_spikes(event, tel_id)
 
             # remove samples at beginning / end of waveform
             start = self.r1_sample_start.tel[tel_id]
@@ -226,7 +255,12 @@ def calibrate(self, event: ArrayEventContainer):
 
         for tel_id in event.r0.tel:
             r1 = event.r1.tel[tel_id]
-            waveform = r1.waveform
+            # if `apply_drs4_corrections` is False, we did not fill in the
+            # waveform yet.
+            if r1.waveform is None:
+                waveform = event.r0.tel[tel_id].waveform.astype(np.float32)
+            else:
+                waveform = r1.waveform
 
             # do gain selection before converting to pe
             # like eventbuilder will do
@@ -378,6 +412,12 @@ def _get_drs4_pedestal_data(path, offset=0):
         The result is cached so we can repeatedly call this method
         using the configured path without reading it each time.
         """
+        if path is None:
+            raise ValueError(
+                "DRS4 pedestal correction requested"
+                " but no file provided for telescope"
+            )
+
         pedestal_data = np.empty(
             (N_GAINS, N_PIXELS_MODULE * N_MODULES, N_CAPACITORS_PIXEL + N_SAMPLES),
             dtype=np.int16

ctapipe_io_lst/tests/test_calib.py

@@ -29,7 +29,7 @@ def test_get_first_capacitor():
     )
 
     tel_id = 1
-    source = LSTEventSource(test_r0_calib_path)
+    source = LSTEventSource(test_r0_calib_path, apply_drs4_corrections=False)
     event = next(iter(source))
 
     first_capacitor_id = event.lst.tel[tel_id].evt.first_capacitor_id

ctapipe_io_lst/tests/test_event_time.py

@@ -56,6 +56,7 @@ def test_ucts_jumps():
     lst = event.lst.tel[tel_id]
     lst.evt.extdevices_presence = 0b1111_1111
     lst.svc.module_ids = np.arange(N_MODULES)
+    lst.evt.module_status = np.ones(N_MODULES)
 
     subarray = LSTEventSource.create_subarray(geometry_version=4, tel_id=1)
 
@@ -173,6 +174,8 @@ def test_extract_reference_values(caplog):
     lst.evt.pps_counter = np.full(N_MODULES, 100)
     lst.evt.tenMHz_counter = np.zeros(N_MODULES)
     lst.evt.tenMHz_counter[CENTRAL_MODULE] = 2
+    lst.evt.module_status = np.ones(N_MODULES)
+
 
     # actually test it
     time = time_calculator(tel_id, event)
@@ -220,6 +223,7 @@ def test_no_reference_values_no_ucts(caplog):
     lst.evt.tenMHz_counter = np.zeros(N_MODULES)
     lst.evt.tenMHz_counter[CENTRAL_MODULE] = 2
     lst.svc.date = run_start.unix
+    lst.evt.module_status = np.ones(N_MODULES)
 
     # actually test it
     time = time_calculator(tel_id, event)

ctapipe_io_lst/tests/test_lsteventsource.py

@@ -26,6 +26,7 @@ def test_loop_over_events():
     source = LSTEventSource(
         input_url=test_r0_path,
         max_events=n_events,
+        apply_drs4_corrections=False,
     )
 
     for i, event in enumerate(source):
@@ -45,7 +46,10 @@ def test_loop_over_events():
 def test_multifile():
     from ctapipe_io_lst import LSTEventSource
 
-    source = LSTEventSource(input_url=test_r0_path_all_streams)
+    source = LSTEventSource(
+        input_url=test_r0_path_all_streams,
+        apply_drs4_corrections=False,
+    )
     assert len(set(source.file_list)) == 4
 
     for i, event in enumerate(source):
@@ -90,7 +94,10 @@ def test_subarray():
 
 def test_missing_modules():
     from ctapipe_io_lst import LSTEventSource
-    source = LSTEventSource(test_missing_module_path)
+    source = LSTEventSource(
+        test_missing_module_path,
+        apply_drs4_corrections=False,
+    )
 
     assert source.lst_service.telescope_id == 1
     assert source.lst_service.num_modules == 264