For a given 2D mask, I would like to draw its outline with no internal lines between adjacent grid cells. Highlighting one cell is straightforward: How to mark cells in matplotlib.pyplot.imshow (drawing cell borders) and highlighting many cells is also straightforward: Python matplotlib - add borders to grid plot based on value. However, I do not want internal boundaries within the mask to be separated by a line, so the second link above is not suitable.
Here is example code to generate a mask:
import matplotlib.pyplot as plt
import numpy as np
N = 50
xbound = np.linspace(0, 10, N + 1)
ybound = np.linspace(0, 10, N + 1)
x = (xbound[:-1] + xbound[1:]) / 2
y = (ybound[:-1] + ybound[1:]) / 2
X, Y = np.meshgrid(x, y)
Z = np.exp(-((X - 2.5)**2 + (Y - 2.5)**2)) + np.exp(-2 * ((X - 7.5)**2 + (Y - 6.5)**2))
mask = Z > 0.2
plt.imshow(mask, origin='lower', extent=(0, 10, 0, 10))
Which generates the image:
I want the boundary to be around each of the pixelated yellow circles, and to follow the same edges as show above (i.e. be parallel to the x/y-axis) - I want to emphasize that the underlying data is gridded. Using:
plt.contour(x, y, mask, levels=[0.5])
comes close, but the contour is at 45° along the staircase edges.
Bonus points if the outline can be filled or shown using cartopy.
This can be done using shapely
shapes and the union operation. The idea is to create a square (or rectangle etc.) for each cell where the mask is true, then combine contiguous squares together. These geometries can then be plotted or used to create matplotlib
polygons. The following shows the outline in red (with the staircase edge as required):
import shapely.geometry
import shapely.ops
geoms = []
for yidx, xidx in zip(*np.where(mask)):
geoms.append(shapely.geometry.box(xbound[xidx], ybound[yidx], xbound[xidx+1], ybound[yidx+1]))
full_geom = shapely.ops.unary_union(geoms)
for geom in full_geom.geoms:
plt.plot(*geom.exterior.xy, linewidth=4, color='r')
plt.imshow(mask, origin='lower', extent=(0, 10, 0, 10))
It would be hard to fill each shape with a colour/pattern, so instead matplotlib
polygons can be used:
from matplotlib.patches import Polygon
polygons = []
for geom in full_geom.geoms:
p = Polygon(list(zip(*geom.exterior.xy)), closed=True, facecolor='none', edgecolor='r', linestyle='-', hatch='//', linewidth=4)
polygons.append(p)
plt.figure()
axes = plt.gca()
for p in polygons:
axes.add_patch(p)
plt.imshow(mask, origin='lower', extent=(0, 10, 0, 10))
Finally, the raw shapely
polygons can be added directly to a cartopy map:
import cartopy.crs as ccrs
fig, ax = plt.subplots(subplot_kw={'projection': ccrs.PlateCarree()})
ax.add_geometries(
full_geom.geoms,
crs=ccrs.PlateCarree(),
facecolor='white',
edgecolor='red',
linewidth=4,
hatch='//',
)
ax.coastlines()
ax.gridlines(draw_labels=True, dms=True, x_inline=False, y_inline=False)
ax.set_xlim((0, 10))
ax.set_ylim((0, 10))