Testing Geospatial CLI Tools

A geospatial command-line tool is testable in seconds rather than minutes once its fixtures stop touching disk and its assertions stop comparing floating-point coordinates for exact equality. This guide covers the fixture patterns, the assertion styles, and the isolation boundaries that make a GIS test suite fast enough to run on every save. It is part of the CLI Architecture & Design Patterns guide.

Prerequisites

• Python 3.10 or later, for tmp_path typing and structural pattern matching in assertion helpers.
• pip install pytest pytest-cov rasterio pyogrio shapely numpy — rasterio and pyogrio both ship manylinux wheels with GDAL bundled, so no system install is needed for the fixtures below.
• A CLI to test. The examples assume the layered structure described in the parent CLI Architecture & Design Patterns reference, where the command function is thin and the geospatial work lives behind it.
• If your tool is built on Click rather than Typer, the harness in Testing Click Commands with CliRunner applies unchanged — typer.testing.CliRunner wraps the same object.

Problem framing

Geospatial test suites go bad in a predictable way. Someone commits a fixture: a 40 MB GeoTIFF checked into the repository because it reproduces a bug. A year later there are nine of them, the clone is slow, and the suite takes four minutes because every test decompresses a real file. Meanwhile the assertions have drifted toward comparing whole output files byte for byte, so a GDAL upgrade that changes a metadata tag turns the whole suite red for reasons unrelated to the code.

Both problems have the same root: the tests are exercising GDAL rather than exercising your tool. What your tool actually does is parse arguments, decide which windows to read, call a transform, and write output with the right georeferencing. All of that can be checked against rasters that exist only in memory and vectors with four features in them.

What each layer of the suite should own

The three layers below correspond to the interface, orchestration and engine split the section describes, and each one has assertions the others cannot make cheaply.

The rule that keeps the layers honest is that each test stubs everything below it. An interface test that reprojects a real raster is really an end-to-end test wearing a disguise: it is slow, and when it fails you do not know whether the parser or the warp is at fault.

Step-by-step implementation

1. Build rasters in memory, never on disk

rasterio.MemoryFile gives a real GDAL dataset backed by an in-memory buffer. It supports the whole reader and writer API, including windowed reads and overviews, so code under test cannot tell the difference.

# tests/conftest.py
from contextlib import contextmanager

import numpy as np
import pytest
import rasterio
from rasterio.io import MemoryFile
from rasterio.transform import from_origin


@contextmanager
def make_raster(width=256, height=256, bands=1, crs="EPSG:4326",
                origin=(-1.0, 51.6), res=0.001, dtype="float32", fill=None):
    """Yield an open, readable dataset with a known transform and CRS."""
    transform = from_origin(origin[0], origin[1], res, res)
    profile = {
        "driver": "GTiff",
        "width": width,
        "height": height,
        "count": bands,
        "dtype": dtype,
        "crs": crs,
        "transform": transform,
        "tiled": True,
        "blockxsize": 128,
        "blockysize": 128,
    }
    data = (np.arange(width * height, dtype=dtype).reshape(height, width)
            if fill is None else np.full((height, width), fill, dtype=dtype))

    with MemoryFile() as memfile:
        with memfile.open(**profile) as dst:
            for band in range(1, bands + 1):
                dst.write(data, band)
        with memfile.open() as src:      # reopen read-only, as callers expect
            yield src


@pytest.fixture
def wgs84_raster():
    with make_raster() as src:
        yield src

Two details in that profile matter more than they look. Setting tiled=True with an explicit block size means tests exercise the same block-aligned code path production does — a striped fixture will happily pass tests that a tiled input would fail. And building the data from arange rather than zeros means a bug that returns the wrong window produces different values rather than the same zeros everywhere.

2. Give the CLI a path when it insists on one

Most commands take a path, not a dataset. Rather than weakening the command’s signature to accept an open dataset, write the fixture to tmp_path once per test that needs it — tmp_path is on a temporary filesystem and the file is small enough that the write is measured in microseconds.

@pytest.fixture
def raster_path(tmp_path):
    """A small, tiled, EPSG:4326 GeoTIFF on disk, for commands that take a path."""
    path = tmp_path / "input.tif"
    with make_raster(width=256, height=256) as src:
        profile = src.profile
        with rasterio.open(path, "w", **profile) as dst:
            dst.write(src.read())
    return path

Note that this is still not a checked-in fixture. The bytes are produced by code the suite can read, so a reviewer can see exactly what shape the input has without opening it in a GIS.

