Test Fixtures for GIS Code: Building Spatial Test Data You Can Trust

Problem statement

You want to test a cleaning function, so you point it at real data:

def test_removes_slivers():
    gdf = gpd.read_file("/data/parcels_2026.gpkg")
    result = remove_slivers(gdf, min_area=10)
    assert len(result) == 8_412_991

Six months later the test fails. Not because the code broke β€” because the supplier reissued the file with 4,000 more parcels. The assertion encoded a fact about a dataset, not a fact about the function.

The other failure is slower to notice. The test suite takes eleven minutes because every test reads a 6 GB GeoPackage. It cannot run in CI because the file is not in the repository. It cannot run on a colleague's laptop because the path is absolute. And when it does fail, the message is assert 8412987 == 8412991, which tells you nothing about which parcel was wrong or why.

A fixture is test data you control. Getting them right is most of what makes a GIS test suite useful.

Quick answer

Build geometry in code, at a scale you can reason about:

import pytest
import geopandas as gpd
from shapely.geometry import Polygon, Point

@pytest.fixture
def parcels():
    """Three parcels: one normal, one a sliver, one invalid."""
    return gpd.GeoDataFrame(
        {"id": [1, 2, 3], "class": ["residential", "sliver", "broken"]},
        geometry=[
            Polygon([(0, 0), (100, 0), (100, 100), (0, 100)]),      # 10,000 mΒ²
            Polygon([(100, 0), (100.5, 0), (100.5, 100), (100, 100)]),  # 50 mΒ²
            Polygon([(0, 0), (10, 10), (10, 0), (0, 10)]),          # bow-tie
        ],
        crs=27700,
    )

def test_removes_slivers(parcels):
    result = remove_slivers(parcels, min_area=100)
    assert set(result["id"]) == {1, 3}
    assert 2 not in set(result["id"])
Checklist of the five properties a good spatial test fixture has.
Five properties. Real data satisfies none of them.
Property Why it matters
built in code no external file, no path, works in CI
small you can hold the whole thing in your head
deliberate every feature is there to exercise something
deterministic the same input every run, forever
named the assertion says why, not just what

Step-by-step solution

1. Make every feature exist for a reason

A fixture is not a sample. Each feature should test a specific behaviour, and the fixture's docstring should say so:

@pytest.fixture
def geometry_cases():
    """One feature per condition the cleaner must handle."""
    from shapely.geometry import Polygon, MultiPolygon
    from shapely import wkt

    return gpd.GeoDataFrame(
        {"case": ["valid", "bowtie", "empty", "null", "sliver",
                  "multipart", "duplicate", "duplicate"]},
        geometry=[
            Polygon([(0, 0), (10, 0), (10, 10), (0, 10)]),          # ordinary
            Polygon([(0, 0), (10, 10), (10, 0), (0, 10)]),          # self-intersecting
            wkt.loads("POLYGON EMPTY"),                              # empty, not null
            None,                                                    # null, not empty
            Polygon([(0, 0), (0.01, 0), (0.01, 10), (0, 10)]),      # 0.1 mΒ²
            MultiPolygon([Polygon([(0, 0), (1, 0), (1, 1), (0, 1)]),
                          Polygon([(5, 5), (6, 5), (6, 6), (5, 6)])]),
            Polygon([(20, 20), (30, 20), (30, 30), (20, 30)]),      # same as next
            Polygon([(20, 20), (30, 20), (30, 30), (20, 30)]),
        ],
        crs=27700,
    )

Eight features, eight distinct conditions. Compare with a 10,000-row sample of real data, in which you do not know whether any of these conditions is present β€” so a passing test proves nothing about them.

The empty-versus-null distinction matters and is easy to conflate; see null, empty, missing and invalid.

2. Use coordinates you can do arithmetic on

Panels contrasting real-world coordinates with round test coordinates and the assertions each allows.
An assertion you can verify by hand is an assertion you can trust.
# ❌ real coordinates β€” what should the area be? nobody knows
Polygon([(351204.117, 381009.882), (351298.443, 381012.004), ...])

# βœ… round numbers β€” the area is 10,000 mΒ², obviously
Polygon([(0, 0), (100, 0), (100, 100), (0, 100)])

This turns opaque assertions into transparent ones:

def test_buffer_area(parcels):
    buffered = parcels.geometry.iloc[0].buffer(10)
    # 100Γ—100 square buffered by 10: 120Γ—120 minus four corners plus a circle
    expected = 120 * 120 - 4 * (10 * 10) + 3.14159 * 100
    assert buffered.area == pytest.approx(expected, rel=0.01)

A reviewer can check that arithmetic. They cannot check assert area == pytest.approx(9412.88).

Use a projected CRS in fixtures unless the test is about geographic coordinates. EPSG:27700 with coordinates near the origin is not geographically valid, but it makes areas come out in square metres and the test does not care where on Earth the shapes are.

3. Pick the smallest fixture that exercises the behaviour

# testing a spatial join needs exactly this much
@pytest.fixture
def join_inputs():
    left = gpd.GeoDataFrame(
        {"id": ["inside", "straddling", "outside"]},
        geometry=[Point(5, 5), Point(10, 5), Point(25, 5)],
        crs=27700,
    )
    right = gpd.GeoDataFrame(
        {"zone": ["west", "east"]},
        geometry=[Polygon([(0, 0), (10, 0), (10, 10), (0, 10)]),
                  Polygon([(10, 0), (20, 0), (20, 10), (10, 10)])],
        crs=27700,
    )
    return left, right

def test_boundary_point_matches_both_zones(join_inputs):
    left, right = join_inputs
    joined = gpd.sjoin(left, right, predicate="intersects")
    straddling = joined[joined["id"] == "straddling"]
    assert len(straddling) == 2, "a point on a shared edge intersects both zones"
    assert set(straddling["zone"]) == {"west", "east"}

Three points and two squares. The test documents a real, surprising behaviour β€” a point on a shared boundary joins twice, which is the root of duplicate rows from a spatial join β€” and it does so in a form anyone can verify by looking at the numbers.

4. Make fixtures composable

pytest fixtures can depend on other fixtures, which lets you build a small vocabulary rather than one big blob:

@pytest.fixture
def crs():
    return 27700

@pytest.fixture
def unit_square(crs):
    return gpd.GeoDataFrame({"id": [1]},
                            geometry=[Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])],
                            crs=crs)

@pytest.fixture
def grid_10x10(crs):
    """A clean 10Γ—10 coverage: no gaps, no overlaps, 100 cells of 1 mΒ²."""
    from shapely.geometry import box
    cells = [box(x, y, x + 1, y + 1) for y in range(10) for x in range(10)]
    return gpd.GeoDataFrame({"id": range(len(cells))}, geometry=cells, crs=crs)

@pytest.fixture
def grid_with_gap(grid_10x10):
    """The same coverage with cell 44 removed β€” one 1 mΒ² hole."""
    return grid_10x10[grid_10x10["id"] != 44].reset_index(drop=True)

@pytest.fixture
def grid_with_overlap(grid_10x10):
    """The same coverage with one cell enlarged so it overlaps its neighbour."""
    from shapely.geometry import box
    gdf = grid_10x10.copy()
    gdf.loc[gdf["id"] == 44, "geometry"] = box(4, 4, 5.5, 5)
    return gdf

Now a coverage test reads like a specification:

def test_clean_coverage_has_no_gaps(grid_10x10):
    assert find_gaps(grid_10x10).empty

def test_gap_is_found_and_measured(grid_with_gap):
    gaps = find_gaps(grid_with_gap)
    assert len(gaps) == 1
    assert gaps.geometry.area.sum() == pytest.approx(1.0)

def test_overlap_is_found_and_measured(grid_with_overlap):
    overlaps = find_overlaps(grid_with_overlap)
    assert len(overlaps) == 1
    assert overlaps.geometry.area.sum() == pytest.approx(0.5)

1.0 and 0.5 are checkable by hand, which is what makes these assertions worth trusting. The behaviour is covered in how to fix gaps and overlaps in a polygon coverage.

5. Reserve real data for a small number of slow tests

Fixtures cannot cover everything. Real data has encoding problems, mixed geometry types, coordinates near the limits of floating point, and combinations nobody would think to construct.

import pytest
from pathlib import Path

SAMPLES = Path(__file__).parent / "data"

@pytest.fixture(scope="session")
def real_sample():
    """A small committed extract of real data. Regenerate with tools/make_sample.py."""
    path = SAMPLES / "parcels_sample.gpkg"
    if not path.exists():
        pytest.skip(f"{path} not present β€” run tools/make_sample.py")
    return gpd.read_file(path)

