Quick guide to correct a single S1 SLC burst with the S1-ETAD product

This guide should provide the basic steps which need to be performed to correct a single Sentinel-1 SLC burst with the timings provided in the S1-ETAD product.

For detailed explanations of the product format and grid definition, please consult the “Product Format Specification” Document (ETAD-DLR-PS-0014) v1.2.

All required information for timing correction is contained in the NetCDF product. The XML annotation can be consulted for informative purposes.

Notebook setup

[1]:

%matplotlib inline

[2]:

%load_ext autoreload
%autoreload 2

[3]:

import sys

sys.path.append("../..")

[4]:

import datetime

import numpy as np
from scipy.constants import c
from matplotlib import pyplot as plt

1. Select a burst of a S1-SLC product and extract relevant parameters

The following parameters have been extracted by the S1 product:

S1B_IW_SLC__1SDV_20200124T221416_20200124T221444_019964_025C43_78B2.SAFE

In particular it has been used the annotation file:

annotations/s1b-iw1-slc-vh-20200124t221416-20200124t221444-019964-025c43-001.xml

[5]:

swath_name = "IW1"
burst_t0 = datetime.datetime.fromisoformat("2020-01-24T22:14:22.029900")
burst_tau0 = 5.330221013350372e-03
burst_dt = 2.055556299999998e-03
burst_range_pixel_spacing = 2.329562
burst_dtau = burst_range_pixel_spacing / (c / 2)
burst_lines = 1488
burst_samples = 20808

burst_duration = burst_lines * burst_dt
burst_center = burst_t0 + datetime.timedelta(
    seconds=0.5 * burst_lines * burst_dt
)
t_ref = burst_center

[6]:

print("Burst parameters")
print("t0:            ", burst_t0)
print("tau0:          ", burst_tau0)
print("dt:            ", burst_dt)
print("dtau:          ", burst_dtau)
print("lines:         ", burst_lines)
print("sammples:      ", burst_samples)
print("burst duration:", burst_duration)
print("referenc time: ", t_ref, "(burst center, arbitrary choice)")

Burst parameters
t0:             2020-01-24 22:14:22.029900
tau0:           0.005330221013350372
dt:             0.002055556299999998
dtau:           1.554116481475995e-08
lines:          1488
sammples:       20808
burst duration: 3.058667774399997
referenc time:  2020-01-24 22:14:23.559234 (burst center, arbitrary choice)

2. Extract the relevant timing information from the ETAD product

The NetCDF file is structured in hierarchical groups. The first level contains groups for each sub-swath (i.e. IW1-IW3 for IW mode). Each swath group is further divided into groups for the individual Bursts which are named and ordered according to their acquisition start time.

The overall organization of the product is reflected by the Python API provided by the s1etad package, but internal details of XML and NetCDF files are hidden.

Load S1-ETAD product

Please refer to s1etad quick start guide for more information about the installation procedure of the s1etad Python package.

[7]:

import s1etad

[8]:

filename = (
    "data/"
    "S1B_IW_ETA__AXDV_20200124T221416_20200124T221444_019964_025C43_0A63.SAFE"
)

[9]:

eta = s1etad.Sentinel1Etad(filename)

[10]:

eta

[10]:

Sentinel1Etad("data/S1B_IW_ETA__AXDV_20200124T221416_20200124T221444_019964_025C43_0A63.SAFE")  # 0x711d3f50cfe0
Number of Sentinel-1 slices: 1
Sentinel-1 products list:
  S1B_IW_SLC__1ADV_20200124T221416_20200124T221444_019964_025C43_95FB.SAFE
Number of swaths: 3
Swath list: IW1, IW2, IW3
Azimuth time:
  min: 2020-01-24 22:14:16.480938
  max: 2020-01-24 22:14:44.428152
Range time:
  min: 0.005328684957372668
  max: 0.006383362874313361
Grid sampling:
  x: 8.131672451354599e-07
  y: 0.02932551319648094
  unit: s
Grid spacing:
  x: 200.0
  y: 200.0
  unit: m
Processing settings:
  troposphericDelayCorrection: True
  ionosphericDelayCorrection: True
  solidEarthTideCorrection: True
  bistaticAzimuthCorrection: True
  dopplerShiftRangeCorrection: True
  FMMismatchAzimuthCorrection: True

Query ETAD product

[11]:

margin = 0.5
t0_query = burst_center - datetime.timedelta(
    seconds=(burst_duration + margin) / 2
)
t1_query = burst_center + datetime.timedelta(
    seconds=(burst_duration + margin) / 2
)

[12]:

selection = eta.query_burst(
    swath="IW1", first_time=t0_query, last_time=t1_query
)

[13]:

selection

[13]:

	bIndex	pIndex	sIndex	productID	swathID	azimuthTimeMin	azimuthTimeMax
2	7	1	1	S1B_IW_SLC__1ADV_20200124T221416_20200124T2214...	IW1	2020-01-24 22:14:21.994134480	2020-01-24 22:14:25.131964392

Get the burst object

[14]:

swath = eta[swath_name]
burst = swath[selection.bIndex.values[0]]

[15]:

burst

[15]:

Sentinel1EtadBurst("/IW1/Burst0007")  0x711daad9b200
Swaths ID: IW1
Burst index: 7
Shape: (108, 402)
Sampling start:
  x: 0.0
  y: 5.513196480938404
  units: s
Sampling:
  x: 8.131672451354599e-07
  y: 0.02932551319648094
  units: s

Compute ETAD time grids for the selected burst

Start time computation

The burst group contains range and azimuth start times which need to be calculated:

Az = rootGroup:azimuthTimeMin + burst:gridStartAzimuthTime
Rg = rootGroup:rangeTimeMin + burst:gridStartRangeTime

The above formula corresponds to the following expression when one uses the s1etad Python API:

[16]:

eta_az_start = (
    eta.min_azimuth_time - t_ref
).total_seconds() + burst.sampling_start["y"]
eta_rg_start = eta.min_range_time + burst.sampling_start["x"]

NOTE: please note that the t_ref is used to be able to express eta_az_ax_start as a float.

[17]:

print("eta_az_start:", eta_az_start, 'sec (relative to "t_ref")')
print("eta_rg_start:", eta_rg_start, "sec")

eta_az_start: -1.5650995190615955 sec (relative to "t_ref")
eta_rg_start: 0.005328684957372668 sec

Grid axis computation

[18]:

eta_az_ax = eta_az_start + np.arange(burst.lines) * burst.sampling["y"]
eta_rg_ax = eta_rg_start + np.arange(burst.samples) * burst.sampling["x"]

A simpler alternative is to use the burst.get_burst_grid method:

[19]:

eta_az_ax_rel, eta_rg_ax_rel = burst.get_burst_grid()

[20]:

eta_az_ax2 = (eta.min_azimuth_time - t_ref).total_seconds() + eta_az_ax_rel
eta_rg_ax2 = eta.min_range_time + eta_rg_ax_rel

NOTE: in this case the burst.sampling_start is not used because it is already taken into account in burst.get_burst_grid.

Extraction of correction data

The relevant timing corrections are named (in the NetCDF file):

sumOfCorrectionsRg
sumOfCorrectionsAz

They can be accessed with the Python API as follows:

[21]:

corrections = burst.get_correction(s1etad.ECorrectionType.SUM)
rg_correction = corrections["x"]
az_correction = corrections["y"]

[22]:

fig, ax = plt.subplots(1, 2, figsize=[13, 5], sharex=True, sharey=True)

img = ax[0].imshow(rg_correction * 1e9, aspect="auto", origin="lower")
ax[0].set_title("Sum of Range Corrections")
ax[0].set_xlabel("Range [samples]")
ax[0].set_ylabel("Azimuth [lines]")
fig.colorbar(img, ax=ax[0]).set_label(r"$[ns]$")

img = ax[1].imshow(az_correction * 1e3, aspect="auto", origin="lower")
ax[1].set_title("Sum of Azimuth Corrections")
ax[1].set_xlabel("Range [samples]")
# ax[1].set_ylabel('Azimuth [lines]')
fig.colorbar(img, ax=ax[1]).set_label(r"$[ms]$")

../_images/notebooks_S1-ETAD-Product_Quick_User_Guide_44_0.png

3. Resample the regularly sampled ETAD burst grids onto the regularly sampled S1-SLC grid with a bi-linear interpolation

Azimuth times \(t\) and range times \(\tau\) of the burst grid are defined by the following NetCDF parameters:

Az = rootGroup:azimuthTimeMin + burst:gridStartAzimuthTime + burst:azimuth
Rg = rootGroup:rangeTimeMin   + burst:gridStartRangeTime   + burst:range

The relevant ETAD NetCDF corrections which need to be applied are:

sumOfCorrectionsRg: \(\Delta_{\tau}(t, \tau)\)
sumOfCorrectionsAz: \(\Delta_{t}(t, \tau)\)

Computation of the regularly sampled S1-SLC grid

[23]:

burst_t0_sec = (burst_t0 - t_ref).total_seconds()
slc_az_ax = burst_t0_sec + np.arange(burst_lines) * burst_dt
slc_rg_ax = burst_tau0 + np.arange(burst_samples) * burst_dtau

slc_rg, slc_az = np.meshgrid(slc_rg_ax, slc_az_ax)

[24]:

eta_box = np.asarray(
    [
        (eta_rg_ax[0], eta_az_ax[0]),
        (eta_rg_ax[0], eta_az_ax[-1]),
        (eta_rg_ax[-1], eta_az_ax[-1]),
        (eta_rg_ax[-1], eta_az_ax[0]),
        (eta_rg_ax[0], eta_az_ax[0]),
    ]
)
slc_box = np.asarray(
    [
        (slc_rg_ax[0], slc_az_ax[0]),
        (slc_rg_ax[0], slc_az_ax[-1]),
        (slc_rg_ax[-1], slc_az_ax[-1]),
        (slc_rg_ax[-1], slc_az_ax[0]),
        (slc_rg_ax[0], slc_az_ax[0]),
    ]
)

fix, ax = plt.subplots()
ax.plot(eta_box[:, 0] * 1e3, eta_box[:, 1], "r", label="eta")
ax.plot(slc_box[:, 0] * 1e3, slc_box[:, 1], "b", label="slc")
ax.set_xlabel("Range [ms]")
ax.set_ylabel("Azimuth [s]")
ax.legend()
ax.grid()

../_images/notebooks_S1-ETAD-Product_Quick_User_Guide_49_0.png

Grid resampling

[25]:

from scipy.interpolate import RectBivariateSpline

rg_corr_interpolator = RectBivariateSpline(
    eta_az_ax, eta_rg_ax, rg_correction, kx=1, ky=1
)  # bi-linear
az_corr_interpolator = RectBivariateSpline(
    eta_az_ax, eta_rg_ax, az_correction, kx=1, ky=1
)  # bi-linear

delta_tau = rg_corr_interpolator(slc_az_ax, slc_rg_ax)
delta_t = az_corr_interpolator(slc_az_ax, slc_rg_ax)

4. Subtract the resampled range corrections \(\Delta_{\tau}(t, \tau)\) from the S1-SLC range times

[26]:

slc_rg -= delta_tau

5. Subtract the resampled azimuth corrections \(\Delta_{t}(t, \tau)\) from the S1-SLC azimuth times

[27]:

slc_az -= delta_t

6. The result is an irregularly spaced burst grid with corrected timings for the S1 SLC data

NOTE: dummy synthetic (float) data are used just to demonstrate the procedure.

[28]:

burst_slc = np.linspace(0, 1, burst_lines).reshape(
    burst_lines, 1
) + np.linspace(0, 1, burst_samples)

NOTE: in this example it is used a standard scipy interpolator. In a real use case a proper interpolator for complex interferometric data shall be used.

NOTE: the interpolation step is very slow. To reduce the computation time data are sub-sampled in the range direction.

[29]:

from scipy.interpolate import RectBivariateSpline

subsamp = 10
slc_interpolator = RectBivariateSpline(
    slc_az_ax, slc_rg_ax[::subsamp], burst_slc[:, ::subsamp]
)

[30]:

interpolated_slc = slc_interpolator(
    slc_az[:, ::subsamp], slc_rg[:, ::subsamp], grid=False
)

[31]:

extent = (0, burst_samples, burst_lines, 0)

fig, ax = plt.subplots(1, 2, figsize=[13, 5], sharex=True, sharey=True)

img = ax[0].imshow(
    burst_slc[:, ::subsamp], aspect="auto", origin="lower", extent=extent
)
ax[0].set_title("Oridinal SLC (dummy) data")
ax[0].set_xlabel("Range [samles]")
ax[0].set_ylabel("Azimuth [lines]")
ax[0].grid()
fig.colorbar(img, ax=ax[0]).set_label("Amplitude")

img = ax[1].imshow(
    burst_slc[:, ::subsamp] - interpolated_slc,
    aspect="auto",
    origin="lower",
    extent=extent,
)
ax[1].set_title(r"$\Delta$ = data - corrected_data")
ax[1].set_xlabel("Range [samles]")
# ax[1].set_ylabel('Azimuth [lines]')
ax[1].grid()
cb = fig.colorbar(img, ax=ax[1])
cb.set_label(r"$\Delta$")
cb.formatter.set_scientific(True)
cb.formatter.set_powerlimits((0, 0))

../_images/notebooks_S1-ETAD-Product_Quick_User_Guide_63_0.png