Rasterio Returned Zero. Was the Point Inside the Raster?

A raster sampler returns zero for two requests. One coordinate selects a valid pixel whose stored value is genuinely zero. The other coordinate falls outside the grid, where the sampling implementation supplies a fallback zero. A consumer that sees only the scalar cannot distinguish an observation from that out-of-range representation.

The repair is not to reject every zero. That would destroy the legitimate observation while leaving the underlying admission problem unresolved. Nor does a more distinctive nodata sentinel prove that a coordinate addressed a pixel. The useful question is whether each original request passed four separate checks: coordinate validity, established coordinate reference system identity, membership in the finite pixel grid, and selected-band validity.

This article qualifies that decision using three owned, in-memory GeoTIFFs and eight identified points. The raster values are deliberately small enough to inspect without a satellite product, download, or scientific interpretation. Two different point IDs select the same valid zero pixel; three other IDs lie outside. A third dataset additionally invalidates one in-bounds, nonzero pixel through its mask.

The complete Python fixture below was executed under CPython 3.13.5, Rasterio 1.5.0, NumPy 2.3.5, and GDAL 3.12.1. Its results are synthetic engine evidence, not field measurements. The operational objective is to preserve valid zeros, explain every excluded or held request, and prevent a returned number from becoming stronger evidence than the lookup actually produced.

Define a valid observation before inspecting its value

For this workflow, an admissible observation requires a correctly formed, finite coordinate pair; trusted evidence that those coordinates use the raster's CRS; row and column indices within the declared raster dimensions; and a valid mask for band 1. All four conditions must hold. The sample value itself is not an admission predicate.

The Rasterio sampling API defines input pairs in the dataset's reference system and provides a masked sampling option. It does not establish the provenance of an arbitrary caller's CRS declaration. In this fixture, that provenance comes from the authored manifest: every main coordinate is explicitly defined in EPSG:3857. An equivalent-looking label supplied without supporting provenance would not establish the same fact in a production pipeline.

Keep coordinate meaning separate from coordinate formatting. The adjacent discussion of coordinate labels versus reprojection concerns that upstream distinction. Here, no coordinates are transformed or relabeled. An unknown or mismatched declared CRS leads to HOLD_CRS before a pixel index is calculated.

“Valid” also has a deliberately limited meaning. IN_BOUNDS_VALID says that the coordinate selects a pixel admitted by this source's selected-band mask. It does not certify calibration, physical interpretation, geolocation accuracy, or suitability for a scientific study. Those are different acceptance decisions with different evidence.

Define the downstream interface accordingly. A record should carry its request identity, spatial evidence, mask flag, and classification, rather than returning only an integer. A zero may then proceed through the same valid-data path as 11 or 33. Outside and masked records remain visible in the ledger, but they are not ordinary observations. This separation makes rejection explainable without quietly changing the requested population.

Pin the engine and construct the two-by-two raster

The executed environment is CPython 3.13.5, Rasterio 1.5.0, NumPy 2.3.5, and GDAL runtime 3.12.1. These are test coordinates, not a latest-version recommendation. The fixture captures the installed versions, the Python implementation, and the platform; it stops on a mismatch instead of printing an inherited pass.

Documentation identity is separate. The stable sampling, transform, and mask pages consulted on 16 October 2026 identify themselves as Rasterio 1.4.4 documentation. The Rasterio 1.5.0 release record identifies another implementation version and lists GDAL 3.12.1 in its wheel packaging notes. That packaging statement does not replace inspection of the actual runtime. A reference-page access date is not a feature launch date.

Each dataset has width 2, height 2, one int16 band, CRS EPSG:3857, and transform from_origin(100, 200, 10, 10). The stored array is [[0, 11], [22, 33]]. Only validity configuration changes:

The constructors create three separate MemoryFile objects. Each writer closes before its MemoryFile reopens for reading. Rasterio's in-memory file documentation describes this GDAL-backed storage abstraction. The executable fixture explicitly keeps the mask internal to the owned GeoTIFF; it does not copy or test transfer assets.

Before sampling, the reader's actual dimensions, dtype, CRS, affine coefficients, stored values, nodata setting, and band mask are recorded and checked. That prevents a mistaken constructor from becoming an apparently successful sampling test. The affine coefficients are recorded in a,b,c,d,e,f order, not silently presented as a GDAL geotransform tuple. Dataset IDs identify the three logical fixtures, while reader names provide supplementary diagnostics, not durable source identity.

Publish the eight-point pixel-index oracle

Write the expected indices before consulting the sampling engine. For this particular north-up transform, the independent formulas are row = floor((200 - y) / 10) and col = floor((x - 100) / 10). Grid membership requires 0 <= row < 2 and 0 <= col < 2.

The GDAL Geotransform Tutorial explains the image-coordinate frame: the upper-left corner is at pixel/line coordinate (0,0), and the upper-left pixel center is at (0.5,0.5). With this fixture's origin and pixel dimensions, the formulas above follow directly. They are authored from the declared geometry, not copied from Rasterio's returned indices.

