Bulk Downloading Satellite Imagery with Python

# requirements.txt
pystac-client==0.8.*
requests==2.32.*
rasterio==1.4.*
tqdm==4.67.*
tenacity==9.*
boto3==1.35.*          # only if archiving to S3

pip install pystac-client requests rasterio tqdm tenacity boto3

import rasterio
from rasterio.env import Env

with Env(GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR", CPL_VSIL_CURL_CACHE_SIZE="200000000"):
    with rasterio.open("https://sentinel-cogs.s3.amazonaws.com/sentinel-s2-l2a-cogs/36/Q/YD/2023/6/S2A_36QYD_20230601_0_L2A/B04.tif") as ds:
        print(ds.profile)  # should print without downloading the full ~100 MB

from __future__ import annotations

import json
import logging
from pathlib import Path
from typing import Iterator

import pystac_client

logger = logging.getLogger(__name__)

STAC_ENDPOINT = "https://earth-search.aws.element84.com/v1"
COLLECTION = "sentinel-2-l2a"


def discover_scenes(
    bbox: list[float],          # [west, south, east, north] in EPSG:4326
    date_range: str,            # e.g. "2023-01-01/2023-12-31"
    max_cloud_cover: float = 20.0,
    cache_path: Path | None = None,
) -> list[dict]:
    """
    Query a STAC catalog and return a deduplicated list of scene asset dicts.
    Results are cached as newline-delimited JSON to avoid re-querying on restart.
    """
    if cache_path and cache_path.exists():
        logger.info("Loading scene list from cache: %s", cache_path)
        return [json.loads(line) for line in cache_path.read_text().splitlines() if line]

    client = pystac_client.Client.open(STAC_ENDPOINT)
    search = client.search(
        collections=[COLLECTION],
        bbox=bbox,
        datetime=date_range,
        query={"eo:cloud_cover": {"lt": max_cloud_cover}},
        sortby=["+datetime"],
        max_items=None,          # let pystac-client paginate
    )

    scenes: list[dict] = []
    for page in search.pages():
        for item in page:
            for band, asset in item.assets.items():
                if band in ("B04", "B03", "B02", "B08", "SCL"):
                    scenes.append({
                        "scene_id": item.id,
                        "band": band,
                        "href": asset.href,
                        "datetime": item.datetime.isoformat() if item.datetime else None,
                    })
        logger.info("Discovered %d assets so far…", len(scenes))

    if cache_path:
        cache_path.parent.mkdir(parents=True, exist_ok=True)
        cache_path.write_text("\n".join(json.dumps(s) for s in scenes))
        logger.info("Cached %d assets to %s", len(scenes), cache_path)

    return scenes

import time
import threading
import requests as req

class TokenManager:
    """Thread-safe OAuth2 access token manager with proactive refresh."""

    def __init__(self, token_url: str, client_id: str, client_secret: str) -> None:
        self._url = token_url
        self._id = client_id
        self._secret = client_secret
        self._token: str = ""
        self._expires_at: float = 0.0
        self._lock = threading.Lock()

    def get_token(self) -> str:
        with self._lock:
            # Refresh 60 s before actual expiry to avoid race conditions
            if time.monotonic() >= self._expires_at - 60:
                self._refresh()
            return self._token

    def _refresh(self) -> None:
        resp = req.post(
            self._url,
            data={
                "grant_type": "client_credentials",
                "client_id": self._id,
                "client_secret": self._secret,
            },
            timeout=15,
        )
        resp.raise_for_status()
        payload = resp.json()
        self._token = payload["access_token"]
        self._expires_at = time.monotonic() + payload["expires_in"]
        logger.debug("OAuth2 token refreshed; expires in %ds", payload["expires_in"])

import hashlib
import os
from concurrent.futures import ThreadPoolExecutor, as_completed
from pathlib import Path
from typing import NamedTuple

import requests as req
from tenacity import retry, stop_after_attempt, wait_exponential, retry_if_exception_type
from tqdm import tqdm

CHUNK_SIZE = 1 << 20   # 1 MiB
MAX_WORKERS = 10
MAX_RETRIES = 5


class DownloadResult(NamedTuple):
    scene_id: str
    band: str
    dest: Path
    sha256: str
    success: bool
    error: str | None


@retry(
    stop=stop_after_attempt(MAX_RETRIES),
    wait=wait_exponential(multiplier=1, min=2, max=60),
    retry=retry_if_exception_type((req.exceptions.RequestException, OSError)),
    reraise=True,
)
def _download_asset(
    href: str,
    dest: Path,
    session: req.Session,
) -> str:
    """
    Download a single asset with HTTP Range resumption.
    Returns the hex SHA-256 digest of the complete file.
    """
    dest.parent.mkdir(parents=True, exist_ok=True)
    existing_bytes = dest.stat().st_size if dest.exists() else 0
    headers = {}
    if existing_bytes:
        headers["Range"] = f"bytes={existing_bytes}-"
        logger.debug("Resuming %s from byte %d", dest.name, existing_bytes)

    with session.get(href, headers=headers, stream=True, timeout=30) as resp:
        if resp.status_code == 416:
            # Server says the file is already complete
            pass
        elif resp.status_code in (200, 206):
            mode = "ab" if existing_bytes and resp.status_code == 206 else "wb"
            with dest.open(mode) as fh:
                for chunk in resp.iter_content(CHUNK_SIZE):
                    if chunk:
                        fh.write(chunk)
        else:
            resp.raise_for_status()

    # Compute digest over the complete file
    sha256 = hashlib.sha256()
    with dest.open("rb") as fh:
        for chunk in iter(lambda: fh.read(CHUNK_SIZE), b""):
            sha256.update(chunk)
    return sha256.hexdigest()


def download_scenes_concurrent(
    scenes: list[dict],
    output_dir: Path,
    token_manager: TokenManager | None = None,
) -> list[DownloadResult]:
    """
    Download all assets in parallel. Returns one DownloadResult per asset.
    """
    session = req.Session()
    adapter = req.adapters.HTTPAdapter(
        pool_connections=MAX_WORKERS,
        pool_maxsize=MAX_WORKERS,
        max_retries=0,   # tenacity handles retries
    )
    session.mount("https://", adapter)

    results: list[DownloadResult] = []

    def _worker(scene: dict) -> DownloadResult:
        if token_manager:
            session.headers.update({"Authorization": f"Bearer {token_manager.get_token()}"})
        dest = output_dir / scene["scene_id"] / f"{scene['band']}.tif"
        try:
            digest = _download_asset(scene["href"], dest, session)
            return DownloadResult(scene["scene_id"], scene["band"], dest, digest, True, None)
        except Exception as exc:
            logger.error("Failed %s/%s: %s", scene["scene_id"], scene["band"], exc)
            return DownloadResult(scene["scene_id"], scene["band"], dest, "", False, str(exc))

    with ThreadPoolExecutor(max_workers=MAX_WORKERS) as pool:
        futures = {pool.submit(_worker, s): s for s in scenes}
        with tqdm(total=len(scenes), unit="asset", desc="Downloading") as pbar:
            for fut in as_completed(futures):
                results.append(fut.result())
                pbar.update(1)

    return results

import rasterio
from rasterio.env import Env

EXPECTED_BAND_COUNT = {"B04": 1, "B03": 1, "B02": 1, "B08": 1, "SCL": 1}
EXPECTED_DTYPE = {"B04": "uint16", "B03": "uint16", "B02": "uint16", "B08": "uint16", "SCL": "uint8"}


def validate_raster(result: DownloadResult, catalog_checksum: str | None = None) -> bool:
    """
    Returns True only when both checksum and structural checks pass.
    Moves corrupt files to a quarantine directory rather than deleting them.
    """
    if catalog_checksum and result.sha256 != catalog_checksum:
        logger.error(
            "Checksum mismatch for %s/%s: expected %s got %s",
            result.scene_id, result.band, catalog_checksum, result.sha256,
        )
        _quarantine(result.dest)
        return False

    try:
        with Env(GDAL_DISABLE_READDIR_ON_OPEN="EMPTY_DIR"):
            with rasterio.open(result.dest) as ds:
                assert ds.count == EXPECTED_BAND_COUNT[result.band], \
                    f"Band count mismatch: {ds.count}"
                assert ds.dtypes[0] == EXPECTED_DTYPE[result.band], \
                    f"Dtype mismatch: {ds.dtypes[0]}"
                assert ds.crs is not None, "CRS is None"
                assert ds.transform != rasterio.transform.IDENTITY, "Transform is identity"
    except Exception as exc:
        logger.error("Structural validation failed for %s: %s", result.dest, exc)
        _quarantine(result.dest)
        return False

    return True


def _quarantine(path: Path) -> None:
    quarantine_dir = path.parent.parent / "_quarantine"
    quarantine_dir.mkdir(exist_ok=True)
    path.rename(quarantine_dir / path.name)
    logger.warning("Quarantined corrupt file: %s", path.name)

import sqlite3
from datetime import datetime, timezone

def open_manifest(manifest_path: Path) -> sqlite3.Connection:
    conn = sqlite3.connect(manifest_path)
    conn.execute("""
        CREATE TABLE IF NOT EXISTS completed_assets (
            scene_id  TEXT NOT NULL,
            band      TEXT NOT NULL,
            sha256    TEXT NOT NULL,
            dest      TEXT NOT NULL,
            validated_at TEXT NOT NULL,
            PRIMARY KEY (scene_id, band)
        )
    """)
    conn.commit()
    return conn


def is_already_validated(conn: sqlite3.Connection, scene_id: str, band: str) -> bool:
    row = conn.execute(
        "SELECT 1 FROM completed_assets WHERE scene_id=? AND band=?",
        (scene_id, band),
    ).fetchone()
    return row is not None


def record_validated(conn: sqlite3.Connection, result: DownloadResult) -> None:
    conn.execute(
        "INSERT OR REPLACE INTO completed_assets VALUES (?,?,?,?,?)",
        (result.scene_id, result.band, result.sha256,
         str(result.dest), datetime.now(timezone.utc).isoformat()),
    )
    conn.commit()

def preflight_check(href: str, session: req.Session) -> tuple[bool, int]:
    """Returns (is_available, content_length_bytes)."""
    try:
        resp = session.head(href, timeout=10, allow_redirects=True)
        if resp.status_code == 200:
            return True, int(resp.headers.get("Content-Length", 0))
        return False, 0
    except req.exceptions.RequestException:
        return False, 0

from shapely.geometry import shape
from shapely.ops import unary_union

def deduplicate_by_coverage(scenes: list[dict], aoi_geom) -> list[dict]:
    """Remove scenes whose footprints are already covered by earlier acquisitions."""
    covered = None
    selected = []
    for scene in sorted(scenes, key=lambda s: s["datetime"]):
        footprint = shape(scene["footprint"])
        if covered is None or not covered.contains(footprint):
            selected.append(scene)
            covered = unary_union([covered, footprint]) if covered else footprint
    return selected

from osgeo import gdal

def validate_cog(path: Path) -> bool:
    """Returns True if the file is a valid, fully-formed COG."""
    ds = gdal.Open(str(path))
    if ds is None:
        return False
    is_cog, _, errors = gdal.GetDriverByName("GTiff").ValidateCOG(ds, {})
    if errors:
        logger.warning("COG validation errors for %s: %s", path.name, errors)
    return bool(is_cog) and not errors

import urllib3

# Increase socket buffer sizes to fill high-BDP links
urllib3.util.connection.TCP_NODELAY = True

session = req.Session()
adapter = req.adapters.HTTPAdapter(
    pool_connections=MAX_WORKERS,
    pool_maxsize=MAX_WORKERS * 2,   # headroom for connection churn
)
session.mount("https://", adapter)

# Request HTTP/2 when the provider supports it (requires httpx or h2 extras)
# For standard requests, ensure Keep-Alive is honoured:
session.headers.update({"Connection": "keep-alive"})