@pytest.mark.slow
def test_pipeline_on_real_extract(real_sample):
    result = clean(real_sample)
    assert result.geometry.is_valid.all()
    assert result.crs == real_sample.crs
    assert len(result) <= len(real_sample)

Three things make this workable. The sample is committed and small β€” a few hundred features, under a megabyte β€” so CI can run it. scope="session" reads it once for the whole run rather than per test. And the assertions are properties rather than exact counts: valid geometry, preserved CRS, no rows invented. Those stay true when the sample is regenerated; assert len(result) == 412 does not.

Mark them and let developers skip them:

# pytest.ini
[pytest]
markers =
    slow: reads real data or takes over a second
pytest -m "not slow"          # the fast loop while developing
pytest                        # everything, in CI

Code examples

Example 1: a fixture library for geometry conditions

# tests/conftest.py
import pytest
import geopandas as gpd
import numpy as np
from shapely import wkt
from shapely.geometry import Point, LineString, Polygon, MultiPolygon, box

CRS = 27700

def gdf(geoms, **cols):
    n = len(geoms)
    data = {k: (v if isinstance(v, list) else [v] * n) for k, v in cols.items()}
    return gpd.GeoDataFrame(data, geometry=geoms, crs=CRS)

@pytest.fixture
def valid_polygons():
    """Three well-formed squares of 100, 400 and 900 mΒ²."""
    return gdf([box(0, 0, 10, 10), box(20, 0, 40, 20), box(50, 0, 80, 30)],
               id=[1, 2, 3], area_m2=[100, 400, 900])

@pytest.fixture
def invalid_polygons():
    """One of each way a polygon can be invalid under the OGC rules."""
    return gdf(
        [
            Polygon([(0, 0), (10, 10), (10, 0), (0, 10)]),            # self-intersection
            Polygon([(0, 0), (10, 0), (10, 10), (0, 10)],
                    [[(2, 2), (2, 8), (8, 8), (8, 2)],
                     [(4, 4), (4, 6), (6, 6), (6, 4)]]),              # nested holes
            Polygon([(0, 0), (10, 0), (0, 0)]),                       # zero area
        ],
        reason=["self-intersection", "nested holes", "degenerate"],
    )

@pytest.fixture
def missing_geometries():
    """Null and empty are different states and must be tested separately."""
    return gdf([box(0, 0, 1, 1), None, wkt.loads("POLYGON EMPTY")],
               kind=["present", "null", "empty"])

@pytest.fixture
def duplicate_geometries():
    """Exact duplicates, and a near-duplicate 1 mm away."""
    return gdf([box(0, 0, 10, 10), box(0, 0, 10, 10), box(0, 0, 10.001, 10)],
               kind=["original", "exact_duplicate", "near_duplicate"])

@pytest.fixture
def mixed_crs():
    """The same square in two CRS β€” for testing that code reprojects."""
    a = gdf([box(351000, 381000, 351100, 381100)], id=[1])
    b = a.to_crs(4326)
    return a, b

@pytest.fixture
def line_network():
    """A network with one dangle, one overshoot and one clean junction."""
    return gdf(
        [
            LineString([(0, 0), (10, 0)]),          # main
            LineString([(10, 0), (20, 0)]),         # continues β€” clean junction
            LineString([(10, 0), (10, 5)]),         # branch β€” clean junction
            LineString([(30, 0), (35, 0)]),         # disconnected dangle
            LineString([(20, 0), (20, -1)]),        # overshoot past the junction
        ],
        kind=["main", "continuation", "branch", "dangle", "overshoot"],
    )

The gdf helper removes the boilerplate that otherwise makes a fixture module tedious to read, and broadcasting a scalar column across every row keeps the definitions to one line each.

Every fixture's docstring states what it is for, so a failing test in six months explains itself. duplicate_geometries having a near-duplicate 1 mm away is the kind of detail that makes a fixture catch a real bug β€” exact-match deduplication passes on it, and it should not.

Example 2: parametrised fixtures for whole classes of input

import pytest
import geopandas as gpd
from shapely.geometry import box

@pytest.fixture(params=[27700, 4326, 3857, 32630])
def any_crs(request):
    """Run the test once per CRS β€” catches code that assumes metres."""
    return request.param

def test_cleaning_preserves_crs(any_crs):
    gdf = gpd.GeoDataFrame({"id": [1]}, geometry=[box(0, 0, 1, 1)], crs=any_crs)
    assert clean(gdf).crs == gdf.crs

@pytest.fixture(params=["Polygon", "MultiPolygon", "LineString", "Point"])
def any_geom_type(request):
    from shapely.geometry import (Point, LineString, Polygon, MultiPolygon)
    shapes = {
        "Polygon": Polygon([(0, 0), (1, 0), (1, 1), (0, 1)]),
        "MultiPolygon": MultiPolygon([Polygon([(0, 0), (1, 0), (1, 1), (0, 1)]),
                                      Polygon([(2, 2), (3, 2), (3, 3), (2, 3)])]),
        "LineString": LineString([(0, 0), (1, 1)]),
        "Point": Point(0.5, 0.5),
    }
    return gpd.GeoDataFrame({"id": [1]}, geometry=[shapes[request.param]], crs=27700)

def test_handles_every_geometry_type(any_geom_type):
    result = clean(any_geom_type)
    assert len(result) == 1
    assert result.geometry.iloc[0].geom_type == any_geom_type.geometry.iloc[0].geom_type
tests/test_clean.py::test_cleaning_preserves_crs[27700] PASSED
tests/test_clean.py::test_cleaning_preserves_crs[4326] PASSED
tests/test_clean.py::test_cleaning_preserves_crs[3857] PASSED
tests/test_clean.py::test_cleaning_preserves_crs[32630] PASSED
tests/test_clean.py::test_handles_every_geometry_type[Polygon] PASSED
tests/test_clean.py::test_handles_every_geometry_type[MultiPolygon] PASSED
tests/test_clean.py::test_handles_every_geometry_type[LineString] FAILED

Sixteen tests from two functions, and one real bug: the cleaner works on polygons and breaks on lines. The params id appears in the test name, so the failure names the case without any extra reporting.

Parametrising over CRS is particularly valuable in GIS code, because the commonest hidden assumption is that coordinates are in metres. A function computing area > 100 is correct in EPSG:27700 and meaningless in EPSG:4326, and only a test that runs under both will say so.

Example 3: generated fixtures with a fixed seed

For property-based testing, generate data β€” but pin the seed so failures reproduce:

import numpy as np
import geopandas as gpd
from shapely.geometry import Point, box
import pytest

def random_points(n, seed, bounds=(0, 0, 1000, 1000), crs=27700):
    """Deterministic random points β€” the same seed gives the same points forever."""
    rng = np.random.default_rng(seed)
    minx, miny, maxx, maxy = bounds
    xs = rng.uniform(minx, maxx, n)
    ys = rng.uniform(miny, maxy, n)
    return gpd.GeoDataFrame(
        {"id": np.arange(n), "value": rng.normal(100, 25, n)},
        geometry=[Point(x, y) for x, y in zip(xs, ys)], crs=crs)

@pytest.fixture(scope="session")
def points_1000():
    return random_points(1000, seed=20260821)

def test_count_in_polygons_totals_correctly(points_1000):
    grid = gpd.GeoDataFrame(
        {"cell": range(100)},
        geometry=[box(x * 100, y * 100, (x + 1) * 100, (y + 1) * 100)
                  for y in range(10) for x in range(10)],
        crs=27700)
    counted = count_points_in_polygons(points_1000, grid)
    # the grid covers the full generation bounds, so every point lands somewhere
    assert counted["count"].sum() == len(points_1000)
    assert (counted["count"] >= 0).all()

np.random.default_rng(seed) rather than np.random.seed() is the important detail: the modern generator is isolated, so another library calling into NumPy's global random state cannot change your fixture. A test whose data depends on global state is a flaky test waiting to happen.

The assertion is a conservation property β€” every point is counted exactly once β€” which holds for any seed and any count. That is the right shape for a generated fixture: assert invariants, not specific values, because the specific values are an artefact of the seed.

Explanation

Grid comparing constructed fixtures with real data on coverage, speed, portability and realism.
The two are complementary, not alternatives. Many of the first, a few of the second.

A fixture is a controlled experiment, and the value of an experiment comes from controlling everything except the thing under test.

