Reason or Problem
The planar CPU slope kernel (_cpu in xrspatial/slope.py) starts every interior cell with
if np.isnan(data[y, x]):
continueOn DEMs with scattered nodata (masked water, cloud holes, void-filled SRTM edges) the branch outcome depends on the data and the CPU can't predict it. What the skip saves is a few float32 multiplies and one arctan, and that is cheaper than the mispredicts.
A spike removed the branch and timed a 2000x4000 float32 DEM (Gaussian bump plus default_rng(71942).normal(0, 2) noise). With 30% random NaN it went from 69 ms to 33 ms (0.47x). With the left half NaN, a contiguous block the predictor handles well, it was 0.87x. With no NaN at all it was within noise, 0.95x to 1.0x. Results were bitwise identical in every case.
Proposal
Drop the early-out. out is pre-filled with NaN, and the 3x3 Horn stencil's arithmetic already carries any neighbour NaN into p. The only case the branch actually protects is a NaN centre with eight valid neighbours (the #2761 regression), and a select after the arithmetic handles that:
r = np.arctan(p) * 57.29578
ctr = data[y, x]
out[y, x] = r if ctr == ctr else np.nan
No fastmath, no dtype change, no reordering of the neighbour loads. The GPU and geodesic kernels are untouched.
Usage: nothing changes for callers. Same slope() signature, same results.
Value: about 2x on speckled-nodata DEMs for the numpy and dask+numpy backends.
Stakeholders and Impacts
Anyone calling slope() on CPU. The dask+numpy path calls the same kernel per chunk, so it gets the same win.
This should also come with a NaN-bearing asv benchmark (SlopeNaN) so the nodata path has a timing from now on. One catch: get_xr_dataarray(include_nan=True) in benchmarks/benchmarks/common.py only sets the single corner cell [0, 0] to NaN, which sits on the border the kernel never visits. The benchmark has to add its own speckle on top or it measures the NaN-free path again.
Drawbacks
Cells with a NaN centre now run the stencil and an arctan before being masked. That is the cost being traded away, and the numbers above say it is smaller than the branch.
Alternatives
Leave the branch in and accept the mispredict cost. A per-row "any NaN" pre-check was also considered; it adds code and I could not see it beating the plain select.
Unresolved Questions
None.
Additional Notes or Context
Gates before merge: pytest xrspatial/tests/test_slope.py, a direct A/B against the current kernel with np.array_equal(..., equal_nan=True) on the three regimes above plus 3x3, 1xN and int16 inputs, and a test that pins NaN to exactly the centre-NaN cells and their 8-neighbours on a speckled raster.
Reason or Problem
The planar CPU slope kernel (
_cpuinxrspatial/slope.py) starts every interior cell withOn DEMs with scattered nodata (masked water, cloud holes, void-filled SRTM edges) the branch outcome depends on the data and the CPU can't predict it. What the skip saves is a few float32 multiplies and one
arctan, and that is cheaper than the mispredicts.A spike removed the branch and timed a 2000x4000 float32 DEM (Gaussian bump plus
default_rng(71942).normal(0, 2)noise). With 30% random NaN it went from 69 ms to 33 ms (0.47x). With the left half NaN, a contiguous block the predictor handles well, it was 0.87x. With no NaN at all it was within noise, 0.95x to 1.0x. Results were bitwise identical in every case.Proposal
Drop the early-out.
outis pre-filled with NaN, and the 3x3 Horn stencil's arithmetic already carries any neighbour NaN intop. The only case the branch actually protects is a NaN centre with eight valid neighbours (the #2761 regression), and a select after the arithmetic handles that:No fastmath, no dtype change, no reordering of the neighbour loads. The GPU and geodesic kernels are untouched.
Usage: nothing changes for callers. Same
slope()signature, same results.Value: about 2x on speckled-nodata DEMs for the numpy and dask+numpy backends.
Stakeholders and Impacts
Anyone calling
slope()on CPU. The dask+numpy path calls the same kernel per chunk, so it gets the same win.This should also come with a NaN-bearing asv benchmark (
SlopeNaN) so the nodata path has a timing from now on. One catch:get_xr_dataarray(include_nan=True)inbenchmarks/benchmarks/common.pyonly sets the single corner cell[0, 0]to NaN, which sits on the border the kernel never visits. The benchmark has to add its own speckle on top or it measures the NaN-free path again.Drawbacks
Cells with a NaN centre now run the stencil and an
arctanbefore being masked. That is the cost being traded away, and the numbers above say it is smaller than the branch.Alternatives
Leave the branch in and accept the mispredict cost. A per-row "any NaN" pre-check was also considered; it adds code and I could not see it beating the plain select.
Unresolved Questions
None.
Additional Notes or Context
Gates before merge:
pytest xrspatial/tests/test_slope.py, a direct A/B against the current kernel withnp.array_equal(..., equal_nan=True)on the three regimes above plus 3x3, 1xN and int16 inputs, and a test that pins NaN to exactly the centre-NaN cells and their 8-neighbours on a speckled raster.