Skip to content

Commit 112cd66

Browse files
authored
Stop surface_distance outputs adopting the dask graph token as .name (#3708) (#3716)
surface_distance, surface_allocation and surface_direction built their result with xr.DataArray(..., coords=, dims=, attrs=) and left name unset. DataArray.__init__ then falls back to getattr(data, "name"), which on a dask array is the graph key, so the dask backends returned '_trim-<hash>' (bounded map_overlap route), 'xrspatial.surface_*-<hash>' (unbounded iterative route) or 'asarray-<hash>' (dask+cupy unbounded) while numpy and cupy returned None. The divergence is user-visible through .to_dataset(): the numpy result raises "unable to convert unnamed DataArray", the dask result silently creates a variable named after the hash, and the hash moves when chunking or max_distance changes. Route the three functions through a shared _wrap_result() helper that resets .name to None after construction, matching the proximity/allocation/direction trio this module mirrors. Same bug class as cost_distance #3344 and pathfinding #3652. Tests cover .name parity over 4 backends x 3 functions x bounded/unbounded max_distance, plus an attrs/coords/dims/dtype preservation guard. The 12 dask cases fail without the fix. Also records the metadata sweep result for this module in .claude/sweep-metadata-state.csv.
1 parent 737f9b4 commit 112cd66

3 files changed

Lines changed: 88 additions & 18 deletions

File tree

.claude/sweep-metadata-state.csv

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -17,6 +17,7 @@ proximity,2026-05-29,2723,MEDIUM,4;5,"Audited 2026-05-29 (agent-a61dbadc2452a200
1717
rasterize,2026-06-09,3087,MEDIUM,1,GeoDataFrame .crs dropped on no-like path (Cat 1); fixed via #3087 emitting attrs crs/crs_wkt when output has no CRS. like-path attrs/coords/dims/nodata verified live on all 4 backends (CUDA available); Cats 2-5 clean.
1818
reproject,2026-06-12,3262,MEDIUM,4,"Re-audited 2026-06-12 (agent-ae420c90e50a23c5c worktree, branch deep-sweep-metadata-reproject-2026-06-12). CUDA available; all 4 backends (numpy/cupy/dask+numpy/dask+cupy) run live end-to-end for reproject() and merge(). Cat 1 attrs (crs/nodata/res/transform/_FillValue/nodatavals refreshed or carried, crs_wkt dropped), Cat 2 coords (pixel-center verified numerically, scalar time + band coord carry, float64), Cat 3 dims (lat/lon names, band-first (band,y,x) round-trip), Cat 4 int16 sentinel parity, and Cat 5 cross-backend attr parity all identical across the 4 backends for reproject(); vertical_crs=4979/vertical_datum verified on numpy + dask; geoid_height_raster carries input attrs per its documented contract. NEW MEDIUM finding #3262 (Cat 4): merge() hardcoded float64 output on every path (_merge_inmemory, _merge_dask template/meta, empty-chunk fills) while reproject() round-trips integer dtypes on all 5 paths (#2185/#2505/#3093/#3096); undocumented and unpinned by tests, so an int16/uint8 mosaic silently promoted (8x memory for uint8, GeoTIFF round-trip changes file dtype). Fix on this branch: shared-integer-dtype inputs now cast back via the reproject round/clip/cast convention (_cast_merged_dtype), output nodata resolved with _detect_nodata dtype hint (NaN->sentinel swap per #2185, explicit out-of-range raises per #2572), dask template/meta + empty-chunk fills use the output dtype (#3096 trap), docstring documents the rule; mixed/float inputs keep float64. 13 new tests in TestMergeIntegerDtype (eager/dask/dask-empty-chunk/cupy, sentinel defaults, mean rounding, out-of-range raise); full reproject suite 514 passed. LOW (documented, not fixed): reproject() docstring says dask inputs are fully lazy but the dask+cupy VRAM-fitting fast path returns an eager cupy array (codified in tests; doc nit). Prior LOW from 2026-06-09 (geoid_height ndarray return for DataArray input) unchanged."
1919
resample,2026-05-27,2542,MEDIUM,2;4;5,"Audited 2026-05-27 (agent-a8135a6a246ecb93c worktree, branch deep-sweep-metadata-resample-2026-05-27). Cat 2 MEDIUM + Cat 4 MEDIUM + Cat 5 MEDIUM all rolled into issue #2542. (a) 2D non-identity path dropped scalar non-dim coords like rioxarrays spatial_ref and squeezed time/band selectors; identity path (scale==1.0, agg.copy()) and 3D path (per-band xr.concat) preserved them, so the bug was path-inconsistent (Cat 5). (b) _resolve_nodata reads attrs[nodata] as a fallback sentinel but the output post-processing only refreshed _FillValue and nodatavals, leaving attrs[nodata]=-9999 alongside data that was now NaN. Fix in resample(): refresh attrs[nodata] to NaN whenever the input had it, and carry across zero-dim non-dim coords on the 2D non-identity path. 7 new tests in TestMetadataPropagation cover nodata-attr refresh, spatial_ref/scalar coord carry, identity-vs-downsample coord parity, and the explicit choice to drop spatially-shaped extra coords. 4-backend (numpy/cupy/dask+numpy/dask+cupy) parity verified for spatial_ref carry; nodata-attr refresh verified on numpy/cupy/dask+numpy (dask+cupy non-NaN nodata masking hits a pre-existing xarray xr.where + cupy.astype quirk unrelated to this audit). Full resample test suite (175 passed) clean."
20+
surface_distance,2026-08-16,3708,MEDIUM,5,"Audited 2026-08-16 (agent-ab3afc98fdd524ae1 worktree, branch deep-sweep-metadata-surface_distance-2026-08-16). CUDA available; all 4 backends (numpy/cupy/dask+numpy/dask+cupy) run live end-to-end for surface_distance/surface_allocation/surface_direction, across both dask routes (bounded map_overlap and unbounded iterative tile Dijkstra). Cat 1 attrs, Cat 2 coords, Cat 3 dims all clean: the three public functions re-emit coords=raster.coords, dims=raster.dims, attrs=raster.attrs, so res/crs/transform/nodatavals/_FillValue, extra scalar coords (spatial_ref), coord values and dtypes, and custom dim names via x=/y= (lat/lon) all survive identically on every backend. Output dtype is float32 on all four backends, matching the docstring. NEW MEDIUM finding #3708 (Cat 5): the three functions left name unset at the xr.DataArray constructor, so DataArray.__init__ fell back to getattr(data, 'name') and the dask backends adopted the graph key as .name -- '_trim-<hash>' from the bounded map_overlap route, 'xrspatial.surface_*-<hash>' from the unbounded iterative route, and 'asarray-<hash>' from the dask+cupy unbounded route -- while numpy and cupy returned None. User-visible via .to_dataset(): numpy raises 'unable to convert unnamed DataArray', dask silently creates a variable named '_trim-fe8156...' whose hash moves with chunking and max_distance. Recurring bug class: slope #2837, aspect #2841, focal #2733, viewshed #2743, zonal #2611, cost_distance #3344, pathfinding #3652. Fix in PR #3709: a shared _wrap_result() helper that sets result.name = None after construction, matching the proximity/allocation/direction trio this module mirrors. 36 new tests (name parity across 4 backends x 3 functions x bounded/unbounded, plus an attrs/coords/dims/dtype preservation guard); the 12 dask cases fail without the fix. Full suite 73 passed. LOW, documented not fixed: attrs are copied verbatim so a user-supplied nodatavals=(-9999,) / _FillValue=-9999 stays on an output that actually uses NaN as its nodata sentinel, and a 'units' attr carries over onto surface_direction output measured in degrees -- both are the library-wide attrs=raster.attrs convention (proximity and cost_distance behave identically), not surface_distance-specific. Also noted out of scope: the module's _dask_task_name_kwargs task names are overwritten by map_overlap's '_trim-' prefix on the bounded route, so the #3256 task-naming work does not reach that path. No CRITICAL or HIGH findings."
2021
viewshed,2026-05-29,2743,MEDIUM,4;5,output .name differed across backends (None/viewshed/dask-token) and dtype float32 on GPU vs float64 on CPU; added name= param and forced float64 on all backends; attrs/coords/dims already preserved
2122
visibility,2026-06-10,3193,HIGH,5,"cupy backend crash in cumulative_viewshed/visibility_frequency (count np vs cupy add) -> no result/metadata emitted; fixed by cupy count branch + cupy tests. numpy/dask preserve coords/dims/attrs incl crs; visibility_frequency keeps attrs through astype/divide. line_of_sight Dataset drops crs/transform (LOW, transect not raster, documented only)."
2223
zonal,2026-05-29,2611,MEDIUM,5,"Audited 2026-05-29 (agent-ae8d8b65cc3a5c40a worktree, branch deep-sweep-metadata-zonal-2026-05-29). CUDA available; all 4 backends (numpy/cupy/dask+numpy/dask+cupy) run live. 5 DataArray-returning functions checked end-to-end: apply, regions, hypsometric_integral, trim, crop. attrs (res/crs/transform/nodatavals), dims, and coords preserved correctly on all 4 backends for every function; trim/crop slice coords with no half-pixel drift. stats() and crosstab() return DataFrames by design so Cat 1-3 DataArray checks N/A. NEW MEDIUM finding #2611 (Cat 5): apply() never set output .name, so numpy/cupy returned None while dask+numpy/dask+cupy inherited a non-deterministic internal dask task name (e.g. _chunk_fn-<hash>). regions/hypsometric_integral/trim/crop all set deterministic names; apply was the outlier. Fix in PR #2611/#2622: add name param (default None) and assign result.name after DataArray construction (setting name= at construction does not override the dask graph name). New parametrized test test_apply_name_consistent_across_backends covers default-None and explicit-name on all 4 backends. Full zonal suite 213 passed. No other CRITICAL/HIGH/MEDIUM findings; no LOW findings to document."

xrspatial/surface_distance.py

Lines changed: 22 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -1435,6 +1435,25 @@ def _compute(raster, elevation, x, y, target_values, max_distance,
14351435
# ---------------------------------------------------------------------------
14361436

14371437

1438+
def _wrap_result(result_data, raster):
1439+
"""Wrap raw output in a DataArray carrying the input's spatial metadata.
1440+
1441+
The name is reset after construction: on dask backends
1442+
``xr.DataArray(..., name=None)`` adopts the dask array's graph key
1443+
(``_trim-<hash>`` from map_overlap, ``xrspatial.surface_*-<hash>`` from
1444+
the iterative path) as ``.name``, while numpy and cupy return None.
1445+
Issue #3708; same fix as cost_distance #3344 and pathfinding #3652.
1446+
"""
1447+
result = xr.DataArray(
1448+
result_data,
1449+
coords=raster.coords,
1450+
dims=raster.dims,
1451+
attrs=raster.attrs,
1452+
)
1453+
result.name = None
1454+
return result
1455+
1456+
14381457
@supports_dataset
14391458
def surface_distance(
14401459
raster: xr.DataArray,
@@ -1488,12 +1507,7 @@ def surface_distance(
14881507
raster, elevation, x, y, target_values, max_distance,
14891508
connectivity, method, DISTANCE,
14901509
)
1491-
return xr.DataArray(
1492-
result_data,
1493-
coords=raster.coords,
1494-
dims=raster.dims,
1495-
attrs=raster.attrs,
1496-
)
1510+
return _wrap_result(result_data, raster)
14971511

14981512

14991513
@supports_dataset
@@ -1530,12 +1544,7 @@ def surface_allocation(
15301544
raster, elevation, x, y, target_values, max_distance,
15311545
connectivity, method, ALLOCATION,
15321546
)
1533-
return xr.DataArray(
1534-
result_data,
1535-
coords=raster.coords,
1536-
dims=raster.dims,
1537-
attrs=raster.attrs,
1538-
)
1547+
return _wrap_result(result_data, raster)
15391548

15401549

15411550
@supports_dataset
@@ -1573,9 +1582,4 @@ def surface_direction(
15731582
raster, elevation, x, y, target_values, max_distance,
15741583
connectivity, method, DIRECTION,
15751584
)
1576-
return xr.DataArray(
1577-
result_data,
1578-
coords=raster.coords,
1579-
dims=raster.dims,
1580-
attrs=raster.attrs,
1581-
)
1585+
return _wrap_result(result_data, raster)

xrspatial/tests/test_surface_distance.py

Lines changed: 65 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -705,6 +705,71 @@ def test_geodesic_basic():
705705
assert 100000 < sd[pos] < 130000 # roughly 100-130 km
706706

707707

708+
# ---------------------------------------------------------------------------
709+
# Metadata propagation (issue #3708)
710+
# ---------------------------------------------------------------------------
711+
712+
713+
def _metadata_backends():
714+
backends = ['numpy']
715+
if da is not None:
716+
backends.append('dask+numpy')
717+
if has_cuda_and_cupy():
718+
backends.append('cupy')
719+
if da is not None:
720+
backends.append('dask+cupy')
721+
return backends
722+
723+
724+
@pytest.mark.parametrize("func", [surface_distance, surface_allocation,
725+
surface_direction])
726+
@pytest.mark.parametrize("max_distance", [3.0, np.inf],
727+
ids=['bounded', 'unbounded'])
728+
@pytest.mark.parametrize("backend", _metadata_backends())
729+
def test_output_name_consistent_across_backends(backend, max_distance, func):
730+
"""Outputs must not adopt the dask graph token as .name.
731+
732+
Without the post-construction reset the dask backends returned
733+
'_trim-<hash>' (bounded map_overlap route),
734+
'xrspatial.surface_*-<hash>' (unbounded iterative route) or
735+
'asarray-<hash>' (dask+cupy unbounded), while numpy and cupy returned
736+
None. Same bug class as cost_distance #3344 and pathfinding #3652.
737+
"""
738+
source = np.zeros((6, 6), dtype=np.float64)
739+
source[0, 0] = 1.0
740+
elev = np.arange(36, dtype=np.float64).reshape(6, 6) * 0.1
741+
raster = _make_raster(source, backend=backend)
742+
elevation = _make_raster(elev, backend=backend)
743+
744+
result = func(raster, elevation, max_distance=max_distance)
745+
assert result.name is None
746+
747+
748+
@pytest.mark.parametrize("func", [surface_distance, surface_allocation,
749+
surface_direction])
750+
@pytest.mark.parametrize("backend", _metadata_backends())
751+
def test_output_preserves_attrs_coords_dims(backend, func):
752+
"""attrs, coords and dims come through unchanged on every backend."""
753+
source = np.zeros((6, 6), dtype=np.float64)
754+
source[0, 0] = 1.0
755+
elev = np.arange(36, dtype=np.float64).reshape(6, 6) * 0.1
756+
raster = _make_raster(source, backend=backend)
757+
elevation = _make_raster(elev, backend=backend)
758+
raster.attrs.update({'crs': 3857, 'nodatavals': (-9999.0,),
759+
'transform': (1.0, 0.0, 0.0, 0.0, -1.0, 0.0)})
760+
raster = raster.assign_coords(spatial_ref=0)
761+
762+
result = func(raster, elevation, max_distance=3.0)
763+
764+
assert result.dims == raster.dims
765+
assert result.attrs == raster.attrs
766+
assert set(result.coords) == set(raster.coords)
767+
for name in raster.coords:
768+
np.testing.assert_array_equal(result.coords[name].values,
769+
raster.coords[name].values)
770+
assert result.dtype == np.float32
771+
772+
708773
# ---------------------------------------------------------------------------
709774
# Memory guard
710775
# ---------------------------------------------------------------------------

0 commit comments

Comments
 (0)