Skip to content

Commit 62d3929

Browse files
brendancolbrendancol
authored
Add morphological raster operators (#949) (#958)
* Add morphological raster operators (#949) New module xrspatial/morphology.py with four functions: morph_erode, morph_dilate, morph_opening, morph_closing. All four support numpy, cupy, dask+numpy, and dask+cupy backends via ArrayTypeFunctionMapping. * Add tests for morphological operators (#949) 32 tests covering correctness, NaN propagation, edge cases, boundary modes, idempotence, Dataset support, and all four backends (numpy, cupy, dask+numpy, dask+cupy). * Add API reference docs for morphological operators (#949) * Add user guide notebook for morphological operators (#949) * Add morphological operators to README feature matrix (#949)
1 parent 699bca0 commit 62d3929

File tree

8 files changed

+1246
-0
lines changed

8 files changed

+1246
-0
lines changed

‎README.md‎

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -171,6 +171,17 @@ In the GIS world, rasters are used for representing continuous phenomena (e.g. e
171171

172172
-------
173173

174+
### **Morphological**
175+
176+
| Name | Description | NumPy xr.DataArray | Dask xr.DataArray | CuPy GPU xr.DataArray | Dask GPU xr.DataArray |
177+
|:----------:|:------------|:----------------------:|:--------------------:|:-------------------:|:------:|
178+
| [Erode](xrspatial/morphology.py) | Morphological erosion (local minimum over structuring element) | ✅️ | ✅️ | ✅️ | ✅️ |
179+
| [Dilate](xrspatial/morphology.py) | Morphological dilation (local maximum over structuring element) | ✅️ | ✅️ | ✅️ | ✅️ |
180+
| [Opening](xrspatial/morphology.py) | Erosion then dilation (removes small bright features) | ✅️ | ✅️ | ✅️ | ✅️ |
181+
| [Closing](xrspatial/morphology.py) | Dilation then erosion (fills small dark gaps) | ✅️ | ✅️ | ✅️ | ✅️ |
182+
183+
-------
184+
174185
### **Fire**
175186

176187
| Name | Description | NumPy xr.DataArray | Dask xr.DataArray | CuPy GPU xr.DataArray | Dask GPU xr.DataArray |

‎docs/source/reference/index.rst‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -15,6 +15,7 @@ Reference
1515
focal
1616
hydrology
1717
interpolation
18+
morphology
1819
multispectral
1920
pathfinding
2021
proximity

‎docs/source/reference/morphology.rst‎

Lines changed: 40 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,40 @@
1+
.. _reference.morphology:
2+
3+
**********
4+
Morphology
5+
**********
6+
7+
Erode
8+
=====
9+
.. autosummary::
10+
:toctree: _autosummary
11+
12+
xrspatial.morphology.morph_erode
13+
14+
Dilate
15+
======
16+
.. autosummary::
17+
:toctree: _autosummary
18+
19+
xrspatial.morphology.morph_dilate
20+
21+
Opening
22+
=======
23+
.. autosummary::
24+
:toctree: _autosummary
25+
26+
xrspatial.morphology.morph_opening
27+
28+
Closing
29+
=======
30+
.. autosummary::
31+
:toctree: _autosummary
32+
33+
xrspatial.morphology.morph_closing
34+
35+
Kernel Construction
36+
===================
37+
.. autosummary::
38+
:toctree: _autosummary
39+
40+
xrspatial.morphology._circle_kernel

‎examples/user_guide/17_Morphological_Operators.ipynb‎

Lines changed: 217 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,217 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "markdown",
5+
"metadata": {},
6+
"source": [
7+
"# Morphological Operators\n",
8+
"\n",
9+
"This notebook demonstrates the four morphological raster operators\n",
10+
"in `xrspatial.morphology`:\n",
11+
"\n",
12+
"- **`morph_erode`** -- local minimum (shrinks bright regions)\n",
13+
"- **`morph_dilate`** -- local maximum (expands bright regions)\n",
14+
"- **`morph_opening`** -- erosion then dilation (removes small bright features)\n",
15+
"- **`morph_closing`** -- dilation then erosion (fills small dark gaps)\n",
16+
"\n",
17+
"These operations work on any numeric raster, including grayscale\n",
18+
"elevation data and binary classification masks. All four backends\n",
19+
"(numpy, cupy, dask+numpy, dask+cupy) are supported."
20+
]
21+
},
22+
{
23+
"cell_type": "code",
24+
"execution_count": null,
25+
"metadata": {},
26+
"outputs": [],
27+
"source": [
28+
"import numpy as np\n",
29+
"import xarray as xr\n",
30+
"import matplotlib.pyplot as plt\n",
31+
"\n",
32+
"from xrspatial.morphology import (\n",
33+
" _circle_kernel,\n",
34+
" morph_closing,\n",
35+
" morph_dilate,\n",
36+
" morph_erode,\n",
37+
" morph_opening,\n",
38+
")"
39+
]
40+
},
41+
{
42+
"cell_type": "markdown",
43+
"metadata": {},
44+
"source": [
45+
"## 1. Synthetic elevation data\n",
46+
"\n",
47+
"Create a small raster with a peak and a pit to show how erosion and\n",
48+
"dilation behave."
49+
]
50+
},
51+
{
52+
"cell_type": "code",
53+
"execution_count": null,
54+
"metadata": {},
55+
"outputs": [],
56+
"source": [
57+
"rng = np.random.default_rng(42)\n",
58+
"y, x = np.mgrid[-3:3:0.05, -3:3:0.05]\n",
59+
"z = np.exp(-(x**2 + y**2)) - 0.5 * np.exp(-((x - 1.5)**2 + (y - 1)**2) / 0.3)\n",
60+
"z += 0.05 * rng.standard_normal(z.shape) # add noise\n",
61+
"\n",
62+
"raster = xr.DataArray(z, dims=['y', 'x'], name='elevation')\n",
63+
"\n",
64+
"fig, ax = plt.subplots(figsize=(6, 5))\n",
65+
"raster.plot(ax=ax, cmap='terrain')\n",
66+
"ax.set_title('Original raster')\n",
67+
"plt.tight_layout()\n",
68+
"plt.show()"
69+
]
70+
},
71+
{
72+
"cell_type": "markdown",
73+
"metadata": {},
74+
"source": [
75+
"## 2. Erosion and dilation"
76+
]
77+
},
78+
{
79+
"cell_type": "code",
80+
"execution_count": null,
81+
"metadata": {},
82+
"outputs": [],
83+
"source": [
84+
"kernel = np.ones((5, 5), dtype=np.uint8)\n",
85+
"\n",
86+
"eroded = morph_erode(raster, kernel=kernel, boundary='nearest')\n",
87+
"dilated = morph_dilate(raster, kernel=kernel, boundary='nearest')\n",
88+
"\n",
89+
"fig, axes = plt.subplots(1, 3, figsize=(15, 4))\n",
90+
"for ax, data, title in zip(\n",
91+
" axes,\n",
92+
" [raster, eroded, dilated],\n",
93+
" ['Original', 'Eroded (local min)', 'Dilated (local max)'],\n",
94+
"):\n",
95+
" data.plot(ax=ax, cmap='terrain', vmin=z.min(), vmax=z.max())\n",
96+
" ax.set_title(title)\n",
97+
"plt.tight_layout()\n",
98+
"plt.show()"
99+
]
100+
},
101+
{
102+
"cell_type": "markdown",
103+
"metadata": {},
104+
"source": [
105+
"## 3. Opening and closing\n",
106+
"\n",
107+
"Opening removes small bright spikes; closing fills small dark pits."
108+
]
109+
},
110+
{
111+
"cell_type": "code",
112+
"execution_count": null,
113+
"metadata": {},
114+
"outputs": [],
115+
"source": [
116+
"opened = morph_opening(raster, kernel=kernel, boundary='nearest')\n",
117+
"closed = morph_closing(raster, kernel=kernel, boundary='nearest')\n",
118+
"\n",
119+
"fig, axes = plt.subplots(1, 3, figsize=(15, 4))\n",
120+
"for ax, data, title in zip(\n",
121+
" axes,\n",
122+
" [raster, opened, closed],\n",
123+
" ['Original', 'Opening', 'Closing'],\n",
124+
"):\n",
125+
" data.plot(ax=ax, cmap='terrain', vmin=z.min(), vmax=z.max())\n",
126+
" ax.set_title(title)\n",
127+
"plt.tight_layout()\n",
128+
"plt.show()"
129+
]
130+
},
131+
{
132+
"cell_type": "markdown",
133+
"metadata": {},
134+
"source": [
135+
"## 4. Binary mask cleanup\n",
136+
"\n",
137+
"A common use case: clean up noisy classification results. Opening\n",
138+
"removes salt noise; closing removes pepper noise."
139+
]
140+
},
141+
{
142+
"cell_type": "code",
143+
"execution_count": null,
144+
"metadata": {},
145+
"outputs": [],
146+
"source": [
147+
"# Generate a binary mask with noise\n",
148+
"mask = np.zeros((100, 100), dtype=np.float64)\n",
149+
"mask[20:80, 20:80] = 1.0 # large square\n",
150+
"\n",
151+
"# Sprinkle salt-and-pepper noise\n",
152+
"noise = rng.random(mask.shape)\n",
153+
"mask[noise < 0.02] = 1.0 # salt\n",
154+
"mask[noise > 0.98] = 0.0 # pepper\n",
155+
"\n",
156+
"noisy = xr.DataArray(mask, dims=['y', 'x'], name='mask')\n",
157+
"cleaned = morph_opening(\n",
158+
" morph_closing(noisy, kernel=kernel, boundary='nearest'),\n",
159+
" kernel=kernel,\n",
160+
" boundary='nearest',\n",
161+
")\n",
162+
"\n",
163+
"fig, axes = plt.subplots(1, 2, figsize=(10, 4))\n",
164+
"noisy.plot(ax=axes[0], cmap='gray')\n",
165+
"axes[0].set_title('Noisy mask')\n",
166+
"cleaned.plot(ax=axes[1], cmap='gray')\n",
167+
"axes[1].set_title('After closing + opening')\n",
168+
"plt.tight_layout()\n",
169+
"plt.show()"
170+
]
171+
},
172+
{
173+
"cell_type": "markdown",
174+
"metadata": {},
175+
"source": [
176+
"## 5. Circular structuring element\n",
177+
"\n",
178+
"Use `_circle_kernel(radius)` to create a disk-shaped kernel\n",
179+
"instead of a square."
180+
]
181+
},
182+
{
183+
"cell_type": "code",
184+
"execution_count": null,
185+
"metadata": {},
186+
"outputs": [],
187+
"source": [
188+
"disk = _circle_kernel(3)\n",
189+
"print('Circular kernel (radius=3):')\n",
190+
"print(disk)\n",
191+
"\n",
192+
"eroded_disk = morph_erode(raster, kernel=disk, boundary='nearest')\n",
193+
"\n",
194+
"fig, axes = plt.subplots(1, 2, figsize=(10, 4))\n",
195+
"raster.plot(ax=axes[0], cmap='terrain')\n",
196+
"axes[0].set_title('Original')\n",
197+
"eroded_disk.plot(ax=axes[1], cmap='terrain')\n",
198+
"axes[1].set_title('Eroded with circular kernel (r=3)')\n",
199+
"plt.tight_layout()\n",
200+
"plt.show()"
201+
]
202+
}
203+
],
204+
"metadata": {
205+
"kernelspec": {
206+
"display_name": "Python 3",
207+
"language": "python",
208+
"name": "python3"
209+
},
210+
"language_info": {
211+
"name": "python",
212+
"version": "3.11.0"
213+
}
214+
},
215+
"nbformat": 4,
216+
"nbformat_minor": 4
217+
}

