边缘检测算法总结及其python实现--一阶检测算子

边缘检测的实质是采用某种算法来提取出图像中对象与背景间的交界线。边缘为图像中灰度发生急剧变化的区域边界。图像灰度的变化情况可以用图像灰度分布的梯度来反映，因此可以用局部图像微分技术来获得边缘检测算子。经典的边缘检测方法，是通过对原始图像中像素的某小邻域构造边缘检测算子来达到检测边缘的目的。
一旦计算出导数之后，下一步要做的就是给出一个阈值来确定哪里是边缘位置。阈值越低，能够检测出的边线越多，结果也就越容易受到图片噪声的影响，并且越容易从图像中挑出不相关的特性。与此相反，一个高的阈值将会遗失细的或者短的线段。一个常用的这种方法是带有滞后作用的阈值选择。这个方法使用不同的阈值去寻找边缘。首先使用一个阈值上限去寻找边线开始的地方。一旦找到了一个开始点，我们在图像上逐点跟踪边缘路径，当大于门槛下限时一直纪录边缘位置，直到数值小于下限之后才停止纪录。这种方法假设边缘是连续的界线，并且我们能够跟踪前面所看到的边缘的模糊部分，而不会将图像中的噪声点标记为边缘。
一阶检测算子
Roberts Cross算子，Prewitt算子，Sobel算子， Kirsch算子，Canny算子，罗盘算子
1.Roberts算子
Roberts算子是对图像的对角方向的边缘感兴趣的梯度算子，是一种最简单的算子，利用局部差分算子寻找边缘。它采用对角线方向相邻两像素之差近似梯度幅值检测边缘。

下面分别用python自带的函数和自己按照Roberts算法编的代码，采用的方法是在原图像的第一行和第一列填充0，然后对其进行卷积计算，采用的公式为：

代码：

from skimage import data,color,filters
import  matplotlib.pyplot as plt
import numpy as np
original_img=data.chelsea()
gray_img=color.rgb2gray(original_img)

'''
#using system function
edge_img=filters.roberts(gray_img)
figure=plt.figure()
plt.subplot(131).set_title('original_img')
plt.imshow(original_img)
plt.subplot(132).set_title('gray_img')
plt.imshow(gray_img)
plt.subplot(133).set_title('roberts_img')
plt.imshow(edge_img)
plt.show()
'''

#self code
x_roberts=np.array([[1,0],
                    [0,-1]])
y_roberts=np.array([[0,1],
                    [-1,0]])
h,w=gray_img.shape
img=np.zeros([h+1,w+1])
img[1:h+1,1:w+1]=gray_img[0:h,0:w]
def robert_cal(img,filter):
    h,w=img.shape
    img_filter=np.zeros([h,w])
    for i in range(h-1):
        for j in range(w-1):
            img_filter[i][j]=img[i][j]*filter[0][0]+img[i][j+1]*filter[0][1]+img[i+1][j]*filter[1][0]+img[i+1][j+1]*filter[1][1]
    return img_filter

x_edge_img=robert_cal(img,x_roberts)
y_edge_img=robert_cal(img,y_roberts)
edge_img=np.zeros([h,w])
for i in range(h):
    for j in range(w):
      edge_img[i][j]=np.sqrt(x_edge_img[i][j]**2+y_edge_img[i][j]**2)/(np.sqrt(2))
plt.figure('imgs')
plt.subplot(321).set_title('original_img')
plt.imshow(original_img)
plt.subplot(322).set_title('gray_img')
plt.imshow(gray_img)
plt.subplot(323).set_title('x_edge_img')
plt.imshow(x_edge_img)
plt.subplot(324).set_title('y_edge_img')
plt.imshow(y_edge_img)
plt.subplot(325).set_title('edge_img')
plt.imshow(edge_img)

plt.show()

结果如图所示：
用python自带的函数：

自己写的代码：

