Property-Based Testing for Round-Trip Transformations

A fixture tests the coordinates someone thought of; a property test tests the ones they did not. This guide, part of testing and CI for coordinate transformation code, applies property-based testing to coordinate transformations — generating positions across the valid domain and asserting invariants that must hold for every one of them, which is how boundary cases get found without anyone having to imagine them first.

The Invariants Worth Asserting #

Four properties hold for any correct transformation implementation, and each fails in a characteristic way.

Figure — four invariants and the specific defect each one fails on.

Round-trip closure. Forward then inverse returns the input to within the numerical floor. It fails on sign errors, transposed matrices, and an inverse that stops iterating too early. The floor is not zero — it is set by the convergence threshold of the inverse and by double precision — so the assertion is a small tolerance, not equality.

Determinism. The same input transformed twice gives bitwise identical output. It fails on threading, on iteration over an unordered collection, and on any dependence on wall-clock time.

Batch invariance. Transforming an array gives the same result as transforming its elements one at a time, and transforming it in two halves gives the same result as transforming it whole. It fails on rounding applied per chunk and on any accidental dependence on array length.

Order independence. Shuffling the input and un-shuffling the output gives the same result. It fails when results are collected in completion order rather than by index — the classic concurrency defect described in concurrent pyproj transformation pipelines.

Complete Runnable Implementation #

from __future__ import annotations

import numpy as np
from hypothesis import given, settings, strategies as st
from pyproj import Transformer

# Bounds of the operation's area of use, not of the coordinate system: generating
# outside the area of use tests extrapolation, which is a refusal, not a round trip.
LAT_MIN, LAT_MAX = 41.0, 49.0
LON_MIN, LON_MAX = -125.0, -116.0

FWD = Transformer.from_crs("EPSG:6318", "EPSG:6339", always_xy=True)
INV = Transformer.from_crs("EPSG:6339", "EPSG:6318", always_xy=True)

positions = st.tuples(
    st.floats(LON_MIN, LON_MAX, allow_nan=False, allow_infinity=False),
    st.floats(LAT_MIN, LAT_MAX, allow_nan=False, allow_infinity=False),
)


@given(positions)
@settings(max_examples=500, deadline=None)
def test_round_trip_closes(pos: tuple[float, float]) -> None:
    """Forward then inverse returns the input to the numerical floor."""
    lon, lat = pos
    e, n = FWD.transform(lon, lat)
    lon2, lat2 = INV.transform(e, n)
    # Degrees, at ~1e-9 -> about 0.1 mm on the ground.
    assert abs(lon2 - lon) < 1e-9, f"longitude closed to {abs(lon2 - lon):.3e} deg"
    assert abs(lat2 - lat) < 1e-9, f"latitude closed to {abs(lat2 - lat):.3e} deg"


@given(positions)
@settings(max_examples=200, deadline=None)
def test_transform_is_deterministic(pos: tuple[float, float]) -> None:
    """Two identical calls give bitwise identical results."""
    first = FWD.transform(*pos)
    second = FWD.transform(*pos)
    assert first == second, "transformation is not deterministic"


@given(st.lists(positions, min_size=1, max_size=200))
@settings(max_examples=100, deadline=None)
def test_batch_matches_scalar(points: list[tuple[float, float]]) -> None:
    """An array call equals the same points transformed one at a time."""
    lons = np.array([p[0] for p in points], dtype=np.float64)
    lats = np.array([p[1] for p in points], dtype=np.float64)
    be, bn = FWD.transform(lons, lats)
    for i, (lon, lat) in enumerate(points):
        se, sn = FWD.transform(lon, lat)
        assert be[i] == se and bn[i] == sn, (
            f"batch and scalar disagree at index {i}: "
            f"{(be[i], bn[i])} vs {(se, sn)}"
        )


@given(st.lists(positions, min_size=2, max_size=200), st.integers(0, 2**32 - 1))
@settings(max_examples=100, deadline=None)
def test_order_independence(points: list[tuple[float, float]], seed: int) -> None:
    """Shuffling the input and un-shuffling the output changes nothing."""
    lons = np.array([p[0] for p in points], dtype=np.float64)
    lats = np.array([p[1] for p in points], dtype=np.float64)
    straight = np.column_stack(FWD.transform(lons, lats))

    rng = np.random.default_rng(seed)
    order = rng.permutation(len(points))
    shuffled = np.column_stack(FWD.transform(lons[order], lats[order]))
    restored = np.empty_like(shuffled)
    restored[order] = shuffled
    assert np.array_equal(restored, straight), "result depends on input order"