‎xrspatial/__init__.py‎

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -43,6 +43,10 @@
4343
from xrspatial.flow_length import flow_length # noqa
4444
from xrspatial.flow_path import flow_path # noqa
4545
from xrspatial.focal import mean # noqa
46+
from xrspatial.morphology import morph_closing # noqa
47+
from xrspatial.morphology import morph_dilate # noqa
48+
from xrspatial.morphology import morph_erode # noqa
49+
from xrspatial.morphology import morph_opening # noqa
4650
from xrspatial.hand import hand # noqa
4751
from xrspatial.hillshade import hillshade # noqa
4852
from xrspatial.mahalanobis import mahalanobis # noqa

‎xrspatial/accessor.py‎

Lines changed: 36 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -181,6 +181,24 @@ def focal_mean(self, **kwargs):
181181
from .focal import mean
182182
return mean(self._obj, **kwargs)
183183

184+
# ---- Morphological ----
185+
186+
def morph_erode(self, **kwargs):
187+
from .morphology import morph_erode
188+
return morph_erode(self._obj, **kwargs)
189+
190+
def morph_dilate(self, **kwargs):
191+
from .morphology import morph_dilate
192+
return morph_dilate(self._obj, **kwargs)
193+
194+
def morph_opening(self, **kwargs):
195+
from .morphology import morph_opening
196+
return morph_opening(self._obj, **kwargs)
197+
198+
def morph_closing(self, **kwargs):
199+
from .morphology import morph_closing
200+
return morph_closing(self._obj, **kwargs)
201+
184202
# ---- Proximity / Distance ----
185203

186204
def proximity(self, **kwargs):
@@ -517,6 +535,24 @@ def focal_mean(self, **kwargs):
517535
from .focal import mean
518536
return mean(self._obj, **kwargs)
519537

538+
# ---- Morphological ----
539+
540+
def morph_erode(self, **kwargs):
541+
from .morphology import morph_erode
542+
return morph_erode(self._obj, **kwargs)
543+
544+
def morph_dilate(self, **kwargs):
545+
from .morphology import morph_dilate
546+
return morph_dilate(self._obj, **kwargs)
547+
548+
def morph_opening(self, **kwargs):
549+
from .morphology import morph_opening
550+
return morph_opening(self._obj, **kwargs)
551+
552+
def morph_closing(self, **kwargs):
553+
from .morphology import morph_closing
554+
return morph_closing(self._obj, **kwargs)
555+
520556
# ---- Diffusion ----
521557

522558
def diffuse(self, **kwargs):

0 commit comments

Comments
 (0)