Skip to content
metnoPublic

About

Sentinel-1/Sentinel-2 SAFE to GeoTIFF quicklook processing (warp, stretch, overviews)

Resources

Stars

4 stars

Watchers

0 watching

Forks

Latest commit

 

History

50 Commits

Folders and files

Repository files navigation

pysent

Sentinel-1 and Sentinel-2 SAFE → GeoTIFF quicklook processing, packaged as a library.

This is the processing core that was previously embedded in the FastAPI-mapserver ingestion service (core/sentinel_s1.py, core/sentinel_s2.py) and duplicated as a vendor/ copy inside the QA notebook bundle. It is now one installable package that both consume, so the service, the notebooks and any new tooling all run the same code path — no vendored copies to re-sync.

from pathlib import Path
from pysent.s2 import S2_DEFAULT_PRODUCTS, process_sentinel_s2_safe

product = "true_color_vegetation"          # ("B4", "B3", "B2")

results = process_sentinel_s2_safe(
    input_dataset="/archive/S2A/2022/03/19/S2A_MSIL1C_20220319T110701.zip",
    output_dir=Path("/out"),
    product_bands={product: S2_DEFAULT_PRODUCTS[product]},
    output_names={product: "scene_true_color.tif"},
    processing_options={"histogram_stretch": True, "compression": "DEFLATE"},
)

What's in it

Module Contents
pysent.s1 Sentinel-1 amplitude quicklooks from SAFE or NetCDF: GCP warp (TPS or polynomial), dB or linear percentile stretch + alpha, tiled/compressed GeoTIFF with overviews. Handles VV/VH, HH/HV and single-pol products.
pysent.s2 Sentinel-2 RGB band combinations from SAFE: one stacked-VRT warp per scene, per-band percentile stretch with gamma (or min/max), tiled/compressed 8-bit RGB with overviews.
pysent.profiles Platform detection (S1 vs S2) from any textual hint, plus per-platform profile defaults.
pysent.archive Catalogue download URL → local archive path; UUID → local SAFE resolution.
pysent.csw CSW record lookup (needs the csw extra).
pysent.qa Benchmark/quality harness: step timing with peak RSS, raster stats, over/under-stretch reports, plotting (needs the qa extra).

Top-level names resolve lazily, so import pysent.profiles works in an environment without the GDAL bindings.

Install

pip install pysent                # processing core
pip install "pysent[csw]"         # + CSW lookup (OWSLib)
pip install "pysent[qa]"          # + benchmark/quality harness (pandas, matplotlib)
pip install "pysent[all]"         # everything except the GDAL bindings

Installing the geo stack

pysent.s1 and pysent.s2 need the GDAL Python bindings (osgeo), which are deliberately not a hard dependency — there is no reliable wheel for them, and every deployment already gets them from the platform. Provide them one of these ways:

# Debian/Ubuntu (what the FastAPI-mapserver image does)
apt-get install python3-gdal python3-rasterio
pip install --no-deps pysent

# conda (recommended for a workstation)
conda install -c conda-forge gdal rasterio
pip install pysent

pip install "pysent[gdal]" exists if you really want pip to build the bindings, but prefer the system or conda package.

Verified against GDAL 3.8, rasterio 1.5, numpy 1.26 on Python 3.12.

Documentation

The docs are executable notebooks, all of which CI runs on every push:

Notebook Covers
01_quickstart the three processing stages, first output
02_sentinel1 polarisations, the GCP warp, stretch percentiles
03_sentinel2 band combinations, percentile vs min/max stretch
04_benchmarks time and memory per step, stretch quality
05_bulk_benchmarks throughput: scenes/hour and memory across worker and thread counts

They run with no Sentinel data at all (using the small committed fixtures), or against a mounted archive for the full pipeline:

docker run --rm -p 8888:8888 -v /path/to/nbsArchive:/data/nbsArchive:ro \
    ghcr.io/metno/pysent-docs:main

See docs/README.md for details, and docs/tuning-and-roadmap.md for the measured baseline and what is worth changing next.

Development

git clone https://github.com/metno/pysent && cd pysent
pip install -e ".[test]"
pytest

The suite runs without any SAFE product: processing tests build synthetic rasters, real-scene tests use the committed fixtures under tests/data/, and the GDAL-dependent ones skip cleanly when osgeo is absent.

To refresh the fixtures from a newer scene:

python scripts/make_test_data.py --days 7

This queries the NBS catalogue for recent Sentinel-1 GRD and Sentinel-2 products over Norway and cuts a window from each directly out of the remote archive using HTTP range requests, so a few hundred KB crosses the network rather than the full 1–8 GB product. tests/data/manifest.json records the provenance.

Bulk processing

For converting many products at once, examples/ holds a reference runner and the measurements behind its defaults:

python examples/bulk_convert.py --input-dir /archive/S2C/2026/09 --output-dir /out

One worker process per scene, resumable (a scene is skipped once its sidecar JSON and outputs are in place), failure-isolated (one bad scene or a worker killed by the OOM killer does not stop the run) and sized from the CPU and memory budget it is actually given. bulk_from_catalogue.py starts from catalogue UUIDs or a search instead of a directory, and slurm/bulk_array.sbatch runs it as a Slurm array job.

