We have created a TileDB-SAR library to perform coherence change detection in a temporal collection of registered SAR images. Our library uses Dask to perform the change detection operation in parallel on a collection of array data. The library also provides a set of command line tools to perform the operations of stacking and processing the input SAR data.
The TileDB-SAR library has been designed to support products from different vendors. Currently, the library supports SLC products from the UAVSAR mission.
Installation
You can run the TileDB SAR code as follows;
docker run -it --rm -u 0 -v /local/path:/data tiledb/tiledb-geospatial /bin/bash
CLI Usage
The TileDB SAR command line tools focus on:
Creating a temporal stack of SAR data
Parallel processing to create an image that highlights change across a scene over time
Stack
rio stack-sar --help
Create TileDB SAR stack.
Options:
--output TEXT Output array.
-t, --type [uavsar] SAR sensor type. [default: 0]
-f, --function [ccd] InSAR function type.
--bands INTEGER... InSAR bands
--window_size INTEGER
--config FILENAME TileDB config.
--tile_x_size INTEGER
--tile_y_size INTEGER
--bbox INTEGER... subset box, minx,miny,maxx,maxy
--n_workers INTEGER number of dask workers
--threads_per_worker INTEGER dask threads per worker
--help Show this message and exit.
An example of stacking two registered scenes from UAVSAR and extracting a bounding box into a TileDB array is:
rio process-stack --help
Process TileDB SAR stack.
Options:
--output TEXT Output array.
-f, --function [ccd] InSAR function type.
--bands INTEGER... InSAR bands
--config FILENAME TileDB config.
--n_workers INTEGER number of dask workers
--threads_per_worker INTEGER dask threads per worker
--help Show this message and exit.
An example of processing a stack of SAR data is:
rio process-stack --function ccd --bands 0 1 --output result stack_arr
In the above command, stack_arr is the input stack above and result is the output array .
The result array contains values between zero and one and indicate the likelihood of their being a change between the two collected SAR scenes. A result value of 1 indicates no change between the two images. This algorithm is most easily demonstrated by creating a complex test image and adding a known shift to the original complex image data:
The code to generate these images is as follows:
import math
import dask.array as da
import matplotlib.pyplot as plt
import numpy as np
from insar import uavsar
from insar.sar import local_ccd
mu, sigma = 0.5, 0.24
def add_change(x):
# use zero for thermal noise in this synthetic dataset
r = np.random.rand()
# weight shift to determine alpha
return x * r + math.sqrt(1 - r**2) * x
np.random.seed(0)
s = np.random.normal(mu, sigma, 10000)
data = s.reshape((100, 100))
arr1 = data + 1.j*data
arr2 = np.copy(arr1)
# treat additive thermal noise as zero in this dataset
# alpha is the change metric we wish to estimate in the interval [0, 1]
# 0 indicates complete change, 1 no change
vfunc = np.vectorize(add_change)
for c in range(0, 100):
if c in range(48, 52):
arr2[:, c] = vfunc(arr2[:, c])
arr2[c, :] = vfunc(arr2[c, :])
da1 = da.from_array(arr1).rechunk((window, window))
da2 = da.from_array(arr2).rechunk((window, window))
change_result = da.map_blocks(local_ccd, da1, da2).compute()
Flight path
It is useful to plot the flight path of a SAR data collection to debug missing output values when necessary.
rio flight-path --help
Usage: rio flight-path [OPTIONS] INPUT_ [INTERVAL]
Create GeoJSON flight path.
Options:
-t, --type [uavsar] SAR sensor type. [default: 0]
--help Show this message and exit.
An example usage;
rio flight-path -t uavsar lathrop/SDelta_15503_04_BC_s3_2x8.llh | fio collect > path.geojson