Both rowcol() and the dataset's index() must match these literal pairs. Agreement between two helpers from the same library would be useful consistency evidence, but it would not be an independent oracle. The Rasterio transform API documents floor as the default conversion from fractional pixel coordinates to whole indices.

P5 is a useful arithmetic check: its fractional column is -0.5, whose floor is -1. Truncating toward zero would incorrectly put it in column 0. Moreover, obtaining an integer index does not imply that the integer is allowed. The range predicate must run before an array lookup can treat that index as an address.

The oracle also supplies a useful separation of responsibilities during review. A wrong index pair points toward coordinate interpretation, transform use, or index conversion. A correct pair with the wrong inside flag points toward the range predicate. Correct geometry followed by the wrong validity classification points toward mask handling or precedence. Those are different defects, even when they eventually produce the same incorrect accepted value. Keeping the intermediate pair in the report makes that distinction inspectable instead of forcing a reviewer to infer it from the final scalar.

State exactly which boundaries the fixture includes

P8 sits precisely at (100,200). Both independent expressions evaluate to zero, so it selects row 0, column 0. P6 produces column 2, which is outside a width-2 array. P7 produces row 2, which is outside a height-2 array. No tolerance, inward movement, or nearest-edge substitution is applied.

These outcomes follow from this affine and its floor-index convention. They are not a general closed geographic bounding-box rule. Nor should the adjacent workflow on vector boundary assignment with Shapely supply a replacement predicate: this experiment is about membership in a finite pixel grid, not ownership by vector polygons. Other raster geometries require their own qualification rather than reuse of these simplified formulas.

Observe the zero collision with nodata unset

For M0, the observed sample(..., indexes=1, masked=False) vector, in P1 through P8 order, was [0, 11, 22, 33, 0, 0, 0, 0]. The same integer appears for two valid requests and three requests that address no pixel.

Inspection of the installed sample.py in the executed Rasterio 1.5.0 environment found the fallback construction based on dataset.nodata or 0. The implementation tests row and column bounds before reading a one-pixel window; an outside index receives the fallback instead. The versioned Rasterio 1.5.0 source is the comparison point for this implementation detail, while the observations here come from the recorded installed runtime.

A consumer that accepts every returned integer therefore admits P5, P6, and P7 incorrectly. A consumer that discards every zero makes the opposite error: it removes P1 and P8 even though both selected a valid source pixel. Neither shortcut can be repaired by rearranging the scalar comparisons, because the distinguishing information is absent from the scalar alone.

The absence of a nodata declaration is similarly insufficient. M0 contains four valid stored pixels, but its eight requested coordinates do not all fall within those pixels. Source validity and request membership are different populations. Validating the former does not automatically validate the latter.

Keep input identity when points select the same cell

P1 and P8 have different coordinates and different IDs. Both map to (0,0), and both legitimately return zero. They must survive as two records because the unit being reconciled is an original point request, not a distinct raster cell.

Deduplicating by value would also merge unrelated nonzero observations in a less distinctive array. Deduplicating by row and column would change the request population even when the raster lookup was correct. A later analysis may deliberately aggregate by pixel, but that needs a separate, explicit policy. The admission stage should preserve the many-to-one relationship between requests and cells, not erase it before the analyst decides how those requests should contribute.

Compare an explicit nodata sentinel

M1 changes only the nodata metadata to -9999. None of its four stored values equals that sentinel. Its observed raw vector was [0, 11, 22, 33, -9999, -9999, -9999, 0].

The three outside positions now look conspicuously different from the in-bounds values. That is useful diagnostic information, but the geometry has not changed. P5 is still column -1, P6 is still column 2, and P7 is still row 2. The valid request set remains P1, P2, P3, P4, and P8.

A hard-coded value != -9999 filter happens to select that set for this particular M1 array. It is not the acceptance contract demonstrated across all three datasets. Apply the same filter to M0 and it admits every outside zero. M2 will additionally expose an in-bounds invalid value that is neither zero nor -9999. The fixture therefore tests a consumer's reasoning across representations, not merely its ability to recognize one magic number.

The sentinel should be retained as source metadata, alongside the raw result, rather than promoted into the only validity rule. Reviewers need to know whether a number came from an addressed pixel or from an outside representation. The row and column evidence answers that question without depending on whether a fallback is visually surprising.

Do not change the raster's nodata declaration simply to make a downstream table look cleaner. Such a change would alter the source contract and would require its own authorization. Likewise, replacing -9999 with zero after sampling does not make an outside request an observation; it merely recreates M0's ambiguity in another representation.

This test does not introduce an in-bounds pixel containing the declared sentinel. Its conclusion is narrower: a distinctive outside fill makes inspection easier, but an explicit membership check remains necessary. The valid zero controls are retained unchanged in both datasets.

Inspect masked sampling without flattening its evidence

For both M0 and M1, the observed masked=True flags in input order were [False, False, False, False, True, True, True, False]. Only P5, P6, and P7 were masked. The two valid zeros, P1 and P8, remained unmasked.

This result gives the consumer information that the M0 raw vector lacks. Nevertheless, the experiment still records indices separately. A mask flag says that this returned element is masked; it does not, by itself, explain whether the coordinate was outside or whether an addressed pixel was invalid. That distinction becomes visible in M2.