Parameter Reference #

Worked Example #

Running the suite against an implementation with a subtly wrong inverse convergence threshold:

Figure — the shrunk counterexample is a diagnosis, not just a failure.

Falsifying example: test_round_trip_closes(
    pos=(-124.99999999999999, 48.99999999999999),
)
AssertionError: latitude closed to 3.412e-08 deg

The shrinker has done the useful work: it reduced the failure to the corner of the domain, which immediately says the problem is at the edge of the area of use rather than everywhere. A fixture set would have had to contain that exact corner to find it.

Validation Check #

def test_out_of_area_is_refused() -> None:
    """The complement of the property: outside the area of use, expect a refusal."""
    far_away = (10.0, 50.0)          # Europe, for a North American operation
    e, n = FWD.transform(*far_away, errcheck=False)
    assert not (np.isfinite(e) and np.isfinite(n)) or abs(e) > 1e7, (
        "an operation applied far outside its area of use returned a plausible "
        "coordinate; the pipeline must range-check before transforming"
    )

Common Mistakes #

Generating over the whole coordinate system instead of the area of use. A property test that produces latitudes in Antarctica for a North American operation will fail, correctly, and the failure says nothing about the code. Bound the generators by the operation’s area of use, and test the outside separately as a refusal.

A round-trip tolerance tight enough to fail on noise. Requiring bitwise closure asserts something that is not true: the inverse of a projection is iterative, and its result differs from the input in the last bits. Set the tolerance from the ground precision you care about — 1e-9 degrees is a tenth of a millimetre — rather than from the floating-point representation.

Ignoring the shrunk example. Hypothesis reports the smallest input it could find that still fails, and that input is a diagnosis: a corner means a domain-edge problem, a zero means an initialisation problem, and a value adjacent to a grid boundary means an indexing problem. Reading it saves the debugging session.

Choosing Generators That Test Something #

The quality of a property test is almost entirely in its generators, and three choices separate a test that finds real defects from one that exercises the middle of the domain forever.

Figure — where a uniform generator spends its budget, and where the defects are.

Bound by the area of use, then sample its edges deliberately. A uniform generator over a rectangle spends most of its budget in the interior, where nothing goes wrong. Mixing a uniform component with an explicit edge component — values at and just inside the boundary — puts the budget where the defects are.

Include the degenerate values that are still valid. Zero latitude, zero longitude, exactly the central meridian, exactly a grid node, exactly a cell boundary. Each of these is a legitimate coordinate and each is where an off-by-one or a division lives.

Generate the shape, not just the values. For batch properties, the interesting variable is the array length: one element, two, a length that is not a multiple of any internal chunk size, and a length large enough to cross whatever threshold the implementation uses to switch strategies. Those thresholds are exactly where batch and scalar paths diverge.

A generator that produces a thousand comfortable coordinates proves that the comfortable case works, which was never in doubt. The value of the technique is entirely in the awkward inputs nobody would have written down, and the generators are what decide whether it ever produces one.

Frequently Asked Questions #

Do property tests replace golden vectors? #

No — they test different things. A property test asserts internal consistency: the code agrees with itself under round trips, reordering and batching. A golden vector asserts external correctness: the code agrees with an independent source about a specific number. An implementation can be perfectly self-consistent and consistently wrong, which is exactly what a changed operation produces.

How many examples should a property test run? #

Enough that a run takes a second or two, which for a fast transformation is several hundred. The value is in the shrinking rather than the volume: a property that fails will usually fail within the first few dozen examples, and the remaining budget is spent confirming the ones that pass.

Can I use property tests on the grid readers? #

Yes, and it is one of the better applications. Generate positions within a synthesised grid’s extent and assert node-exactness, cell-boundary agreement and monotonicity where the surface is monotone; generate positions outside it and assert refusal. A synthetic grid makes all of that reproducible without a large data dependency.

What about float edge cases like NaN and infinity? #

Exclude them from the generators and test them explicitly, because their correct handling is a policy decision rather than a numerical property. A pipeline should mask non-finite inputs before transforming and restore them afterwards; a property test that generates NaN will simply rediscover that the library propagates it, which is not the question.

Related #

• Testing and CI for coordinate transformation code — the parent topic
• Writing golden vector tests for coordinate transforms
• Concurrent pyproj transformation pipelines
• Map projection forward and inverse in Python
• Benchmarking NumPy vectorized vs looped transforms