Fix finding ids of stereo events for MC data.

ctapipe_io_magic/__init__.py

@@ -1370,13 +1370,26 @@ class MarsRun:
         else:
             stereo_m1_data = self.event_data['M1']['stereo_event_number'][np.where(self.event_data['M1']['trigger_pattern'] == MC_TRIGGER_PATTERN)]
             stereo_m2_data = self.event_data['M2']['stereo_event_number'][np.where(self.event_data['M2']['trigger_pattern'] == MC_TRIGGER_PATTERN)]
+            # remove events with 0 stereo number, which are mono events
             stereo_m1_data = stereo_m1_data[np.where(stereo_m1_data != 0)]
             stereo_m2_data = stereo_m2_data[np.where(stereo_m2_data != 0)]
 
             stereo_numbers = np.intersect1d(stereo_m1_data, stereo_m2_data)
 
-            m1_ids = np.searchsorted(self.event_data['M1']['stereo_event_number'], stereo_numbers)
-            m2_ids = np.searchsorted(self.event_data['M2']['stereo_event_number'], stereo_numbers)
+            # because of IDs equal to 0, we must find indices in a slight different way
+            # see https://stackoverflow.com/questions/8251541/numpy-for-every-element-in-one-array-find-the-index-in-another-array
+
+            index_m1 = np.argsort(self.event_data['M1']['stereo_event_number'])
+            index_m2 = np.argsort(self.event_data['M2']['stereo_event_number'])
+
+            sort_stereo_events_m1 = self.event_data['M1']['stereo_event_number'][index_m1]
+            sort_stereo_events_m2 = self.event_data['M2']['stereo_event_number'][index_m2]
+
+            sort_index_m1 = np.searchsorted(sort_stereo_events_m1, stereo_numbers)
+            sort_index_m2 = np.searchsorted(sort_stereo_events_m2, stereo_numbers)
+
+            m1_ids = np.take(index_m1, sort_index_m1)
+            m2_ids = np.take(index_m2, sort_index_m2)
 
             stereo_ids = list(zip(m1_ids, m2_ids))