python - 找到 map 中每个像素的最近邻居

Question

我试图在 python 中找到一个简单的方法，其中对于 2dim 掩码中的每个像素，我可以获得最近的非零邻居的索引。在 Matlab 中有一个 bwdist 可以准确返回。 例如:如果我的输入是:

array [[0 0 0 0 0 0 0]
       [0 1 0 0 0 0 0]
       [0 0 0 0 0 1 0]
       [0 0 0 0 0 0 0]]

我的输出应该是:

array [[(1,1) (1,1) (1,1) (1,1) (2,5) (2,5) (2,5)]
       [(1,1) (1,1) (1,1) (1,1) (2,5) (2,5) (2,5)]
       [(1,1) (1,1) (1,1) (2,5) (2,5) (2,5) (2,5)]
       [(1,1) (1,1) (1,1) (2,5) (2,5) (2,5) (2,5)]]

该函数还可以像 Matlab 中的 bwdist 一样返回绝对索引(对于 1dim 数组)。

谢谢!

编辑:到目前为止，我已经尝试了一些与 scipy 相关的潜在解决方案，例如 distance_transform_edt，但它只能找到到最近像素的距离，而不是像素本身。 如果相关的话，我还在代码的其他地方使用 OpenCV 和 VLfeat。

Answer 1

来自文档:

OpenCV 有 distanceTransform() 和 distanceTransformWithLabels() 函数，它们的工作方式类似，但与这个 Matlab 函数有一些不同。来自Matlab docs for bwdist :

D = bwdist(BW) computes the Euclidean distance transform of the binary image BW. For each pixel in BW, the distance transform assigns a number that is the distance between that pixel and the nearest nonzero pixel of BW.

将此与 OpenCV docs for distanceTransformWithLabels() 进行比较:

Calculates the distance to the closest zero pixel for each pixel of the source image.

所以Matlab给出了到最近的非零像素的距离，而OpenCV给出了到最近的零像素的距离。所以你需要为 OpenCV 反转图像。此外，带有标签的 Matlab 的可选输出给出了对应于最近像素的线性索引:

[D,idx] = bwdist(BW) also computes the closest-pixel map in the form of an index array, idx. Each element of idx contains the linear index of the nearest nonzero pixel of BW. The closest-pixel map is also called the feature map, feature transform, or nearest-neighbor transform.

使用 OpenCV，获得输出的标签不是图像的坐标，也不是索引。相反，它只是一个数字标签，类似于连通分量标签，与像素位置/索引完全无关。

This variant of the function does not only compute the minimum distance for each pixel (x,y) but also identifies the nearest connected component consisting of zero pixels (labelType==DIST_LABEL_CCOMP) or the nearest zero pixel (labelType==DIST_LABEL_PIXEL).

这意味着您必须使用这个标记图像来屏蔽您的输入并找到与该标签对应的像素(据我所知，至少这是最好的方法)。

解决方案:

因此，为了让我们了解如何到达我们想要的位置，让我们看看这个函数将我们带到哪里(如前所述，使用倒置图像作为输入):

In [138]: img
Out[138]:
array([[  0,   0,   0,   0,   0,   0,   0],
       [  0, 255,   0,   0,   0,   0,   0],
       [  0,   0,   0,   0,   0, 255,   0],
       [  0,   0,   0,   0,   0,   0,   0]], dtype=uint8)

In [139]: dist, labels = cv2.distanceTransformWithLabels(~a, distanceType=cv2.DIST_L2, maskSize=3)

In [140]: print(dist)
[[1.3999939 1.        1.3999939 2.1968994 2.1968994 2.        2.1968994]
 [1.        0.        1.        2.        1.3999939 1.        1.3999939]
 [1.3999939 1.        1.3999939 2.        1.        0.        1.       ]
 [2.1968994 2.        2.1968994 2.1968994 1.3999939 1.        1.3999939]]

In [141]: print(labels)
[[1 1 1 1 2 2 2]
 [1 1 1 1 2 2 2]
 [1 1 1 2 2 2 2]
 [1 1 1 2 2 2 2]]

好吧，如果我们只是遍历标签中的唯一值，为它们中的每一个创建一个 mask ， mask 原始图像......然后找到该标记区域内的白色像素，我们将有索引:

