LUT retrieval scheme changed; updates in accordance with the comments…

… / suggestions from @Liang Yu

src/s1reader/s1_annotation.py

@@ -11,6 +11,7 @@
 
 import numpy as np
 
+from scipy.interpolate import InterpolatedUnivariateSpline, interp1d
 from packaging import version
 
 #IPF version from which the NADS has azimuth noise vector annotation
@@ -470,6 +471,49 @@ def from_noise_annotation(cls, noise_annotation: NoiseAnnotation, azimuth_time:
                    line_from, line_to)
 
 
+    def compute_thermal_noise_lut(self, shape_lut):
+        '''Calculate thermal noise LUT whose shape is `shape_lut`
+
+        Parameter:
+        ----------
+        shape_lut: tuple or list
+            Shape of the output LUT
+
+        Returns
+        -------
+        arr_lut_total: np.ndarray
+            2d array containing thermal noise correction look up table values
+        '''
+
+        nrows, ncols = shape_lut
+
+        # Interpolate the range noise vector
+        rg_lut_interp_obj = InterpolatedUnivariateSpline(self.range_pixel,
+                                                         self.range_lut,
+                                                         k=1)
+        if self.azimuth_last_range_sample is not None:
+            vec_rg = np.arange(self.azimuth_last_range_sample + 1)
+        else:
+            vec_rg = np.arange(ncols)
+        rg_lut_interpolated = rg_lut_interp_obj(vec_rg)
+
+        # Interpolate the azimuth noise vector
+        if (self.azimuth_line is None) or (self.azimuth_lut is None):
+            az_lut_interpolated = np.ones(nrows)
+        else:  # IPF >= 2.90
+            az_lut_interp_obj = InterpolatedUnivariateSpline(self.azimuth_line,
+                                                             self.azimuth_lut,
+                                                             k=1)
+            vec_az = np.arange(self.line_from, self.line_to + 1)
+            az_lut_interpolated = az_lut_interp_obj(vec_az)
+
+        arr_lut_total = np.matmul(az_lut_interpolated[..., np.newaxis],
+                                  rg_lut_interpolated[np.newaxis, ...])
+
+        return arr_lut_total
+
+
+
 @dataclass
 class BurstCalibration:
     '''Calibration information for Sentinel-1 IW SLC burst
@@ -531,7 +575,6 @@ class BurstEAP:
     Sentinel-1 IW SLC burst
     '''
     # from LADS
-    num_sample: int  # number of samples
     freq_sampling: float  # range sampling rate
     eta_start: datetime.datetime
     tau_0: float  # imageInformation/slantRangeTime