3. Assert on the contract, not on the bytes

Output comparison is where geospatial suites acquire their fragility. Three assertions carry almost all the value, and none of them compares whole files.

# tests/test_reproject.py
import numpy as np
import rasterio
from pyproj import CRS

from gistools.cli import app
from typer.testing import CliRunner

runner = CliRunner()


def test_reproject_writes_target_crs(raster_path, tmp_path):
    out = tmp_path / "out.tif"
    result = runner.invoke(app, ["raster", "warp", str(raster_path), str(out),
                                 "--dst-crs", "EPSG:3857"])

    assert result.exit_code == 0
    with rasterio.open(out) as dst:
        # 1. the georeferencing contract
        assert CRS(dst.crs) == CRS.from_epsg(3857)
        # 2. the shape contract — non-degenerate, and roughly square as the input was
        assert dst.width > 0 and dst.height > 0
        assert 0.8 < (dst.width / dst.height) < 1.25
        # 3. the data contract — values survived, within resampling tolerance
        band = dst.read(1, masked=True)
        assert band.count() > 0
        assert np.isfinite(band.compressed()).all()

Each assertion answers a question a user would ask. Is it in the CRS I requested? Did it produce a real raster rather than a one-pixel stub? Did the pixel values survive, or is it full of NaN? A byte comparison answers none of those and breaks whenever GDAL changes a default.

4. Compare coordinates with a tolerance derived from the CRS

Floating-point equality is never the right assertion for a transformed coordinate. The tolerance should come from the units of the target CRS, not from a constant someone picked.

from pyproj import CRS


def coordinate_tolerance(crs) -> float:
    """A sensible absolute tolerance for coordinates in `crs`.

    Projected CRSs are usually metres, where a millimetre is far below any
    meaningful precision. Geographic CRSs are degrees, where the same physical
    distance is about 1e-8.
    """
    axis = CRS(crs).axis_info[0]
    unit = (axis.unit_name or "").lower()
    if unit.startswith("degree"):
        return 1e-8
    if "foot" in unit:
        return 1e-2
    return 1e-3      # metres


def assert_bounds_close(actual, expected, crs):
    tol = coordinate_tolerance(crs)
    for a, e, name in zip(actual, expected, ("left", "bottom", "right", "top")):
        assert abs(a - e) < tol, f"{name}: {a} != {e} (tol {tol})"

Deriving the tolerance rather than hardcoding it means the same helper works for a British National Grid test and a WGS 84 test without anyone remembering to change the constant.

5. Isolate the environment, explicitly

GDAL reads a great deal of configuration from the process environment, and a test that passes because of a variable set on the developer’s machine will fail in CI. The fixture below pins the environment for the whole session.

@pytest.fixture(autouse=True)
def clean_gdal_env(monkeypatch):
    """Remove ambient GDAL tuning so tests see identical conditions everywhere."""
    for var in ("GDAL_CACHEMAX", "GDAL_NUM_THREADS", "CPL_DEBUG",
                "GDAL_DISABLE_READDIR_ON_OPEN", "AWS_PROFILE"):
        monkeypatch.delenv(var, raising=False)
    monkeypatch.setenv("GDAL_CACHEMAX", "64")   # small, deterministic

Marking it autouse is deliberate. An opt-in fixture protects only the tests that remember to ask for it, and the tests that forget are exactly the ones that will fail mysteriously later. For the broader picture of which variables matter and where they come from, see Environment Variable Sync.

Configuration integration

A CLI that layers defaults, a config file, environment variables and flags — the arrangement described in Configuration File Management — needs tests at each layer boundary rather than one test per setting. Three tests cover the whole precedence chain: one asserting that a config-file value is used when nothing overrides it, one asserting that an environment variable beats the file, and one asserting that a flag beats the environment.

def test_flag_beats_environment_and_file(tmp_path, monkeypatch, raster_path):
    (tmp_path / "gis.toml").write_text('[raster]\ndst_crs = "EPSG:27700"\n')
    monkeypatch.chdir(tmp_path)
    monkeypatch.setenv("GISTOOL_DST_CRS", "EPSG:3035")

    out = tmp_path / "out.tif"
    result = runner.invoke(app, ["raster", "warp", str(raster_path), str(out),
                                 "--dst-crs", "EPSG:3857"])

    assert result.exit_code == 0
    with rasterio.open(out) as dst:
        assert CRS(dst.crs).to_epsg() == 3857

Testing every setting through every layer is combinatorially large and adds nothing: the merge code does not know which key it is merging. What can go wrong is the ordering, and three tests pin the ordering.