The Rasterio mask documentation distinguishes two opposite conventions. A GDAL-style valid-data mask uses nonzero entries, typically 255, for valid pixels. A NumPy masked-array flag uses True for invalid elements. Reading those conventions backward would reverse admission even if every coordinate and raw sample were correct.

The fixture therefore extracts each mask using np.ma.getmaskarray(sample). The NumPy 2.3 getmaskarray reference specifies a full boolean array, including a full False array when no mask is present. For a one-band sample, indexing that returned array at position 0 is explicit and shape-consistent; assuming that sample.mask always has an indexable array representation is unnecessary.

Order is part of this comparison. The unmasked and masked calls consume the same eight coordinates in the same sequence, and each returned element is reattached to the original ID. Comparing two correctly shaped vectors is insufficient if one has been reordered. The fixture checks lengths and uses strict pairing when assembling records, so a missing element cannot quietly shorten the ledger. Its final key comparison then checks the dataset/point sequence as well as uniqueness. A valid mask attached to the wrong request would still be invalid evidence for that request.

Preserve the mask separately from the scalar

The ledger's raw_value comes from the unmasked sampling call. Its mask_flag comes from the separate masked call. The code does not cast a masked scalar to an integer or fill a masked result with zero and then call the replacement source data.

Those fields serve different purposes. Raw values help explain what an unsafe consumer might have ingested. Mask flags support admission. Retaining a raw value beside a True mask is acceptable diagnostic bookkeeping only while the record's invalid classification remains attached and downstream arithmetic respects it.

Review serialization as part of the consumer repair. A table export that preserves raw_value but discards mask_flag and classification recreates the original ambiguity. The same is true of a helper that returns only the underlying array after a mask-aware call. The stronger interface is the complete classified record, with any accepted numerical view derived from it rather than used as a substitute for it.

Distinguish an inside masked pixel from outside extent

M2 has nodata=None and the explicit valid-data mask [[255, 255], [0, 255]]. The zero in this mask invalidates row 1, column 0. It does not invalidate the stored zero at row 0, column 0; that pixel's mask entry remains 255.

The observed raw M2 vector was [0, 11, 22, 33, 0, 0, 0, 0], unchanged from M0. Its observed masked flags were [False, False, True, False, True, True, True, False]. P3 is the additional masked request, even though its raw value is the nonzero integer 22.

P3 therefore defeats another shortcut: “inside plus a plausible nonzero value means valid.” Its independently derived indices are (1,0), within the grid, but the selected band's mask rejects the addressed pixel. Geometry establishes that a pixel exists at the address; the mask establishes whether this source admits that pixel.

The fixture reads and checks the actual selected-band mask before interpreting the sample flags. For in-bounds requests, it also checks agreement between the flag and the corresponding band-mask entry. The adjacent workflow on GeoTIFF validity-mask preservation concerns preserving that source meaning across delivery. Here the source remains intact; only sampling admission is under test.

Apply classification precedence explicitly

Classify an out-of-range row or column as OUTSIDE_EXTENT first. Otherwise, classify an invalid selected-band mask as IN_BOUNDS_MASKED. Only the remaining case is IN_BOUNDS_VALID.

Under that order, P3 in M2 is inside but masked. P5, P6, and P7 are outside even though their masked sampling flags are also True. A single generic “missing” category would lose the distinction between an invalid existing pixel and the absence of any addressed cell.

The distinction changes the question sent to the responsible owner. P3 calls for an explanation of source validity, not a larger numeric threshold. P5 through P7 call for a coverage or coordinate-policy decision, not a mask edit. P1 and P8 need neither intervention: they are legitimate valid-zero inputs and stay admitted. The classifier preserves all three states without changing coordinates or guessing why the source owner marked a pixel invalid.

Reject malformed coordinates and unproved CRS before sampling

Four preflight controls sit outside the 24 main dataset/point combinations. Their observed outcomes were:

No control reaches rowcol, index, or sample. The same preflight function also guards every main request before those calls, rather than existing only as an isolated demonstration. A failed preflight in the supposedly fixed main manifest stops the harness and retains the offending record as diagnostics.

The input policy is intentionally explicit. Accept a tuple or list of length two whose members satisfy numbers.Real, excluding Python and NumPy booleans. Convert those already-numeric values to floats with conversion errors handled, then require both results to be finite. Python's math.isfinite documentation includes ordinary zero as finite and excludes NaN and infinities.

A numeric-looking string such as "105" is not admitted. Parsing text coordinates may be a legitimate upstream responsibility, but doing it implicitly here would mix schema repair with spatial qualification. The same applies to a one-element tuple, a complex value, or a boolean: a plausible conversion is not proof of the intended coordinate schema.

For CRS identity, missing declarations or missing raster CRS produce HOLD_CRS. Unparseable declarations and mismatches also remain held. Rasterio's CRS construction API supplies the parser used for comparison; it does not verify the historical origin of the coordinates. Trusted upstream provenance is still a precondition.

