GainSelection should be applied to R0 waveform, not R1

[kosack]

Right now the CameraCalibrator applies the GainSelector (in calibrated_dl0()) and it applies it to the R1 waveforms (which are in "rough photoelectrons" after the R1 calibration). This means that the threshold is calibration dependent (e.g. it's something like 7.0 photoelectrons/ns for an LST). In reality, this should be applied in the step from R0 → R1, where the threshold is in Digital Counts and thus the threshold for switching should be (2^12-1) = 4095 for a 12-bit ADC and assuming the front-end amplifiers are linear in the whole digitization range.

I suggest we either:

• remove the gain selector from CameraCalibrator and move it to SimTelEventSource (not ideal, since then you'd have to configure it, but that is currently where the R1 calibration happens)
• just change CameraCalibrator to apply it to the r0 waveform (but remember in the future, this would not be available for real data)

That way we don't have to "calibrate" the gain threshold for each camera each time the DC/PE or pedestals change.

[kosack]

Also, I suggest we change the threshold config attribute to have a unit in the name (e.g threshold_dc, whereas in the current implementation, it should be threshold_pe). This explains why in #1163 , I never got any low-gain pixels (I've now updated the config file example to set the threshold in PE, but so far it's not configurable per camera)

[HealthyPear]

This is for NectarCam, using 50 events from R0 waveforms.

This is for NectarCam, using 50 events from R1 waveforms.

[watsonjj]

Should the gain selection be applied after the mc pedestal subtraction (but before the mc.dc_to_pe is applied)?

ctapipe/ctapipe/io/simteleventsource.py

Line 282 in 5b81bb6

ped = mc.pedestal[..., np.newaxis] / n_samples

[vuillaut]

Right now the CameraCalibrator applies the GainSelector (in calibrated_dl0()) and it applies it to the R1 waveforms (which are in "rough photoelectrons" after the R1 calibration). This means that the threshold is calibration dependent (e.g. it's something like 7.0 photoelectrons/ns for an LST). In reality, this should be applied in the step from R0 → R1, where the threshold is in Digital Counts and thus the threshold for switching should be (2^12-1) = 4095 for a 12-bit ADC and assuming the front-end amplifiers are linear in the whole digitization range.

I suggest we either:

• remove the gain selector from CameraCalibrator and move it to SimTelEventSource (not ideal, since then you'd have to configure it, but that is currently where the R1 calibration happens)
• just change CameraCalibrator to apply it to the r0 waveform (but remember in the future, this would not be available for real data)

That way we don't have to "calibrate" the gain threshold for each camera each time the DC/PE or pedestals change.

I did not realise this.
Thanks for the info.

[watsonjj]

@HealthyPear Is this issue solved since #1167?

[HealthyPear]

@HealthyPear Is this issue solved since #1167?

Since protopipe relies only on frozen releases of ctapipe, for the moment I hardcoded what I needed.
The changes you are referring to are currently hosted in the development (master) version of ctapipe, so I cannot use them now anyway.

What I did is the following,

[...]
from ctapipe.calib.camera.gainselection import GainSelector
[...]
class MyCameraCalibrator(CameraCalibrator):
    """Create a child class of CameraCalibrator."""

    def _calibrate_dl0(self, event, telid):
        """Override class method to perform gain selection at R0 instead of R1.

        Then select R1 waveforms using R0-selected gain channels.
        """
        waveforms_r0 = event.r0.tel[telid].waveform
        _, selected_gain_channel = self.gain_selector(waveforms_r0)

        waveforms_r1 = event.r1.tel[telid].waveform
        if self._check_r1_empty(waveforms_r1):
            return

        # Use R0-selected gain channels to select R1 waveforms
        _, n_pixels, _ = waveforms_r1.shape
        waveforms_gs = waveforms_r1[selected_gain_channel, np.arange(n_pixels)]
        if selected_gain_channel is not None:
            event.r1.tel[telid].selected_gain_channel = selected_gain_channel
        else:
            if event.r1.tel[telid].selected_gain_channel is None:
                raise ValueError(
                    "EventSource is loading pre-gainselected waveforms "
                    "without filling the selected_gain_channel container"
                )

        reduced_waveforms = self.data_volume_reducer(waveforms_gs)
        event.dl0.tel[telid].waveform = reduced_waveforms
[...]
        cfg = Config()
        cfg["ChargeExtractorFactory"]["window_width"] = 5
        cfg["ChargeExtractorFactory"]["window_shift"] = 2
        cfg["ThresholdGainSelector"]["threshold"] = 4000.0
        extractor = LocalPeakWindowSum(config=cfg)
        gain_selector = GainSelector.from_name("ThresholdGainSelector", config=cfg)

        self.calib = MyCameraCalibrator(
            config=cfg, gain_selector=gain_selector, image_extractor=extractor
        )
[...]

What you decided to do is anyway better for the future, so as for any new feature in ctapipe, I will simply import it and get rid of my hardcoded part.

[kosack]

fixed in #1167