Enable Sentinel-1 Extra Wide (EW) SLC compatibility

src/s1reader/constants.py

@@ -0,0 +1,3 @@
+SUPPORTED_SENSOR_MODES = ["iw", "ew"]
+SENSOR_MODE_SUBSWATHS = {"iw": [1, 2, 3], "ew": [1, 2, 3, 4, 5]}
+SENSOR_MODE_MID_SWATH = {"iw": 2, "ew": 3}

src/s1reader/s1_annotation.py

@@ -1,5 +1,5 @@
 """
-A module to load annotation files for Sentinel-1 IW SLC SAFE data
+A module to load annotation files for Sentinel-1 IW/EW SLC SAFE data
 To be used for the class "Sentinel1BurstSlc"
 """
 
@@ -273,7 +273,7 @@ def from_et(cls, et_in: ET, path_annotation: str):
 @dataclass
 class NoiseAnnotation(AnnotationBase):
     """
-    Reader for Noise Annotation Data Set (NADS) for IW SLC
+    Reader for Noise Annotation Data Set (NADS) for IW/EW SLC
     Based on ESA documentation: "Thermal Denoising of Products Generated by the S-1 IPF"
     """
 
@@ -367,7 +367,7 @@ def from_et(cls, et_in: ET, ipf_version: version.Version, path_annotation: str):
 @dataclass
 class ProductAnnotation(AnnotationBase):
     """
-    Reader for L1 Product annotation for IW SLC
+    Reader for L1 Product annotation for IW/EW SLC
     For Elevation Antenna Pattern (EAP) correction
     """
 
@@ -477,8 +477,9 @@ def load_from_zip_file(cls, path_aux_cal_zip: str, pol: str, str_swath: str):
             Path to the AUX_CAL .zip file
         pol: str {'vv','vh','hh','hv'}
             Polarization of interest
-        str_swath: {'iw1','iw2','iw3'}
-            IW subswath of interest
+        str_swath: IW -> {'iw1','iw2','iw3'},
+            EW -> {'ew1','ew2','ew3','ew4','ew5'}
+            subswath of interest
 
         Returns
         -------
@@ -928,7 +929,7 @@ def compute_thermal_noise_lut(self, shape_lut):
 
 @dataclass
 class BurstCalibration:
-    """Calibration information for Sentinel-1 IW SLC burst"""
+    """Calibration information for Sentinel-1 IW/EW SLC burst"""
 
     basename_cads: str
     azimuth_time: datetime.datetime = None
