python - 从光栅图像中提取街道网络

Question

我有一张 512x512 的街道网格图像：

我想为此图像中的每条街道提取折线（大蓝点 = 交叉点，小蓝点 = 沿折线的点）：

我尝试了一些技巧！一个想法是首先skeletonize将街道压缩到 1px 宽的线条：

from skimage import morphology
morphology.skeletonize(streets_data))

不幸的是，这有一些漏洞会破坏街道网络的连通性；我不完全确定为什么，但我猜这是因为某些街道在某些地方窄了 1 像素，而在其他地方则宽了 1 像素。（更新：间隙不是真实的；它们完全是我如何展示骨架的人工制品。有关悲伤的故事，请参阅此评论。骨架连接良好。）

我可以使用 a 来修补这些binary_dilation，但代价是让街道的宽度再次有所变化：

out = morphology.skeletonize(streets_data)
out = morphology.binary_dilation(out, morphology.selem.disk(1))

使用重新连接的网格，我可以运行霍夫变换来查找线段：

import cv2
rho = 1  # distance resolution in pixels of the Hough grid
theta = np.pi / 180  # angular resolution in radians of the Hough grid
threshold = 8  # minimum number of votes (intersections in Hough grid cell)
min_line_length = 10  # minimum number of pixels making up a line
max_line_gap = 2  # maximum gap in pixels between connectable line segments

# Run Hough on edge detected image
# Output "lines" is an array containing endpoints of detected line segments
lines = cv2.HoughLinesP(
    out, rho, theta, threshold, np.array([]),
    min_line_length, max_line_gap
)

line_image = streets_data.copy()
for i, line in enumerate(lines):
    for x1,y1,x2,y2 in line:
        cv2.line(line_image,(x1,y1),(x2,y2), 2, 1)

这会产生一整串重叠的线段，以及一些间隙（查看右侧的 T 交叉点）：

在这一点上，我可以尝试对重叠的线段进行重复数据删除，但我并不清楚这是一条通往解决方案的道路，尤其是考虑到这一差距。

是否有更直接的方法可以访问我正在寻找的折线网络？具体来说，有哪些方法可以：

1. 寻找交叉口（四向交叉口和 T 交叉口）。
2. 将街道缩小到 1 像素宽，允许可能有一些可变宽度。
3. 寻找交叉点之间的折线。

标签： pythonopencvscikit-imagehough-transform

Answer 1

如果你想提高你的“骨架化”，你可以尝试下面的算法来获得“1-px 宽的街道”：

import imageio
import numpy as np
from matplotlib import pyplot as plt
from scipy.ndimage import distance_transform_edt
from skimage.segmentation import watershed

# read image
image_rgb = imageio.imread('1mYBD.png')

# convert to binary
image_bin = np.max(image_rgb, axis=2) > 0

# compute the distance transform (only > 0)
distance = distance_transform_edt(image_bin)

# segment the image into "cells" (i.e. the reciprocal of the network)
cells = watershed(distance)

# compute the image gradients
grad_v = np.pad(cells[1:, :] - cells[:-1, :], ((0, 1), (0, 0)))
grad_h = np.pad(cells[:, 1:] - cells[:, :-1], ((0, 0), (0, 1)))

# given that the cells have a constant value,
# only the edges will have non-zero gradient
edges = (abs(grad_v) > 0) + (abs(grad_h) > 0)

# extract points into (x, y) coordinate pairs
pos_v, pos_h = np.nonzero(edges)

# display points on top of image
plt.imshow(image_bin, cmap='gray_r')
plt.scatter(pos_h, pos_v, 1, np.arange(pos_h.size), cmap='Spectral')

该算法适用于“街区”而不是“街道”，看一下cells图像：