Skip to content

Commit 4436f14

Browse files
committed
Apply circular distance mask to max_distance output
The window extraction uses a square (Chebyshev) region but max_distance is a Euclidean radius. Add a circular mask after computing the window result so cells at the corners beyond max_distance are set to INVISIBLE. Update tests to account for the circular boundary.
1 parent 4091b6a commit 4436f14

File tree

2 files changed

+27
-4
lines changed

2 files changed

+27
-4
lines changed

xrspatial/tests/test_viewshed.py

Lines changed: 14 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -281,9 +281,14 @@ def test_viewshed_dask_max_distance_lazy_output():
281281
assert isinstance(v.data, da.Array)
282282
assert v.shape == (ny, nx)
283283
# Only compute a small slice to verify correctness
284-
center = v.isel(y=slice(49990, 50011), x=slice(49990, 50011)).values
285-
assert center[10, 10] == 180.0 # observer cell
286-
assert (center > INVISIBLE).all() # flat terrain, all visible
284+
# max_distance=10 → radius_cells=10 → window is obs ± 10
285+
center = v.isel(y=slice(49989, 50012), x=slice(49989, 50012)).values
286+
# Observer is at index 11 within this 23-cell slice
287+
assert center[11, 11] == 180.0 # observer cell
288+
# Cells within the circle should be visible (flat terrain, observer up)
289+
assert center[11, 13] > INVISIBLE # 2 cells away
290+
# Corner (49989,49989) is sqrt(11^2+11^2) ≈ 15.6 from observer → outside
291+
assert center[0, 0] == INVISIBLE
287292

288293

289294
def test_viewshed_numpy_max_distance():
@@ -342,10 +347,15 @@ def test_viewshed_max_distance_matches_full(backend):
342347
dist_vals = v_dist.values if isinstance(v_dist.data, np.ndarray) \
343348
else v_dist.data.get()
344349

345-
# Within the window, results should match
350+
# Within the circular radius, results should match
346351
obs_r, obs_c = 5, 4
352+
max_d = 3.5
347353
for r in range(max(0, obs_r - 3), min(ny, obs_r + 4)):
348354
for c in range(max(0, obs_c - 3), min(nx, obs_c + 4)):
355+
dr = (r - obs_r) * 1.0 # ns_res = 1
356+
dc = (c - obs_c) * 1.0 # ew_res = 1
357+
if np.sqrt(dr**2 + dc**2) > max_d:
358+
continue # outside circle — correctly INVISIBLE
349359
np.testing.assert_allclose(
350360
dist_vals[r, c], full_vals[r, c],
351361
atol=0.03,

xrspatial/viewshed.py

Lines changed: 13 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2090,6 +2090,19 @@ def _viewshed_windowed(raster, x, y, observer_elev, target_elev,
20902090
local_result = _viewshed_cpu(
20912091
window, x, y, observer_elev, target_elev)
20922092

2093+
# Mask cells beyond max_distance (the window is a square, not a circle)
2094+
win_y = local_result.coords['y'].values
2095+
win_x = local_result.coords['x'].values
2096+
wx, wy = np.meshgrid(win_x, win_y)
2097+
dist_sq = (wx - x) ** 2 + (wy - y) ** 2
2098+
outside = dist_sq > max_distance ** 2
2099+
if isinstance(local_result.data, np.ndarray):
2100+
local_result.values[outside] = INVISIBLE
2101+
else:
2102+
# cupy path
2103+
import cupy as cp
2104+
local_result.data[cp.asarray(outside)] = INVISIBLE
2105+
20932106
# Embed in full-size INVISIBLE output, preserving array type
20942107
is_dask = has_dask_array() and isinstance(raster.data, da.Array)
20952108

0 commit comments

Comments
 (0)