From 1e95eeb1d1d8cfec0a15b58d4d8c624afff33de2 Mon Sep 17 00:00:00 2001 From: Brendan Collins Date: Thu, 25 Jun 2026 12:24:29 -0700 Subject: [PATCH 1/2] Fill all-NaN templates in generate_terrain numpy/cupy paths (#3525) The numpy and cupy backends zeroed their work buffer with data * 0, but NaN * 0 == NaN, so an all-NaN template (what from_template() returns) was never cleared and the octave loop accumulated noise onto NaN, leaving the whole output NaN. Saved GeoTIFFs read as 100% nodata in QGIS. Use np.zeros_like / cupy.zeros_like instead. The dask and dask+cupy paths already allocated fresh zeros per chunk and were unaffected. Adds tests covering all four backends: an all-NaN template now yields finite terrain identical to an all-zeros template. --- xrspatial/terrain.py | 10 +++++-- xrspatial/tests/test_terrain.py | 53 +++++++++++++++++++++++++++++++++ 2 files changed, 61 insertions(+), 2 deletions(-) diff --git a/xrspatial/terrain.py b/xrspatial/terrain.py index 234762c93..a0937e034 100644 --- a/xrspatial/terrain.py +++ b/xrspatial/terrain.py @@ -181,7 +181,10 @@ def _terrain_numpy(data, seed, x_range_scaled, y_range_scaled, zfactor, warp_strength, warp_octaves, worley_blend, worley_seed): _check_memory(data.shape[0], data.shape[1]) - data = data * 0 + # zeros_like (not data * 0) so an all-NaN template — e.g. from + # from_template() — is overwritten instead of propagating NaN + # everywhere (NaN * 0 == NaN). + data = np.zeros_like(data) data[:] = _gen_terrain( data, seed, x_range=x_range_scaled, y_range=y_range_scaled, octaves=octaves, lacunarity=lacunarity, persistence=persistence, @@ -437,7 +440,10 @@ def _terrain_cupy(data, seed, x_range_scaled, y_range_scaled, zfactor, warp_strength, warp_octaves, worley_blend, worley_seed): _check_gpu_memory(data.shape[0], data.shape[1]) - data = data * 0 + # zeros_like (not data * 0) so an all-NaN template — e.g. from + # from_template() — is overwritten instead of propagating NaN + # everywhere (NaN * 0 == NaN). + data = cupy.zeros_like(data) data[:] = _terrain_gpu( data, seed, x_range=x_range_scaled, y_range=y_range_scaled, octaves=octaves, lacunarity=lacunarity, persistence=persistence, diff --git a/xrspatial/tests/test_terrain.py b/xrspatial/tests/test_terrain.py index 74d99d710..34befc89b 100644 --- a/xrspatial/tests/test_terrain.py +++ b/xrspatial/tests/test_terrain.py @@ -470,3 +470,56 @@ def test_generate_terrain_dask_cupy_preserves_coords_attrs(): np.testing.assert_array_equal(terrain.x.data, src.x.data) assert terrain.attrs['crs'] == 'EPSG:5070' assert terrain.attrs['res'] == (30, 30) + + +# ===================================================================== +# Issue #3525: an all-NaN template (e.g. from_template) must be filled, +# not propagated. data * 0 left NaN * 0 == NaN on numpy / cupy. +# ===================================================================== + +def _nan_template(backend='numpy'): + """An empty all-NaN grid, like from_template() returns.""" + data = np.full((50, 50), np.nan, dtype=np.float32) + raster = xr.DataArray(data, dims=['y', 'x']) + if has_cuda_and_cupy() and 'cupy' in backend: + import cupy + raster.data = cupy.asarray(raster.data) + if 'dask' in backend and da is not None: + raster.data = da.from_array(raster.data, chunks=(10, 10)) + return raster + + +def test_terrain_all_nan_template_is_finite(): + """All-NaN input must produce all-finite terrain, not all-NaN.""" + terrain = generate_terrain(_nan_template()) + assert np.isfinite(terrain.data).all() + + +def test_terrain_all_nan_template_matches_zeros_template(): + """generate_terrain ignores input values: an all-NaN template and an + all-zeros template must yield identical terrain.""" + from_nan = generate_terrain(_nan_template()) + from_zeros = generate_terrain(create_test_arr()) + np.testing.assert_array_equal(from_nan.data, from_zeros.data) + + +@dask_array_available +def test_terrain_all_nan_template_dask_is_finite(): + terrain = generate_terrain(_nan_template(backend='dask')) + assert np.isfinite(terrain.compute().data).all() + + +@cuda_and_cupy_available +def test_terrain_all_nan_template_cupy_matches_numpy(): + t_np = generate_terrain(_nan_template()) + t_cupy = generate_terrain(_nan_template(backend='cupy')) + assert np.isfinite(t_cupy.data.get()).all() + np.testing.assert_allclose(t_np.data, t_cupy.data.get(), + rtol=1e-5, atol=1e-7) + + +@cuda_and_cupy_available +@dask_array_available +def test_terrain_all_nan_template_dask_cupy_is_finite(): + terrain = generate_terrain(_nan_template(backend='dask+cupy')) + assert np.isfinite(terrain.compute().data.get()).all() From ac34ba9f1c7f4466e038b39fcc9f19e0ac8618a4 Mon Sep 17 00:00:00 2001 From: Brendan Collins Date: Thu, 25 Jun 2026 12:26:44 -0700 Subject: [PATCH 2/2] Address review: NaN-template cross-backend parity + drop em dashes (#3525) - dask / dask+cupy NaN-template tests now assert allclose against the numpy result instead of only checking finiteness. - Replace em dashes in the two new code comments with parentheses. --- xrspatial/terrain.py | 8 ++++---- xrspatial/tests/test_terrain.py | 17 +++++++++++------ 2 files changed, 15 insertions(+), 10 deletions(-) diff --git a/xrspatial/terrain.py b/xrspatial/terrain.py index a0937e034..c66b32f8a 100644 --- a/xrspatial/terrain.py +++ b/xrspatial/terrain.py @@ -181,8 +181,8 @@ def _terrain_numpy(data, seed, x_range_scaled, y_range_scaled, zfactor, warp_strength, warp_octaves, worley_blend, worley_seed): _check_memory(data.shape[0], data.shape[1]) - # zeros_like (not data * 0) so an all-NaN template — e.g. from - # from_template() — is overwritten instead of propagating NaN + # zeros_like (not data * 0) so an all-NaN template (e.g. from + # from_template()) is overwritten instead of propagating NaN # everywhere (NaN * 0 == NaN). data = np.zeros_like(data) data[:] = _gen_terrain( @@ -440,8 +440,8 @@ def _terrain_cupy(data, seed, x_range_scaled, y_range_scaled, zfactor, warp_strength, warp_octaves, worley_blend, worley_seed): _check_gpu_memory(data.shape[0], data.shape[1]) - # zeros_like (not data * 0) so an all-NaN template — e.g. from - # from_template() — is overwritten instead of propagating NaN + # zeros_like (not data * 0) so an all-NaN template (e.g. from + # from_template()) is overwritten instead of propagating NaN # everywhere (NaN * 0 == NaN). data = cupy.zeros_like(data) data[:] = _terrain_gpu( diff --git a/xrspatial/tests/test_terrain.py b/xrspatial/tests/test_terrain.py index 34befc89b..797faa4d5 100644 --- a/xrspatial/tests/test_terrain.py +++ b/xrspatial/tests/test_terrain.py @@ -504,9 +504,11 @@ def test_terrain_all_nan_template_matches_zeros_template(): @dask_array_available -def test_terrain_all_nan_template_dask_is_finite(): - terrain = generate_terrain(_nan_template(backend='dask')) - assert np.isfinite(terrain.compute().data).all() +def test_terrain_all_nan_template_dask_matches_numpy(): + t_np = generate_terrain(_nan_template()) + t_dask = generate_terrain(_nan_template(backend='dask')).compute() + assert np.isfinite(t_dask.data).all() + np.testing.assert_allclose(t_np.data, t_dask.data, rtol=1e-5, atol=1e-7) @cuda_and_cupy_available @@ -520,6 +522,9 @@ def test_terrain_all_nan_template_cupy_matches_numpy(): @cuda_and_cupy_available @dask_array_available -def test_terrain_all_nan_template_dask_cupy_is_finite(): - terrain = generate_terrain(_nan_template(backend='dask+cupy')) - assert np.isfinite(terrain.compute().data.get()).all() +def test_terrain_all_nan_template_dask_cupy_matches_numpy(): + t_np = generate_terrain(_nan_template()) + t_dc = generate_terrain(_nan_template(backend='dask+cupy')).compute() + assert np.isfinite(t_dc.data.get()).all() + np.testing.assert_allclose(t_np.data, t_dc.data.get(), + rtol=1e-4, atol=1e-4)