These statuses are policy results, not successful raster observations. INVALID_INPUT and HOLD_CRS are expected control outcomes, so they do not make a correct test run fail. Conversely, an unexpected exception, missing main record, or oracle mismatch is a harness failure, not a newly invented input category. Keep that distinction visible in reporting so a failed test cannot be renamed into a successful exclusion.

Assemble the complete point-level admission ledger

The observed compact ledger below contains one row for each of the 24 main combinations. The full report additionally retains input coordinates, actual transform, CRS, band, expected classification, source metadata, and coordinate-helper comparisons.

Reconcile exact point sets, not only summary counts

M0 and M1 each admit exactly {P1, P2, P3, P4, P8}; M2 admits exactly {P1, P2, P4, P8} and marks {P3} inside but masked. Every dataset has exactly {P5, P6, P7} outside. Counts of five/three or four/one/three are secondary checks, because swapped identities could preserve those counts.

A report must also distinguish completion from partial progress. The writer begins with NOT_COMPLETED and changes to PASS only after source checks, sampling, status reconciliation, and valid-zero assertions succeed. Ordinary exceptions retain their type, message, and traceback under HARNESS_FAILURE and are re-raised. Recorded vectors or a partially populated ledger remain diagnostic evidence, not a completed qualification. In particular, an empty accepted set in a failed run must not be interpreted as proof that every requested coordinate was outside.

Save the following complete fixture as rasterio_admission.py in the provisioned environment. Run it without optimization flags and use a new report filename:

python rasterio_admission.py sample-run-001.json

The report path is opened exclusively to avoid overwriting prior evidence. Diagnostics are serialized and flushed inside the owning ExitStack before its readers and MemoryFiles close, including ordinary exception paths. This is not a promise of recovery from abrupt process termination or storage failure. The printed PASS means that the declared fixture reconciled; it does not mean that all 24 samples were accepted.

from future import annotations

import argparse
import json
import math
import numbers
import os
import platform
import traceback
from contextlib import ExitStack
from datetime import datetime, timezone
from pathlib import Path
from typing import Any

import numpy as np
import rasterio
from rasterio.crs import CRS
from rasterio.errors import CRSError
from rasterio.io import MemoryFile
from rasterio.transform import from_origin, rowcol

FIXTURE = "rasterio-admission-r1"
WIDTH, HEIGHT, BAND = 2, 2, 1
CRS_LABEL = "EPSG:3857"
# Literal (ID, x, y, expected row, expected column), before engine calls.
POINTS = (
    ("P1", 105, 195, 0, 0),
    ("P2", 115, 195, 0, 1),
    ("P3", 105, 185, 1, 0),
    ("P4", 115, 185, 1, 1),
    ("P5", 95, 195, 0, -1),
    ("P6", 120, 195, 0, 2),
    ("P7", 105, 180, 2, 0),
    ("P8", 100, 200, 0, 0),
)
VALUES = np.array([[0, 11], [22, 33]], dtype="int16")
VALID_MASK = np.array([[255, 255], [0, 255]], dtype="uint8")
EXPECTED_VERSIONS = {
    "implementation": "CPython", "python": "3.13.5",
    "rasterio": "1.5.0", "numpy": "2.3.5", "gdal": "3.12.1",
}
EXPECTED_RAW = {
    "M0": [0, 11, 22, 33, 0, 0, 0, 0],
    "M1": [0, 11, 22, 33, -9999, -9999, -9999, 0],
    "M2": [0, 11, 22, 33, 0, 0, 0, 0],
}
EXPECTED_MASK = {
    "M0": [False, False, False, False, True, True, True, False],
    "M1": [False, False, False, False, True, True, True, False],
    "M2": [False, False, True, False, True, True, True, False],
}
VALID, MASKED, OUTSIDE = (
    "IN_BOUNDS_VALID", "IN_BOUNDS_MASKED", "OUTSIDE_EXTENT"
)
EXPECTED_SETS = {
    "M0": {VALID: {"P1", "P2", "P3", "P4", "P8"},
           MASKED: set(), OUTSIDE: {"P5", "P6", "P7"}},
    "M1": {VALID: {"P1", "P2", "P3", "P4", "P8"},
           MASKED: set(), OUTSIDE: {"P5", "P6", "P7"}},
    "M2": {VALID: {"P1", "P2", "P4", "P8"},
           MASKED: {"P3"}, OUTSIDE: {"P5", "P6", "P7"}},
}


def independent_index(x: float, y: float) -> tuple[int, int]:
    # Valid only for this fixture's declared north-up transform.
    return math.floor((200 - y) / 10), math.floor((x - 100) / 10)