@@ -995,7 +996,7 @@ def from_calibration_annotation(
 
 @dataclass
 class BurstEAP:
-    """EAP correction information for Sentinel-1 IW SLC burst"""
+    """EAP correction information for Sentinel-1 IW/EW SLC burst"""
 
     # from LADS
     freq_sampling: float  # range sampling rate

src/s1reader/s1_burst_id.py

@@ -2,14 +2,19 @@
 from dataclasses import dataclass
 from typing import ClassVar
 from s1reader.s1_orbit import T_ORBIT
+from s1reader.constants import SENSOR_MODE_MID_SWATH
 import numpy as np
 
 
 @dataclass(frozen=True)
 class S1BurstId:
     # Constants in Table 9-7 of Sentinel-1 SLC Detailed Algorithm Definition
-    T_beam: ClassVar[float] = 2.758273  # interval of one burst [s]
-    T_pre: ClassVar[float] = 2.299849  # Preamble time interval [s]
+    # https://sentiwiki.copernicus.eu/__attachments/1673968/DI-MPC-IPFDPM%20-%20Sentinel-1%20Level%201%20Detailed%20Algorithm%20Definition%202022%20-%202.5.pdf?inst-v=4318c067-be91-4544-be2e-16af66246c9f
+    # Leave IW values as T_beam, T_pre for backwards compatibility
+    T_beam: ClassVar[float] = 2.758273  # interval of one IW burst [s]
+    T_pre: ClassVar[float] = 2.299849  # Preamble time interval IW [s]
+    T_beam_ew: ClassVar[float] = 3.038376  # interval of one EW burst [s]
+    T_pre_ew: ClassVar[float] = 2.299970  # Preamble time interval EW [s]
     T_orb: ClassVar[float] = T_ORBIT  # Nominal orbit period [s]
     track_number: int
     esa_burst_id: int
@@ -41,7 +46,7 @@ def from_burst_params(
             Relative orbit number at the start of the acquisition, from 1-175.
         end_track : int
             Relative orbit number at the end of the acquisition.
-        subswath : str, {'IW1', 'IW2', 'IW3'}
+        subswath : str, {'IW1', 'IW2', 'IW3'} or {'EW1', 'EW2', 'EW3', 'EW4', 'EW5'}
             Name of the subswath of the burst (not case sensitive).
 
         Returns
@@ -60,34 +65,48 @@ def from_burst_params(
         ESA Sentinel-1 Level 1 Detailed Algorithm Definition
         https://sentinels.copernicus.eu/documents/247904/1877131/S1-TN-MDA-52-7445_Sentinel-1+Level+1+Detailed+Algorithm+Definition_v2-4.pdf/83624863-6429-cfb8-2371-5c5ca82907b8
         """
+        # map the subswath onto the sensor mode
+        sensor_mode = {"I": "iw", "E": "ew"}[subswath[0].upper()]
         swath_num = int(subswath[-1])
         # Since we only have access to the current subswath, we need to use the
         # burst-to-burst times to figure out
-        #   1. if IW1 crossed the equator, and
-        #   2. The mid-burst sensing time for IW2
+        #   1. if IW1/EW1 crossed the equator, and
+        #   2. The mid-burst sensing time for IW2/EW3
+        # for IW Mode
         # IW1 -> IW2 takes ~0.83220 seconds
         # IW2 -> IW3 takes ~1.07803 seconds
         # IW3 -> IW1 takes ~0.84803 seconds
-        burst_times = np.array([0.832, 1.078, 0.848])
-        iw1_start_offsets = [
-            0,
-            -burst_times[0],
-            -burst_times[0] - burst_times[1],
-        ]
-        offset = iw1_start_offsets[swath_num - 1]
-        start_iw1 = sensing_time + datetime.timedelta(seconds=offset)
-
-        start_iw1_to_mid_iw2 = burst_times[0] + burst_times[1] / 2
-        mid_iw2 = start_iw1 + datetime.timedelta(seconds=start_iw1_to_mid_iw2)
+        # for EW mode
+        # EW1 -> EW2 takes ~ 0.68268 seconds
+        # EW2 -> EW3 takes ~ 0.55873 seconds
+        # EW3 -> EW4 takes ~ 0.61234 seconds
+        # EW4 -> EW5 takes ~ 0.56538 seconds
+        # EW5 -> EW1 takes ~ 0.61925 seconds
+        EW_BURST_TIMES = np.array([0.68268, 0.55873, 0.61234, 0.56538, 0.61925])
+        IW_BURST_TIMES = np.array([0.832, 1.078, 0.848])
+        burst_times = {"iw": IW_BURST_TIMES, "ew": EW_BURST_TIMES}[sensor_mode]
+
+        # generalise offset for swath 1
+        s1_start_offsets = np.concatenate([[0], -np.cumsum(burst_times[:-1])])
+        s1_start_offset = s1_start_offsets[swath_num - 1]
+        start_s_t = sensing_time + datetime.timedelta(seconds=s1_start_offset)
+
+        # generalise offset to mid swath, array indexed at zero (e.g. swath 1 = idx 0)
+        mid_swath = SENSOR_MODE_MID_SWATH[sensor_mode]  # e.g. 2 for IW, 3 for EW
+        ref_idx = mid_swath - 1  # 0-based index of reference swath
+        start_s_to_mid_s_offset = (
+            float(np.sum(burst_times[:ref_idx])) + burst_times[ref_idx] / 2
+        )
+        mid_s_t = start_s_t + datetime.timedelta(seconds=start_s_to_mid_s_offset)
 
         has_anx_crossing = (end_track == start_track + 1) or (
             end_track == 1 and start_track == 175
         )
 
-        time_since_anx_iw1 = (start_iw1 - ascending_node_dt).total_seconds()
-        time_since_anx = (mid_iw2 - ascending_node_dt).total_seconds()
+        time_since_anx_start_s = (start_s_t - ascending_node_dt).total_seconds()
+        time_since_anx = (mid_s_t - ascending_node_dt).total_seconds()
 
-        if (time_since_anx_iw1 - cls.T_orb) < 0:
+        if (time_since_anx_start_s - cls.T_orb) < 0:
             # Less than a full orbit has passed
             track_number = start_track
         else:
@@ -104,7 +123,9 @@ def from_burst_params(
         dt_b = time_since_anx + (start_track - 1) * cls.T_orb
 
         # Eq. 9-91 :   1 + floor((∆tb − T_pre) / T_beam )
-        esa_burst_id = 1 + int(np.floor((dt_b - cls.T_pre) / cls.T_beam))
+        T_pre = {"iw": cls.T_pre, "ew": cls.T_pre_ew}[sensor_mode]
+        T_beam = {"iw": cls.T_beam, "ew": cls.T_beam_ew}[sensor_mode]
+        esa_burst_id = 1 + int(np.floor((dt_b - T_pre) / T_beam))
 
         return cls(track_number, esa_burst_id, subswath)
 

src/s1reader/s1_burst_slc.py

@@ -236,7 +236,7 @@ class Sentinel1BurstSlc:
     azimuth_time_interval: float
     slant_range_time: float
     starting_range: float
-    iw2_mid_range: float
+    iw2_mid_range: Optional[float]
     range_sampling_rate: float
     range_pixel_spacing: float
     shape: tuple()
@@ -279,6 +279,8 @@ class Sentinel1BurstSlc:
 
     burst_misc_metadata: SimpleNamespace
 
+    ew3_mid_range: Optional[float]
+
     def __str__(self):
         return f"Sentinel1BurstSlc: {self.burst_id} at {self.sensing_start}"
 
@@ -611,7 +613,11 @@ def bistatic_delay(self, range_step=1, az_step=1):
 
         pri = 1.0 / self.prf_raw_data
         tau0 = self.rank * pri
-        tau_mid = self.iw2_mid_range * 2.0 / isce3.core.speed_of_light
+
+        if self.iw2_mid_range is not None:
+            tau_mid = self.iw2_mid_range * 2.0 / isce3.core.speed_of_light
+        else:
+            tau_mid = self.ew3_mid_range * 2.0 / isce3.core.speed_of_light
 
         slant_vec, az_vec = self._steps_to_vecs(range_step, az_step)
 
@@ -624,6 +630,8 @@ def bistatic_delay(self, range_step=1, az_step=1):
         # currently we have not been able to verify this from ESA documents.
         # This implementation follows the Gisinger et al. (2021) for now, we
         # can revise when we hear back from ESA folks.
+        # For EW, the mid of third subswath is used and will similarly need to
+        # be verified
         bistatic_correction_vec = tau_mid / 2 + tau / 2 - tau0
         ny = az_vec.size
         bistatic_correction = np.tile(bistatic_correction_vec.reshape(1, -1), (ny, 1))
@@ -1142,7 +1150,7 @@ def width(self):
 
     @property
     def swath_name(self):
-        """Swath name in iw1, iw2, iw3."""
+        """Swath name in iw1, iw2, iw3, ew1, ew2, ew3, ew4, ew5"""
         return self.burst_id.subswath.lower()
 
     @property

src/s1reader/s1_info.py

@@ -15,18 +15,30 @@
 import shapely.ops
 
 import s1reader
+from s1reader.constants import SENSOR_MODE_SUBSWATHS
 
 
 def get_bursts(
-    filename: Union[Path, str], pol: str = "vv", iw: Optional[int] = None
+    filename: Union[Path, str],
+    pol: str = "vv",
+    iw: Optional[int] = None,
+    ew: Optional[int] = None,
 ) -> list[s1reader.Sentinel1BurstSlc]:
-    if iw is not None:
-        iws = [iw]
+
+    # get the sensory acquisition mode
+    sensor_mode = str(Path(filename).name).split("_")[1].lower()
+
+    if (iw is not None) and (sensor_mode == "iw"):
+        swath_nums = [iw]
+    elif (ew is not None) and (sensor_mode == "ew"):
+        swath_nums = [ew]
     else:
-        iws = [1, 2, 3]
+        # set to all slc swaths
+        swath_nums = SENSOR_MODE_SUBSWATHS[sensor_mode]
+
     burst_nested_list = [
-        s1reader.load_bursts(filename, None, iw, pol, flag_apply_eap=False)
-        for iw in iws
+        s1reader.load_bursts(filename, None, swath_num, pol, flag_apply_eap=False)
+        for swath_num in swath_nums
     ]
     return list(chain.from_iterable(burst_nested_list))
 
@@ -196,9 +208,16 @@ def get_cli_args():
         "-i",
         "--iw",
         type=int,
-        choices=[1, 2, 3],
+        choices=SENSOR_MODE_SUBSWATHS["iw"],
         help="Print only the burst IDs for the given IW.",
     )
+    parser.add_argument(
+        "-e",
+        "--ew",
+        type=int,
+        choices=SENSOR_MODE_SUBSWATHS["ew"],
+        help="Print only the burst IDs for the given EW.",
+    )
     parser.add_argument(
         "-b",
         "--burst-id",
@@ -245,7 +264,7 @@ def main():
             all_files.append(path)
         elif path.is_dir():
             # Get all matching files within the directory
-            files = path.glob("S1[ABCD]_IW*")
+            files = path.glob("S1[ABCD]_[IE]W*")
             all_files.extend(list(sorted(files)))
         else:
             warnings.warn(f"{path} is not a file or directory. Skipping.")
@@ -263,7 +282,7 @@ def main():
         print(f"Bursts in {path}:")
         print("-" * 80)
         # Do we want to pretty-print this with rich?
-        for burst in get_bursts(path, args.pol, args.iw):
+        for burst in get_bursts(path, args.pol, args.iw, args.ew):
             if args.burst_id:
                 print(burst.burst_id, end=" ")
             else: