diff --git a/xrspatial/geotiff/_backends/vrt.py b/xrspatial/geotiff/_backends/vrt.py index 405ee17d5..344fe0720 100644 --- a/xrspatial/geotiff/_backends/vrt.py +++ b/xrspatial/geotiff/_backends/vrt.py @@ -571,15 +571,15 @@ def read_vrt(source: str, *, _vrt_is_rotated = ( gt is not None and (gt[2] != 0.0 or gt[4] != 0.0) ) + height, width = arr.shape[:2] + if window is not None: + r0 = max(0, window[0]) + c0 = max(0, window[1]) + coord_window = (r0, c0, r0 + height, c0 + width) + else: + coord_window = None if gt is not None: origin_x, res_x, _, origin_y, _, res_y = gt - height, width = arr.shape[:2] - if window is not None: - r0 = max(0, window[0]) - c0 = max(0, window[1]) - coord_window = (r0, c0, r0 + height, c0 + width) - else: - coord_window = None # Rotated VRTs emit int64 pixel coords to match the eager # non-VRT rotated path. Without this gate # the VRT branch handed back float projected coords while @@ -594,7 +594,18 @@ def read_vrt(source: str, *, has_georef=not _vrt_is_rotated, ) else: - coords = {} + # No ```` at all: synthesise integer pixel coords + # ``0..N-1`` so a no-georef VRT carries the same x/y arrays the + # non-VRT no-georef read path emits (see ``_coords.py`` + # ``has_georef=False`` fallback). Previously this branch handed + # back an empty coord dict, so downstream code that assumes x/y + # exist broke on no-georef VRTs but worked on the equivalent + # plain GeoTIFF. + coords = _coords_from_pixel_geometry( + 0.0, 0.0, 1.0, 1.0, height, width, + window=coord_window, + has_georef=False, + ) # Select the per-band nodata sentinel. When a specific band is # selected, source its nodata from that band's ```` @@ -1099,13 +1110,12 @@ def _read_vrt_chunked(source, *, window, band, name, chunks, gpu, dtype, # extent. Mirrors the eager branch in ``read_vrt`` so chunked and # eager reads share the same x/y arrays. gt = vrt.geo_transform - coords = {} _vrt_is_rotated = ( gt is not None and (gt[2] != 0.0 or gt[4] != 0.0) ) + coord_window = (win_r0, win_c0, win_r0 + full_h, win_c0 + full_w) if gt is not None: origin_x, res_x, _, origin_y, _, res_y = gt - coord_window = (win_r0, win_c0, win_r0 + full_h, win_c0 + full_w) # Rotated VRTs emit int64 pixel coords to match the eager # non-VRT rotated path. Without this gate the # chunked VRT branch handed back float projected coords on a @@ -1117,6 +1127,17 @@ def _read_vrt_chunked(source, *, window, band, name, chunks, gpu, dtype, window=coord_window, has_georef=not _vrt_is_rotated, ) + else: + # No ```` at all: synthesise integer pixel coords + # so the chunked no-georef VRT read carries the same x/y arrays + # the eager and non-VRT no-georef paths emit (see ``_coords.py`` + # ``has_georef=False`` fallback). Previously this branch left + # ``coords = {}``, dropping x/y on no-georef VRTs. + coords = _coords_from_pixel_geometry( + 0.0, 0.0, 1.0, 1.0, full_h, full_w, + window=coord_window, + has_georef=False, + ) # Surface the nodata sentinel for the selected band. The chunked # path declares ``float64`` up front whenever any selected band has diff --git a/xrspatial/geotiff/tests/vrt/test_parity.py b/xrspatial/geotiff/tests/vrt/test_parity.py index 84ed7770a..a210bc6d5 100644 --- a/xrspatial/geotiff/tests/vrt/test_parity.py +++ b/xrspatial/geotiff/tests/vrt/test_parity.py @@ -743,6 +743,61 @@ def test_vrt_chunked_none_matches_dask(tmp_path): assert tiff_attrs == vrt_attrs +def test_vrt_eager_none_synthesizes_pixel_coords(tmp_path): + """Issue #2818: a no-```` VRT read eagerly must + synthesise integer x/y pixel coords, matching the non-VRT + no-georef read instead of dropping coords entirely.""" + tiff, vrt = _make_none_pair(tmp_path, 'none_coords_2818') + ref = open_geotiff(tiff) + vrt_da = read_vrt(vrt) + assert ref.attrs['georef_status'] == GEOREF_STATUS_NONE + assert vrt_da.attrs['georef_status'] == GEOREF_STATUS_NONE + assert 'x' in vrt_da.coords and 'y' in vrt_da.coords + np.testing.assert_array_equal( + _coord_view(vrt_da, 'x'), _coord_view(ref, 'x')) + np.testing.assert_array_equal( + _coord_view(vrt_da, 'y'), _coord_view(ref, 'y')) + assert vrt_da.coords['x'].dtype == ref.coords['x'].dtype + assert vrt_da.coords['y'].dtype == ref.coords['y'].dtype + + +def test_vrt_chunked_none_synthesizes_pixel_coords(tmp_path): + """Issue #2818: a no-```` VRT read chunked must + synthesise the same integer x/y pixel coords as the non-VRT dask + no-georef read rather than dropping coords entirely.""" + tiff, vrt = _make_none_pair(tmp_path, 'none_coords_chunked_2818') + ref = read_geotiff_dask(tiff, chunks=2) + vrt_da = read_vrt(vrt, chunks=2) + assert ref.attrs['georef_status'] == GEOREF_STATUS_NONE + assert vrt_da.attrs['georef_status'] == GEOREF_STATUS_NONE + assert 'x' in vrt_da.coords and 'y' in vrt_da.coords + np.testing.assert_array_equal( + _coord_view(vrt_da, 'x'), _coord_view(ref, 'x')) + np.testing.assert_array_equal( + _coord_view(vrt_da, 'y'), _coord_view(ref, 'y')) + assert vrt_da.coords['x'].dtype == ref.coords['x'].dtype + assert vrt_da.coords['y'].dtype == ref.coords['y'].dtype + + +def test_vrt_none_windowed_synthesizes_offset_pixel_coords(tmp_path): + """Issue #2818: a windowed no-georef VRT read shifts the synthesised + integer coords to the window offset, matching the non-VRT windowed + read, on both the eager and chunked paths.""" + tiff, vrt = _make_none_pair(tmp_path, 'none_coords_win_2818') + window = (1, 2, 4, 4) + ref = open_geotiff(tiff, window=window) + eager = read_vrt(vrt, window=window) + chunked = read_vrt(vrt, window=window, chunks=2) + for label, actual in (('eager', eager), ('chunked', chunked)): + assert 'x' in actual.coords and 'y' in actual.coords, label + np.testing.assert_array_equal( + _coord_view(actual, 'x'), _coord_view(ref, 'x'), + err_msg=label) + np.testing.assert_array_equal( + _coord_view(actual, 'y'), _coord_view(ref, 'y'), + err_msg=label) + + def test_vrt_chunked_rotated_dropped(tmp_path): _, vrt = _make_rotated_pair(tmp_path, 'rot_chunked_2180') attrs = dict(read_vrt(vrt, allow_rotated=True, chunks=2).attrs)