Update the reading of spectral response, i.e. central wavelength (cw) and fwhm for HCO coregistered products

Author: telespazio-tim

prisma4sen2like/prisma/adapter.py

@@ -95,7 +95,6 @@ def interpolate_gps_position(v):
     for rec in u[1:]:
         i = np.where(v == rec)[0][0]
         indice.append(i)
-    print(indice)
     x_0 = np.linspace(1, v.shape[0], v.shape[0])
     # linear interpolation :
     x = np.array(indice)
@@ -357,9 +356,9 @@ def _extract_mask(self, grid, values: list[int], output_file_path: str) -> NDArr
             dest_mask[mask == val] = 1
             logger.debug(dest_mask[dest_mask == val])
 
-        logger.info("NB maked px %s", len(dest_mask[dest_mask >= 1]))
+        logger.debug("NB maked px %s", len(dest_mask[dest_mask >= 1]))
 
-        logger.info("maked px %s", dest_mask[dest_mask >= 1])
+        logger.debug("maked px %s", dest_mask[dest_mask >= 1])
 
         image = Image.fromarray(dest_mask)
         image.save(output_file_path)
@@ -472,8 +471,6 @@ def _compute_orbital_model(self):
         pos_y = interpolate_gps_position(np.array(list(self._product.product["Info/Ancillary/PVSdata/Wgs84_pos_y"])))
         pos_z = interpolate_gps_position(np.array(list(self._product.product["Info/Ancillary/PVSdata/Wgs84_pos_z"])))
 
-        print(pos_x[0], pos_y[0], pos_z[0])
-
         lat_lon_grids = self._product.lat_lon_grids
 
         lat = np.pi / 180.0 * np.rot90(lat_lon_grids[0], k=-1)
@@ -511,7 +508,6 @@ def _compute_orbital_model(self):
         r = res / res2
         zenith = (np.arccos(res / res2)) * 180 / np.pi
         # zenith -> off nadir
-        print(f"zenith {np.mean(zenith)}")
 
         # compute azimuth:
         # consider directionnal vector (d)
@@ -528,14 +524,9 @@ def _compute_orbital_model(self):
 
         azimuth[azimuth < 0] = azimuth[azimuth < 0] + 360
 
-        print(f"azimuth {np.mean(azimuth)}")
-
         self._viewing_zenith_angle = np.mean(zenith)
         self._viewing_azimuth_angle = np.mean(azimuth)
 
-        print(f"Sun_azimuth_angle {self._product.product.attrs['Sun_azimuth_angle']}")
-        print(f"Sun_zenith_angle {self._product.product.attrs['Sun_zenith_angle']}")
-
     @property
     def viewing_angle_grid(self) -> AngleGrid:
         # _viewing_zenith_angle and _viewing_azimuth_angle

prisma4sen2like/prisma/spectral_aggregation_functions.py

@@ -36,14 +36,31 @@ def generate_aggregation_coefficients_prisma_s2(product_file):
     n_detectors = 1000
 
     # Compute spectral aggregation coefficients (PRISMA => S2-MSI-A)
-    index_start_vnir = 96
-    index_stop_vnir = 162
-    index_start_swir = 81
-    index_stop_swir = 254
-    prisma_vnir_wl = product_file["KDP_AUX/Cw_Vnir_Matrix"][:, index_start_vnir:index_stop_vnir]
-    prisma_vnir_fwhm = product_file["KDP_AUX/Fwhm_Vnir_Matrix"][:, index_start_vnir:index_stop_vnir]
-    prisma_swir_wl = product_file["KDP_AUX/Cw_Swir_Matrix"][:, index_start_swir:index_stop_swir]
-    prisma_swir_fwhm = product_file["KDP_AUX/Fwhm_Swir_Matrix"][:, index_start_swir:index_stop_swir]
+
+    # Read matrix of central wavelength and fwhm for each detector (case of not coregistered products)
+    # This part is commented as not used when working with coregistered products (HCO)
+    # This part could be uncommented in future if working with not coregistered products.
+    # index_start_vnir = 96
+    # index_stop_vnir = 162
+    # index_start_swir = 81
+    # index_stop_swir = 254
+    # prisma_vnir_wl = product_file["KDP_AUX/Cw_Vnir_Matrix"][:, index_start_vnir:index_stop_vnir]
+    # prisma_vnir_fwhm = product_file["KDP_AUX/Fwhm_Vnir_Matrix"][:, index_start_vnir:index_stop_vnir]
+    # prisma_swir_wl = product_file["KDP_AUX/Cw_Swir_Matrix"][:, index_start_swir:index_stop_swir]
+    # prisma_swir_fwhm = product_file["KDP_AUX/Fwhm_Swir_Matrix"][:, index_start_swir:index_stop_swir]
+
+    # For coregistered cubes (HCO), read central wavelength and fwhm (VNIR & SWIR) using
+    # the GLOBAL attributes: List_Cw_Vnir, List_Fwhm_Vnir, List_Cw_Swir, List_Fwhm_Swir
+    prisma_vnir_wl_list = np.float32(product_file.attrs.get("List_Cw_Vnir"))
+    prisma_vnir_fwhm_list = np.float32(product_file.attrs.get("List_Fwhm_Vnir"))
+    prisma_swir_wl_list = np.float32(product_file.attrs.get("List_Cw_Swir"))
+    prisma_swir_fwhm_list = np.float32(product_file.attrs.get("List_Fwhm_Swir"))
+
+    # Replicate prisma_vnir_wl values per column with single List_cw_vnir value
+    prisma_vnir_wl = np.tile(prisma_vnir_wl_list, (n_detectors, 1))
+    prisma_vnir_fwhm = np.tile(prisma_vnir_fwhm_list, (n_detectors, 1))
+    prisma_swir_wl = np.tile(prisma_swir_wl_list, (n_detectors, 1))
+    prisma_swir_fwhm = np.tile(prisma_swir_fwhm_list, (n_detectors, 1))
 
     n_bands_prisma_vnir = prisma_vnir_wl.shape[1]
     n_bands_prisma_swir = prisma_swir_wl.shape[1]

prisma4sen2like/prisma/version.py

@@ -18,4 +18,4 @@
 """version info module."""
 
 
-__version__ = "00.01.00"
+__version__ = "00.01.01"