2.Prewitt算子
Prewitt算子利用像素点上下、左右邻点的灰度差，在边缘处达到极值检测边缘，去掉部分伪边缘，对噪声具有平滑作用 。其原理是在图像空间利用两个方向模板与图像进行邻域卷积来完成的，这两个方向模板一个检测水平边缘，一个检测垂直边缘。

下面分别用python自带的函数和自己按照prewitt算法编的代码（采用的是在原图像的前两行和前两列填充0，然后对其进行卷积计算）进行了图像的边缘提取，其中自己写代码时，考虑了四种方法：1.用某点的x向梯度和y向梯度的最大值代表该点梯度；2.用x向梯度和y向梯度的和代替该点梯度；3.用x向梯度和y向梯度的平方根代替该点的梯度（同python自带的运算）；4.用x向梯度和y向梯度绝对值的和代替该点的梯度。
代码：

from skimage import data,color,filters
import matplotlib.pyplot as plt
import numpy as np
np.set_printoptions(threshold=100000000)

#input image
original_img=data.chelsea()
gray_img=color.rgb2gray(original_img)

'''
#using prewitt in package
edge_img=filters.prewitt(gray_img)
plt.figure('imgs')
plt.subplot(131).set_title('original_img')
plt.imshow(original_img)
plt.subplot(132).set_title('gray_img')
plt.imshow(gray_img)
plt.subplot(133).set_title('edge_img')
plt.imshow(edge_img)
plt.show()
'''


#self code
x_prewitt=np.array([[1,0,-1],
                   [1,0,-1],
                   [1,0,-1]])
y_prewitt=np.array([[1,1,1],
                   [0,0,0],
                   [-1,-1,-1]])
def conv_calculate(img,filter):
    h,w=img.shape
    conv_img=np.zeros([h-2,w-2])
    for i in range(h-2):
        for j in range(w-2):
            conv_img[i][j]=img[i][j]*filter[0][0]+img[i][j+1]*filter[0][1]+img[i][j+2]*filter[0][2]+\
                    img[i+1][j]*filter[1][0]+img[i+1][j+1]*filter[1][1]+img[i+1][j+2]*filter[1][2]+\
                    img[i+2][j]*filter[2][0]+img[i+2][j+1]*filter[2][1]+img[i+2][j+2]*filter[2][2]

    return conv_img
h,w=gray_img.shape
img=np.zeros([h+2,w+2])
img[2:h+2,2:w+2]=gray_img[0:h]
edge_x_img=conv_calculate(img,x_prewitt)
edge_y_img=conv_calculate(img,y_prewitt)

#p(i,j)=max[edge_x_img,edge_y_img]
edge_img_max=np.zeros([h,w])
for i in range(h):
    for j in range(w):
        if edge_x_img[i][j]>edge_y_img[i][j]:
            edge_img_max=edge_x_img[i][j]
        else:
            edge_img_max=edge_y_img

#p(i,j)=edge_x_img+edge_y_img
edge_img_sum=np.zeros([h,w])
for i in range(h):
    for j in range(w):
        edge_img_sum[i][j]=edge_x_img[i][j]+edge_y_img[i][j]

# p(i,j)=|edge_x_img|+|edge_y_img|
edge_img_abs = np.zeros([h, w])
for i in range(h):
    for j in range(w):
        edge_img_abs[i][j] = abs(edge_x_img[i][j]) + abs(edge_y_img[i][j])


#p(i,j)=sqrt(edge_x_img**2+edge_y_img**2)
edge_img_sqrt=np.zeros([h,w])
for i in range(h):
    for j in range(w):
        edge_img_sqrt[i][j]=np.sqrt((edge_x_img[i][j])**2+(edge_y_img[i][j])**2)

