Add slope algorithm (#1088)

* Add slope algorithm

* Add test for slope algorithm

src/titiler/core/tests/test_algorithms.py

@@ -168,6 +168,29 @@ def test_hillshade():
     assert out.array[0, 0, 0] is numpy.ma.masked
 
 
+def test_slope():
+    """test slope."""
+    algo = default_algorithms.get("slope")()
+
+    arr = numpy.random.randint(0, 5000, (1, 262, 262), dtype="uint16")
+    img = ImageData(arr)
+    out = algo(img)
+    assert out.array.shape == (1, 256, 256)
+    assert out.array.dtype == "float32"
+
+    arr = numpy.ma.MaskedArray(
+        numpy.random.randint(0, 5000, (1, 262, 262), dtype="uint16"),
+        mask=numpy.zeros((1, 262, 262), dtype="bool"),
+    )
+    arr.mask[0, 0:100, 0:100] = True
+
+    img = ImageData(arr)
+    out = algo(img)
+    assert out.array.shape == (1, 256, 256)
+    assert out.array.dtype == "float32"
+    assert out.array[0, 0, 0] is numpy.ma.masked
+
+
 def test_contours():
     """test contours."""
     algo = default_algorithms.get("contours")()

src/titiler/core/titiler/core/algorithm/__init__.py

@@ -11,12 +11,13 @@
 
 from titiler.core.algorithm.base import AlgorithmMetadata  # noqa
 from titiler.core.algorithm.base import BaseAlgorithm
-from titiler.core.algorithm.dem import Contours, HillShade, TerrainRGB, Terrarium
+from titiler.core.algorithm.dem import Contours, HillShade, Slope, TerrainRGB, Terrarium
 from titiler.core.algorithm.index import NormalizedIndex
 from titiler.core.algorithm.ops import CastToInt, Ceil, Floor
 
 default_algorithms: Dict[str, Type[BaseAlgorithm]] = {
     "hillshade": HillShade,
+    "slope": Slope,
     "contours": Contours,
     "normalizedIndex": NormalizedIndex,
     "terrarium": Terrarium,

src/titiler/core/titiler/core/algorithm/dem.py

@@ -9,7 +9,7 @@
 
 from titiler.core.algorithm.base import BaseAlgorithm
 
-__all__ = ["HillShade", "Contours", "Terrarium", "TerrainRGB"]
+__all__ = ["HillShade", "Slope", "Contours", "Terrarium", "TerrainRGB"]
 
 
 class HillShade(BaseAlgorithm):
@@ -63,6 +63,53 @@ def __call__(self, img: ImageData) -> ImageData:
         )
 
 
+class Slope(BaseAlgorithm):
+    """Slope calculation."""
+
+    title: str = "Slope"
+    description: str = "Calculate degrees of slope from DEM dataset."
+
+    # parameters
+    buffer: int = Field(3, ge=0, le=99, description="Buffer size for edge effects")
+
+    # metadata
+    input_nbands: int = 1
+    output_nbands: int = 1
+    output_dtype: str = "float32"
+
+    def __call__(self, img: ImageData) -> ImageData:
+        """Calculate degrees slope from DEM dataset."""
+        # Get the pixel size from the transform
+        pixel_size_x = abs(img.transform[0])
+        pixel_size_y = abs(img.transform[4])
+
+        x, y = numpy.gradient(img.array[0])
+        dx = x / pixel_size_x
+        dy = y / pixel_size_y
+
+        slope = numpy.arctan(numpy.sqrt(dx * dx + dy * dy)) * (180 / numpy.pi)
+
+        bounds = img.bounds
+        if self.buffer:
+            slope = slope[self.buffer : -self.buffer, self.buffer : -self.buffer]
+
+            window = windows.Window(
+                col_off=self.buffer,
+                row_off=self.buffer,
+                width=slope.shape[1],
+                height=slope.shape[0],
+            )
+            bounds = windows.bounds(window, img.transform)
+
+        return ImageData(
+            slope.astype(self.output_dtype),
+            assets=img.assets,
+            crs=img.crs,
+            bounds=bounds,
+            band_names=["slope"],
+        )
+
+
 class Contours(BaseAlgorithm):
     """Contours.