From 1b25c0517ade7369fc189256f46146f73d68d7ae Mon Sep 17 00:00:00 2001 From: Seongsu Jeong Date: Thu, 21 Sep 2023 10:00:28 -0700 Subject: [PATCH] Patch for "no polynomial" case for azimuth FM rate mitigation & handling the burst without burst polygon (#134) * Take care of the case when the az. time of the polynomial list does not cover sensing start / stop * Let `InterpolatedUnivariateSpline` to return boundary values when the input parameter is out of the interpolator's range * extrapolate the burst center point the border polygon when the burst polygon data is missing in the metadata * Interface to turn on / off the correction * ditch out burst border polygon extrapolation, and fill that with `None` * Place empty geometry instead of `None` * whitespace removal * Apply suggestions from code review Co-authored-by: Liang Yu * add description about "no polynomials for burst" case * get rid of unnecessary module in `s1-reader.py` * bump the version * more clarification about the missing polynomial case --------- Co-authored-by: Seongsu Jeong Co-authored-by: Liang Yu --- src/s1reader/s1_annotation.py | 23 ++++++++++++++++++++--- src/s1reader/s1_burst_slc.py | 14 ++++++++++++++ src/s1reader/s1_reader.py | 33 ++++++++++++++++++++++++--------- src/s1reader/version.py | 1 + 4 files changed, 59 insertions(+), 12 deletions(-) diff --git a/src/s1reader/s1_annotation.py b/src/s1reader/s1_annotation.py index 220e790..5da874d 100644 --- a/src/s1reader/s1_annotation.py +++ b/src/s1reader/s1_annotation.py @@ -1180,13 +1180,15 @@ def from_polynomial_lists(cls, fm_rate_coeff_burst_arr, fm_rate_tau0_burst_vec) = cls.extract_polynomial_sequence(az_fm_rate_list, sensing_start, - sensing_end) + sensing_end, + handle_out_of_range=True) (dc_aztime_burst_vec, dc_coeff_burst_arr, dc_tau0_burst_vec) = cls.extract_polynomial_sequence(doppler_centroid_list, sensing_start, - sensing_end) + sensing_end, + handle_out_of_range=True) return cls(fm_rate_aztime_burst_vec, fm_rate_coeff_burst_arr, fm_rate_tau0_burst_vec, dc_aztime_burst_vec, dc_coeff_burst_arr, dc_tau0_burst_vec) @@ -1195,7 +1197,8 @@ def from_polynomial_lists(cls, @classmethod def extract_polynomial_sequence(cls, polynomial_list: list, datetime_start: datetime.datetime, - datetime_end: datetime.datetime): + datetime_end: datetime.datetime, + handle_out_of_range=True): ''' Scan `vec_azimuth_time` end find indices of the vector that covers the period defined with @@ -1246,6 +1249,20 @@ def extract_polynomial_sequence(cls, polynomial_list: list, # Scale factor to convert range (in meters) to seconds (tau) range_to_tau = 2.0 / speed_of_light + + # Take care of the case when the az. time of the polynomial list does not cover + # the sensing start/stop + if (index_end == index_start) and handle_out_of_range: + # 0--1--2--3--4--5 <- az. time of polynomial list (index shown on the left) + #|--| <- sensing start / stop + if index_start == 0: + index_end += 1 + + # 0--1--2--3--4--5 <- az. time of polynomial list (index shown on the left) + # |--| <- sensing start / stop + else: + index_start -= 1 + for poly in polynomial_list[index_start:index_end+1]: vec_aztime_sequence.append(poly[0]) arr_coeff_sequence.append(poly[1].coeffs) diff --git a/src/s1reader/s1_burst_slc.py b/src/s1reader/s1_burst_slc.py index 5b8db5b..8736244 100644 --- a/src/s1reader/s1_burst_slc.py +++ b/src/s1reader/s1_burst_slc.py @@ -911,15 +911,29 @@ def az_fm_rate_mismatch_from_llh(self, for datetime_vec in self.extended_coeffs.dc_aztime_vec] # calculate splined interpolation of the coeffs. and tau_0s + if len(fm_rate_aztime_sec_vec) <= 1: + # Interpolator object cannot be created with only one set of polynomials. + # Such case happens when there is no polygon that falls between a burst's start / stop + # which was found from very few Sentinel-1 IW SLCs processed by IPF 002.36 + # Return an empty LUT2d in that case + return isce3.core.LUT2d() interpolator_tau0_ka = InterpolatedUnivariateSpline( fm_rate_aztime_sec_vec, self.extended_coeffs.fm_rate_tau0_vec, + ext=3, k=1) tau0_ka_interp = interpolator_tau0_ka(vec_t)[..., np.newaxis] + if len(dc_aztime_sec_vec) <=1: + # Interpolator object cannot be created with only one set of polynomials. + # Such case happens when there is no polygon that falls between a burst's start / stop + # which was found from very few Sentinel-1 IW SLCs processed by IPF 002.36 + # Return an empty LUT2d in that case + return isce3.core.LUT2d() interpolator_tau0_fdc_interp = InterpolatedUnivariateSpline( dc_aztime_sec_vec, self.extended_coeffs.dc_tau0_vec, + ext=3, k=1) tau0_fdc_interp = interpolator_tau0_fdc_interp(vec_t)[..., np.newaxis] diff --git a/src/s1reader/s1_reader.py b/src/s1reader/s1_reader.py index 46bf361..7d063ee 100644 --- a/src/s1reader/s1_reader.py +++ b/src/s1reader/s1_reader.py @@ -19,7 +19,7 @@ from nisar.workflows.stage_dem import check_dateline from s1reader import s1_annotation # to access __file__ -from s1reader.s1_annotation import (RFI_INFO_AVAILABLE_FROM, +from s1reader.s1_annotation import (RFI_INFO_AVAILABLE_FROM, CalibrationAnnotation, AuxCal, BurstCalibration, BurstEAP, BurstNoise, BurstExtendedCoeffs, @@ -221,16 +221,20 @@ def calculate_centroid(lons, lats): return shapely.geometry.Point(llh_centroid[:2]) -def get_burst_centers_and_boundaries(tree): +def get_burst_centers_and_boundaries(tree, num_bursts=None): '''Parse grid points list and calculate burst center lat and lon - Parameters: - ----------- + Parameters + ---------- tree : Element Element containing geolocation grid points. - - Returns: - -------- + num_bursts: int or None + Expected number of bursts in the subswath. + To check if the # of polygon is the same as the parsed polygons. + When it's None, the number of burst is the same as # polygons found in this function. + + Returns + ------- center_pts : list List of burst centroids ass shapely Points boundary_pts : list @@ -252,6 +256,9 @@ def get_burst_centers_and_boundaries(tree): lons[i] = float(grid_pt[5].text) unique_line_indices = np.unique(lines) + + if num_bursts is None: + num_bursts = len(unique_line_indices) - 1 n_bursts = len(unique_line_indices) - 1 center_pts = [[]] * n_bursts boundary_pts = [[]] * n_bursts @@ -273,6 +280,14 @@ def get_burst_centers_and_boundaries(tree): poly = shapely.geometry.Polygon(zip(burst_lons, burst_lats)) boundary_pts[i] = check_dateline(poly) + num_border_polygon = len(unique_line_indices) - 1 + if num_bursts > num_border_polygon: + warnings.warn('Inconsistency between # bursts in subswath, and the # polygons. ') + num_missing_polygons = num_bursts - num_border_polygon + + center_pts += [shapely.Point()] * num_missing_polygons + boundary_pts += [[shapely.Polygon()]] * num_missing_polygons + return center_pts, boundary_pts @@ -853,8 +868,6 @@ def burst_from_xml(annotation_path: str, orbit_path: str, tiff_path: str, orbit_number = int(tree.find('adsHeader/absoluteOrbitNumber').text) - center_pts, boundary_pts = get_burst_centers_and_boundaries(tree) - wavelength = isce3.core.speed_of_light / radar_freq starting_range = slant_range_time * isce3.core.speed_of_light / 2 range_pxl_spacing = isce3.core.speed_of_light / (2 * range_sampling_rate) @@ -912,6 +925,8 @@ def burst_from_xml(annotation_path: str, orbit_path: str, tiff_path: str, n_bursts = int(burst_list_elements.attrib['count']) bursts = [[]] * n_bursts + center_pts, boundary_pts = get_burst_centers_and_boundaries(tree, num_bursts=n_bursts) + for i, burst_list_element in enumerate(burst_list_elements): # Zero Doppler azimuth time of the first line of this burst sensing_start = as_datetime(burst_list_element.find('azimuthTime').text) diff --git a/src/s1reader/version.py b/src/s1reader/version.py index 1c914ee..15132fa 100644 --- a/src/s1reader/version.py +++ b/src/s1reader/version.py @@ -4,6 +4,7 @@ # release history Tag = collections.namedtuple('Tag', 'version date') release_history = ( + Tag('0.2.3', '2023-09-21'), Tag('0.2.2', '2023-09-08'), Tag('0.2.1', '2023-08-23'), Tag('0.2.0', '2023-07-25'),