@@ -569,7 +612,6 @@ def from_product_annotation_and_aux_cal(cls, product_annotation: ProductAnnotati
         '''
         id_closest = closest_block_to_azimuth_time(product_annotation.antenna_pattern_azimuth_time,
                                                    azimuth_time)
-        num_sample = product_annotation.number_of_samples
         freq_sampling = product_annotation.range_sampling_rate
         eta_start = azimuth_time
         tau_0 = product_annotation.slant_range_time
@@ -580,11 +622,49 @@ def from_product_annotation_and_aux_cal(cls, product_annotation: ProductAnnotati
 
         ascending_node_time = product_annotation.ascending_node_time
 
-        return cls(num_sample, freq_sampling, eta_start, tau_0, tau_sub, theta_sub,
+        return cls(freq_sampling, eta_start, tau_0, tau_sub, theta_sub,
                    azimuth_time, ascending_node_time,
                    gain_eap, delta_theta)
 
 
+    def compute_eap_compensation_lut(self, num_sample):
+        '''Returns LUT for EAP compensation whose size is `num_sample`.
+        Based on ESA docuemnt :
+        "Impact of the Elevation Antenna Pattern Phase Compensation
+         on the Interferometric Phase Preservation"
+
+        Document URL:
+        https://sentinel.esa.int/documents/247904/1653440/Sentinel-1-IPF_EAP_Phase_correction
+
+        Parameter:
+        num_sample: int
+            Size of the output LUT
+
+        Return:
+        -------
+            gain_eap_interpolatd: EAP phase for the burst to be compensated
+
+        '''
+
+        n_elt = len(self.gain_eap)
+
+        theta_am = (np.arange(n_elt) - (n_elt - 1) / 2) * self.delta_theta
+
+        delta_anx = self.eta_start - self.ascending_node_time
+        theta_offnadir = self._anx2roll(delta_anx.seconds + delta_anx.microseconds * 1.0e-6)
+
+        theta_eap = theta_am + theta_offnadir
+
+        tau = self.tau_0 + np.arange(num_sample) / self.freq_sampling
+
+        theta = np.interp(tau, self.tau_sub, self.theta_sub)
+
+        interpolator_gain = interp1d(theta_eap, self.gain_eap)
+        gain_eap_interpolated = interpolator_gain(theta)
+
+        return gain_eap_interpolated
+
+
     def _anx2roll(self, delta_anx):
         '''
         Returns the Platform nominal roll as function of elapsed time from
@@ -599,7 +679,7 @@ def _anx2roll(self, delta_anx):
 
         Returns
         -------
-        _: float
+        nominal_roll: float
             Estimated nominal roll (degrees)
         '''
         # Estimate altitude based on time elapsed since ANX
@@ -637,7 +717,7 @@ def _anx2height(self, delta_anx):
 
         Returns:
         --------
-        _: float
+        h_t: float
             nominal height of the platform
 
         '''

src/s1reader/s1_burst_slc.py

@@ -8,7 +8,6 @@
 import isce3
 import numpy as np
 from osgeo import gdal
-from scipy.interpolate import InterpolatedUnivariateSpline, interp1d
 
 from s1reader import s1_annotation
 
@@ -625,77 +624,25 @@ def swath_name(self):
 
     @property
     def thermal_noise_lut(self):
-        '''Returns the burst-sized LUT for thermal noise correction
-
-        Returns
-        -------
-        arr_lut_total: np.array
-            2d array containing thermal noise correction look up table values
         '''
-
+        Returns the LUT for thermal noise correction for the burst
+        '''
         if self.burst_noise is None:
             raise ValueError('burst_noise is not defined for this burst.')
 
-        nrows, ncols = self.shape
-
-        # Interpolate the range noise vector
-        rg_lut_interp_obj = InterpolatedUnivariateSpline(self.burst_noise.range_pixel,
-                                                         self.burst_noise.range_lut,
-                                                         k=1)
-        if self.burst_noise.azimuth_last_range_sample is not None:
-            vec_rg = np.arange(self.burst_noise.azimuth_last_range_sample + 1)
-        else:
-            vec_rg = np.arange(ncols)
-        rg_lut_interpolated = rg_lut_interp_obj(vec_rg).reshape((1, ncols))
-
-
-        # Interpolate the azimuth noise vector
-        if (self.burst_noise.azimuth_line is None) or (self.burst_noise.azimuth_lut is None):
-            az_lut_interpolated = np.ones(nrows).reshape((nrows, 1))
-        else:  # IPF >= 2.90
-            az_lut_interp_obj = InterpolatedUnivariateSpline(self.burst_noise.azimuth_line,
-                                                             self.burst_noise.azimuth_lut,
-                                                             k=1)
-            vec_az = np.arange(self.burst_noise.line_from, self.burst_noise.line_to + 1)
-            az_lut_interpolated = az_lut_interp_obj(vec_az).reshape((nrows, 1))
-
-        arr_lut_total = np.matmul(az_lut_interpolated, rg_lut_interpolated)
-
-        return arr_lut_total
+        return self.burst_noise.compute_thermal_noise_lut(self.shape)
 
     @property
     def eap_compensation_lut(self):
         '''Returns LUT for EAP compensation.
-        Based on ESA docuemnt :
-        "Impact of the Elevation Antenna Pattern Phase Compensation
-         on the Interferometric Phase Preservation"
-
-        Document URL:
-        https://sentinel.esa.int/documents/247904/1653440/Sentinel-1-IPF_EAP_Phase_correction
-
+
         Returns:
         -------
-            gain_eap: EAP phase for the burst to be compensated
+            _: Interpolated EAP gain for the burst's lines
 
         '''
+        if self.burst_eap is None:
+            raise ValueError('burst_eap is not defined for this burst.'
+                            f' IPF version = {self.ipf_version}')
 
-        n_elt = len(self.burst_eap.gain_eap)
-
-        theta_am = (np.arange(n_elt) - (n_elt - 1) / 2) * self.burst_eap.delta_theta
-
-        delta_anx = self.burst_eap.eta_start-self.burst_eap.ascending_node_time
-        theta_offnadir = self.burst_eap._anx2roll(delta_anx.seconds + delta_anx.microseconds * 1.0e-6)
-
-        theta_eap = theta_am + theta_offnadir
-
-        tau = self.burst_eap.tau_0 + np.arange(self.burst_eap.num_sample) / self.burst_eap.freq_sampling
-
-        theta = np.interp(tau, self.burst_eap.tau_sub, self.burst_eap.theta_sub)
-
-        interpolator_gain = interp1d(theta_eap, self.burst_eap.gain_eap)
-        gain_eap_interpolated = interpolator_gain(theta)
-        phi_eap = np.angle(gain_eap_interpolated)
-        c_j = np.complex64(1.0j)
-        gain_eap = np.exp(c_j * phi_eap)
-
-        return gain_eap
+        return self.burst_eap.compute_eap_compensation_lut(self.width)

src/s1reader/s1_reader.py

@@ -301,10 +301,11 @@ def get_path_aux_cal(directory_aux_cal: str, str_platform: str, instrument_cfg_i
     Returns:
     --------
     path_aux_cal: AUX_CAL file that corresponds to the criteria provided
+
     '''
     list_aux_cal = glob.glob(f'{directory_aux_cal}/{str_platform}_AUX_CAL_V*.SAFE.zip')
     list_aux_cal.sort()
-    path_aux_cal = list_aux_cal[inst_config_id-1]
+    path_aux_cal = list_aux_cal[instrument_cfg_id-1]
 
     return path_aux_cal