Multi-temporal Sentinel-2 super-resolution by optimization

Jump to: Quick start | How it works | Scripts | Repo layout | Citation

Reconstruct a high-resolution Sentinel-2 reflectance scene from ~30 low-resolution observations of the same area by exploiting the sub-pixel jitter between repeat passes. The scene is parameterized as a continuous field of 2D Gaussian splats; each S2 observation is the analytic integral of that field over shifted pixel footprints convolved with the sensor PSF. We jointly optimize the splat weights and the per-observation sub-pixel shifts in PyTorch. Pure optimization, with no training data and no neural network.

This repo provides:

• A PyTorch SR solver() that fits a Gaussian splat field to a stack of Sentinel-2 observations using one of four optimization strategies (sr.py

adam

,lbfgs

,adam_lbfgs

,lbfgs_adam

), with PSF / TV / scale / coarse-to-fine knobs exposed on the command line. - A parameterized STAC down() that pulls a cloud-free Sentinel-2 L2A time series over an arbitrary AOI / TOI from the Microsoft Planetary Computer and saves the cropped scenes as Cloud-Optimized GeoTIFFs.download_data.py

Figure 1. A scene near Redmond, WA reconstructed at 8× from 8 Sentinel-2 observations. (Left) One raw S2 observation at the original 10 m grid, nearest-neighbor displayed at the SR output size. (** Middle**) Bicubic upsample of the temporal mean across the 8 observations — clean but soft. (** Right**) The Gaussian-splat fit at 1.25 m, recovering edge structure that no single observation contains. Reproduce with python experiments/psf_sweep.py

; the per-σ outputs are committed under experiments/psf_sweep/.

The achievable resolution gain is bounded by the ratio of the per-observation jitter to the sensor PSF width. For Sentinel-2 the jitter is roughly 1.5 m std and the PSF is roughly 3–5 m wide, so the ~30 shifted observations mostly contain redundant spatial information. Multi-temporal optimization produces clean denoised renders at 2–5 m, but it cannot recover genuine 1 m structure from the data alone. A sharper sensor, larger jitter, or an external structural prior (such as an aerial basemap) would be needed to push further.

git clone https://github.com/calebrob6/s2-superres.git
cd s2-superres
pip install -r requirements.txt

python download_data.py

python sr.py --crop 256 --save-png

ls output/sr_*/

The default sr.py

invocation runs LBFGS coarse-to-fine to recover the splat weights and per-observation shifts, then Adam sharpens the result with a TV anneal followed by an unregularized phase. On a single GPU a 256 × 256 LR crop takes a few minutes; full ~1000 × 1000 scenes take longer.

To run on a different area, pass an AOI bbox and a date range to download_data.py

:

python download_data.py \
    --bbox -122.450,37.700,-122.350,37.800 \
    --start 2024-04-01 --end 2024-10-31 \
    --output-dir data_sf/

python sr.py --data-dir data_sf/ --crop 256 --save-png

Each LR observation is modeled as

y[t,c,i,j] = sum_k w_k[c] * erf_box(mu_k, sigma_eff, pixel(i,j) - delta_t)

where w_k[c]

is the per-channel weight of splat k

, erf_box(...)

is the analytic Gaussian integral over a square LR pixel footprint, and delta_t

is a learnable per-observation 2D shift. Both the splat basis and the PSF are Gaussian, so their convolution collapses to sigma_eff = sqrt(sigma_splat^2 + sigma_psf^2)

evaluated at the splat-pixel offset. There is no HR pixel grid in the forward model and no numerical integration. The gradient flows analytically through erf

to both the weights and the shifts.

The optimizer picks one of four strategies:

Strategy	Description
adam
Joint Adam on weights and shifts at one splat scale. Simple baseline.
lbfgs
Block-coordinate LBFGS, alternating between weights and shifts, run coarse-to-fine across --c2f-levels .
adam_lbfgs
Short Adam warmup at the coarsest level (escapes bad basins), then LBFGS C2F.
lbfgs_adam
LBFGS C2F to convergence, then a long Adam phase that anneals TV to zero and then runs unregularized for sharpening. Default; highest measured edge contrast.

The shared core is in src/:

—src/splat_field.py

SplatField

nn.Module

with the analyticerf

forward model.—src/strategies.py

run_strategy(...)

dispatches to the four strategies above.—src/sharpness.py

edge_contrast

,laplacian_var

,gradient_energy

for evaluating real-data outputs.— phase-correlation shift initialization with parabolic sub-pixel refinement.src/shift_estimation.py

— STAC GeoTIFF stack with cloud masking and COG writer.src/data_.py

Searches the Microsoft Planetary Computer STAC catalog and saves cropped Sentinel-2 L2A scenes as COGs.

python download_data.py \
    --bbox -122.186,47.629,-122.056,47.719 \
    --start 2025-04-01 --end 2025-10-31 \
    --max-cloud-cover 1.0 \
    --bands B02,B03,B04,B08 \
    --min-scenes 32 \
    --output-dir data/

Pass --basemap

to also fetch an ESRI World Imagery aerial mosaic over the same AOI for visual comparison. Run python download_data.py --help

for the full flag list.

python sr.py \
    --data-dir data/ \
    --output-dir output/my_run/ \
    --strategy lbfgs_adam \
    --psf 0.5 \
    --tv 1e-3 \
    --c2f-levels 2,4,8 \
    --n-obs 8 \
    --crop 256 \
    --save-png

Useful flags:

Flag	Description
--strategy
One of {adam, lbfgs, adam_lbfgs, lbfgs_adam} .
--psf
Sensor PSF Gaussian σ in LR pixels (Sentinel-2 is approximately 0.4–0.6).
--tv
Huber-TV regularization weight on the splat weights.
--c2f-levels
Comma-separated splat-scale ladder for the LBFGS strategies.
--n-obs
Number of lowest-cloud observations to use.
--crop
Center-crop size in LR pixels (0 = full scene).
--save-cog
Write a georeferenced Cloud-Optimized GeoTIFF (only when --crop 0 ).

Run python sr.py --help

for everything.

Sweeps sigma_psf ∈ {0.3, 0.4, 0.5, 0.6, 0.7, 0.8}

over the central 256 × 256 crop with lbfgs_adam

and otherwise-default settings. Writes per-σ RGB PNGs, a 2 × 4 comparison grid (input + bicubic + 6 σ values), an edge_contrast

vs sigma_psf

plot, and a metrics CSV to experiments/psf_sweep/. The PNGs in this README come from this script.

python experiments/psf_sweep.py

See experiments/psf_sweep/README.md for the most recent results.

sr.py                       SR entry point
download_data.py            Parameterized STAC down
src/                        Reusable library (SplatField, strategies, sharpness, ...)
experiments/                Experiment scripts and their committed PNG/CSV outputs
images/                     Figures embedded in this README

data/

and output/

are gitignored. Fill them locally with download_data.py

and sr.py

runs.

If you use this repo, please cite it:

@misc{robinson2026s2superres,
  author       = {Robinson, Caleb},
  title        = {{s2-superres}: multi-temporal {Sentinel-2} super-resolution by optimization},
  year         = {2026},
  howpublished = {\url{https://github.com/calebrob6/s2-superres}}
}

MIT. See LICENSE.

source & further reading

github.com — original article