In [146]: for l in np.unique(labels):
     ...:     mask = label == l
     ...:     i = np.where(img * mask)
     ...:     print(i)
     ...:
(array([1]), array([1]))
(array([2]), array([5]))

这不是您要求的确切输出，但它是一个索引列表，并且您有标签。所以现在我们只需要映射这些。我要做的是创建一个空的双 channel 矩阵来保存索引值，然后根据标签中的掩码填充它:

In [177]: index_img = np.zeros((*img.shape, 2), dtype=np.intp)

In [178]: for l in np.unique(labels):
     ...:     mask = label == l
     ...:     index_img[mask] = np.dstack(np.where(img * mask))

这是一个包含您想要的信息的双 channel 数组。结构略有不同(不为每个条目使用元组)，但通常是其他 OpenCV 函数(双 channel 数组)所需的结构:

In [204]: index_img[:, :, 0]
Out[204]:
array([[1, 1, 1, 1, 2, 2, 2],
       [1, 1, 1, 1, 2, 2, 2],
       [1, 1, 1, 2, 2, 2, 2],
       [1, 1, 1, 2, 2, 2, 2]])

In [205]: index_img[:, :, 1]
Out[205]:
array([[1, 1, 1, 1, 5, 5, 5],
       [1, 1, 1, 1, 5, 5, 5],
       [1, 1, 1, 5, 5, 5, 5],
       [1, 1, 1, 5, 5, 5, 5]])

综合考虑

这是一个执行此操作的函数，并且可以选择输出这两个 channel 输出或像 Matlab 那样输出线性输出:

def bwdist(img, metric=cv2.DIST_L2, dist_mask=cv2.DIST_MASK_5, label_type=cv2.DIST_LABEL_CCOMP, ravel=True):
    """Mimics Matlab's bwdist function.

    Available metrics:
        https://docs.opencv.org/3.4/d7/d1b/group__imgproc__misc.html#gaa2bfbebbc5c320526897996aafa1d8eb
    Available distance masks:
        https://docs.opencv.org/3.4/d7/d1b/group__imgproc__misc.html#gaaa68392323ccf7fad87570e41259b497
    Available label types:
        https://docs.opencv.org/3.4/d7/d1b/group__imgproc__misc.html#ga3fe343d63844c40318ee627bd1c1c42f
    """
    flip = cv2.threshold(img, 0, 255, cv2.THRESH_BINARY_INV)[1]
    dist, labeled = cv2.distanceTransformWithLabels(flip, metric, dist_mask)

    # return linear indices if ravel == True (default)
    if ravel:  
        idx = np.zeros(img.shape, dtype=np.intp)  # np.intp type is for indices
        for l in np.unique(labeled):
            mask = labeled == l
            idx[mask] = np.flatnonzero(img * mask)
        return dist, idx

    # return two-channel indices if ravel == False
    idx = np.zeros((*img.shape, 2), dtype=np.intp)  
    for l in np.unique(labeled):
        mask = labeled == l
        idx[mask] = np.dstack(np.where(img * mask))
    return dist, idx

以及 Matlab 在文档中给出的示例:

In [241]: bw = np.zeros((5, 5), dtype=np.uint8)
     ...: bw[1, 1] = 1
     ...: bw[3, 3] = 1
     ...: print(bw)
     ...:
[[0 0 0 0 0]
 [0 1 0 0 0]
 [0 0 0 0 0]
 [0 0 0 1 0]
 [0 0 0 0 0]]

In [244]: d, idx = bwdist(bw)

In [245]: print(d)
[[1.3999939 1.        1.3999939 2.1968994 3.1968994]
 [1.        0.        1.        2.        2.1968994]
 [1.3999939 1.        1.3999939 1.        1.3999939]
 [2.1968994 2.        1.        0.        1.       ]
 [3.1968994 2.1968994 1.3999939 1.        1.3999939]]

In [246]: print(idx)
[[ 6  6  6  6 18]
 [ 6  6  6  6 18]
 [ 6  6  6 18 18]
 [ 6  6 18 18 18]
 [ 6 18 18 18 18]]