如何从二进制骨架化图像中找到分支点

Question

我使用 Python OpenCV 像这样对图像进行骨架化：

我想找到骨架的分支点

我不知道该怎么做。有什么想法吗？

Answer 1

这个问题已经很老了，但是如果其他人偶然发现这个问题并希望得到一个不依赖于其他包并使用简单形态学操作的答案，您可能会发现以下内容很有帮助。

这个想法只是应用命中或未命中变换来搜索满足分支点条件的像素。骨架中的分支点是连接到三个或四个其他像素的像素。给定一个适当的结构元素列表selems，您可以在一个输出图像中优雅地组合几个命中或未命中转换，如下所示。

import numpy as np
import scipy.ndimage as ndi


branches = np.zeros_like(skeleton, dtype=bool)
for selem in selems:
    branches |= ndi.binary_hit_or_miss(skeleton, selem)

这非常节省空间，因为您直接将每个转换的结果添加到同一个结果数组中。现在的问题是如何创建结构元素列表。一种解决方案如下。

selems = list()
selems.append(np.array([[0, 1, 0], [1, 1, 1], [0, 0, 0]]))
selems.append(np.array([[1, 0, 1], [0, 1, 0], [1, 0, 0]]))
selems.append(np.array([[1, 0, 1], [0, 1, 0], [0, 1, 0]]))
selems.append(np.array([[0, 1, 0], [1, 1, 0], [0, 0, 1]]))
selems.append(np.array([[0, 0, 1], [1, 1, 1], [0, 1, 0]]))
selems = [np.rot90(selems[i], k=j) for i in range(5) for j in range(4)]

selems.append(np.array([[0, 1, 0], [1, 1, 1], [0, 1, 0]]))
selems.append(np.array([[1, 0, 1], [0, 1, 0], [1, 0, 1]]))

仅当您还想检测具有四个分支的分支点时才需要最后两行。如果您只对三个分支感兴趣，可以省略最后两行。

完整的解决方案将是

import numpy as np
import scipy.ndimage as ndi


selems = list()
selems.append(np.array([[0, 1, 0], [1, 1, 1], [0, 0, 0]]))
selems.append(np.array([[1, 0, 1], [0, 1, 0], [1, 0, 0]]))
selems.append(np.array([[1, 0, 1], [0, 1, 0], [0, 1, 0]]))
selems.append(np.array([[0, 1, 0], [1, 1, 0], [0, 0, 1]]))
selems.append(np.array([[0, 0, 1], [1, 1, 1], [0, 1, 0]]))
selems = [np.rot90(selems[i], k=j) for i in range(5) for j in range(4)]

selems.append(np.array([[0, 1, 0], [1, 1, 1], [0, 1, 0]]))
selems.append(np.array([[1, 0, 1], [0, 1, 0], [1, 0, 1]]))

branches = np.zeros_like(skeleton, dtype=bool)
for selem in selems:
    branches |= ndi.binary_hit_or_miss(skeleton, selem)

同样，您可以使用以下结构元素列表搜索端点。

selems = list()
selems.append(np.array([[0, 1, 0], [0, 1, 0], [0, 0, 0]]))
selems.append(np.array([[1, 0, 0], [0, 1, 0], [0, 0, 0]]))
selems = [np.rot90(selems[i], k=j) for i in range(2) for j in range(4)]

Answer 2

这里的主要想法是在附近寻找。您可以为每个像素使用一个 8 连接邻域（[1 1 1; 1 1 1; 1 1 1]，中心是正在探索的邻域的像素！）。

在每个分支点，像素的度数将 > 2，而常规像素的度数为 2（即，连接到其邻域中的 2 个像素）。

在您的情况下，您的目标是找到度数 > 3 的“交叉点”。

Answer 3

您可以使用库Skan 来查找分支点。

from skan import skeleton_to_csgraph
from skimage import io, morphology

# loading image
image = io.imread('https://i.stack.imgur.com/FpSt9.jpg')

# make image binary
image_binary = image >= 200

# skeletonize
sk = morphology.skeletonize(image_binary).astype(bool)

# receive a degree matrix
_, _, degrees = skeleton_to_csgraph(sk)

# consider all values larger than two as intersection
intersection_matrix = degrees > 2

这里的结果不是一个像素，因为多个分支在一个点相遇。可以用regionprops 确定这个中心区域，然后如果需要一个像素，就可以找到质心。

Answer 4

输出（覆盖在输入图像上）：

代码：

#!/usr/bin/env python
# -*- coding: utf-8 -*-

"""
Find branch point in example image.
"""

import numpy as np
import matplotlib.pyplot as plt
from mahotas.morph import hitmiss as hit_or_miss
from skimage.morphology import medial_axis as skeletonize
# from scipy.ndimage import binary_hit_or_miss as hit_or_miss