Configuration

Processing parameters are passed explicitly via processing_options; the environment variables below only supply defaults when an option is omitted.

Variable Effect
S1_PARALLEL_MODE, S2_PARALLEL_MODE serial (default) / threads / processes (S1 only) fan-out across products
S1_PRODUCT_WORKERS, S2_PRODUCT_WORKERS Worker count for that fan-out
S1_GDAL_NUM_THREADS, S2_GDAL_NUM_THREADS GDAL threads per product: warp NUM_THREADS, plus GDAL_NUM_THREADS (decoding, compression) when products run one at a time
S1_WARP_MEMORY_LIMIT_MB, S2_WARP_MEMORY_LIMIT_MB GDAL warpMemoryLimit
S1_USE_NUMBA Enable the numba fast path for the S1 stretch
GDAL_CACHEMAX GDAL block cache size, applied even if set after GDAL started
GDAL_NUM_THREADS Kept unless gdal_num_threads is given. Don't combine it with parallel_mode="threads": GDAL 3.8 can segfault (a warning is raised)
NBS_ARCHIVE_ROOT Local mount of the archive (pysent.archive)
NBS_SENTINEL_CSW_ENDPOINT, CSW_ENDPOINT Catalogue endpoints (pysent.csw)
NBS_SENTINEL_PLATFORM_PROFILES_JSON Profile overrides when the caller supplies none

Unattended and bulk runs

  • work_dir (processing option): where intermediates go, default output_dir. Each call uses its own temporary directory there, removed whether it succeeds or fails. Outputs are written under a hidden temporary name and renamed when complete, so a killed worker never leaves a truncated file under the final name. For Sentinel-2, budget about 2 bytes per pixel per band the products need: roughly 1.2 GB for the three default products of a 10 m scene. intermediate_compression (default none) trades time for space.
  • gdal_cachemax_mb (processing option): GDAL block cache during the call. GDAL's default of 5 % of RAM applies per process.
  • preset (processing option): full (default) or quicklook, which renders Sentinel-2 at 60 m and Sentinel-1 at 160 m. A three-product S2 scene takes 5 s instead of 39 s, because GDAL then reads the JP2 overviews.
  • CPU budget: defaults come from the CPUs the process may use (affinity mask and cgroup quota, e.g. docker --cpus or Slurm), not the host total. Products run one at a time per process by default, so each gets every thread for its compression; parallel_mode="threads" overlaps them instead. Each product of a call may use all of them, so when running several calls side by side, set gdal_num_threads to each call's share. parallel_mode="processes" runs serially inside a worker process. Its pool uses forkserver or spawn, so guard scripts with if __name__ == "__main__":.
  • Errors (pysent.errors): every product of a call is attempted. If some of several fail, PartialFailure carries .results (the finished products) and .errors. An S2 scene without valid pixels raises EmptySceneError. Both subclass RuntimeError, and GDAL's own message is included.
  • A ready-made runner with these settings, resume and failure isolation is in examples/bulk_convert.py.

Reducing output size

Three processing options decide how large a product is. The platforms want different settings, because Sentinel-1 already defaults to its smallest codec and Sentinel-2 defaults to a lossless one.

Option Default Effect on size
target_resolution 10 m (S2, the finest band used), 40 m (S1) The biggest lever: doubling it quarters the pixel count. An S2 true-colour product is 175 MB at 10 m and 6 MB at 60 m.
preset="quicklook" full Shorthand for target_resolution 60 m (S2) / 160 m (S1), nothing else. An explicit target_resolution wins over it.
compression DEFLATE (S2), JPEG (S1) S2: JPEG or WEBP are several times smaller but lossy. S1: leave it alone.
overview_factors (2, 4, 8, 16) Overviews add about a third; (4, 8, 16) drops most of that, at the cost of a slower first zoom-out. An empty list falls back to the default, so overviews cannot be turned off.
# Sentinel-2
processing_options={
    "target_resolution": 20,
    "compression": "JPEG",              # or "WEBP" for the smallest file
    "overview_factors": (4, 8, 16),
}

# Sentinel-1: keep the default JPEG
processing_options={
    "target_resolution": 80,
    "overview_factors": (4, 8, 16),
}
  • Lossy Sentinel-2 output can fringe at the data edge. The product has no alpha band; 0 is the NoData value (see below), and JPEG and WEBP do not keep exact zeros next to valid pixels. Where transparency at a swath edge matters, stay on DEFLATE and reduce the resolution instead.
  • Another lossless codec gains nothing. On the Sentinel-2 test window ZSTD is the size of DEFLATE and LZW is 20 % larger, against JPEG at about a fifth and WEBP smaller still.
  • Sentinel-1 is gray + alpha, two bands. WEBP refuses that layout, and DEFLATE or ZSTD come out about 2.5 times the size of the default JPEG.
  • There is no quality setting: compression is passed straight to GDAL, so lossy codecs run at GDAL's default quality. block_size is not a size lever.

