Decouple datashader from xarray-spatial

[Melissari1997]

Closes #3640. Implements the consolidated datashader decoupling across #3450, #3451, and #3453–#3464 (#3452 was already closed).

What changed

• Phase 1 — source cleanup (Remove canvas_like and datashader references from source code #3451): removed canvas_like() from xrspatial/utils.py (the only library-source datashader import; no callers, never exported). Rewrote the suggest_zonal_canvas docstring example to use template.xrs.rasterize(gdf, column='id', merge='max'). Reworded the pixel-center comment in xrspatial/terrain.py.
• Phase 2 — notebook/script migrations (Migrate surface user guide off datashader #3453–Migrate animated_hillshade example script off datashader #3464): all 9 user-guide notebooks and 3 example notebooks/script migrated. Vector rasterization (Canvas().points/line/polygons + ds.mean/max/sum/min) → template.xrs.rasterize(gdf, column=, merge=) on GeoDataFrames. Rendering (shade/stack/dynspread/set_background + Elevation/Set1/Set3/inferno) → matplotlib _shade/_stack/_show helpers.
  • Migrate multispectral user guide off datashader #3458 multispectral: Canvas().raster(agg='mean') resampling → xarray .coarsen/.interp; orient_array/Images → the helpers.
  • Migrate pharmacy-deserts example off datashader #3462 pharmacy-deserts: merge='mean' is not a supported reducer, so non-overlapping census polygons use merge='last' (equivalent for a single covering polygon per pixel); counties use merge='max' to preserve the ds.max intent.
• Phase 3 — packaging (Drop datashader from packaging and install docs #3450, strictly last): replaced datashader >= 0.15.0 in the examples extra with matplotlib, geopandas, shapely. Updated installation.rst to match.

Backend coverage

No backend code changed — this is a packaging/docs/internal-cleanup effort with no change to any public xrspatial.* function. numpy / cupy / dask+numpy / dask+cupy paths are untouched.

Test plan

• test_canvas_like_removed_3451_3640 regression test added — asserts canvas_like stays absent and import xrspatial.utils does not load datashader.
• pytest xrspatial/tests/test_utils.py test_terrain.py test_zonal.py -k "not dask" — 139 passed, 42 skipped (dask backends skipped because dask is not installed in this env; unrelated to the change).
• xrspatial/ source has zero datashader imports (only an explanatory docstring comment in zonal.py).
• All 11 targeted notebooks are valid JSON with zero datashader imports; animated_hillshade.py compiles.
• Notebook output cells were cleared (matplotlib/nbconvert could not be installed in this env to re-execute). Source cells are run-ready; CI or a maintainer can re-execute to regenerate preview images.

Deliberate datashader references kept

• The A/B benchmark in benchmarks/rasterizer_benchmarks.py (import-guarded via HAS_DATASHADER).
• Historical-credit prose in docs/source/getting_started/why_xarray_spatial.rst and docs/source/index.rst.

Notes

• CHANGELOG.md was not touched (per rockout rules; maintained separately at release time).
• No new public API; no behavior changes.

[Melissari1997]

PR Review: Decouple datashader from xarray-spatial

Reviewed the diff against main (18 files). The change removes datashader from library source, all 9 user-guide notebooks, 3 example notebooks/script, and the examples packaging extra, with no public API changes. The approach (substitute template.xrs.rasterize(...) for vector rasterization, matplotlib _shade/_stack/_show helpers for rendering, xarray .coarsen/.interp for raster resampling) is sound and matches the consolidated plan.

Blockers (must fix before merge)

None.

Suggestions (should fix, not blocking)

• examples/animated_hillshade.py:115-124 — _composite has a dead alpha branch. Line 120 reassigns base to an RGB array (top[..., :3] * a + base[..., :3] * (1 - a), shape (H, W, 3)), so the guard on line 121 (base.shape[-1] == 4) is always false and the output alpha is hard-coded to 255. The visual output happens to look right because the terrain/water base layers are fully opaque, but the compositing is incorrect for any layer stack with a transparent base. The notebook _stack helper (e.g. focal.ipynb, classification.ipynb) does this correctly by keeping rgb separate and rebuilding base with np.concatenate([rgb, out_a[..., None]], axis=-1). Recommend aligning _composite with that implementation so the alpha math is right and the dead branch goes away.

• examples/pharmacy-deserts.ipynb cell that rasterizes block groups — the comment says merge='mean' is "equivalent to datashader's ds.mean('PCT_ABOVE_65') here" because block groups tile the plane. That equivalence holds only if polygons do not overlap and each pixel is covered by exactly one polygon. Census block groups are designed to tile, so this is correct in practice, but worth a one-line note that the equivalence relies on non-overlap (the current comment is close but could state the assumption more sharply). No code change needed if the comment is tightened.

Nits (optional improvements)

• xrspatial/zonal.py docstring example — ys = np.linspace(y_range[0], y_range[1], height) uses ascending y. The original ds.Canvas y-orientation differs, but the actual_min_pixels count assertion is invariant to a y-flip so the doctest still passes if it were run. No action required; flagging only so a future reader doesn't assume the y-direction was carefully matched.

• The _shade helper is duplicated verbatim across 8 notebooks. This is acceptable for self-contained notebooks, but if a shared notebook helper module ever exists, these would be good candidates to centralize. Not blocking.

• examples/animated_hillshade.py:7 docstring says "Requires: numpy, xarray, matplotlib, pillow." — also needs xrspatial (imported on line 16). Trivial docstring fix.

What looks good

• canvas_like() removal is clean — it had no callers in xrspatial/ and was never exported, and the regression test (test_canvas_like_removed_3451_3640) guards against re-introduction, including checking that import xrspatial.utils does not load datashader via sys.modules.
• The packaging change correctly adds matplotlib, geopandas, shapely to the examples extra (the deps the migrated notebooks actually use) rather than just dropping datashader.
• Deliberate datashader references (the A/B benchmark in benchmarks/rasterizer_benchmarks.py and the historical-credit prose) were correctly preserved.
• The merge= mapping for the multi-reducer cases is correct: ds.min→merge='min', ds.sum→merge='sum', ds.max→merge='max'. The ds.mean case is correctly flagged as unsupported and handled with merge='last' on non-overlapping polygons.
• All 11 notebooks are valid JSON with zero datashader imports; all Python files compile.

Checklist

• N/A — no algorithm changes (packaging/docs/internal cleanup only)
• N/A — no backend code changed; numpy/cupy/dask paths untouched
• N/A — no neighborhood operations changed
• Regression test added for the one source removal (canvas_like)
• N/A — no dask chunk changes
• No premature materialization introduced (notebook helpers operate on numpy arrays)
• N/A — no new benchmarkable function
• N/A — no new public function (README matrix unchanged, correctly)
• Docstring rewritten for suggest_zonal_canvas; comment reworded in terrain.py

[Melissari1997]

Follow-up review (after review-driven fixes)

Commit c2a85c85 addresses the two Suggestions and the actionable Nit from the first review.

Disposition of original findings

• Suggestion 1 — animated_hillshade.py _composite dead alpha branch: fixed. Rewrote _composite to keep rgb separate and rebuild base with np.concatenate([rgb, out_a], axis=-1), matching the notebook _stack helper. The output alpha is now tracked correctly through the layer stack instead of being hard-coded to 255.

• Suggestion 2 — pharmacy-deserts.ipynb merge='last' comment: fixed. Sharpened the comment to state the non-overlap assumption ("each pixel is covered by exactly one polygon") that makes merge='last' equivalent to ds.mean for the census block groups.

• Nit — animated_hillshade.py docstring Requires: line: fixed. Added xrspatial to the requires list.

• Nit — duplicated _shade helper across 8 notebooks: dismissed. Self-contained notebooks are the project convention (each user-guide notebook runs standalone without importing a shared helper module). Centralizing would require a new support module shipped with the docs, which is out of scope for this decoupling effort.

• Nit — ys linspace direction in the zonal.py docstring: dismissed. The actual_min_pixels count assertion is invariant to a y-flip, and the accessor handles orientation internally. No correctness impact.

What looks good now

• All review findings are either fixed or dismissed with a recorded reason.
• animated_hillshade.py compiles; pharmacy-deserts.ipynb is valid JSON.
• No new findings from the follow-up diff.