def find_branch_points(skel):
    X=[]
    #cross X
    X0 = np.array([[0, 1, 0],
                   [1, 1, 1],
                   [0, 1, 0]])
    X1 = np.array([[1, 0, 1],
                   [0, 1, 0],
                   [1, 0, 1]])
    X.append(X0)
    X.append(X1)
    #T like
    T=[]
    #T0 contains X0
    T0=np.array([[2, 1, 2],
                 [1, 1, 1],
                 [2, 2, 2]])

    T1=np.array([[1, 2, 1],
                 [2, 1, 2],
                 [1, 2, 2]])  # contains X1

    T2=np.array([[2, 1, 2],
                 [1, 1, 2],
                 [2, 1, 2]])

    T3=np.array([[1, 2, 2],
                 [2, 1, 2],
                 [1, 2, 1]])

    T4=np.array([[2, 2, 2],
                 [1, 1, 1],
                 [2, 1, 2]])

    T5=np.array([[2, 2, 1],
                 [2, 1, 2],
                 [1, 2, 1]])

    T6=np.array([[2, 1, 2],
                 [2, 1, 1],
                 [2, 1, 2]])

    T7=np.array([[1, 2, 1],
                 [2, 1, 2],
                 [2, 2, 1]])
    T.append(T0)
    T.append(T1)
    T.append(T2)
    T.append(T3)
    T.append(T4)
    T.append(T5)
    T.append(T6)
    T.append(T7)
    #Y like
    Y=[]
    Y0=np.array([[1, 0, 1],
                 [0, 1, 0],
                 [2, 1, 2]])

    Y1=np.array([[0, 1, 0],
                 [1, 1, 2],
                 [0, 2, 1]])

    Y2=np.array([[1, 0, 2],
                 [0, 1, 1],
                 [1, 0, 2]])

    Y2=np.array([[1, 0, 2],
                 [0, 1, 1],
                 [1, 0, 2]])

    Y3=np.array([[0, 2, 1],
                 [1, 1, 2],
                 [0, 1, 0]])

    Y4=np.array([[2, 1, 2],
                 [0, 1, 0],
                 [1, 0, 1]])
    Y5=np.rot90(Y3)
    Y6 = np.rot90(Y4)
    Y7 = np.rot90(Y5)
    Y.append(Y0)
    Y.append(Y1)
    Y.append(Y2)
    Y.append(Y3)
    Y.append(Y4)
    Y.append(Y5)
    Y.append(Y6)
    Y.append(Y7)

    bp = np.zeros(skel.shape, dtype=int)
    for x in X:
        bp = bp + hit_or_miss(skel,x)
    for y in Y:
        bp = bp + hit_or_miss(skel,y)
    for t in T:
        bp = bp + hit_or_miss(skel,t)

    return bp

def find_end_points(skel):
    endpoint1=np.array([[0, 0, 0],
                        [0, 1, 0],
                        [2, 1, 2]])

    endpoint2=np.array([[0, 0, 0],
                        [0, 1, 2],
                        [0, 2, 1]])

    endpoint3=np.array([[0, 0, 2],
                        [0, 1, 1],
                        [0, 0, 2]])

    endpoint4=np.array([[0, 2, 1],
                        [0, 1, 2],
                        [0, 0, 0]])

    endpoint5=np.array([[2, 1, 2],
                        [0, 1, 0],
                        [0, 0, 0]])

    endpoint6=np.array([[1, 2, 0],
                        [2, 1, 0],
                        [0, 0, 0]])

    endpoint7=np.array([[2, 0, 0],
                        [1, 1, 0],
                        [2, 0, 0]])

    endpoint8=np.array([[0, 0, 0],
                        [2, 1, 0],
                        [1, 2, 0]])

    ep1=hit_or_miss(skel,endpoint1)
    ep2=hit_or_miss(skel,endpoint2)
    ep3=hit_or_miss(skel,endpoint3)
    ep4=hit_or_miss(skel,endpoint4)
    ep5=hit_or_miss(skel,endpoint5)
    ep6=hit_or_miss(skel,endpoint6)
    ep7=hit_or_miss(skel,endpoint7)
    ep8=hit_or_miss(skel,endpoint8)
    ep = ep1+ep2+ep3+ep4+ep5+ep6+ep7+ep8
    return ep

# --------------------------------------------------------------------------------
# script
# --------------------------------------------------------------------------------

img = plt.imread("test.jpg")

# for some reason (screenshot?), example image is not binary
binary = img >= 250

skeleton = skeletonize(binary)

branch_pts = find_branch_points(skeleton)

plt.imshow(binary + branch_pts, cmap='gray', interpolation='nearest'); plt.show()