diff --git a/src/ctapipe_io_lst/__init__.py b/src/ctapipe_io_lst/__init__.py index 597cf1cb..6c33ab34 100644 --- a/src/ctapipe_io_lst/__init__.py +++ b/src/ctapipe_io_lst/__init__.py @@ -479,27 +479,30 @@ def fill_from_cta_r1(self, array_event, zfits_event): offset = self.data_stream.waveform_offset pixel_id_map = self.camera_config.pixel_id_map - # FIXME: missing modules / pixels - # FIXME: DVR? should not happen in r1 but dl0, but our own converter uses the old R1 - # reorder to nominal pixel order pixel_status = _reorder_pixel_status( - zfits_event.pixel_status, pixel_id_map, set_dvr_bits=True + zfits_event.pixel_status, pixel_id_map, set_dvr_bits=not self.dvr_applied, ) n_channels = zfits_event.num_channels - n_pixels = zfits_event.num_pixels n_samples = zfits_event.num_samples + if self.dvr_applied: + stored_pixels = (zfits_event.pixel_status & PixelStatus.DVR_STATUS) > 0 + n_pixels = np.count_nonzero(stored_pixels) + else: + stored_pixels = slice(None) # all pixels stored + n_pixels = zfits_event.num_pixels + readout_shape = (n_channels, n_pixels, n_samples) raw_waveform = zfits_event.waveform.reshape(readout_shape) waveform = raw_waveform.astype(np.float32) / scale - offset reordered_waveform = np.full((n_channels, N_PIXELS, n_samples), 0.0, dtype=np.float32) - reordered_waveform[:, pixel_id_map] = waveform + reordered_waveform[:, pixel_id_map[stored_pixels]] = waveform waveform = reordered_waveform - # FIXME, check using evb_preprocessing and make ctapipe support 2 gains + if zfits_event.num_channels == 2: selected_gain_channel = None else: @@ -529,22 +532,23 @@ def fill_from_cta_r1(self, array_event, zfits_event): if DataLevel.R0 in self.datalevels: reordered_raw_waveform = np.full((n_channels, N_PIXELS, n_samples), 0, dtype=np.uint16) - reordered_raw_waveform[:, pixel_id_map] = raw_waveform + reordered_raw_waveform[:, pixel_id_map[stored_pixels]] = raw_waveform array_event.r0.tel[self.tel_id] = R0CameraContainer( waveform=reordered_raw_waveform, ) def fill_lst_from_ctar1(self, zfits_event): + pixel_status = _reorder_pixel_status( + zfits_event.pixel_status, + self.pixel_id_map, + set_dvr_bits=not self.dvr_applied, + ) evt = LSTEventContainer( - pixel_status=zfits_event.pixel_status.copy(), + pixel_status=pixel_status, first_capacitor_id=zfits_event.first_cell_id, calibration_monitoring_id=zfits_event.calibration_monitoring_id, local_clock_counter=zfits_event.module_hires_local_clock_counter, ) - # set bits for dvr if not already set - if not self.dvr_applied: - not_broken = get_channel_info(evt.pixel_status) != 0 - evt.pixel_status[not_broken] |= PixelStatus.DVR_STATUS_0 if zfits_event.debug is not None: debug = zfits_event.debug @@ -971,7 +975,7 @@ def tag_flatfield_events(self, array_event): ''' tel_id = self.tel_id waveform = array_event.r1.tel[tel_id].waveform - + if waveform.ndim == 3: image = waveform[HIGH_GAIN].sum(axis=1) else: @@ -981,7 +985,7 @@ def tag_flatfield_events(self, array_event): n_in_range = np.count_nonzero(in_range) looks_like_ff = n_in_range >= self.min_flatfield_pixel_fraction * image.size - + if looks_like_ff: # Tag as FF only events with 2-gains waveforms: both gains are needed for calibration if waveform.ndim == 3: @@ -1091,7 +1095,7 @@ def fill_r0r1_camera_container(self, zfits_event): if CTAPIPE_0_20: # reorder to nominal pixel order pixel_status = _reorder_pixel_status( - zfits_event.pixel_status, pixel_id_map, set_dvr_bits=True + zfits_event.pixel_status, pixel_id_map, set_dvr_bits=not self.dvr_applied ) r1.pixel_status = pixel_status r1.event_time = cta_high_res_to_time( diff --git a/src/ctapipe_io_lst/tests/test_lsteventsource.py b/src/ctapipe_io_lst/tests/test_lsteventsource.py index 424077f4..90bea5f1 100644 --- a/src/ctapipe_io_lst/tests/test_lsteventsource.py +++ b/src/ctapipe_io_lst/tests/test_lsteventsource.py @@ -21,6 +21,7 @@ test_r0_path_all_streams = test_r0_dir / 'LST-1.1.Run02008.0000_first50.fits.fz' test_missing_module_path = test_data / 'real/R0/20210215/LST-1.1.Run03669.0000_first50.fits.fz' +test_missing_module_r1v1_path = test_data / 'real/R0/20240310/LST-1.1.Run17016.0000_first10.fits.fz' test_drive_report = test_data / 'real/monitoring/DrivePositioning/DrivePosition_log_20200218.txt' @@ -142,6 +143,35 @@ def test_missing_modules(): assert np.all(event.r0.tel[1].waveform[:, 514] == fill) +def test_missing_modules_r1v1(): + from ctapipe_io_lst import LSTEventSource + source = LSTEventSource( + test_missing_module_r1v1_path, + apply_drs4_corrections=False, + pointing_information=False, + ) + + assert source.lst_service.telescope_id == 1 + assert source.lst_service.num_modules == 264 + + n_events = 0 + for event in source: + n_events += 1 + # one module missing, so 7 pixels + assert np.count_nonzero(event.mon.tel[1].pixel_status.hardware_failing_pixels) == N_PIXELS_MODULE * N_GAINS + assert np.count_nonzero(event.r0.tel[1].waveform == 0.0) == N_PIXELS_MODULE * N_SAMPLES * N_GAINS + + missing_gain, missing_pixel = np.nonzero(event.mon.tel[1].pixel_status.hardware_failing_pixels) + # 514 is one of the missing pixels + for gain, pixel in zip(missing_gain, missing_pixel): + np.testing.assert_equal(event.r0.tel[1].waveform[gain, pixel], 0.0) + + if CTAPIPE_0_20: + np.testing.assert_equal(event.lst.tel[1].evt.pixel_status, event.r1.tel[1].pixel_status) + + assert n_events == 40 + + def test_gain_selected(): from ctapipe_io_lst import LSTEventSource