The codec ratios are indicative: they come from the 512 × 512 fixtures under tests/data/, not from full scenes. From the bulk runner the same options are --preset quicklook and --option target_resolution=20, --option overview_factors=4,8,16 (examples/README.md).

How Sentinel-2 products are rendered

Each band is stretched from its 0.5–99.5 percentile range through a gamma of 0.7 into [1, 255], with 0 reserved for NoData — so real fill stays transparent and no valid pixel ever is. On a hazy scene that lifts the mean from 17/255 (plain min/max, which one bright cloud is enough to flatten) to 44/255, clipping 0.6 % of pixels, nearly all of them cloud tops.

processing_options={
    "histogram_stretch": True,
    "stretch_method": "percentile",     # or "minmax" for the full range
    "stretch_percentiles": (0.5, 99.5),
    "stretch_gamma": 0.7,               # 1.0 for a straight linear ramp
}

How Sentinel-1 products are rendered

Each polarisation is stretched in dB (20*log10(amplitude)) between the 1st and 99th percentile, which is how SAR is normally read: the bright end is compressed so that water, radar shadow and smooth ground show detail instead of going flat black. Validity lives in an explicit alpha band, so S1 NoData was never ambiguous.

processing_options={
    "stretch_method": "db",             # or "linear" for the amplitude itself
    "stretch_percentiles": (1.0, 99.0),
    "speckle_filter": "lee",            # off by default; "none" to keep it off
    "speckle_window": 5,
}

stretch_method="linear" with stretch_percentiles=(2.0, 98.0) reproduces the rendering from before September 2026 exactly.

Speckle filtering is available and off by default. A Lee filter makes a browse image much easier to read — equivalent looks rise from 16 to 139 on a real 40 m scene — but it alters the pixel values (point targets are dimmed, fine texture is averaged away) and adds about 9 s per polarisation:

S1 scene, VV+VH, 8 cores time peak RAM
default (no filter) 20.2 s 1.35 GB
speckle_filter="lee" 38.1 s 1.55 GB
"lee" + preset="quicklook" 23.5 s 0.92 GB

The number of looks is estimated from the warped raster rather than taken from the product type: a GRDH scene has about 4.4 looks natively, but warping to 40 m already averages it to about 16, and a filter told to expect 4.4 would smooth away real texture.

What a product records about itself

Every product carries the settings that determined its pixels as PYSENT_* items in the default GDAL metadata domain, so plain gdalinfo shows them and they travel with the file (a byte copy keeps them):

Metadata:
  PYSENT_PRESET=full
  PYSENT_SPECKLE_FILTER=none
  PYSENT_STRETCH_METHOD=db
  PYSENT_STRETCH_PERCENTILES=1.0,99.0
  PYSENT_TARGET_EPSG=EPSG:32661
  PYSENT_TARGET_RESOLUTION=40.0
  PYSENT_VERSION=0.1.0
  PYSENT_WARP_TRANSFORM=thin-plate-spline
  TIFFTAG_SOFTWARE=pysent 0.1.0

PYSENT_WARP_TRANSFORM (thin-plate-spline or polynomial-gcp, from use_tps) is written for Sentinel-1 SAFE only; PYSENT_SPECKLE_WINDOW and PYSENT_SPECKLE_LOOKS (the looks actually used) only when the Lee filter ran; PYSENT_STRETCH_GAMMA for Sentinel-2, whose PYSENT_STRETCH_METHOD is none when histogram_stretch is off. Each result dict returns the same mapping as provenance, so a caller can record it without reopening the file. The pixels are unchanged: a real S1 and S2 scene give the same pixel digest before and after.

Performance

One 16-core machine, real scenes end to end, Sentinel-2 counting three products per scene:

scenes/hour per worker
Sentinel-2, full resolution 312 ~1.0 GB
Sentinel-2, preset="quicklook" 2090 ~0.4 GB
Sentinel-1, VV+VH 598 ~1.4 GB

Two GDAL threads per worker measured fastest nearly everywhere; the bulk runner sizes itself that way. examples/README.md has the full sweep and the memory budget, and 05_bulk_benchmarks.ipynb re-runs it on your own machine.

Known tuning work

Each is measurable with pysent.qa:

  • The S1 warp dominates its cost (~27 s of 31 s for an S1 GRD). use_tps=False swaps the thin-plate spline for the polynomial GCP transform and saves about 23 %. It shifts features by 62 m RMS (194 m at worst) against the product's own GCP grid — 1.5 px at the default 40 m, which is why it is not the default, but only 0.4 px at the 160 m quicklook grid, where it is close to free. Full distribution in docs/tuning-and-roadmap.md.
  • S1 warps to EPSG:32661 (UPS North) at 40 m by default, which suits the Nordic archive; scenes further south want their own UTM zone via target_epsg.
  • COG output was evaluated and not adopted: slower, for the same size, and the GeoTIFF already carries internal overviews.

Licence

GNU General Public License v3.0 or later — see LICENSE.

About

Sentinel-1/Sentinel-2 SAFE to GeoTIFF quicklook processing (warp, stretch, overviews)

Resources

Stars

4 stars

Watchers

0 watching

Forks

Releases

Packages

Used by

Contributors

Languages