plt.figure('imgs')
plt.subplot(331).set_title('original_img')
plt.imshow(original_img)
plt.subplot(332).set_title('gray_img')
plt.imshow(gray_img)
plt.subplot(333).set_title('x_edge_img')
plt.imshow(edge_x_img)
plt.subplot(334).set_title('y_edge_img')
plt.imshow(edge_y_img)
plt.subplot(335).set_title('edge_img_max')
plt.imshow(edge_img_max)
plt.subplot(336).set_title('edge_img_sum')
plt.imshow(edge_img_max)
plt.subplot(337).set_title('edge_img_sqrt')
plt.imshow(edge_img_sqrt)
plt.subplot(338).set_title('edge_img_abs')
plt.imshow(edge_img_abs)
plt.show()

结果如图所示：
用python自带的函数计算：

自己写的代码：

3.Sobel算子

下面分别用python自带的函数和自己按照Sobel算法编的代码，采用的方法是在原图的首行和首列前填充两行/两列0，在用Sobel算子对其进行卷积计算，采用的公式为：

代码：

from skimage import data,color,filters
import  matplotlib.pyplot as plt
import numpy as np
original_img=data.chelsea()
gray_img=color.rgb2gray(original_img)

'''
#using system function
edge_img=filters.sobel(gray_img)
figure=plt.figure()
plt.subplot(131).set_title('original_img')
plt.imshow(original_img)
plt.subplot(132).set_title('gray_img')
plt.imshow(gray_img)
plt.subplot(133).set_title('sobel_img')
plt.imshow(edge_img)
plt.show()
'''

#self code
x_sobel=np.array([[-1,0,1],
                    [-2,0,2],
                    [-1,0,1]])
y_sobel=np.array([[-1,-2,-1],
                    [0,0,0],
                    [1,2,1]])
h,w=gray_img.shape
img=np.zeros([h+2,w+2])
img[2:h+2,2:w+2]=gray_img[0:h,0:w]
def sobel_cal(img,filter):
    h,w=img.shape
    img_filter=np.zeros([h,w])
    for i in range(h-2):
        for j in range(w-2):
            img_filter[i][j]=img[i][j]*filter[0][0]+img[i][j+1]*filter[0][1]+img[i][j+2]*filter[0][2]+\
                    img[i+1][j]*filter[1][0]+img[i+1][j+1]*filter[1][1]+img[i+1][j+2]*filter[1][2]+\
                    img[i+2][j]*filter[2][0]+img[i+2][j+1]*filter[2][1]+img[i+2][j+2]*filter[2][2]
    return img_filter

x_edge_img=sobel_cal(img,x_sobel)
y_edge_img=sobel_cal(img,y_sobel)
edge_img=np.zeros([h,w])
for i in range(h):
    for j in range(w):
      edge_img[i][j]=np.sqrt(x_edge_img[i][j]**2+y_edge_img[i][j]**2)/(np.sqrt(2))
plt.figure('imgs')
plt.subplot(321).set_title('original_img')
plt.imshow(original_img)
plt.subplot(322).set_title('gray_img')
plt.imshow(gray_img)
plt.subplot(323).set_title('x_edge_img')
plt.imshow(x_edge_img)
plt.subplot(324).set_title('y_edge_img')
plt.imshow(y_edge_img)
plt.subplot(325).set_title('edge_img')
plt.imshow(edge_img)

plt.show()

结果如图所示：
python自带的函数：

自己写的代码：

4.Krisch算子
Kirsch边缘算子由八个方向的卷积核构成，这8个模板代表8个方向，对图像上的8个特定边缘方向作出最大响应，运算中取最大值作为图像的边缘输出：

代码：

from skimage import data,color
import matplotlib.pyplot as plt
import numpy as np
original_img=data.chelsea()
gray_img=color.rgb2gray(original_img)