Speed, and where it goes

The reason to care about suite speed is not tidiness — it is that a suite people run before pushing catches more than a suite people run in CI. The chart below shows where the time actually goes in a typical geospatial suite before and after the fixture changes above.

That leftover import cost is worth attacking too, but through session-scoped fixtures and -p no:cacheprovider rather than through anything clever. Beyond a point the honest answer is that a Python process that uses GDAL takes about a second to start, and the suite should be structured so it starts once.

Fixture scope is the other lever, and getting it wrong is how a fast suite becomes slow again. The three scopes below cover every case a geospatial suite needs.

When in doubt, start at function scope and promote only what a --durations run shows to be expensive. Promoting early is how a suite acquires the kind of order-dependent failure that only reproduces on someone else’s machine.

Error handling and gotchas

A masked array’s == is not a boolean. Comparing two masked arrays with == yields another masked array, and assert a == b then raises a truth-value error or, worse, passes on an empty comparison. Use numpy.testing.assert_allclose on .filled(fill_value) or compare .compressed() explicitly.

tmp_path is not shared across parametrised cases. Each parametrised run gets its own directory, which is usually what you want, but a test that builds an expensive fixture and parametrises over ten assertions will build it ten times. Move the build into a scope="module" fixture and keep only the assertion parametrised.

A CliRunner result hides exceptions by default. result.exit_code == 1 with an unhelpful message usually means the command raised something the runner caught. result.exception holds it, and asserting result.exception is None before asserting on the exit code turns a mystery into a traceback.

Nodata does not survive every driver. Writing a nodata value and asserting it round-trips works for GeoTIFF and fails for formats that do not carry the tag. If the test is about nodata handling, pin the driver in the fixture rather than letting it default.

Geometry equality is not geometric equality. Two shapely polygons describing the same area can differ in vertex order or start point, so a == b is False while a.equals(b) is True. For transformed geometry, a.equals_exact(b, tolerance) is usually what you mean.

Verification

A suite is doing its job when three things are true, and each is cheap to check.

# 1. It is fast enough to run on save.
pytest -q --durations=10

# 2. It fails when the code is wrong — mutate one line and confirm red.
#    (a deliberate typo in the CRS argument should break exactly one test)

# 3. It covers the paths that matter, not just the lines.
pytest --cov=gistools --cov-report=term-missing --cov-fail-under=85

The middle one is the check people skip and the one that matters most. Coverage says a line ran; it does not say an assertion would have noticed if the line were wrong. Deliberately breaking one behaviour and confirming that exactly one test turns red is the only cheap evidence that the suite has teeth.

Performance notes

Session-scoped fixtures are the single largest lever. A MemoryFile raster built once per session and reopened per test costs a few microseconds per reopen against a few milliseconds per build. Because MemoryFile datasets are read-only once written, sharing one across tests is safe as long as no test mutates it.

Parallelism with pytest-xdist helps less than it does for pure-Python suites, because each worker pays the GDAL import cost separately. On a four-core machine, expect roughly a two-and-a-half-times speed-up rather than four. It is still worth having, but fixing a slow fixture beats adding workers.

Finally, keep any test that touches the network behind a marker and out of the default run. A suite whose speed depends on someone’s connection is a suite that will be skipped. Guidance on faking that boundary is in Mocking S3 and Network I/O in Geospatial CLI Tests.

FAQ

Keeping the suite trustworthy as it grows

Two failure modes appear once a geospatial suite passes a few hundred tests, and both are worth naming because the fix is structural rather than a matter of discipline.

The first is the accumulating skip. A test that only runs on a particular GDAL version, or only when a network fixture is available, quietly stops running everywhere and nobody notices. Printing the skip reasons in the CI summary — pytest -ra does it — turns an invisible gap into a line someone reads on every run. A skip that appears in every environment is not a conditional test; it is a deleted test that still costs collection time.

The second is the test that passes for the wrong reason. A raster assertion against an output that was never written will fail on the open, which is fine; an assertion against an output written by a previous test in the same directory will pass, which is not. Scoping every output path to tmp_path prevents it structurally, and it is the strongest argument for never letting a test write to a fixed relative path even when doing so would be convenient.

Related

• Testing Click Commands with CliRunner for GIS Tools — the harness these fixtures plug into.
• Matrix Testing a Geospatial CLI Across GDAL Versions — running this suite against the versions your users have.
• CLI Architecture & Design Patterns — the layering that makes each test layer possible.