From 8a492555fdbfa46c14289f23577cac032c13cf1e Mon Sep 17 00:00:00 2001 From: Maximilian Linhoff Date: Fri, 25 Nov 2022 12:43:29 +0100 Subject: [PATCH] Implement software stereo trigger This trigger allows it to remove lone LST telescope events from array events (as required by the hardware stereo trigger in real life), and to reject events with only one telescope after a subarray has been selected. --- ctapipe/instrument/__init__.py | 2 + ctapipe/instrument/tests/test_trigger.py | 52 +++++++++++ ctapipe/instrument/trigger.py | 113 +++++++++++++++++++++++ ctapipe/tools/process.py | 26 +++--- 4 files changed, 181 insertions(+), 12 deletions(-) create mode 100644 ctapipe/instrument/tests/test_trigger.py create mode 100644 ctapipe/instrument/trigger.py diff --git a/ctapipe/instrument/__init__.py b/ctapipe/instrument/__init__.py index db941ce9be0..16c96403d9b 100644 --- a/ctapipe/instrument/__init__.py +++ b/ctapipe/instrument/__init__.py @@ -4,6 +4,7 @@ from .optics import FocalLengthKind, OpticsDescription, ReflectorShape, SizeType from .subarray import SubarrayDescription, UnknownTelescopeID from .telescope import TelescopeDescription +from .trigger import SoftwareTrigger __all__ = [ "CameraDescription", @@ -19,4 +20,5 @@ "FocalLengthKind", "ReflectorShape", "SizeType", + "SoftwareTrigger", ] diff --git a/ctapipe/instrument/tests/test_trigger.py b/ctapipe/instrument/tests/test_trigger.py new file mode 100644 index 00000000000..3371cdfa587 --- /dev/null +++ b/ctapipe/instrument/tests/test_trigger.py @@ -0,0 +1,52 @@ +from ctapipe.containers import ArrayEventContainer + + +def test_software_trigger(subarray_prod5_paranal): + from ctapipe.instrument.trigger import SoftwareTrigger + + subarray = subarray_prod5_paranal + trigger = SoftwareTrigger( + subarray=subarray, + min_telescopes=2, + min_telescopes_of_type=[ + ("type", "*", 0), + ("type", "LST*", 2), + ], + ) + + # only one telescope, no SWAT + event = ArrayEventContainer() + event.trigger.tels_with_trigger = [5] + assert trigger(event) == False + assert event.trigger.tels_with_trigger == [] + + # 1 LST + 1 MST, 1 LST would not have triggered LST hardware trigger + # and after LST is removed, we only have 1 telescope, so no SWAT either + event = ArrayEventContainer() + event.trigger.tels_with_trigger = [1, 6] + assert trigger(event) == False + assert event.trigger.tels_with_trigger == [] + + # two MSTs and 1 LST, -> remove single LST + event = ArrayEventContainer() + event.trigger.tels_with_trigger = [1, 5, 6] + assert trigger(event) == True + assert event.trigger.tels_with_trigger == [5, 6] + + # two MSTs, nothing to change + event = ArrayEventContainer() + event.trigger.tels_with_trigger = [5, 6] + assert trigger(event) == True + assert event.trigger.tels_with_trigger == [5, 6] + + # three LSTs, nothing to change + event = ArrayEventContainer() + event.trigger.tels_with_trigger = [1, 2, 3] + assert trigger(event) == True + assert event.trigger.tels_with_trigger == [1, 2, 3] + + # thee LSTs, plus MSTs, nothing to change + event = ArrayEventContainer() + event.trigger.tels_with_trigger = [1, 2, 3, 5, 6, 7] + assert trigger(event) == True + assert event.trigger.tels_with_trigger == [1, 2, 3, 5, 6, 7] diff --git a/ctapipe/instrument/trigger.py b/ctapipe/instrument/trigger.py new file mode 100644 index 00000000000..31120160fdb --- /dev/null +++ b/ctapipe/instrument/trigger.py @@ -0,0 +1,113 @@ +from ctapipe.containers import ArrayEventContainer +from ctapipe.core import TelescopeComponent +from ctapipe.core.traits import Integer, IntTelescopeParameter + + +class SoftwareTrigger(TelescopeComponent): + """ + A stereo trigger that can remove telescope events from subarray events. + + This class is needed to correctly handle super-arrays as simulated for + CTA and still handle the LST hardware stereo trigger and the normal stereo + trigger correctly. + + When selecting subarrays from simulations that contain many more telescopes, + as is done in all major CTA productions to date, the stereo trigger is not + correctly simulated as in that after selecting a realistic subarray, events + are still in the data stream where only one telescope of the selected subarray + triggered, which would in reality not trigger the stereo trigger. + + An additional complexity is the LST hardware stereo trigger, that forces that + an array event has always to contain no or at least two LST telescope events. + + This means that after selectig a subarray, we need to: + - Remove LST telescope events from the subarray if only one LST triggered + - Ignore events with only 1 telescope after this has been applied + + With the default settings, this class is a no-op. To get the correct behavior + for CTA simulations, use the following configuration: + + .. + SoftwareTrigger: + min_telescopes: 2 + min_telescopes_of_type: + - ["type", "*", 0] + - ["type", "LST*", 2] + """ + + min_telescopes = Integer( + default_value=1, + help=( + "Minimum number of telescopes required globally." + " Events with fewer telescopes will be filtered out completely." + ), + ).tag(config=True) + + min_telescopes_of_type = IntTelescopeParameter( + default_value=0, + help=( + "Minimum number of telescopes required for a specific type." + " In events with fewer telescopes of that type" + " , those telescopes will be removed from the array event." + " This might result in the event not fullfilling ``min_telescopes`` anymore" + " and thus being filtered completely." + ), + ).tag(config=True) + + def __init__(self, subarray, *args, **kwargs): + super().__init__(subarray, *args, **kwargs) + + self._ids_by_type = { + str(type): set(self.subarray.get_tel_ids_for_type(type)) + for type in self.subarray.telescope_types + } + + def __call__(self, event: ArrayEventContainer) -> bool: + """ + Remove telescope events that have not the required number of telescopes of + a given type from the subarray event and decide if the event would + have triggered the stereo trigger. + + Data is cleared from events that did not trigger. + + Returns + ------- + triggered : bool + Whether or not this event would have triggered the stereo trigger + """ + + for tel_type in self.subarray.telescope_types: + tel_type_str = str(tel_type) + min_tels = self.min_telescopes_of_type.tel[tel_type_str] + + # no need to check telescopes for which we have no min requirement + if min_tels == 0: + continue + + tels_with_trigger = set(event.trigger.tels_with_trigger) + tel_ids = self._ids_by_type[tel_type_str] + tels_in_event = tels_with_trigger.intersection(tel_ids) + if len(tels_in_event) < min_tels: + for tel_id in tels_in_event: + self.log.debug( + "Removing tel_id %d of type %s from event due to type requirement", + tel_id, + tel_type_str, + ) + + # remove from tels_with_trigger + event.trigger.tels_with_trigger.remove(tel_id) + + # remove any related data + for container in ("trigger", "r0", "r1", "dl0", "dl1", "dl2"): + tel_map = getattr(event, container).tel + if tel_id in tel_map: + del tel_map[tel_id] + + if len(event.trigger.tels_with_trigger) < self.min_telescopes: + event.trigger.tels_with_trigger = [] + for container in ("trigger", "r0", "r1", "dl0", "dl1", "dl2"): + tel_map = getattr(event, container).tel + tel_map.clear() + return False + return True diff --git a/ctapipe/tools/process.py b/ctapipe/tools/process.py index acfdd776e54..f62ff5a89c6 100644 --- a/ctapipe/tools/process.py +++ b/ctapipe/tools/process.py @@ -11,6 +11,7 @@ from ..core.traits import Bool, classes_with_traits, flag from ..image import ImageCleaner, ImageModifier, ImageProcessor from ..image.extractor import ImageExtractor +from ..instrument import SoftwareTrigger from ..io import ( DataLevel, DataWriter, @@ -144,6 +145,7 @@ class ProcessorTool(Tool): ShowerProcessor, metadata.Instrument, metadata.Contact, + SoftwareTrigger, ] + classes_with_traits(EventSource) + classes_with_traits(ImageCleaner) @@ -169,17 +171,11 @@ def setup(self): ) sys.exit(1) - self.calibrate = CameraCalibrator( - parent=self, subarray=self.event_source.subarray - ) - self.process_images = ImageProcessor( - subarray=self.event_source.subarray, parent=self - ) - - self.process_shower = ShowerProcessor( - subarray=self.event_source.subarray, parent=self - ) - + subarray = self.event_source.subarray + self.software_trigger = SoftwareTrigger(parent=self, subarray=subarray) + self.calibrate = CameraCalibrator(parent=self, subarray=subarray) + self.process_images = ImageProcessor(subarray=subarray, parent=self) + self.process_shower = ShowerProcessor(subarray=subarray, parent=self) self.write = DataWriter(event_source=self.event_source, parent=self) # add ml reco classes if model paths were supplied via cli and not already configured @@ -196,7 +192,7 @@ def setup(self): reconstructor = Reconstructor.from_name( name, parent=self.process_shower, - subarray=self.event_source.subarray, + subarray=subarray, ) self.process_shower.reconstructors.append(reconstructor) self.process_shower.reconstructor_types.append(name) @@ -298,6 +294,12 @@ def start(self): if not self.event_type_filter(event): continue + if not self.software_trigger(event): + self.log.debug( + "Skipping event %i due to software trigger", event.index.event_id + ) + continue + if self.should_calibrate: self.calibrate(event)