Dolphin basic walkthrough¶

This notebook demonstrates the basic usage of the dolphin command line tool to execute the stack-based phase linking workflow. In this notebook, we will

• Download geocoded, co-registered single-look complex (CSLC) radar images from ASF
• Prepare a configuration file for a stack of coregistered single-look complex (SLC) radar images with dolphin config
• Run this configuration file with dolphin run
• Inspect the resulting output interferograms
• Show the most common parameters you may want to change or customize

Setup¶

We first need to install dolphin as outlined in the Getting Started section of the documentation.

If you are running this in Colab, you can install dolphin using pip

In [ ]:

!pip install dolphin

!pip install dolphin

We can check that we have the command line tool correctly installed by running:

In [1]:

!dolphin --help

!dolphin --help

usage: dolphin [-h] [--version] {run,config,unwrap,timeseries,filter} ...

options:
  -h, --help            show this help message and exit
  --version             show program's version number and exit

subcommands:
  {run,config,unwrap,timeseries,filter}

or by importing dolphin in python

In [2]:

import dolphin

dolphin.show_versions()

import dolphin dolphin.show_versions()

dolphin version: 0.36.1.post1.dev1+g4835cec

Python deps:
        h5py: 3.12.1
         jax: 0.4.35
       numba: 0.60.0
       numpy: 2.0.2
 opera-utils: 0.16.0
    pydantic: 2.10.5
      pyproj: 3.7.0
    rasterio: 1.4.3
 ruamel_yaml: None
       scipy: 1.15.1
threadpoolctl: 3.5.0
        tqdm: 4.67.1
  osgeo.gdal: 3.10.1

System:
      python: 3.12.8 | packaged by conda-forge | (main, Dec  5 2024, 14:19:53) [Clang 18.1.8 ]
  executable: /Users/staniewi/miniconda3/envs/mapping-312/bin/python
     machine: macOS-15.3.2-arm64-arm-64bit
Unwrapping packages:
      snaphu: 0.4.1.dev18
       spurt: None
       isce3: 0.24.2
       tophu: None
   whirlwind: None
optional GPU info:
         jax: 0.4.35
gpu_is_available: False

If you have a GPU available to you, you can follow the extra installation set up so that the GPU verion of the workflow run. This can be 5-20x faster than the CPU version, depending on the size of your workstation. Here we will be processing a relatively small area, so the CPU will suffice.

Input dataset¶

To find input data, you can use the ASF search UI to explore and get a list of URLs to download; for our purposes, we will use the OPERA Co-regisred Single-look Complex product, which dolphin can directly process.

The helper functions in the opera-utils library provide wrappers over the ASF library to make it easy to download OPERA CSLCs over a certain region.

We will use a stack of Sentinel-1 SLCs from track 78 over West Texas, where wastewater injection lead to a huge jump in surface displacement in 2022:

In [ ]:

!pip install opera-utils asf_search

!pip install opera-utils asf_search

Since ASF requires a login to download data, you must add your username/password here:

In [4]:

import getpass

asf_username = input("Enter your ASF username: ")
asf_password = getpass.getpass("Enter ASF Password: ")

import getpass asf_username = input("Enter your ASF username: ") asf_password = getpass.getpass("Enter ASF Password: ")

Enter your ASF username:  
Enter ASF Password:  ········

In [ ]:

from pathlib import Path

import asf_search as asf
import opera_utils.download

aoi = (
    "POLYGON((-102.8136 31.3039,-102.5927 31.3039,-102.5927 31.532,-102.8136"
    " 31.532,-102.8136 31.3039))"
)
results, options = opera_utils.download.search_cslcs(
    aoi_polygon=aoi,
    # We want to have the same set of dates for each Burst ID (spatial footprint)
    check_missing_data=True,
    track=78,
    start="2021-06-01",
    end="2022-06-01",
)
best_option = options[0]

slc_dir = Path("input_slcs")
slc_dir.mkdir(parents=True, exist_ok=True)
asf.download_urls(best_option.inputs, path=slc_dir)