def preflight(coords: Any, declared_crs: Any,
              dataset_crs: CRS | None) -> dict[str, Any]:
    base = {"original_coordinates": repr(coords),
            "declared_crs": str(declared_crs)}

    def result(status: str, reason: str) -> dict[str, Any]:
        return {**base, "status": status, "reason": reason}

    if not isinstance(coords, (tuple, list)) or len(coords) != 2:
        return result("INVALID_INPUT", "a tuple/list of arity 2 is required")
    if any(not isinstance(v, numbers.Real)
           or isinstance(v, (bool, np.bool_)) for v in coords):
        return result("INVALID_INPUT", "real numbers required; no strings/bools")
    try:
        x, y = (float(v) for v in coords)
    except (TypeError, ValueError, OverflowError):
        return result("INVALID_INPUT", "coordinates cannot become finite floats")
    if not (math.isfinite(x) and math.isfinite(y)):
        return result("INVALID_INPUT", "coordinates must be finite")
    if declared_crs is None or dataset_crs is None:
        return result("HOLD_CRS", "same-CRS provenance is unknown")
    try:
        point_crs = CRS.from_user_input(declared_crs)
    except (CRSError, TypeError, ValueError):
        return result("HOLD_CRS", "declared CRS cannot be interpreted")
    if point_crs != dataset_crs:
        return result("HOLD_CRS", "declared CRS differs from raster CRS")
    return {**result("READY", "finite pair; declared CRS matches"),
            "x": x, "y": y}


def classify(inside: bool, masked: bool) -> str:
    if not isinstance(inside, bool) or not isinstance(masked, bool):
        raise TypeError("classification needs established boolean evidence")
    if not inside:
        return OUTSIDE
    return MASKED if masked else VALID


def make_raster(stack: ExitStack, nodata: int | None,
                mask: np.ndarray | None) -> MemoryFile:
    mem = stack.enter_context(MemoryFile())
    with mem.open(driver="GTiff", width=WIDTH, height=HEIGHT,
                  count=1, dtype="int16", crs=CRS_LABEL,
                  transform=from_origin(100, 200, 10, 10),
                  nodata=nodata) as writer:
        writer.write(VALUES, BAND)
        if mask is not None:
            writer.write_mask(mask)
    return mem  # The writer is closed; the MemoryFile is still owned.


def make_m0(stack: ExitStack) -> MemoryFile:
    return make_raster(stack, None, None)


def make_m1(stack: ExitStack) -> MemoryFile:
    return make_raster(stack, -9999, None)


def make_m2(stack: ExitStack) -> MemoryFile:
    return make_raster(stack, None, VALID_MASK)


def check_source(dataset_id: str, src: Any,
                 report: dict[str, Any]) -> np.ndarray:
    transform = [float(v) for v in tuple(src.transform)[:6]]
    band_mask = src.read_masks(BAND)
    meta = {
        "fixture": FIXTURE, "dataset_id": dataset_id,
        "reader_name": src.name, "driver": src.driver,
        "crs": str(src.crs), "transform": transform,
        "transform_order": "a,b,c,d,e,f",
        "width": src.width, "height": src.height,
        "count": src.count, "dtypes": list(src.dtypes),
        "nodata": src.nodata, "values": src.read(BAND).tolist(),
        "valid_data_mask": band_mask.tolist(),
        "mask_flags": [[flag.name for flag in group]
                       for group in src.mask_flag_enums],
    }
    report["datasets"][dataset_id] = meta  # Retain before assertions.
    assert src.driver == "GTiff"
    assert (src.width, src.height, src.count) == (2, 2, 1)
    assert src.dtypes == ("int16",)
    assert src.crs == CRS.from_epsg(3857)
    assert transform == [10, 0, 100, 0, -10, 200]
    assert meta["values"] == VALUES.tolist()
    assert src.nodata == (-9999 if dataset_id == "M1" else None)
    expected = VALID_MASK if dataset_id == "M2" else np.full(
        (2, 2), 255, dtype="uint8"
    )
    assert np.array_equal(band_mask, expected)
    return band_mask


def sample_dataset(dataset_id: str, src: Any,
                   report: dict[str, Any]) -> None:
    band_mask = check_source(dataset_id, src, report)
    meta = report["datasets"][dataset_id]
    meta["main_preflight"], meta["coordinate_checks"] = [], []
    ready, checks = [], meta["coordinate_checks"]
    for point_id, x, y, expected_row, expected_col in POINTS:
        gate = preflight((x, y), CRS_LABEL, src.crs)
        meta["main_preflight"].append({"point_id": point_id, gate})
        if gate["status"] != "READY":
            raise RuntimeError(f"main input failed preflight: {point_id}")
        x1, y1 = gate["x"], gate["y"]
        row, col = independent_index(x1, y1)
        rr, cc = rowcol(src.transform, x1, y1)
        ir, ic = src.index(x1, y1)
        check = {
            "point_id": point_id, "oracle": [expected_row, expected_col],
            "independent": [row, col],
            "rowcol": [int(rr), int(cc)], "index": [int(ir), int(ic)],
        }
        checks.append(check)
        assert check["independent"] == check["oracle"]
        assert check["rowcol"] == check["oracle"]
        assert check["index"] == check["oracle"]
        ready.append((x1, y1))

    raw_samples = list(src.sample(ready, indexes=1, masked=False))
    masked_samples = list(src.sample(ready, indexes=1, masked=True))
    raw = [int(np.asarray(s)[0]) for s in raw_samples]
    flags = [bool(np.ma.getmaskarray(s)[0]) for s in masked_samples]
    meta.update(raw_vector=raw, mask_vector=flags)
    assert len(raw) == len(flags) == len(POINTS)
    for point, check, value, mask_flag in zip(
        POINTS, checks, raw, flags, strict=True
    ):
        point_id, x, y, ,  = point
        row, col = check["independent"]
        inside = 0 <= row < src.height and 0 <= col < src.width
        expected = next(status for status, ids in
                        EXPECTED_SETS[dataset_id].items() if point_id in ids)
        record = {
            "fixture": FIXTURE, "dataset_id": dataset_id,
            "transform": meta["transform"], "transform_order": "a,b,c,d,e,f",
            "crs": meta["crs"], "band": BAND, "point_id": point_id,
            "x": x, "y": y, "row": row, "col": col,
            "in_bounds": inside, "raw_value": value, "mask_flag": mask_flag,
            "classification": classify(inside, mask_flag),
            "expected_classification": expected,
        }
        report["ledger"].append(record)
        if inside:
            assert mask_flag == bool(band_mask[row, col] == 0)
    assert raw == EXPECTED_RAW[dataset_id]
    assert flags == EXPECTED_MASK[dataset_id]


def run_controls(report: dict[str, Any], dataset_crs: CRS) -> None:
    controls = (
        ("C_NAN", (float("nan"), 195), CRS_LABEL, "INVALID_INPUT"),
        ("C_INF", (105, float("inf")), CRS_LABEL, "INVALID_INPUT"),
        ("C_MALFORMED", (105,), CRS_LABEL, "INVALID_INPUT"),
        ("C_CRS", (105, 195), "EPSG:4326", "HOLD_CRS"),
    )
    for control_id, coords, crs, expected in controls:
        result = {"control_id": control_id,
                  preflight(coords, crs, dataset_crs),
                  "expected_status": expected}
        report["preflight"].append(result)
        assert result["status"] == expected
    assert [r["status"] for r in report["preflight"]] == [
        "INVALID_INPUT", "INVALID_INPUT", "INVALID_INPUT", "HOLD_CRS"
    ]


def reconcile(report: dict[str, Any]) -> None:
    ledger = report["ledger"]
    keys = [(r["dataset_id"], r["point_id"]) for r in ledger]
    expected_keys = [(d, p[0]) for d in EXPECTED_SETS for p in POINTS]
    assert keys == expected_keys and len(set(keys)) == len(keys) == 24
    report["status_sets"], report["summaries"] = {}, {}
    for dataset_id, expected_sets in EXPECTED_SETS.items():
        rows = [r for r in ledger if r["dataset_id"] == dataset_id]
        assert all(r["classification"] == r["expected_classification"]
                   for r in rows)
        observed = {status: {r["point_id"] for r in rows
                            if r["classification"] == status}
                    for status in (VALID, MASKED, OUTSIDE)}
        assert observed == expected_sets
        report["status_sets"][dataset_id] = {
            status: sorted(ids) for status, ids in observed.items()
        }
        zeros = [r for r in rows if r["point_id"] in {"P1", "P8"}]
        assert len(zeros) == 2
        assert all((r["row"], r["col"]) == (0, 0)
                   and r["raw_value"] == 0 and not r["mask_flag"]
                   and r["classification"] == VALID for r in zeros)
        accepted = [r for r in rows if r["classification"] == VALID]
        raw_sum = sum(r["raw_value"] for r in rows)
        valid_sum = sum(r["raw_value"] for r in accepted)
        report["summaries"][dataset_id] = {
            "raw_count": len(rows), "raw_sum": raw_sum,
            "raw_mean": raw_sum / len(rows),
            "accepted_ids": [r["point_id"] for r in accepted],
            "accepted_count": len(accepted), "accepted_sum": valid_sum,
            "accepted_mean": valid_sum / len(accepted),
        }
    assert len(report["preflight"]) == 4
    assert len({r["control_id"] for r in report["preflight"]}) == 4
    assert [report["summaries"][d]["accepted_sum"]
            for d in EXPECTED_SETS] == [66, 66, 44]
    assert [report["summaries"][d]["accepted_count"]
            for d in EXPECTED_SETS] == [5, 5, 4]


def run(output: Path) -> dict[str, Any]:
    versions = {
        "implementation": platform.python_implementation(),
        "python": platform.python_version(), "rasterio": rasterio.__version__,
        "numpy": np.__version__, "gdal": rasterio.__gdal_version__,
    }
    report: dict[str, Any] = {
        "fixture": FIXTURE,
        "started_utc": datetime.now(timezone.utc).isoformat(),
        "versions": versions, "expected_versions": EXPECTED_VERSIONS,
        "platform": platform.platform(), "outcome": "NOT_COMPLETED",
        "oracle": [dict(point_id=p[0], x=p[1], y=p[2], row=p[3], col=p[4],
                        in_bounds=0 <= p[3] < 2 and 0 <= p[4] < 2)
                   for p in POINTS],
        "datasets": {}, "ledger": [], "preflight": [],
    }
    # Exclusive creation refuses to overwrite any previous evidence file.
    with output.open("x", encoding="utf-8") as out, ExitStack() as stack:
        try:
            if not debug:
                raise RuntimeError("Run without -O; fixture assertions required")
            if versions != EXPECTED_VERSIONS:
                raise RuntimeError("version mismatch; investigate before rerun")
            # Publish the literal oracle before CRS/transform/sampling calls.
            print(json.dumps({"literal_oracle": report["oracle"]}), flush=True)
            stack.enter_context(rasterio.Env(GDAL_TIFF_INTERNAL_MASK=True))
            run_controls(report, CRS.from_epsg(3857))
            constructors = {"M0": make_m0, "M1": make_m1, "M2": make_m2}
            memfiles = {d: build(stack) for d, build in constructors.items()}
            for dataset_id, mem in memfiles.items():
                report["phase"] = f"sampling {dataset_id}"
                src = stack.enter_context(mem.open())
                sample_dataset(dataset_id, src, report)
            report["phase"] = "reconciliation"
            reconcile(report)
            report["outcome"] = "PASS"
        except Exception as exc:
            report["outcome"] = "HARNESS_FAILURE"
            report["error"] = {
                "type": type(exc).__name__, "message": str(exc),
                "traceback": traceback.format_exc(),
            }
            raise
        finally:
            # Runs inside ExitStack, before readers/MemoryFiles are cleaned up.
            json.dump(report, out, indent=2, allow_nan=False)
            out.write("\n")
            out.flush()
            os.fsync(out.fileno())
    return report


def main() -> None:
    parser = argparse.ArgumentParser()
    parser.add_argument("report", type=Path, help="new JSON evidence filename")
    args = parser.parse_args()
    result = run(args.report)
    print(f"{result['outcome']}: 24 main records; 4 preflight controls")


if name == "__main__":
    main()

Repair the consumer while preserving source requests

The authorized repair replaces a raw-value heuristic with the explicit admission gate. It does not alter the raster, move coordinates, rewrite the point manifest, or treat an outside fallback as an alternative measurement. Start from retained original coordinates and the intended source, then rerun qualification with the corrected consumer.

Preserve the previous report and any contaminated output before changing the implementation. Record which consumer version admitted values incorrectly, which dataset identity was used, and which request population was affected. That evidence determines whether the problem can be repaired by a replay or whether missing provenance requires a hold.

Expose policy decisions as data. Preflight returns INVALID_INPUT, HOLD_CRS, or READY with a reason. For qualified coordinates with established mask evidence, the classifier returns one of the three main statuses. Only IN_BOUNDS_VALID records supply values to the accepted numerical view. Test assertions check the fixed experiment; they are not the production-facing policy interface.

Do not silently clamp P5's column -1 to 0, snap P6 to the last column, or move P7 upward. Those changes would sample different positions. Replacing masked values with zero is equally unauthorized, because zero is already a legitimate source value in the fixture. A domain owner may later approve another policy, but that policy needs its own specification and evidence.

Keep excluded IDs available to downstream users. A report with only the five accepted M0 records cannot show whether three requests were intentionally excluded, lost in transport, or never submitted. Retaining eight classified records makes the accounting explicit while allowing analysis to consume a filtered view.

Review the replay at the point where values actually enter the downstream calculation. A corrected helper does not repair a caller that continues to use its old raw-value list. Compare the accepted IDs with the ledger at that boundary and confirm that the retained zero records contribute normally. This also establishes a practical rollback limit: restoring an earlier consumer version would restore its admission behavior, not recover a scientifically justified population. Preserve the corrected and previous artifacts separately rather than overwriting the evidence needed to explain their difference.

For a deployed interface, also specify what happens when expected evidence is unavailable. A missing mask flag must not be coerced to False; an incomplete source identity must not inherit the last known raster's metadata. Hold the affected admission path and preserve the unresolved request. The repair is complete only when a replay demonstrates correct point membership, valid-zero retention, and an auditable reason for every non-admitted input, not merely when a new conditional appears in the code.

Rebuild affected summaries from qualified samples

Correcting the sampler does not retroactively correct a table or summary already computed from inadmissible values. The executable fixture calculates a deliberately naive mean over all eight raw returns and a qualified mean over admitted point records. The observed synthetic bookkeeping was:

M0 is particularly instructive. Removing three outside zeros leaves the sum unchanged at 66, but the admitted denominator is five rather than eight. An unchanged sum therefore does not demonstrate an unaffected result. M1's sentinel makes the contaminated mean conspicuous, while M2 shows that rejecting the inside invalid 22 changes the numerator as well.

These are point-request summaries. P1 and P8 each contribute their valid zero because the input population contains both IDs. They do not constitute two independent raster cells. Whether a scientific analysis should weight requests, distinct cells, or another unit is a separate analytical policy; this admission test does not choose that policy on the analyst's behalf.

Rebuild from the classified records, not from an old vector whose point identities have been lost. Retain the excluded populations alongside the new summary: P5, P6, and P7 for every dataset, plus P3 as inside but masked for M2. Record both the old and replacement artifact identities so consumers can tell which output supersedes which.

The adjacent workflow on rasterization overlap priority addresses a different source of raster-result changes. It is not part of this experiment. Here the stored pixels remain fixed; the correction concerns which requested samples may contribute to a downstream result.

