geotiff: validate VRT writer source compatibility (#1733)

write_vrt previously read metadata from every source but used only
the first source's pixel size, sample format + bps (dtype), band
count, and CRS to write the VRT. A mismatched source would silently
produce a VRT that misplaced pixels or mis-typed data downstream.

Enforce the docstring contract: every source must agree with the
first on pixel size (within 1e-6 relative tolerance), sample format +
bps, band count, and CRS WKT (when both sides set it). Mismatches
raise ValueError naming the offending source and field.

Document nodata precedence: caller arg wins over the first source's
per-band nodata.

Adds regression coverage in test_vrt_writer_source_compat_1733.py.

Author: brendancol

xrspatial/geotiff/_vrt.py

@@ -783,7 +783,10 @@ def write_vrt(vrt_path: str, source_files: list[str], *,
     """Generate a VRT file that mosaics multiple GeoTIFF tiles.
 
     Each source file is placed in the virtual raster based on its
-    geo transform. Files must share the same CRS and pixel size.
+    geo transform. All sources must share the same pixel size, dtype
+    (sample format + bits-per-sample), band count, and CRS. Mismatches
+    raise ``ValueError`` rather than producing a misplaced or mis-typed
+    mosaic.
 
     Parameters
     ----------
@@ -796,10 +799,12 @@ def write_vrt(vrt_path: str, source_files: list[str], *,
     crs_wkt : str or None
         CRS as WKT string. If None, taken from the first source.
     nodata : float, int, or None
-        NoData value. If None, taken from the first source. Integer
-        sentinels (e.g. ``65535`` for uint16, ``-9999`` for int32) are
-        accepted so the surface lines up with the ``nodata`` kwarg on
-        ``to_geotiff`` and ``write_geotiff_gpu``.
+        NoData value applied to every band of the mosaic. Caller-supplied
+        value takes precedence; when ``None``, the first source's
+        per-band nodata is used. Integer sentinels (e.g. ``65535`` for
+        uint16, ``-9999`` for int32) are accepted so the surface lines up
+        with the ``nodata`` kwarg on ``to_geotiff`` and
+        ``write_geotiff_gpu``.
 
     Returns
     -------
@@ -849,6 +854,65 @@ def write_vrt(vrt_path: str, source_files: list[str], *,
     res_x = first['transform'].pixel_width
     res_y = first['transform'].pixel_height
 
+    # Enforce the docstring contract: every source must agree with the
+    # first on pixel size, sample format + bits-per-sample (i.e. dtype),
+    # band count, and CRS WKT (when both sides set it). Without this,
+    # build_vrt would silently produce a syntactically valid VRT that
+    # misplaces or mis-types data downstream (issue #1733).
+    #
+    # Pixel size is compared with a small relative tolerance: TIFFs
+    # written by different tools occasionally round the GeoTransform
+    # slightly, and the existing mosaic-extent math rounds the final
+    # raster size anyway. Sample format + bps must match exactly because
+    # the VRT dtype is taken from ``first`` and applied to every band.
+    _PIXEL_SIZE_RTOL = 1e-6
+
+    def _pixel_size_mismatch(a: float, b: float) -> bool:
+        # Compare magnitudes; pixel_height is negative for the common
+        # north-up case so a direct equality test would mis-flag a
+        # source with the opposite sign.
+        if a == b:
+            return False
+        denom = max(abs(a), abs(b))
+        if denom == 0.0:
+            return abs(a - b) > 0.0
+        return abs(a - b) / denom > _PIXEL_SIZE_RTOL
+
+    first_crs = first.get('crs_wkt')
+    for m in sources_meta[1:]:
+        t = m['transform']
+        if _pixel_size_mismatch(t.pixel_width, res_x) \
+                or _pixel_size_mismatch(t.pixel_height, res_y):
+            raise ValueError(
+                f"VRT source {m['path']!r} has pixel size "
+                f"({t.pixel_width}, {t.pixel_height}) which does not "
+                f"match the first source ({res_x}, {res_y}). All sources "
+                f"in a VRT must share the same pixel size."
+            )
+        if m['sample_format'] != first['sample_format'] \
+                or m['bps'] != first['bps']:
+            raise ValueError(
+                f"VRT source {m['path']!r} has sample_format="
+                f"{m['sample_format']}, bps={m['bps']} which does not "
+                f"match the first source (sample_format="
+                f"{first['sample_format']}, bps={first['bps']}). All "
+                f"sources in a VRT must share the same dtype."
+            )
+        if m['bands'] != first['bands']:
+            raise ValueError(
+                f"VRT source {m['path']!r} has {m['bands']} band(s) "
+                f"which does not match the first source "
+                f"({first['bands']} band(s)). All sources in a VRT "
+                f"must share the same band count."
+            )
+        m_crs = m.get('crs_wkt')
+        if first_crs and m_crs and m_crs != first_crs:
+            raise ValueError(
+                f"VRT source {m['path']!r} has CRS WKT that does not "
+                f"match the first source. All sources in a VRT must "
+                f"share the same CRS."
+            )
+
     # Compute the bounding box of all sources
     all_x0, all_y0, all_x1, all_y1 = [], [], [], []
     for m in sources_meta:

xrspatial/geotiff/tests/test_vrt_writer_source_compat_1733.py

@@ -0,0 +1,118 @@
+"""Regression tests for issue #1733.
+
+``write_vrt`` previously trusted the first source for resolution,
+sample format + bps (dtype), band count, and CRS. A mismatched source
+would silently produce a VRT that placed pixels incorrectly or
+re-interpreted bytes as the wrong dtype downstream.
+
+These tests assert that ``write_vrt`` now rejects mismatched sources
+with a clear ``ValueError`` covering each of those properties, and
+still accepts sources that match within a small float tolerance on
+pixel size.
+"""
+from __future__ import annotations
+
+import os
+import uuid
+
+import numpy as np
+import pytest
+import xarray as xr
+
+from xrspatial.geotiff import to_geotiff
+from xrspatial.geotiff._vrt import write_vrt
+
+
+def _unique_dir(tmp_path, label: str) -> str:
+    d = tmp_path / f"vrt_1733_{label}_{uuid.uuid4().hex[:8]}"
+    d.mkdir()
+    return str(d)
+
+
+def _write_tif(path: str, *, h: int, w: int, dtype, bands: int = 1,
+               px: float = 1.0, py: float = -1.0,
+               origin_x: float = 0.0, origin_y: float = 100.0,
+               crs: int | None = 4326) -> None:
+    if bands == 1:
+        arr = np.arange(h * w, dtype=dtype).reshape(h, w)
+        dims = ['y', 'x']
+    else:
+        arr = np.arange(h * w * bands, dtype=dtype).reshape(h, w, bands)
+        dims = ['y', 'x', 'band']
+    y = origin_y + (np.arange(h) + 0.5) * py
+    x = origin_x + (np.arange(w) + 0.5) * px
+    coords = {'y': y, 'x': x}
+    attrs = {}
+    if crs is not None:
+        attrs['crs'] = crs
+    da = xr.DataArray(arr, dims=dims, coords=coords, attrs=attrs)
+    to_geotiff(da, path, compression='none')
+
+
+def test_mismatched_pixel_size_raises(tmp_path):
+    d = _unique_dir(tmp_path, "px")
+    a = os.path.join(d, "a.tif")
+    b = os.path.join(d, "b.tif")
+    _write_tif(a, h=4, w=4, dtype=np.float32, px=1.0, py=-1.0)
+    # Place b adjacent so the geometry would otherwise work, but the
+    # pixel size disagrees.
+    _write_tif(b, h=4, w=4, dtype=np.float32, px=2.0, py=-2.0,
+               origin_x=4.0)
+    vrt = os.path.join(d, "out.vrt")
+    with pytest.raises(ValueError, match="pixel size"):
+        write_vrt(vrt, [a, b])
+
+
+def test_mismatched_dtype_raises(tmp_path):
+    d = _unique_dir(tmp_path, "dtype")
+    a = os.path.join(d, "a.tif")
+    b = os.path.join(d, "b.tif")
+    _write_tif(a, h=4, w=4, dtype=np.float32)
+    _write_tif(b, h=4, w=4, dtype=np.int16, origin_x=4.0)
+    vrt = os.path.join(d, "out.vrt")
+    with pytest.raises(ValueError, match="dtype|sample_format|bps"):
+        write_vrt(vrt, [a, b])
+
+
+def test_mismatched_band_count_raises(tmp_path):
+    d = _unique_dir(tmp_path, "bands")
+    a = os.path.join(d, "a.tif")
+    b = os.path.join(d, "b.tif")
+    _write_tif(a, h=4, w=4, dtype=np.float32, bands=1)
+    _write_tif(b, h=4, w=4, dtype=np.float32, bands=3, origin_x=4.0)
+    vrt = os.path.join(d, "out.vrt")
+    with pytest.raises(ValueError, match="band count"):
+        write_vrt(vrt, [a, b])
+
+
+def test_compatible_sources_succeed(tmp_path):
+    d = _unique_dir(tmp_path, "ok")
+    a = os.path.join(d, "a.tif")
+    b = os.path.join(d, "b.tif")
+    _write_tif(a, h=4, w=4, dtype=np.float32)
+    _write_tif(b, h=4, w=4, dtype=np.float32, origin_x=4.0)
+    vrt = os.path.join(d, "out.vrt")
+    write_vrt(vrt, [a, b])
+    assert os.path.exists(vrt)
+
+
+def test_pixel_size_within_tolerance_accepted(tmp_path):
+    d = _unique_dir(tmp_path, "tol")
+    a = os.path.join(d, "a.tif")
+    b = os.path.join(d, "b.tif")
+    _write_tif(a, h=4, w=4, dtype=np.float32, px=1.0, py=-1.0)
+    # Drift well below the 1e-6 relative tolerance.
+    _write_tif(b, h=4, w=4, dtype=np.float32,
+               px=1.0 + 1e-10, py=-1.0, origin_x=4.0)
+    vrt = os.path.join(d, "out.vrt")
+    write_vrt(vrt, [a, b])
+    assert os.path.exists(vrt)
+
+
+def test_single_source_still_works(tmp_path):
+    d = _unique_dir(tmp_path, "one")
+    a = os.path.join(d, "a.tif")
+    _write_tif(a, h=4, w=4, dtype=np.float32)
+    vrt = os.path.join(d, "out.vrt")
+    write_vrt(vrt, [a])
+    assert os.path.exists(vrt)