Instrumenting gdalwarp with Prometheus Counters

To instrument a reprojection, wrap the gdal.Warp call (or a subprocess gdalwarp invocation) so it increments a Counter for pixels reprojected and bytes read and observes a Histogram of warp seconds, tagging each sample with low-cardinality source and target EPSG labels. Read the output dimensions from the warped dataset to compute pixel counts, keep the EPSG labels drawn from a known projection set, and confirm the numbers land by curling the /metrics endpoint and running a PromQL rate() query. This page shows the full wrapper for both the Python-API and subprocess paths.

When to Use This Pattern

You run reprojections with gdal.Warp or shell out to gdalwarp and need throughput visibility per CRS pair.
You want pixels-per-second and warp-duration percentiles broken down by source and target projection.
Your warp step is a discrete task whose regressions currently only surface as slow downstream tiles.
You need the counters to feed the same Prometheus stack described in Prometheus metrics for raster throughput.

Complete Working Example

The wrapper below runs a warp through the GDAL Python API, reads the resulting dataset to count reprojected pixels, and records three metrics. EPSG labels are stringified integer codes from a bounded set, so cardinality stays flat no matter how many rasters flow through.

import os
import time
import subprocess
from osgeo import gdal
from prometheus_client import Counter, Histogram

gdal.UseExceptions()

PIXELS_REPROJECTED = Counter(
    "gdalwarp_pixels_reprojected_total",
    "Pixels written by a warp, by source and target EPSG",
    labelnames=("src_epsg", "dst_epsg"),
)
BYTES_READ = Counter(
    "gdalwarp_source_bytes_read_total",
    "Source file bytes read by a warp",
    labelnames=("src_epsg",),
)
WARP_SECONDS = Histogram(
    "gdalwarp_seconds",
    "Wall-clock duration of a gdalwarp reprojection",
    labelnames=("src_epsg", "dst_epsg"),
    buckets=(0.5, 1, 2, 4, 8, 16, 32, 64),
)

def warp_with_metrics(src_path: str, dst_path: str,
                      dst_epsg: int, resample: str = "bilinear") -> int:
    # Read source EPSG up front so the label is known before work starts.
    src_ds = gdal.Open(src_path)
    src_srs = src_ds.GetSpatialRef()
    src_epsg: str = src_srs.GetAuthorityCode(None) or "unknown"
    d_epsg: str = str(dst_epsg)
    BYTES_READ.labels(src_epsg=src_epsg).inc(os.path.getsize(src_path))

    start = time.perf_counter()
    warp_opts = gdal.WarpOptions(
        dstSRS=f"EPSG:{dst_epsg}",
        resampleAlg=resample,
        dstNodata=src_ds.GetRasterBand(1).GetNoDataValue(),  # preserve nodata
        multithread=True,
    )
    out_ds = gdal.Warp(dst_path, src_ds, options=warp_opts)
    elapsed = time.perf_counter() - start

    pixels: int = out_ds.RasterXSize * out_ds.RasterYSize * out_ds.RasterCount
    WARP_SECONDS.labels(src_epsg=src_epsg, dst_epsg=d_epsg).observe(elapsed)
    PIXELS_REPROJECTED.labels(src_epsg=src_epsg, dst_epsg=d_epsg).inc(pixels)

    out_ds = None  # flush to disk
    src_ds = None
    return pixels

def warp_subprocess_with_metrics(src_path: str, dst_path: str,
                                 src_epsg: int, dst_epsg: int) -> None:
    # Subprocess variant: gdalwarp emits no metrics itself, so we time the call.
    s_epsg, d_epsg = str(src_epsg), str(dst_epsg)
    cmd = [
        "gdalwarp", "-s_srs", f"EPSG:{src_epsg}", "-t_srs", f"EPSG:{dst_epsg}",
        "-r", "bilinear", "-dstnodata", "0", "-overwrite", src_path, dst_path,
    ]
    BYTES_READ.labels(src_epsg=s_epsg).inc(os.path.getsize(src_path))
    start = time.perf_counter()
    subprocess.run(cmd, check=True, capture_output=True)
    WARP_SECONDS.labels(src_epsg=s_epsg, dst_epsg=d_epsg).observe(
        time.perf_counter() - start
    )
    ds = gdal.Open(dst_path)
    PIXELS_REPROJECTED.labels(src_epsg=s_epsg, dst_epsg=d_epsg).inc(
        ds.RasterXSize * ds.RasterYSize * ds.RasterCount
    )
    ds = None

Both paths preserve nodata — through dstNodata in the API options and -dstnodata on the command line — so fill pixels are never mistaken for valid data, and both derive labels from EPSG codes rather than free-form projection strings.

Parameter & Option Reference

Parameter	Type	Default	Spatial notes
`dst_epsg`	`int`	required	Target EPSG; becomes the `dst_epsg` label. Keep to a known set.
`resample`	`str`	`"bilinear"`	GDAL resampling algorithm; `nearest` for categorical rasters.
`dstNodata`	`float`	source nodata	Must be set or reprojected edges corrupt band statistics.
`buckets`	`tuple`	0.5–64 s	Histogram bounds; widen if warps exceed 64 s.
`src_epsg` label	`str`	read from source	Authority code from the source SRS; `"unknown"` if absent.
`multithread`	`bool`	`True`	Uses multiple cores for the warp; does not affect metrics.

Verification & Testing

First confirm the source and target CRS of a test warp with gdalinfo so you know the EPSG labels the wrapper should emit:

gdalinfo output.tif | grep -A2 "Coordinate System"

Then scrape the metrics endpoint the worker exposes and check the counters carry the expected labels and non-zero values:

curl -s localhost:9105/metrics | grep gdalwarp_
# gdalwarp_pixels_reprojected_total{dst_epsg="3857",src_epsg="4326"} 6.7108864e+07
# gdalwarp_seconds_bucket{dst_epsg="3857",le="8.0",src_epsg="4326"} 3

A minimal correctness check asserts the counter moved and the output CRS actually matches the requested target:

from osgeo import gdal

def test_warp_emits_metrics_and_correct_crs(tmp_path):
    before = PIXELS_REPROJECTED.labels(src_epsg="4326", dst_epsg="3857")._value.get()
    warp_with_metrics("fixture_4326.tif", str(tmp_path / "out.tif"), dst_epsg=3857)
    after = PIXELS_REPROJECTED.labels(src_epsg="4326", dst_epsg="3857")._value.get()
    assert after > before
    ds = gdal.Open(str(tmp_path / "out.tif"))
    assert ds.GetSpatialRef().GetAuthorityCode(None) == "3857"

Finally, once several warps have run, verify live throughput in the Prometheus expression browser with a rate() query over the counter:

rate(gdalwarp_pixels_reprojected_total{dst_epsg="3857"}[5m])

Common Pitfalls

Labeling on the full PROJ string instead of the EPSG code. A PROJ pipeline string is effectively unbounded and detonates cardinality. Always reduce to the integer authority code before using it as a label.
Counting pixels from the source instead of the output. A warp changes raster dimensions, so counting source pixels understates or overstates delivered work. Read RasterXSize/RasterYSize from the warped dataset.
Omitting dstNodata. Without it, reprojected fill areas become valid zero-value pixels that skew band statistics and downstream analytics — a silent correctness bug no latency metric will catch.
Forgetting to flush the output dataset. Leaving out_ds open (not setting it to None) can leave the file partially written when you immediately reopen it to count pixels, producing wrong or zero counts.

Related: This wrapper feeds the metric family defined in Prometheus metrics for raster throughput, pairs with measuring tile generation latency percentiles for the tiling stage, and slots into the broader Observability & Monitoring for Geospatial Pipelines design. See the GDAL Warp documentation for the full option set.

← Back to Prometheus Metrics for Raster Throughput

When to Use This Pattern#

Complete Working Example#

Parameter & Option Reference#

Verification & Testing#

Common Pitfalls#

Related in this section

When to Use This Pattern

Complete Working Example

Parameter & Option Reference

Verification & Testing

Common Pitfalls