from pathlib import Path import asf_search as asf import opera_utils.download aoi = ( "POLYGON((-102.8136 31.3039,-102.5927 31.3039,-102.5927 31.532,-102.8136" " 31.532,-102.8136 31.3039))" ) results, options = opera_utils.download.search_cslcs( aoi_polygon=aoi, # We want to have the same set of dates for each Burst ID (spatial footprint) check_missing_data=True, track=78, start="2021-06-01", end="2022-06-01", ) best_option = options[0] slc_dir = Path("input_slcs") slc_dir.mkdir(parents=True, exist_ok=True) asf.download_urls(best_option.inputs, path=slc_dir)

In the input_slcs directory, we have stored the NetCDF-format SLCs:

In [3]:

!ls input_slcs | head -3

!ls input_slcs | head -3

OPERA_L2_CSLC-S1_T078-165573-IW2_20210606T005130Z_20240726T084639Z_S1B_VV_v1.1.h5
OPERA_L2_CSLC-S1_T078-165573-IW2_20210618T005131Z_20240726T114540Z_S1B_VV_v1.1.h5
OPERA_L2_CSLC-S1_T078-165573-IW2_20210630T005131Z_20240726T145129Z_S1B_VV_v1.1.h5

See the OPERA CSLC documentation for the full filename convention, but the main points are

• T078 is Sentinel-1 track (relative orbit) 78
• 165573 the Burst IDs from ESA's Burst database
• IW2 indicates these are from the first subswath (since the "Burst ID" is repeated for subswaths IW1,2,3.)
• 20210606T005130Z is the acquisition datetime

Note that we specified the data we want is in /data/VV. This is not necessary for other SLC formats (e.g. binary files from ISCE2) which have only one raster layer.

You can process one single stack, or multiple geocoded stacks. If you have different spatial regions, as is the case with OPERA CSLCs, dolphin will form burst-wise interferograms and stitch them before unwrapping.

Let's make a configuration file for all of the bursts:

In [4]:

!dolphin config --slc-files input_slcs/*.h5 --subdataset "/data/VV"

!dolphin config --slc-files input_slcs/*.h5 --subdataset "/data/VV"

Saving configuration to dolphin_config.yaml

Common configuration¶

Strides¶

You can create a downsampled version of the output using --output-options.strides.y and --output-options.strides.x (or --sy and --sx) This will save time/space by creation an output with coarser pixel spacing than your SLCs. For COMPASS outputs, adding --sx 6 --sy 3 will convert the inputs at (5m, 10m) in (x, y) (aka (easting, northing)) to a (30 meter, 30 meter) output.

Specify the working directory¶

Use --work-directory to save all rasters to a different directory other than the one you call dolphin run from.

Specify how many bursts, or patches, to process in parallel¶

By adding the --n-parallel-bursts, you can process separate geocoded bursts at the same time (assuming sufficient resources are available). On a large single stack (like the flat binary files from isce2, or NISAR swaths), this option controls the number of parallel blocks to process during phase linking.

Additionally, use --threads-per-worker to limit the number of CPUs you would like the workflow to use, especially on large workstations.

Mini-stack size for phase linking¶

When processing many (dozens to hundreds) of SLCs, dolphin will use the strategy of batching SLCs into "mini-stacks", rather than processing the entire stack at once. The default is to use batch of 15 dates at a time, but this can be changed with --phase-linking.ministack-size or --ms. For example, you may want to try processing all dates at once (likely a slower processing time), or you may want to adjust the number slightly (if, e.g., you have 46 dates and would rather not leave the last batch with one remaining).

Phase unwrapping¶

dolphin supports multiple options for phase unwrapping. Here, we will use the Python wrapper for SNAPHU, one of the most widely used phase unwrapping algorithms.

Bounds subsetting¶

If you are interested in a smaller subset of the area, you can specify the --output-bounds to force the raster boundaries.

In [ ]:

%%bash
dolphin config \
    --slc-files input_slcs/*.h5  --subdataset "/data/VV" \
    --sx 6 --sy 3 --n-parallel-bursts 2 \
    --output-options.bounds -102.8 31.3 -102.6 31.5 \
    --work-directory work-walkthrough

%%bash dolphin config \ --slc-files input_slcs/*.h5 --subdataset "/data/VV" \ --sx 6 --sy 3 --n-parallel-bursts 2 \ --output-options.bounds -102.8 31.3 -102.6 31.5 \ --work-directory work-walkthrough

Saving configuration to dolphin_config.yaml

Running the workflow¶

Now that we have created the dolphin_config.yaml file, we can run it using dolphin run

In [ ]:

# Make sure we have the SNAPHU wrapper installed:
!pip install snaphu

# Make sure we have the SNAPHU wrapper installed: !pip install snaphu

In [6]:

%%capture outputs

# To see all the output logging, remove the line above

import os

# TQDM doesn't play nicely with notebook outputs, not like terminal
os.environ["TQDM_DISABLE"] = "1"

# Running from the command line
!dolphin run dolphin_config.yaml

%%capture outputs # To see all the output logging, remove the line above import os # TQDM doesn't play nicely with notebook outputs, not like terminal os.environ["TQDM_DISABLE"] = "1" # Running from the command line !dolphin run dolphin_config.yaml

Visualizing the displacement¶

The outputs can be plotted using any tool capable of reading GDAL-compatible rasters. You can also use the dolphin.io.load_gdal function for convenience.

In [7]:

file_list = sorted(Path("work-walkthrough/timeseries/").glob("2*.tif"))
print(f"Found {len(file_list)} timeseries files")

file_list = sorted(Path("work-walkthrough/timeseries/").glob("2*.tif")) print(f"Found {len(file_list)} timeseries files")

Found 25 timeseries files

In [8]:

import matplotlib.pyplot as plt
import numpy as np
from opera_utils import get_dates

from dolphin.io import load_gdal

%matplotlib inline

import matplotlib.pyplot as plt import numpy as np from opera_utils import get_dates from dolphin.io import load_gdal %matplotlib inline

In [25]:

fig, axes = plt.subplots(ncols=2, figsize=(11, 5))
ax = axes[0]

last_displacement = load_gdal(file_list[-1], masked=True)

axim = ax.imshow(
    last_displacement - np.ma.median(last_displacement),
    cmap="RdBu_r",
    vmax=0.15,
    vmin=-0.15,
)
fig.colorbar(axim, label="LOS [meters]")

ax = axes[1]
row, col = np.unravel_index(np.argmax(last_displacement), last_displacement.shape)
dates = [get_dates(f)[1] for f in file_list]
timeseries_values = [
    load_gdal(f, rows=slice(row, row + 1), cols=slice(col, col + 1)) for f in file_list
]

ax.plot(dates, np.array(timeseries_values).squeeze(), ".-")
fig.tight_layout()

fig, axes = plt.subplots(ncols=2, figsize=(11, 5)) ax = axes[0] last_displacement = load_gdal(file_list[-1], masked=True) axim = ax.imshow( last_displacement - np.ma.median(last_displacement), cmap="RdBu_r", vmax=0.15, vmin=-0.15, ) fig.colorbar(axim, label="LOS [meters]") ax = axes[1] row, col = np.unravel_index(np.argmax(last_displacement), last_displacement.shape) dates = [get_dates(f)[1] for f in file_list] timeseries_values = [ load_gdal(f, rows=slice(row, row + 1), cols=slice(col, col + 1)) for f in file_list ] ax.plot(dates, np.array(timeseries_values).squeeze(), ".-") fig.tight_layout()

Apppendix: More on the output folders¶

For each stack of SLCs (which may be > 1 when processing OPERA geocoded burst SLCs), the workflow creates a folder for

1. persistent scatter outputs (PS)
2. linked phase optimized SLCs (linked_phase)
3. (virtual) interferograms formed using the optimized SLCs (interferograms)

Here we have two of these subdirectories named t042_088905_iw1 and t042_088906_iw1. At the top level, we have

• The slc_stack.vrt is a VRT file pointing to the input SLCs for that burst stack.
• The nodata_mask.tif has been created from the OPERA CSLC metadata to skip over the nan regions
• The bounds_mask.tif was created to more quickly process a subset of the CSLCS.

Last, there are a top-level directorie for interferograms that have been stitched together, unwrapped folder for the outputs of phase unwrapping, and timeseries for the result of inverting in unwrapped interferograms (if a network of redundant interferograms were created).

$ tree -L 1 work-walkthrough/
work-walkthrough/
├── dolphin.log
├── interferograms
├── t078_165573_iw2
├── t078_165574_iw2
├── timeseries
└── unwrapped

$ tree -L 1 work-walkthrough/t078_165574_iw2/
work-walkthrough/t078_165574_iw2/
├── PS
├── bounds_mask.tif
├── combined_mask.tif
├── interferograms
├── linked_phase
├── nodata_mask.tif
├── slc_stack.vrt
├── timeseries
└── unwrapped
    

In [ ]: