Image Resampling Tutorial

One of the key operations for SLC co-registration is resampling SLC images from one geometry to another. Resampling is the process of moving pixels from one place in an image to another in a new image while accounting for fractional pixel indices via interpolation. Currently, we represent the transformation with a pixel-by-pixel map of “range” offsets (i.e., offsets in the horizontal direction) and “azimuth” offsets (i.e., offsets in the vertical direction). For SLC co-registration, in addition to resampling of the complex pixel data, the user may also need to account for non-zero carrier phases in the azimuth direction (e.g., for native Doppler images) and flattening of the complex phase to account for differences in center frequency between the image and a master image. The isce::image::ResampSlc class contains all relevant operations for SLC co-registration.

Example 1: SLC resampling

For the basic task of image resampling without carrier phase or flattening considerations, let’s look at the following example. Here, we have a crop of an Envisat SLC, azimuth offsets that represent a contraction of the image in the vertical dimension, and range offsets that represent a shearing in the horizontal dimension.

The following is example code to perform the resampling.

import isceextension
import gdal

# Instantiate a ResampSlc object
resamp = isceextension.pyResampSlc()

# Open rasters for input files
inputSlc = isceextension.pyRaster('input.slc')
rgOff = isceextension.pyRaster('range.off')
azOff = isceextension.pyRaster('azimuth.off')

# Create raster for output resampled SLC
outputSlc = isceextension.pyRaster('output.slc', access=1, width=rgOff.width,
                                   length=rgOff.length, numBands=1, driver='ISCE',
                                   dtype=gdal.GDT_CFloat32)

# Run resamp
resamp.resamp(inSlc=inputSlc, outSlc=outputSlc, rgoffRaster=rgOff, azoffRaster=azOff)

First, we created pyRaster objects for all input rasters: the input SLC image, the pixel-by-pixel range offsets, and the pixel-by-pixel azimuth offsets. We then created an output SLC image with the output geometry determined by the either the range or azimuth offset raster. These objects are then passed to a default pyResampSlc instance to perform the resampling. After resampling, we obtain the following image:

Example 2: SLC resampling with carrier phase and flattening

If the user wishes to account for SLCs with native Doppler frequencies in the azimuth direction and flattening with respect to a master SLC, the following code can be used.

import isceextension
import gdal

# Create a polynomial for native Doppler
# Note: 0th order in azimuth, 2nd order in range
doppler = isceextension.pyPoly2d(azimuthOrder=0, rangeOrder=2)
doppler.coeffs = [301.35306906319204, -0.04633312447837377, 2.044436266418998e-06]

# Create an ImageMode for the input image using Envisat parameters
mode = isceextension.pyImageMode()
# Set relevant parameters
mode.wavelength = 0.056
mode.startingRange = 826988.69
mode.rangePixelSpacing = 7.80
mode.prf = 1652.416

# Create an ImageMode for the reference master image
modeRef = isceextension.pyImageMode()
# Set relevant parameters for reference
modeRef.wavelength = 0.057
modeRef.startingRange = 826991.0
modeRef.rangePixelSpacing = 7.80
modeRef.prf = 1652.416

# Instantiate a ResampSlc object
resamp = isceextension.pyResampSlc(doppler=doppler, mode=mode)
resamp.refImageMode = modeRef

# Open rasters for input files
inputSlc = isceextension.pyRaster('input.slc')
rgOff = isceextension.pyRaster('range.off')
azOff = isceextension.pyRaster('azimuth.off')

# Create raster for output resampled SLC
outputSlc = isceextension.pyRaster('output.slc', access=1, width=rgOff.width,
                                   length=rgOff.length, numBands=1, driver='ISCE',
                                   dtype=gdal.GDT_CFloat32)

# Run resamp
resamp.resamp(inSlc=inputSlc, outSlc=outputSlc, rgoffRaster=rgOff, azoffRaster=azOff)

Example 3: SLC resampling with NumPy arrays

As we saw in the raster tutorials, rasters can be created from NumPy arrays via the GDAL gdal_array class. Therefore, range and azimuth offset arrays can be created in memory, dressed as pyRaster objects, and passed to pyResampSlc as before.

import isceextension
import numpy as np
from osgeo import gdal_array
import gdal

# Instantiate a ResampSlc object
resamp = isceextension.pyResampSlc()

# Open rasters for input files
inputSlc = isceextension.pyRaster('input.slc')

# Meshgrids for coordinates for offsets
xvals = np.arange(500)
yvals = np.arange(500)
X, Y = np.meshgrid(xvals, yvals)

# Create offsets
roff = (Y - 250.0) / 5.0
aoff = (Y - 250.0) / 2.5

# Dress offset numpy arrays with gdalarray
roff_ds = gdal_array.OpenNumPyArray(roff)
aoff_ds = gdal_array.OpenNumPyArray(aoff)

# Pass gdal datasets to pyRaster
roff_raster = isceextension.pyRaster('', dataset=roff_ds)
aoff_raster = isceextension.pyRaster('', dataset=aoff_ds)

# Create raster for output resampled SLC
outputSlc = isceextension.pyRaster('output.slc', access=1, width=roff_raster.width,
                                   length=roff_raster.length, numBands=1, driver='ISCE',
                                   dtype=gdal.GDT_CFloat32)

# Run resamp
resamp.resamp(inSlc=inputSlc, outSlc=outputSlc, rgoffRaster=roff_raster,
              azoffRaster=aoff_raster)