Merge pull request #19 from FrancaCassol/new_time_and_geometry

New time and geometry

ctapipe_io_lst/__init__.py

@@ -31,7 +31,7 @@
 def load_camera_geometry(version=3):
     ''' Load camera geometry from bundled resources of this repo '''
     f = resource_filename(
-        'ctapipe_io_lst', 'resources/LSTCam-{:03d}'.format(version)
+        'ctapipe_io_lst', 'resources/LSTCam-{:03d}.camgeom.fits.gz'.format(version)
     )
     return CameraGeometry.from_table(f)
 
@@ -101,6 +101,8 @@ def __init__(self, **kwargs):
 
         self.multi_file = MultiFiles(self.file_list)
 
+        self.geometry_version=3
+
         self.camera_config = self.multi_file.camera_config
         self.log.info(
             "Read {} input files".format(
@@ -125,7 +127,7 @@ def subarray(self):
         optics = OpticsDescription.from_name("LST")
 
         # camera info from LSTCam-[geometry_version].camgeom.fits.gz file
-        camera = load_camera_geometry(version=3)
+        camera = load_camera_geometry(version=self.geometry_version)
 
         tel_descr = TelescopeDescription(
             name='LST', tel_type='LST', optics=optics, camera=camera
@@ -156,14 +158,11 @@ def _generator(self):
         # Instrument information
         for tel_id in self.data.lst.tels_with_data:
 
-            assert (tel_id == 0 or tel_id == 1) # only LST1 (for the moment id = 0)
-
             # optics info from standard optics.fits.gz file
             optics = OpticsDescription.from_name("LST")
 
             # camera info from LSTCam-[geometry_version].camgeom.fits.gz file
-            geometry_version = 2
-            camera = CameraGeometry.from_name("LSTCam", geometry_version)
+            camera = load_camera_geometry(version=self.geometry_version)
 
             tel_descr = TelescopeDescription(
                 name='LST', tel_type='LST', optics=optics, camera=camera
@@ -175,6 +174,7 @@ def _generator(self):
             # LSTs telescope position taken from MC from the moment
             tel_pos = {tel_id: [50., 50., 16] * u.m}
 
+
         subarray = SubarrayDescription("LST1 subarray")
         subarray.tels = tels
         subarray.positions = tel_pos
@@ -248,7 +248,6 @@ def fill_lst_service_container_from_zfile(self):
         svc_container.num_samples = self.camera_config.num_samples
         svc_container.pixel_ids = self.camera_config.expected_pixels_id
         svc_container.data_model_version = self.camera_config.data_model_version
-
         svc_container.num_modules = self.camera_config.lstcam.num_modules
         svc_container.module_ids = self.camera_config.lstcam.expected_modules_id
         svc_container.idaq_version = self.camera_config.lstcam.idaq_version
@@ -283,37 +282,38 @@ def fill_lst_event_container_from_zfile(self, event):
             event_container.tib_tenMHz_counter = unpacked_tib[2]
             event_container.tib_stereo_pattern = unpacked_tib[3]
             event_container.tib_masked_trigger = unpacked_tib[4]
-        else:
-            # missing TIB data
-            event_container.tib_event_counter = -1
-            event_container.tib_pps_counter = -1
-            event_container.tib_tenMHz_counter = -1
-            event_container.tib_stereo_pattern = -1
-            event_container.tib_masked_trigger = -1
 
         # if UCTS data are there
         if event_container.extdevices_presence & 2:
 
-            # unpack UCTS-CDTS data
-            rec_fmt = '=IIIQQBBB'
-            unpacked_cdts =  struct.unpack(rec_fmt, event.lstcam.cdts_data)
-            event_container.ucts_event_counter = unpacked_cdts[0]
-            event_container.ucts_pps_counter = unpacked_cdts[1]
-            event_container.ucts_clock_counter = unpacked_cdts[2]
-            event_container.ucts_timestamp = unpacked_cdts[3]
-            event_container.ucts_camera_timestamp = unpacked_cdts[4]
-            event_container.ucts_trigger_type = unpacked_cdts[5]
-            event_container.ucts_white_rabbit_status = unpacked_cdts[6]
-
-        else:
-            event_container.ucts_event_counter = -1
-            event_container.ucts_pps_counter = -1
-            event_container.ucts_clock_counter = -1
-            event_container.ucts_timestamp = -1
-            event_container.ucts_camera_timestamp = -1
-            event_container.ucts_trigger_type = -1
-            event_container.ucts_white_rabbit_status = -1
+            if int(self.data.lst.tel[self.camera_config.telescope_id].svc.idaq_version) > 37201:
+
+                # unpack UCTS-CDTS data (new version)
+                rec_fmt = '=QIIIIIBBBBI'
+                unpacked_cdts = struct.unpack(rec_fmt, event.lstcam.cdts_data)
+                event_container.ucts_timestamp = unpacked_cdts[0]
+                event_container.ucts_address = unpacked_cdts[1]        # new
+                event_container.ucts_event_counter = unpacked_cdts[2]
+                event_container.ucts_busy_counter = unpacked_cdts[3]   # new
+                event_container.ucts_pps_counter = unpacked_cdts[4]
+                event_container.ucts_clock_counter = unpacked_cdts[5]
+                event_container.ucts_trigger_type = unpacked_cdts[6]
+                event_container.ucts_white_rabbit_status = unpacked_cdts[7]
+                event_container.ucts_stereo_pattern = unpacked_cdts[8] # new
+                event_container.ucts_num_in_bunch = unpacked_cdts[9]   # new
+                event_container.ucts_cdts_version = unpacked_cdts[10]  # new
 
+            else:
+                # unpack UCTS-CDTS data (old version)
+                rec_fmt = '=IIIQQBBB'
+                unpacked_cdts =  struct.unpack(rec_fmt, event.lstcam.cdts_data)
+                event_container.ucts_event_counter = unpacked_cdts[0]
+                event_container.ucts_pps_counter = unpacked_cdts[1]
+                event_container.ucts_clock_counter = unpacked_cdts[2]
+                event_container.ucts_timestamp = unpacked_cdts[3]
+                event_container.ucts_camera_timestamp = unpacked_cdts[4]
+                event_container.ucts_trigger_type = unpacked_cdts[5]
+                event_container.ucts_white_rabbit_status = unpacked_cdts[6]
 
         # if SWAT data are there
         if event_container.extdevices_presence & 4:
@@ -329,16 +329,6 @@ def fill_lst_event_container_from_zfile(self, event):
             event_container.swat_array_flag = unpacked_swat[6]
             event_container.swat_array_event_num = unpacked_swat[7]
 
-        else:
-            event_container.swat_timestamp =-1
-            event_container.swat_counter1 = -1
-            event_container.swat_counter2 = -1
-            event_container.swat_event_type = -1
-            event_container.swat_camera_flag = -1
-            event_container.swat_camera_event_num = -1
-            event_container.swat_array_flag = -1
-            event_container.swat_array_event_num = -1
-
         # unpack Dragon counters
         rec_fmt = '=HIIIQ'
         rec_len = struct.calcsize(rec_fmt)
@@ -374,13 +364,20 @@ def fill_r0_camera_container_from_zfile(self, r0_container, event):
         #    r0_container.trigger_time = self.data.lst.tel[self.camera_config.telescope_id].evt.ucts_timestamp/1e9
 
         # consider for the moment only TIB time since UCTS seems not correct
-        if self.data.lst.tel[self.camera_config.telescope_id].evt.tib_pps_counter > 0:
-            r0_container.trigger_time = (
-                self.data.lst.tel[self.camera_config.telescope_id].svc.date +
-                self.data.lst.tel[self.camera_config.telescope_id].evt.tib_pps_counter +
-                self.data.lst.tel[self.camera_config.telescope_id].evt.tib_tenMHz_counter * 10**(-7))
-        else:
-            r0_container.trigger_time = 0
+        #if self.data.lst.tel[self.camera_config.telescope_id].evt.tib_pps_counter > 0:
+        #    r0_container.trigger_time = (
+        #        self.data.lst.tel[self.camera_config.telescope_id].svc.date +
+        #        self.data.lst.tel[self.camera_config.telescope_id].evt.tib_pps_counter +
+        #        self.data.lst.tel[self.camera_config.telescope_id].evt.tib_tenMHz_counter * 10**(-7))
+        #else:
+        #    r0_container.trigger_time = 0
+
+        #consider for the moment trigger time from central dragon module
+        module_rank = np.where(self.data.lst.tel[self.camera_config.telescope_id].svc.module_ids == 132)
+        r0_container.trigger_time = (
+                    self.data.lst.tel[self.camera_config.telescope_id].svc.date +
+                    self.data.lst.tel[self.camera_config.telescope_id].evt.pps_counter[module_rank] +
+                    self.data.lst.tel[self.camera_config.telescope_id].evt.tenMHz_counter[module_rank] * 10**(-7))
 
         # look for correct trigger type first in UCTS and then in TIB
         #if self.data.lst.tel[self.camera_config.telescope_id].evt.ucts_trigger_type > 0:

ctapipe_io_lst/containers.py

@@ -19,6 +19,7 @@
 ]
 
 
+
 class LSTDriveContainer(Container):
     """
     Drive report container
@@ -74,28 +75,33 @@ class LSTEventContainer(Container):
     module_status = Field([], "status of the modules (n_modules)")
     extdevices_presence = Field(None, "presence of data for external devices")
 
-    tib_event_counter = Field(None, "TIB event counter")
-    tib_pps_counter = Field(None, "TIB pps counter")
-    tib_tenMHz_counter = Field(None, "TIB 10 MHz counter")
-    tib_stereo_pattern = Field(None, "TIB stereo pattern")
-    tib_masked_trigger = Field(None, "TIB trigger mask")
-
-    ucts_event_counter =  Field(None, "UCTS event counter")
-    ucts_pps_counter = Field(None, "UCTS pps counter")
-    ucts_clock_counter = Field(None, "UCTS clock counter")
-    ucts_timestamp = Field(None, "UCTS timestamp")
-    ucts_camera_timestamp = Field(None, "UCTS camera timestamp")
-    ucts_trigger_type = Field(None, "UCTS trigger type")
-    ucts_white_rabbit_status = Field(None, "UCTS whiteRabbit status")
-
-    swat_timestamp = Field(None, "SWAT timestamp")
-    swat_counter1 = Field(None, "SWAT event counter 1")
-    swat_counter2 = Field(None, "SWAT event counter 2")
-    swat_event_type = Field(None, "SWAT event type")
-    swat_camera_flag = Field(None, "SWAT camera flag ")
-    swat_camera_event_num = Field(None, "SWAT camera event number")
-    swat_array_flag = Field(None, "SWAT array negative flag")
-    swat_array_event_num = Field(None, "SWAT array event number")
+    tib_event_counter = Field(-1, "TIB event counter")
+    tib_pps_counter = Field(-1, "TIB pps counter")
+    tib_tenMHz_counter = Field(-1, "TIB 10 MHz counter")
+    tib_stereo_pattern = Field(-1, "TIB stereo pattern")
+    tib_masked_trigger = Field(-1, "TIB trigger mask")
+
+    ucts_event_counter =  Field(-1, "UCTS event counter")
+    ucts_pps_counter = Field(-1, "UCTS pps counter")
+    ucts_clock_counter = Field(-1, "UCTS clock counter")
+    ucts_timestamp = Field(-1, "UCTS timestamp")
+    ucts_camera_timestamp = Field(-1, "UCTS camera timestamp")
+    ucts_trigger_type = Field(-1, "UCTS trigger type")
+    ucts_white_rabbit_status = Field(-1, "UCTS whiteRabbit status")
+    ucts_address = Field(-1,"UCTS address")
+    ucts_busy_counter = Field(-1, "UCTS busy counter")
+    ucts_stereo_pattern = Field(-1, "UCTS stereo pattern")
+    ucts_num_in_bunch = Field(-1, "UCTS num in bunch (for debugging)")
+    ucts_cdts_version = Field(-1, "UCTS CDTS version")
+
+    swat_timestamp = Field(-1, "SWAT timestamp")
+    swat_counter1 = Field(-1, "SWAT event counter 1")
+    swat_counter2 = Field(-1, "SWAT event counter 2")
+    swat_event_type = Field(-1, "SWAT event type")
+    swat_camera_flag = Field(-1, "SWAT camera flag ")
+    swat_camera_event_num = Field(-1, "SWAT camera event number")
+    swat_array_flag = Field(-1, "SWAT array negative flag")
+    swat_array_event_num = Field(-1, "SWAT array event number")
 
     pps_counter= Field([], "Dragon pulse per second counter (n_modules)")
     tenMHz_counter = Field([], "Dragon 10 MHz counter (n_modules)")

setup.py

@@ -17,7 +17,7 @@
     install_requires=[
         'astropy',
         'ctapipe',
-        'protozfits',
+        'protozfits'
     ],
     package_data={
         'ctapipe_io_lst': ['resources/*'],