Real data controls nothing. Its contents change when the supplier reissues it, its size makes tests slow, its location makes them unportable, and β€” most damagingly β€” you do not know what conditions it contains. A test that passes on 10,000 real parcels tells you the function did not crash on those parcels. It does not tell you the function handles null geometry, because you do not know whether any was present.

Constructed fixtures invert this. Eight carefully chosen features cover eight conditions with certainty, and a failure names the condition. The trade is that constructed data is unrealistic β€” it lacks the encoding oddities, coordinate extremes and structural surprises that real data supplies for free. Hence the two-tier approach: many fast constructed fixtures for behaviour, a few slow real-data tests for the unknown unknowns, marked so the fast loop stays fast.

Round coordinates matter more than they look. An assertion's value depends on a reader being able to check it. assert result.area == pytest.approx(10_000) on a 100 Γ— 100 square is self-evidently right; assert result.area == pytest.approx(9412.883) is a number someone once copied from a test run, and it will be copied again the next time the test fails β€” which converts the test from a specification into a record of current behaviour. That is how a suite full of green tests ends up asserting a bug.

Determinism is not optional and is easy to lose. np.random.seed() sets a global that any imported library can also set. Unseeded generation makes failures unreproducible. Even something as small as iterating a set of geometries can vary between runs. A test that fails one time in twenty is worse than no test, because it trains people to re-run rather than investigate.

Parametrised fixtures are unusually valuable in GIS, because the discipline's most common latent bug is a unit assumption. Code that computes area > 100 or distance < 500 is correct in a projected CRS and meaningless in a geographic one, and nothing in the type system distinguishes them. Running every test across several CRS makes those assumptions fail loudly instead of surfacing months later as an inexplicable result β€” the class of problem described in the geographic CRS warning.

Finally, fixtures are documentation that cannot go stale. A fixture called grid_with_gap whose docstring says "cell 44 removed β€” one 1 mΒ² hole" tells a new reader what a gap is, what the code should do about it, and how big the answer should be. Prose documentation drifts from the code; a fixture is executed on every run.

Edge cases or notes

  • Null and empty geometry are different states. None is null; wkt.loads("POLYGON EMPTY") is empty. Test both.
  • scope="session" builds a fixture once per run. Only for immutable data β€” a test that mutates a session fixture corrupts every later test.
  • Never mutate a fixture in place. Return a copy, or use function scope so each test gets a fresh object.
  • pytest.approx is required for coordinate arithmetic. Floating point makes exact equality unreliable β€” see coordinate precision explained.
  • np.random.default_rng(seed) is isolated; np.random.seed() is global and can be changed by any library.
  • Fixtures near the origin in EPSG:27700 are not geographically valid. That is fine unless the test is about location.
  • A committed sample must be small. Under a megabyte, a few hundred features, and regenerable by a script in the repo.
  • conftest.py shares fixtures across a directory without imports; nested directories can override.
  • pytest --fixtures lists everything available, with docstrings β€” which is why the docstrings matter.
  • Assert properties, not counts, on real data. Counts encode a dataset version; properties encode a specification.

FAQ

Should I test against real data?

A little, in a few tests marked slow. Real data supplies surprises no constructed fixture will, but it changes, it is large, and you cannot know which conditions it contains. Most tests should use fixtures you built.

Why do my tests break when the data is updated?

They assert facts about a dataset β€” a row count, a specific area β€” rather than about the function. Assert properties that survive a data refresh: valid geometry, preserved CRS, no rows invented.

How big should a fixture be?

Small enough to reason about entirely. Three to ten features is usually right; a hundred is already too many to hold in your head when the test fails.

Why use round coordinates?

So the expected values are checkable by hand. area == 10_000 for a 100 Γ— 100 square is verifiable; a number copied from a previous test run is not, and it will be re-copied the next time it fails.

How do I make random test data reproducible?

np.random.default_rng(seed) with a fixed seed. Avoid np.random.seed(), which sets a global that any library can change, and assert invariants rather than seed-specific values.

Should fixtures be in conftest.py?

Yes for anything shared across several test files β€” pytest finds them without imports. Keep single-use fixtures next to the test that uses them.

How do I test that my code works in any CRS?

Parametrise a fixture over several CRS. It is the most effective way to catch code that assumes coordinates are in metres, which is the commonest hidden assumption in GIS code.