Implement software trigger

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",
 ]

ctapipe/instrument/tests/test_trigger.py

@@ -1,7 +1,7 @@
 from ctapipe.containers import ArrayEventContainer
 
 
-def test_softare_trigger(subarray_prod5_paranal):
+def test_software_trigger(subarray_prod5_paranal):
     from ctapipe.instrument.trigger import SoftwareTrigger
 
     subarray = subarray_prod5_paranal

ctapipe/instrument/trigger.py

@@ -4,6 +4,18 @@
 
 
 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.
+
+    LSTs have a hardware trigger, that requires at least two LSTs to have a trigger
+    in each event. Thus, after selecting a subarray, it is needed to
+    - Remove single LST telescope events from events
+    - Ignore events with only 1 telescope after this has been applied
+    """
 
     min_telescopes = Integer(
         default_value=1,
@@ -13,16 +25,71 @@ class SoftwareTrigger(TelescopeComponent):
         ),
     ).tag(config=True)
 
-    min_telescopes_by_type = IntTelescopeParameter(
+    min_telescopes_of_type = IntTelescopeParameter(
         default_value=1,
         help=(
             "Minimum number of telescopes required for a specific type."
-            " In events with fewer telescopes of a given 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 == 1:
+                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

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)