OPERA Sentinel-1 RTC backscatter (NASA Earthdata)¶
Download OPERA Radiometric-Terrain-Corrected Sentinel-1 backscatter from the ASF DAAC over a small area and map the VV channel in decibels — the kind of SAR layer used for flood and land-surface monitoring. Cloud-Optimized GeoTIFF raster, read with pyramids. Live query — needs the [earthdata] extra and EDL credentials, and the ASF application authorized for your EDL account (see Authentication); wrapped for nbval-lax safety offline.
ASF auth note: ASF's datapool uses an EDL OAuth redirect that drops a bearer token across hosts (HTTP 401). ASF downloads therefore need username/password (or a
~/.netrcentry) soearthaccesscan hold the session — a bareEARTHDATA_TOKENis not enough — or in-region S3. You must also authorize the Alaska Satellite Facility Data Access application for your EDL account (see Authentication).
Setup¶
The imports and the output directory. earthlens provides the unified EarthLens entry point; pyramids reads the Cloud-Optimized GeoTIFF later, and matplotlib/numpy render the VV channel in decibels.
from datetime import datetime
from pathlib import Path
import numpy as np
from pyramids.dataset import Dataset
from pyramids.feature.bbox import transform as transform_bbox
from pyramids.plot import ColorBar
from earthlens.core import EarthLens
OUT_DIR = Path('earthdata_output')
OUT_DIR.mkdir(exist_ok=True)
AOI = (-121.8, 36.3, -121.5, 36.6) # W, S, E, N -- inland of Monterey Bay
1 · Download the OPERA RTC-S1 scene¶
First describe the request: the earthdata backend, the OPERA RTC-S1 collection, the VV polarization, a short January 2024 window, and a small bounding box near Monterey Bay. Keeping the constructor on its own statement makes the request easy to read and re-run.
scene = EarthLens(
data_source='earthdata',
dataset='OPERA_L2_RTC-S1_V1',
variables=['VV'],
start='2024-01-01',
end='2024-01-13',
aoi=list(AOI),
path=OUT_DIR,
)
download() performs the live ASF query and writes the COG file(s) under OUT_DIR, returning the list of
written paths. A failed query or missing EDL credentials raises here rather than being swallowed, so the cause
is visible.
paths = scene.download(progress_bar=False)
print(len(paths), 'file(s):', [Path(p).name for p in paths][:4])
2 — Pick the granule that covers the AOI, then crop to it¶
The AOI decides which bursts are downloaded, not what they contain: every
OPERA granule is a whole Sentinel-1 burst, roughly 100 x 45 km. Most of the
ones returned here only clip a corner of the box, so the first _VV.TIF in
the list is the wrong scene to map. Keep the granules whose footprint
contains the AOI, then crop one down to it.
vv_paths = [p for p in paths if str(p).upper().endswith('_VV.TIF')]
covering = []
for candidate in vv_paths:
granule = Dataset.read_file(candidate)
west, south, east, north = transform_bbox(granule.bbox, granule.epsg, 4326)
granule.close()
if west <= AOI[0] and south <= AOI[1] and east >= AOI[2] and north >= AOI[3]:
covering.append(candidate)
print(f'{len(covering)} of {len(vv_paths)} VV granules contain the AOI')
if not covering:
raise RuntimeError(
f'no VV granule fully contains {AOI}; widen the date window or relax '
'the containment test above'
)
selected = covering[0]
# OPERA encodes the acquisition instant in the file name, so the plot below can
# name the date it is actually showing rather than repeating a literal.
acquired = datetime.strptime(selected.name.split('_')[4], '%Y%m%dT%H%M%SZ')
print(f'using {selected.name}')
print(f'acquired {acquired:%d %B %Y %H:%M} UTC')
source = Dataset.read_file(selected)
aoi_scene = source.crop(bbox=AOI, epsg=4326)
backscatter = aoi_scene.apply(lambda a: 10.0 * np.log10(np.where(a <= 0, np.nan, a)))
print('cropped to', backscatter.rows, 'x', backscatter.columns, 'cells')
Map the VV channel in decibels¶
Linear power becomes decibels (10 * log10), and the grey ramp is stretched to
the range this scene actually occupies (-16 to -2 dB) rather than a
generic —25 to 0: the median here is about —9 dB, so the wider window
would crush more than half the pixels into black. Dark is smooth ground that
reflects the pulse away — open water, and flooded land; bright is rough,
vegetated or built-up. Greyscale is deliberate: a colour ramp turns SAR speckle
into noise rather than texture.
glyph = backscatter.plot(
cmap='gray',
vmin=-16,
vmax=-2,
colorbar=ColorBar(label='VV backscatter (dB)'),
title=f'OPERA RTC-S1 VV backscatter — {acquired:%d %B %Y}',
)
glyph.ax.xaxis.set_ticks_position('bottom')
glyph.ax.xaxis.set_label_position('bottom')