The rollback boundary is important. Updating a classification column does not undo numbers already exported elsewhere. A derived output without sufficient lineage may need to remain held until its inputs can be reconstructed. No scientific effect size, operational saving, or mission-quality improvement is measured by these arithmetic examples. They demonstrate the accounting needed to rebuild honestly after admission changes.

Choose accept, repair, rebuild, or hold

The decision belongs to the evidence state, not to how reassuring the raw output looks. Use the following matrix for this declared workflow:

ACCEPT is a per-record decision. M0 and M1 have five accepted point IDs; M2 has four. A passing harness can include expected outside and masked records, because the harness tests correct classification of those records rather than universal acceptance.

REPAIR and REBUILD are different actions. A repaired consumer may produce correct new records while an old summary remains wrong. Track completion of both actions where necessary. Likewise, a successful replay against a newly chosen raster does not prove what an earlier, unidentified raster would have returned. That historical gap remains a provenance issue.

HOLD should identify what would change the decision. Examples include a trusted CRS declaration, the intended source revision, an authoritative mask, or a complete replayable point manifest. “The values look reasonable” supplies none of those missing facts. Missing ledger entries are also not evidence of successful exclusion; they are unaccounted requests.

Keep the acceptance claim within the tested raster model

The executed result concerns three local, synthetic, single-band GeoTIFFs, one fixed north-up affine, and known same-CRS inputs. It qualifies direct cell selection and selected-band validity under the recorded engine. It does not establish multiband alpha-mask precedence, rotated-grid behavior, reprojection, interpolation, remote products, or transfer completeness.

A changed runtime or source configuration requires a new qualification record. Do not edit expected outputs until a divergent result looks like a pass. First determine whether the changed behavior, the fixture, or the acceptance assumptions account for the difference. Retain the failed run and its diagnostics. The acceptance claim should remain as narrow as the evidence that actually survived reconciliation.

Assign coordinate, raster, and downstream ownership

The coordinate producer owns original point IDs, submitted coordinates, coordinate order, and CRS provenance. It should provide enough evidence for the consumer to establish what the numbers mean without guessing from their magnitude. A request with an unproved CRS remains a request worth accounting for, but not a coordinate eligible for this sampling path.

The raster owner owns the intended source identity, declared CRS, affine transform, dimensions, band selection, and mask meaning. A reader name is useful for diagnosing which open handle was used, but an ephemeral in-memory pathname is not a durable definition of the source. In this fixture, the authored values, geometry, validity configuration, and fixture revision supply that definition.

The consumer owns the admission sequence and its evidence retention. It must apply preflight to the real inputs, compare index results with the declared expectations during qualification, preserve raw and masked evidence separately, and reconcile the original request population. A helper that is correct in isolation is insufficient if the production call path bypasses it.

The downstream analyst or reviewer owns the interpretation of exclusions and the choice of analytical unit. Software can establish that P5 does not address a cell in this raster. It cannot decide whether the larger study should acquire another source, change its coverage policy, or treat the absence as informative. Those are domain decisions, not automatic consequences of the sample's fallback value.

These responsibilities can sit within a broader SAR data workflow, but the current evidence does not depend on that workflow's mission products or acquisition interfaces. Keeping ownership separate prevents a successful low-level lookup from being mistaken for approval of the entire processing chain.

Agree who may authorize a changed input rather than a rerun of the original one. Replacing a CRS declaration, changing a validity mask, or moving a coordinate changes the evidence being admitted. Such changes should produce a new request or source revision and retain a relationship to the earlier record. Without that boundary, a replay can appear to fix the sampler while actually answering a different spatial question.

Make the handoff concrete: preserve the request manifest, source definition, runtime record, admission ledger, and identities of derived outputs. An owner should be able to answer why a particular ID was admitted or held without reverse-engineering a mean. When that explanation is missing, asking for more provenance is more useful than changing a numerical threshold. The next engineer should inherit a reproducible decision, not an unexplained list of plausible values.

Strengthen remote-sensing quality practice

A valid zero is a demanding regression control precisely because it looks ordinary. In M0, the raw zeros cannot distinguish P1 and P8 from P5 through P7. The distinction becomes defensible only when request identity, independently checked indices, and selected-band validity remain attached to each result.

Keep that control when maintaining the consumer. A change that removes outside fallbacks but also loses the valid zeros has not repaired admission. A change that preserves the counts but swaps point identities has not preserved the requested population. A passing fixture should demonstrate the declared memberships, not merely produce attractive totals.

The same discipline supports practical remote-sensing education: state a contract, build a small counterexample, retain diagnostic evidence, and explain the limits of the result. No enrollment is needed to run this fixture or apply its admission rule.

For structured practice with the tools named in its syllabus, readers may review Refonte Learning's Remote Sensing Scientist/Engineer program. Its published page describes three months at 10–12 hours per week, relevant technical foundations and Python comfort, tools including QGIS/ArcGIS Pro, Google Earth Engine, GDAL/OGR and rasterio, applied capstone work, and a Training Certificate plus a Certificate of Internship on successful completion. These are the provider's stated program features, not a guaranteed external appointment.