def conv_cal(img,filter):
    h,w=img.shape
    img_filter=np.zeros([h,w])
    for i in range(h-2):
        for j in range(w-2):
            img_filter[i][j]=img[i][j]*filter[0][0]+img[i][j+1]*filter[0][1]+img[i][j+2]*filter[0][2]+\
                    img[i+1][j]*filter[1][0]+img[i+1][j+1]*filter[1][1]+img[i+1][j+2]*filter[1][2]+\
                    img[i+2][j]*filter[2][0]+img[i+2][j+1]*filter[2][1]+img[i+2][j+2]*filter[2][2]
    return img_filter

krisch1=np.array([[5,5,5],
                  [-3,0,-3],
                  [-3,-3,-3]])
krisch2=np.array([[-3,-3,-3],
                  [-3,0,-3],
                  [5,5,5]])
krisch3=np.array([[5,-3,-3],
                  [5,0,-3],
                  [5,-3,-3]])
krisch4=np.array([[-3,-3,5],
                  [-3,0,5],
                  [-3,-3,5]])
krisch5=np.array([[-3,-3,-3],
                  [-3,0,5],
                  [-3,5,5]])
krisch6=np.array([[-3,-3,-3],
                  [5,0,-3],
                  [5,5,-3]])
krisch7=np.array([[-3,5,5],
                  [-3,0,5],
                  [-3,-3,-3]])
krisch8=np.array([[5,5,-3],
                  [5,0,-3],
                  [-3,-3,-3]])

w,h=gray_img.shape
img=np.zeros([w+2,h+2])
img[2:w+2,2:h+2]=gray_img[0:w,0:h]
edge1=conv_cal(img,krisch1)
edge2=conv_cal(img,krisch2)
edge3=conv_cal(img,krisch3)
edge4=conv_cal(img,krisch4)
edge5=conv_cal(img,krisch5)
edge6=conv_cal(img,krisch6)
edge7=conv_cal(img,krisch7)
edge8=conv_cal(img,krisch8)
edge_img=np.zeros([w,h])
for i in range(w):
    for j in range(h):
        edge_img[i][j]=max(list([edge1[i][j],edge2[i][j],edge3[i][j],edge4[i][j],edge5[i][j],edge6[i][j],edge7[i][j],edge8[i][j]]))

plt.figure('imgs')
plt.subplot(431).set_title('edge1')
plt.imshow(edge1)
plt.subplot(432).set_title('edge2')
plt.imshow(edge2)
plt.subplot(433).set_title('edge3')
plt.imshow(edge3)
plt.subplot(434).set_title('edge4')
plt.imshow(edge4)
plt.subplot(435).set_title('edge5')
plt.imshow(edge5)
plt.subplot(436).set_title('edge6')
plt.imshow(edge6)
plt.subplot(437).set_title('edge7')
plt.imshow(edge7)
plt.subplot(438).set_title('edge8')
plt.imshow(edge8)
plt.subplot(439).set_title('edge')
plt.imshow(edge_img)
plt.show()

计算结果：

5.Canny算子
Canny边缘检测算法的流程如下：

之所以把Canny算子归为一阶算子，我想是因为在计算梯度值和方向的过程中使用的是一阶检测算子的原因。
下面的代码直接用了opencv中的函数（会单独写一篇Canny检测边缘的文章）：

from skimage import data,color,filters
import  matplotlib.pyplot as plt
import cv2
original_img=data.clock()


#using system function
edge_img=cv2.Canny(original_img,5,30,apertureSize=3)
figure=plt.figure()
plt.subplot(121).set_title('original_img')
plt.imshow(original_img)
plt.subplot(122).set_title('canny_img')
plt.imshow(edge_img)
plt.show()

结果：

6.罗盘算子
罗盘算子性能超过Canny算子，相比其他算子，它以像素颜色分布代替均值，利用EMD距离找到罗盘直径方向，该方向的直径能最大化圆形窗口中被直径分隔开的两半区域的差异，此方向即为边缘方向，两半圆区域之间的差异度即为边缘幅值，两者构成向量，类似于梯度。
关于罗盘算子的介绍比较少，查找资料后会单独写一篇罗盘算子检测边缘的文章。