Load & xarray¶
Read a clipped or decimated SAR scene straight into an xarray.DataArray, or
write a GeoTIFF. Requires the [load] extra.
to_stack is the multi-date companion: it co-registers several acquisitions
onto one shared grid and returns a (time, y, x) datacube — the step
stackstac / odc-stac play elsewhere in the STAC ecosystem, which can't be
pointed at Umbra because successive passes over a site arrive in whatever UTM
zone and extent each acquisition used. The shared grid is lon/lat by default;
crs="utm" (or any CRS) builds it in projected units instead, so its cells are
equal-area and a cell count is a measurement.
stack_stats reduces such a cube to the answer most multi-date searches are
after: one record per pass (its distribution and the signed decibel change
against the pass before it) plus a net first-to-last record, all plain JSON. It
is what umbra stack --stats prints and what the stack_stats agent tool
returns over MCP / LangChain / LlamaIndex. Its blocks=N argument adds the
spatial half of the answer — the scene cut into an N×N grid, each block
reporting its own net change, a compass label and lon/lat centre, and the pair
of passes it moved most between — so a change confined to one corner, which the
scene-wide mean dilutes, reads as where and when. umbra stack --blocks N
prints the same breakdown. Adding block_series=True (umbra stack
--block-series) keeps each block's whole pass-to-pass sequence rather than
only the interval it moved most in, which is what distinguishes a steady drift
from a single step.
A cube costs max_size² × the number of passes in memory, which is what
caps how much series can be stacked sharp. to_stack(lazy=True) (umbra
stack --lazy, the [dask] extra) defers each pass's read into one dask
chunk, and the consumers that reduce a cube — stack_stats and
stack_to_geotiff — walk it a slice at a time, so peak memory follows the
grid rather than the length of the series. The numbers are identical; only
what is resident differs. chunk_size=N (umbra stack --lazy --chunk-size N)
takes the same step within a pass: each slice is cut into N-square windows
read and written independently, so max_size stops being bounded by how much of
one scene fits in memory.
stack_stats(windowed=True) (umbra stack --stats-windowed) measures those
windows rather than whole passes, so a cube stacked sharper than a slice you can
hold is measurable and not only writable. Every count, mean, standard deviation
and change number stays exact — each is a sum, so a window folds in — while each
pass's median / p5 / p95 become histogram estimates good to about 0.05 dB,
because a percentile is the one statistic that needs the whole pass at once. The
summary says which it is (quantile_method / quantile_bin_db, plus a caveat),
so the two kinds of number are never confused.
to_stack(speckle_filter=...) (umbra stack --speckle-filter) averages
speckle down in every pass before the series is assembled — the one
uncertainty in those numbers that no correction in the conversion pipeline
touches, and the largest: a single look's power scatters about its surface's
true backscatter as widely as its own mean, so an unfiltered cell-to-cell
difference is mostly interference rather than change. "boxcar" averages the
window unconditionally (the multilook); "lee" averages only where a window is
no more variable than speckle alone explains, so edges and points survive. It is
the only surface that reaches Umbra's published GEC rasters —
umbra convert --speckle-filter filters complex products in the
radar's own image space, before geocoding. Opt-in, because what it spends is
resolution, and both halves are recorded in the cube's provenance (the same
speckle_filter / speckle_window keys a converted raster carries), so
stack_stats states the trade and a later stack refuses to difference a
filtered cube against an unfiltered pass.
Load an Umbra SAR image as a georeferenced :class:xarray.DataArray.
Reads band 1 of the item's geocoded GeoTIFF (the GEC cloud-optimized
GeoTIFF by default) and returns it as a 2D DataArray with y / x
coordinate axes in the raster's native CRS, ready for the scientific Python
stack. Only the requested window and resolution are streamed via HTTP range
requests -- no full download.
Parameters¶
item:
The acquisition to load.
asset:
Which product to read. Defaults to "GEC"; "CSI" (the
single-band color-sub-aperture GeoTIFF) also works. The complex
SICD / CPHD products are not amplitude rasters and aren't
supported here.
bbox:
Optional (min_lon, min_lat, max_lon, max_lat) in EPSG:4326. When
given, only that geographic window is read (reprojected to the
raster's CRS first). Useful for pulling a small area out of a large
scene without reading the whole thing.
max_size:
Optional cap on the longest output side in pixels. The raster is
decimated to fit (GDAL pulls the matching cloud-optimized GeoTIFF
overview, so this is cheap). None reads full resolution -- which
for a multi-GB scene can be a lot of data; pair large reads with a
bbox or a max_size.
db:
Convert linear amplitude to decibels (20*log10(amplitude)), the
radiometrically-meaningful SAR scale. Implies masked=True for the
non-positive pixels log10 can't represent.
masked:
Replace nodata and non-positive amplitudes with NaN so they don't
contaminate statistics. The array is always returned as float32.
Returns¶
xarray.DataArray
Dimensions ("y", "x") with descending y (north-up) and
ascending x cell-center coordinates. attrs carry the CRS
(crs, a WKT/PROJ string), the affine transform (a 6-tuple),
the geographic bounds, units, and acquisition metadata
(item_id, datetime, platform, product_type), plus the
Umbra license and attribution you must carry with derived
products. The CRS string round-trips through rasterio.crs.CRS /
pyproj and rioxarray (da.rio.write_crs(da.attrs["crs"])).
A raster ``umbra convert`` produced also carries a ``provenance`` dict
-- exactly what :func:`~umbra_py.convert.read_conversion_tags` reads off
the file (``calibration``, ``rtc_model``, ``scale``, ``dem``, ...) -- so
an array knows whether its values are a physical backscatter coefficient
or relative brightness. The key is absent for Umbra's published products,
which carry no such tags.
to_stack(items, *, asset='GEC', bbox=None, max_size=1024, db=False, extent='intersection', crs=None, lazy=False, chunk_size=None, speckle_filter=None, speckle_window=SPECKLE_WINDOW_DEFAULT)
Co-register several acquisitions into one (time, y, x) datacube.
The multi-date companion to :func:to_xarray, and the missing primitive
between search and analysis: hand it a search result and get back a
labelled :class:xarray.DataArray whose slices are pixel-aligned, so
cube.mean("time"), cube.std("time") and cube.diff("time") are
honest per-ground-cell statistics rather than per-scene ones.
Alignment is real work, not a reshape: Umbra's passes over a site are
delivered in whatever UTM zone each acquisition used and at whatever extent
it happened to cover, so every scene is warped to one shared grid derived
from the requested extent -- EPSG:4326 (lon/lat) by default, or the
projected CRS crs names when the cells have to be equal-area. Only
decimated overviews are streamed via HTTP range requests -- no full download.
Parameters¶
items:
The acquisitions to stack. Order doesn't matter (the result is sorted
oldest-first); each must carry a datetime.
asset:
Which product to read, as in :func:to_xarray. Stack one
polarization: mixing VV and VH puts a polarization difference on the
time axis where you'll read it as change. UmbraItem.polarizations
and search(polarizations=...) are how you keep them apart.
bbox:
Optional (min_lon, min_lat, max_lon, max_lat) in EPSG:4326 clipping
the cube to a sub-area of whatever extent selected.
max_size:
Longest side of the output grid in pixels. The grid is shared by every
slice, so this caps the whole cube (bytes fetched grow ~quadratically).
db:
Return the decibel scale (20*log10(amplitude)) instead of linear
amplitude -- the radiometrically meaningful scale for differencing,
where a ratio of backscatter becomes a subtraction.
extent:
One of :data:STACK_EXTENTS. "intersection" (default) keeps only
the ground every acquisition covers, so no cell has a gap; it raises
if the footprints don't all overlap. "union" keeps all ground any
acquisition covers and fills each slice outside its own footprint with
NaN.
crs:
CRS of the shared output grid. None (default) builds it in lon/lat
(EPSG:4326), whose cells are not equal-area -- see the note below.
:data:STACK_AUTO_CRS ("utm") picks the UTM zone containing the
stacked ground, giving square metre-sized cells without your having to
know the zone; any other value is a CRS name ("EPSG:32633", a PROJ
or WKT string) warped to as given. bbox stays lon/lat either way.
lazy:
Defer each pass's read into a dask task instead of streaming the
whole series up front -- one chunk per acquisition. The grid is still
resolved eagerly (the sources' footprints decide it, and a bad extent
or bbox still raises here rather than at compute time), but no pixels
are fetched until something asks for them, and a reduction that walks the
cube one slice at a time -- cube.mean("time"),
:func:stack_stats, :func:stack_to_geotiff -- never holds more than a
few slices at once. That is what lifts the ceiling this function
otherwise has: an eager cube costs max_size² × the number of
acquisitions in RAM, so a long series has to be stacked coarse. The
values are identical either way. Requires the dask extra
(pip install "umbra-py[dask]").
chunk_size:
Cut each pass into chunk_size-square windows instead of reading it
as one slab -- the second half of the ceiling lazy lifts. One
chunk per acquisition makes the unit of work a whole slice, so a single
pass at a large max_size is still read and held whole (a 8192-pixel
grid is 256 MB of float32 per slice); windowing makes the unit
chunk_size² instead, so how sharp a cube can be stacked stops
depending on how much of one scene fits in memory. Requires lazy.
It costs range requests -- each window opens the source and reads its
own bytes, so a pass costs ⌈h/c⌉ × ⌈w/c⌉ reads rather than one -- which
is why it is opt-in and why the window wants to be a decent fraction of
the grid (512–2048), not a tile. The values are unchanged (including
under speckle_filter, which reads each window with a halo -- see it
below for the one number a chunked "lee" estimates rather than reads).
speckle_filter:
Average speckle down before the series is assembled: one of
:data:~umbra_py.convert.SPECKLE_FILTERS ("boxcar", the multilook,
or "lee", which averages only where a window is no more variable than
speckle alone explains), applied to each pass on the shared grid in the
power domain. Speckle is not sensor noise -- it is the interference
pattern coherent illumination makes on a rough surface, so a single
look's power scatters about the surface's true backscatter with a
standard deviation equal to its mean, and averaging is the only thing
that removes it. It is therefore the dominant uncertainty in every number
:func:stack_stats reports from an unfiltered cube, and the reason a
cell-by-cell difference between two passes is mostly interference rather
than change.
Opt-in, because what it spends is resolution: a window that averages N
cells reports ground N cells across. Both halves are recorded in the
cube's ``provenance`` (``speckle_filter`` / ``speckle_window``, the same
keys ``umbra convert`` writes), so :func:`stack_stats` states the trade,
a GeoTIFF written from the cube carries it, and a later stack refuses to
difference this cube against an unfiltered one.
The filtering happens *after* co-registration -- this is the first point
a source exists on the cube's own grid -- so the window averages the
cells the cube reports rather than the source's own pixels. Filtering
earlier, in the radar's image space where speckle is one independent
sample per pixel, is ``umbra convert --speckle-filter``'s job; sources
that already record one are refused here rather than filtered twice.
Composes with ``chunk_size``: each window is read with a half-window halo
and cropped after filtering, so a filtered window holds the cells the
whole-pass filter would have put there, and ``"lee"``'s speckle parameter
is resolved once per pass rather than per window (from a fixed sample of
it, since a chunked build is the case where the pass does not fit in
memory -- a pass small enough to sample whole gives the identical
number). ``boxcar`` needs no such parameter, so a chunked ``boxcar`` cube
is cell-for-cell the unchunked one.
speckle_window:
Edge of the odd, centred window speckle_filter averages over, in
cells of the shared grid. Wider removes more speckle and more detail; it
costs no more to compute (the filters use a summed-area table).
Returns¶
xarray.DataArray
Dimensions ("time", "y", "x"): ascending time, descending y
(north-up) and ascending x cell-center coordinates in the cube's CRS
(degrees by default, projected units under crs), plus an item_id
coordinate along time so every slice keeps its provenance. Nodata and
non-positive pixels are NaN and the dtype is always float32.
attrs mirror :func:to_xarray's (crs, transform,
bounds, units, license, attribution), plus the
provenance the sources agree on when they carry one. Backed by NumPy,
or -- with lazy=True -- by a dask array chunked one slice per
acquisition (or chunk_size-square windows within each slice), which
.compute() / .load() turn into the former.
Notes¶
Sources that disagree about what their pixel values are are refused rather
than stacked. A raster umbra convert produced records its calibration,
RTC model and amplitude scale in UMBRA_* GeoTIFF tags, and stacking a
calibrated pass against an uncalibrated one (or against a published GEC,
which carries no tags at all) would put the difference between the two
conversions on the time axis where you would read it as change on the
ground. The rule is :data:MEASUREMENT_PROVENANCE_KEYS, the refusal names
the disagreement and the acquisitions on each side, and what the sources
do agree on is carried into attrs["provenance"] -- so a measurement
from :func:stack_stats, and any GeoTIFF written from the cube, can say
which conversion produced it. speckle_filter= adds to that record rather
than sitting outside it: a cube that averaged its own slices says so in the
same two keys, and is refused against an unfiltered one for the same reason.
The default lon/lat grid stretches with latitude (cells are not equal-area),
the same quick-look approximation umbra change / umbra timescan make.
That is fine at scene scale and for comparing a cell to itself across
dates, which is what a time series does -- but it makes a cell count a poor
proxy for an area, and it distorts distances. Pass crs="utm" (or a
projected CRS of your own) when the answer is "how many hectares changed":
every cell then covers the same ground, so counting them is measuring.
Summarize a datacube's time axis as a JSON-ready statistics series.
The reporting companion to :func:to_stack: the cube itself is an array, but
the question a multi-date search is usually asking — how did this site
change, and by how much? — has a small numeric answer. This reduces the
(time, y, x) cube to one record per pass (distribution statistics plus
the signed change against the pass before it) and one net baseline → latest
record, all plain JSON so it fits a CLI print, a manifest, or an agent tool
result without carrying pixels around.
Complementary to :func:~umbra_py.narrate.compute_change_stats, which cuts
two passes into spatial blocks to say where change sits. This walks the
whole series to say when it happened and how much ground moved — and
with blocks=N it does both at once, cutting every pass into the same
coarse grid so each block answers where and when.
Parameters¶
cube:
A cube from :func:to_stack — dimensions ("time", "y", "x") with the
item_id coordinate and crs / transform / units attributes
it sets. Slices must already be co-registered; nothing is re-gridded here.
change_threshold_db:
How many decibels a cell has to move between two passes to count as
changed. 3 dB (a doubling of backscatter power) is the same default
umbra change --narrate grounds its narration on.
blocks:
Cut the cube into a blocks × blocks grid and report each block
separately (0, the default, skips the breakdown entirely). A scene-wide
mean hides a change that moved one corner hard, so this is what turns
"the site changed 1.4 dB" into "the northeast corner brightened 9 dB,
between the March and April passes".
block_series:
Keep each block's whole pass-to-pass sequence, not just the interval
it moved most in. The steps are computed either way — this only decides
whether they are reported — so it costs payload, not arithmetic:
blocks × blocks × (count − 1) records at most. Requires
blocks. Ask for it when the question is the shape of a block's
history — did it move once and stay, or drift every pass? — which a
single peak interval cannot answer.
windowed:
Measure the cube one window at a time instead of one pass at a time,
following the cube's own chunks (:func:to_stack's chunk_size). The
default reads a whole slice per pass, so a cube stacked sharper than
memory can be written but not measured; this drops the resident
footprint to three windows and lifts that last ceiling.
The trade is stated rather than hidden: every count, mean, standard
deviation and change number is still exact (they are sums, so a window
folds in), but the per-pass ``median``/``p5``/``p95`` become histogram
estimates, good to about one ``_QUANTILE_BIN_DB`` bin — a quantile needs
the whole distribution, which is the one thing a window-by-window walk
never has.
The summary says so (``quantile_method`` / ``quantile_bin_db`` plus a
caveat), so a consumer can tell the two kinds of number apart. An
unchunked cube is one window, i.e. the default read with estimated
percentiles.
Each pass's ``looks`` is a median over measuring blocks either way, and
the blocks are cut from whatever array is in hand — so a window whose
edge is not a whole number of blocks drops its remainder and the two
modes can differ in the last decimal. That is a property of the
diagnostic rather than of this mode: ``looks`` is a read of the scene,
like ``umbra convert``'s ENL pair, not one of the exact sums beside it.
Returns¶
dict
{count, units, product_type, grid, passes, net_change,
change_threshold_db, license, attribution, caveats}, plus
provenance when the cube carries one (:func:to_stack) — the
conversion its slices were made by, which is also what decides whether
the first caveat calls these decibels relative or calibrated. Each entry in
passes carries item_id, datetime, valid_fraction and the
distribution of that pass (mean/median/std/p5/p95, in
the cube's own units), plus looks — that pass's equivalent number
of looks, read off its own blocks — and change_vs_previous (None
for the first pass). net_change compares the first pass to the last.
A multi-pass cube whose looks could be read also carries ``detection``:
what speckle alone does to a change measured at ``change_threshold_db``.
Speckle is not an error bar on the mean, it is the dominant variation in
a single cell — at one look the pass-to-pass decibel difference of
*unchanged* ground has a 7.9 dB spread — so ``changed_fraction`` is only
evidence to the degree it stands clear of what interference produces by
itself. ``detection`` says by how much: ``looks`` (the representative
equivalent looks, the median of the passes that gave a reading),
``cell_sigma_db`` (that spread), ``false_alarm_fraction`` (the share of
unchanged cells speckle alone pushes past the threshold) and
``target_threshold_db`` (the threshold that would hold that share to
``false_alarm_target``, :data:`DETECTION_FALSE_ALARM_TARGET`). The looks
are read off the *cube's* cells rather than the source products, because
they describe the numbers being quoted — :func:`to_stack` decimates onto
a shared grid, which averages speckle down — and a textured scene reads
low, so the floor is an upper bound on the false alarms. It is absent
rather than null on a single-pass cube (no comparison to weigh) and on
one no block could be read from (nothing measured it).
With ``blocks``, an extra ``spatial`` key carries the grid: one record
per block with its ``row``/``col``, a plain-language ``compass`` label,
``bounds`` in the cube's CRS, a ``center_lonlat`` to map or geocode it
by, its ``net_change`` (first → last, same fields as the top-level one)
and its ``peak_interval`` — the consecutive pair of passes that block
moved most between, named by item id and timestamp. Alongside them,
``peak_block`` names the block that moved most overall and ``grid_text``
renders the net signed change as a north-up ASCII heat-grid. With
``block_series`` each block additionally carries the ``series`` those
peaks were picked from — every consecutive step, oldest first, in the
same shape as ``peak_interval``.
The document is public API, pinned by
``docs/schemas/stack-stats.schema.json`` — the same shape ``umbra stack
--stats --json`` carries under its manifest's ``stats`` key, ``POST
/artifacts/stats`` returns, and the ``stack_stats`` agent tool hands a
model.
Change is **always** reported in decibels — a ratio of backscatter is a
difference on the log scale — whether the cube holds dB or linear
amplitude, so the numbers mean the same thing either way.
``changed_area_km2`` is ``None`` unless the cube's grid is projected
(see :func:`to_stack`'s ``crs``), because counting geographic cells
measures nothing.
Load an Umbra SAR image and write it to dest as a GeoTIFF.
A file-producing companion to :func:to_xarray for users who want a
clipped / decimated raster on disk (for QGIS, GDAL, or any GIS) rather
than an in-memory array. Same windowing and resolution options: bbox
clips to a lon/lat area, max_size decimates via the cloud-optimized
GeoTIFF overviews, db writes the decibel scale. Only the requested
window/resolution is streamed (no full download).
The output is a single-band float32 GeoTIFF in the source raster's
native CRS, with nodata / non-positive pixels written as NaN
(nodata=NaN) so masking survives the round-trip. Deflate-compressed
and tiled.
stack_to_geotiff(items, dest, *, asset='GEC', bbox=None, max_size=1024, db=False, extent='intersection', crs=None, lazy=False, chunk_size=None, speckle_filter=None, speckle_window=SPECKLE_WINDOW_DEFAULT)
Co-register several acquisitions and write the cube to a GeoTIFF.
The file-producing companion to :func:to_stack, mirroring what
:func:to_geotiff is to :func:to_xarray. The output is a multi-band
float32 GeoTIFF in the cube's CRS (EPSG:4326 unless crs names
another, e.g. "utm" for equal-area cells) -- one band per acquisition,
oldest first -- with each band described by its acquisition timestamp and
the item ids carried in the file tags, so the time axis survives the trip
into QGIS, GDAL or any GIS. Nodata is NaN; deflate-compressed and tiled.
lazy (see :func:to_stack) makes this the memory-bounded path to a big
file: bands are written one at a time, so a series long enough to blow up an
in-memory cube still writes, at the resolution it deserves. Add
chunk_size and a band is written one window at a time too, so a grid
too large for one slice to be resident still writes. The file is
byte-identical however it was read.
speckle_filter (see :func:to_stack) averages speckle down in every
band before it is written, and the file records that it did -- so the raster
a GIS opens says what its values are, and re-stacking it against an
unfiltered product is refused rather than measured.