Canny边缘检测

边缘检测一般步骤

1. 使用高斯滤波
2. 在x和y方向上使用Sobel滤波器，在此之上求出边缘的强度和边缘的梯度
3. 对梯度幅值进行非极大值抑制来使边缘变得更细
4. 使用滞后阈值来对阈值进行处理

def Canny_step1(img):#Canny边缘检测第一步
    # Gray scale
    def BGR2GRAY(img):#图像转灰度
        b = img[:, :, 0].copy()
        g = img[:, :, 1].copy()
        r = img[:, :, 2].copy()
        # Gray scale
        out = 0.2126 * r + 0.7152 * g + 0.0722 * b#将b、g、r分量组合成灰度，三通道合一
        out = out.astype(np.uint8)
        return out#将灰度图像返回
    # Gaussian filter for grayscale，对灰度图像进行高斯滤波
    def gaussian_filter(img, K_size=3, sigma=1.3):#设置高斯核大小，标准差sigma的值
        if len(img.shape) == 3:
            H, W, C = img.shape
            gray = False
        else:
            img = np.expand_dims(img, axis=-1)
            H, W, C = img.shape
            gray = True
        ## Zero padding
        pad = K_size // 2   #对边界进行补充
        out = np.zeros([H + pad * 2, W + pad * 2, C], dtype=np.float)
        out[pad: pad + H, pad: pad + W] = img.copy().astype(np.float)
        ## prepare Kernel
        K = np.zeros((K_size, K_size), dtype=np.float) #生成高斯核
        for x in range(-pad, -pad + K_size):
            for y in range(-pad, -pad + K_size):
                K[y + pad, x + pad] = np.exp(- (x ** 2 + y ** 2) / (2 * (sigma ** 2)))
        K /= (2 * np.pi * sigma * sigma)
        K /= K.sum() #将高斯核进行归一化处理
        tmp = out.copy()
        # filtering
        for y in range(H):
            for x in range(W):
                for c in range(C): #对图像进行高斯卷积滤波
                    out[pad + y, pad + x, c] = np.sum(K * tmp[y: y + K_size, x: x + K_size, c])
        out = np.clip(out, 0, 255) #将out值控制在0-255之间
        out = out[pad: pad + H, pad: pad + W] #将out边部切去，生成与原图等大的滤波后的图片
        # out = out.astype(np.uint8)
        if gray:
            out = out[..., 0]
        return out
    # sobel filter sobel滤波
    def sobel_filter(img, K_size=3):
        if len(img.shape) == 3:
            H, W, C = img.shape
        else:
            # img = np.expand_dims(img, axis=-1)
            H, W = img.shape
        # Zero padding
        pad = K_size // 2
        out = np.zeros((H + pad * 2, W + pad * 2), dtype=np.float)
        out[pad: pad + H, pad: pad + W] = img.copy().astype(np.float)
        tmp = out.copy()
        out_v = out.copy()
        out_h = out.copy()
        ## Sobel vertical垂直滤波器
        Kv = [[1., 2., 1.], [0., 0., 0.], [-1., -2., -1.]]
        ## Sobel horizontal水平滤波器
        Kh = [[1., 0., -1.], [2., 0., -2.], [1., 0., -1.]]
        # filtering
        for y in range(H):
            for x in range(W):#分别使用水平sobel滤波器和垂直sobel滤波器对图像进行滤波
                out_v[pad + y, pad + x] = np.sum(Kv * (tmp[y: y + K_size, x: x + K_size]))
                out_h[pad + y, pad + x] = np.sum(Kh * (tmp[y: y + K_size, x: x + K_size]))
        out_v = np.clip(out_v, 0, 255)
        out_h = np.clip(out_h, 0, 255)
        out_v = out_v[pad: pad + H, pad: pad + W].astype(np.uint8)
        out_h = out_h[pad: pad + H, pad: pad + W].astype(np.uint8)
        return out_v, out_h #返回经垂直滤波和水平滤波之后的图像
    def get_edge_angle(fx, fy):#获取角度的函数
        # get edge strength 获取边缘强度
        edge = np.sqrt(np.power(fx, 2) + np.power(fy, 2))
        fx = np.maximum(fx, 1e-5)#np.maximum函数接收两个参数，逐次比较返回较大者，入np.maximum([2,2,3],12)返回[12,12,12]
        #此处为防止fx出现负值或为0，直接到1e-5处将其截断
        # get edge angle 获取边缘角度
        angle = np.arctan(fy / fx)#获取边缘角度
        return edge, angle#返回边缘强度和角度
    def angle_quantization(angle):#角度量化函数，将梯度角度量化
        angle = angle / np.pi * 180#将弧度值转为角度值
        angle[angle < -22.5] = 180 + angle[angle < -22.5]
        _angle = np.zeros_like(angle, dtype=np.uint8)#指定_angle与angle尺寸一样，并将其初始化为0
        _angle[np.where(angle <= 22.5)] = 0#将角度量化至0-45-90-135这四个值上面
        _angle[np.where((angle > 22.5) & (angle <= 67.5))] = 45
        _angle[np.where((angle > 67.5) & (angle <= 112.5))] = 90
        _angle[np.where((angle > 112.5) & (angle <= 157.5))] = 135
        return _angle
    # grayscale
    gray = BGR2GRAY(img)#调用前面的函数进行图像转灰度处理
    # gaussian filtering
    gaussian = gaussian_filter(gray, K_size=5, sigma=1.4)#对灰度图像进行高斯滤波
    # sobel filtering
    fy, fx = sobel_filter(gaussian, K_size=3)#对高斯滤波后的图像进行两个方向上的sobel滤波，得到fy梯度图像及fx梯度图像
    # get edge strength, angle
    edge, angle = get_edge_angle(fx, fy)#将fx、fy传入获取角度的函数，得到每个像素点位置的梯度强度和弧度制角度信息
    # angle quantization梯度方向量化
    angle = angle_quantization(angle)
    return edge, angle#返回梯度强度图像和梯度方向量化之后的图像

第三步 非极大值抑制

• 非极大值抑制是对除去非极大值以外的值的操作的总称（这个术语在其它的任务中也经常出现）。
• 在这里，我们比较我们所关注的地方梯度的法线方向邻接的三个像素点的梯度幅值，如果该点的梯度值不比其它两个像素大，那么这个地方的值设置为0。
• 也就是说，我们在注意梯度幅值edge(x,y)的时候，可以根据下式由梯度方向angle(x,y)来变换edge(x,y)
• canny算子中非最大抑制(Non-maximum suppression)是回答这样一个问题: “当前的梯度值在梯度方向上是一个局部最大值吗?” 所以,要把当前位置的梯度值与梯度方向上两侧的梯度值进行比较.

1. angle(y,x)=0

• 如果在edge(y,x)、edge(y,x-1)、edge(y,x+1)中edge(y,x)不是最大的，那么edge(y,x)=0；

2. angle(y,x)=45

• 如果在edge(y,x)、edge(y+1,x-1)、edge(y-1,x+1)中edge(y,x)不是最大的，那么edge(x,y)=0；

3. angle(y,x)=90

• 如果在edge(y,x)、edge(y−1,x)、edge(y+1,x)中edge(y,x)不是最大的，那么edge(y,x)=0；

4. angle(y,x)=135

• 如果在edge(y,x)、edge(y−1,x-1)、edge(y+1,x+1)中edge(y,x)不是最大的，那么edge(x,y)=0；

#非极大值抑制
    def non_maximum_suppression(angle, edge):#传入梯度强度图像及梯度方向图像
        H, W = angle.shape#获取图像尺寸
        _edge = edge.copy()#复制一个新图像备用
        for y in range(H):
            for x in range(W):#依次检测每个像素点的梯度方向，并赋予dx1 dy1 dx2 dy2不同的值用来做非极大值抑制
                    if angle[y, x] == 0:
                            dx1, dy1, dx2, dy2 = -1, 0, 1, 0
                    elif angle[y, x] == 45:
                            dx1, dy1, dx2, dy2 = -1, 1, 1, -1
                    elif angle[y, x] == 90:
                            dx1, dy1, dx2, dy2 = 0, -1, 0, 1
                    elif angle[y, x] == 135:
                            dx1, dy1, dx2, dy2 = -1, -1, 1, 1
                    if x == 0:#处理图像边界
                            dx1 = max(dx1, 0)
                            dx2 = max(dx2, 0)
                    if x == W-1:#处理图像边界
                            dx1 = min(dx1, 0)
                            dx2 = min(dx2, 0)
                    if y == 0:#处理图像边界
                            dy1 = max(dy1, 0)
                            dy2 = max(dy2, 0)
                    if y == H-1:#处理图像边界
                            dy1 = min(dy1, 0)
                            dy2 = min(dy2, 0)
                    if max(max(edge[y, x], edge[y + dy1, x + dx1]), edge[y + dy2, x + dx2]) != edge[y, x]:#如果所在像素不是参与比较的三个像素中的最大值，就将其像素值置为0
                            _edge[y, x] = 0
        return _edge#返回非极大值抑制后的图像

最后一步，滞后阈值

• 在这里我们将通过设置高阈值（HT：high threshold）和低阈值（LT：low threshold）来将梯度幅值二值化。
• 如果梯度幅值edge(x,y)大于高阈值的话，令edge(x,y)=255；
• 如果梯度幅值edge(x,y)小于低阈值的话，令edge(x,y)=0；
• 如果梯度幅值edge(x,y)介于高阈值和低阈值之间并且周围8邻域内有比高阈 值高的像素点存在，令edge(x,y)=255；
• 在这里，我们使高阈值为100，低阈值为20。顺便说一句，阈值的大小需要边看结果边调整。

#滞后阈值函数
    def hysterisis(edge, HT=100, LT=30):#传入经过非极大值抑制的梯度幅值图像，高阈值和低阈值参数
        H, W = edge.shape
        # Histeresis threshold
        edge[edge >= HT] = 255#对图像进行二值化处理
        edge[edge <= LT] = 0#对图像进行二值化处理
        _edge = np.zeros((H + 2, W + 2), dtype=np.float32)
        _edge[1 : H + 1, 1 : W + 1] = edge
        ## 8 - Nearest neighbor 构建检测八邻域内是否有比中心像素点高的数组
        nn = np.array(((1., 1., 1.), (1., 0., 1.), (1., 1., 1.)), dtype=np.float32)
        for y in range(1, H+2):
                for x in range(1, W+2):
                        if _edge[y, x] < LT or _edge[y, x] > HT:
                                continue#不在范围之内，就跳过
                        if np.max(_edge[y-1:y+2, x-1:x+2] * nn) >= HT:
                                _edge[y, x] = 255#在范围之内，如果被检测到，将像素点置为255，否则置为0
                        else:
                                _edge[y, x] = 0
        edge = _edge[1:H+1, 1:W+1]          
        return edge

总结一下

1. 先进性高斯滤波
2. 在高斯滤波基础上进行x方向和y方向的sobel滤波
3. 将x方向的sobel滤波与y方向上的sobel滤波结合，得到梯度幅度图及梯度方向图
4. 将梯度幅度图和梯度方向图结合起来进行非极大值抑制
5. 对经过了非极大值抑制的图像进行最后一步滞后阈值处理-大于高阈值的置为255，小于低阈值的置为0，在此之间的看像素八邻域之间有没有大于高阈值的，有就置255，否则置0