参考自：2.canny edge detection（canny边缘检测）_哔哩哔哩_bilibili

本篇博客是通过学习上述链接视频，将其内容转成文字化的描述，我强烈推荐观看该视频和他提供的代码。下面是1-4是我思考的一些思路。如果不想看可以直接看5，他是表示canny的完整流程。

1. 图像边缘检测问题和目标

已知条件：

一幅图像。

问题和目标：

找出图像中的边缘。

什么是边缘？

如上图中红线区域，可以理解为是图像变化较大的部分可以认为是边缘。

解决方法：

假设一个连续可导函数在某点的变化情况，是用梯度进行表示。但是图像是一个点阵，是一个离散情况，需要将其转化。

数学上的梯度求法，

其实呢，是可以直接用x和y方向的梯度的，但是这样不够直观，不能直接感受梯度具体是朝哪个方向，变化值是多大。因此就将x和y方向的梯度合并起来，其中梯度的模，也就是幅值是

梯度的方向为

综上所述，只需要根据像素点和如下公式就可以求出梯度幅值和方向。

2. 根据公式转化成实际运算方法

求梯度的方式有很多种，有人使用[-1,0, 1]和[-1, 0, 1]T当做Gx和Gy，进行卷积计算得到相应方向梯度，也就是当前像素的后一个像素减去前一个像素。但常规使用的是一个方阵作为卷积运算模板。例如作为Gy和Gx。

首先实现卷积运算：

def convolve(filter,mat,padding,strides):''':param filter:卷积核，必须为二维(2 x 1也算二维) 否则返回None:param mat:图片:param padding:对齐:param strides:移动步长:return:返回卷积后的图片。(灰度图，彩图都适用)@author:bilibili-会飞的吴克'''result = Nonefilter_size = filter.shapemat_size = mat.shapeif len(filter_size) == 2:if len(mat_size) == 3:channel = []for i in range(mat_size[-1]):pad_mat = np.pad(mat[:,:,i], ((padding[0], padding[1]), (padding[2], padding[3])), 'constant')temp = []for j in range(0,mat_size[0],strides[1]):temp.append([])for k in range(0,mat_size[1],strides[0]):val = (filter*pad_mat[j:j+filter_size[0],k:k+filter_size[1]]).sum()temp[-1].append(val)channel.append(np.array(temp))channel = tuple(channel)result = np.dstack(channel)elif len(mat_size) == 2:channel = []pad_mat = np.pad(mat, ((padding[0], padding[1]), (padding[2], padding[3])), 'constant')for j in range(0, mat_size[0], strides[1]):channel.append([])for k in range(0, mat_size[1], strides[0]):val = (filter * pad_mat[j:j + filter_size[0],k:k + filter_size[1]]).sum()channel[-1].append(val)result = np.array(channel)return result

在使用卷积运算实现横向X轴和竖向Y轴的求导。

img = plt.imread(filename)if i[-4:] == '.png':img = img*255img = img.mean(axis=-1) #this is a way to get a gray image. sobel_kernel_y = np.array([[-1, -2, -1], [0, 0, 0], [1, 2, 1]])sobel_kernel_x = np.array([[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]])img3 = convolve(sobel_kernel_x,img2,[1,1,1,1],[1,1])img4 = convolve(sobel_kernel_y,img2,[1,1,1,1],[1,1])

但这个边缘检测似乎有点粗糙，似乎能达到一定的效果，但是这个边缘是内容太多了，我们需要的就是一条线类似这样的边缘，需要进行再一步处理。

3. 边缘的进一步确定

问题

边缘区域不够明确，产生了很多不合理的区域。不合理区域包含2个部分

（1）红色区域产生了很多的不合理点，我们并不认为他是边缘；

（2）绿色区域是边缘部分，但是他太粗了，我们希望获得一个唯一边缘的点。

解决方法

来看看canny是怎么解决的？

（1）在梯度上找到极大值，具体如下图所示，C点为当前点梯度幅值，theta为梯度角度，选取8领域内的梯度幅值，其中g1，g2是左上角的2个梯度幅值，下图使用线性插值的方式计算出dTmp1=(g2-g1)*(1/tantheta) + g2，dTmp2 = (g4-g3)*(1/tantheta)+g3，公式不是一定这样，而是根据斜率的正负号以及涉及到的8领域中的点有关。

（2）设定2个阈值，如下图所示，当小于等于minval值即为抛弃的梯度值将该点设为0，当大于maxval值时认为该点设为255，即表示该点为确定边缘点。当该点的幅值处于maxval和minval中间时，判断该点周边时候存在确定边缘点，如果存在则认为该点是边缘点的一部分，否则该点被剔除设为0，这是为了让边缘点能够相对连续。

该函数即实现了上述功能

def judgeConnect(m2,threshold):e = 0.01s = []cood = []for i in range(m2.shape[0]):cood.append([])for j in range(m2.shape[1]):cood[-1].append([i,j])if abs(m2[i,j] - 255) < e:s.append([i,j])cood = np.array(cood)while not len(s) == 0:index = s.pop()jud = m2[max(0, index[0] - 1):min(index[0] + 2, m2.shape[1]), max(0, index[1] - 1):min(index[1] + 2, m2.shape[0])]jud_i = cood[max(0, index[0] - 1):min(index[0] + 2, cood.shape[1]), max(0, index[1] - 1):min(index[1] + 2, cood.shape[0])]jud = (jud > threshold[0])&(jud < threshold[1])jud_i = jud_i[jud]for i in range(jud_i.shape[0]):s.append(list(jud_i[i]))m2[jud_i[i][0],jud_i[i][1]] = 255return m2def DecideAndConnectEdge(g_l,g_t,threshold = None):if threshold == None:lower_boundary = g_l.mean()*0.5threshold = [lower_boundary,lower_boundary*3]result = np.zeros(g_l.shape)for i in range(g_l.shape[0]):for j in range(g_l.shape[1]):isLocalExtreme = Trueeight_neiborhood = g_l[max(0,i-1):min(i+2,g_l.shape[0]),max(0,j-1):min(j+2,g_l.shape[1])]if eight_neiborhood.shape == (3,3):if g_t[i,j] <= -1:x = 1/g_t[i,j]first = eight_neiborhood[0,1] + (eight_neiborhood[0,1] - eight_neiborhood[0,0])*xx = -xsecond = eight_neiborhood[2,1] + (eight_neiborhood[2,2] - eight_neiborhood[2,1])*xif not (g_l[i,j] > first and g_l[i,j] > second):isLocalExtreme = Falseelif g_t[i,j] >= 1:x = 1 / g_t[i, j]first = eight_neiborhood[0, 1] + (eight_neiborhood[0, 2] - eight_neiborhood[0, 1]) * xx = -xsecond = eight_neiborhood[2, 1] + (eight_neiborhood[2, 1] - eight_neiborhood[2, 0]) * xif not (g_l[i, j] > first and g_l[i, j] > second):isLocalExtreme = Falseelif g_t[i,j] >= 0 and g_t[i,j] < 1:y = g_t[i, j]first = eight_neiborhood[1, 2] + (eight_neiborhood[0, 2] - eight_neiborhood[1, 2]) * yy = -ysecond = eight_neiborhood[1, 0] + (eight_neiborhood[1, 0] - eight_neiborhood[2, 0]) * yif not (g_l[i, j] > first and g_l[i, j] > second):isLocalExtreme = Falseelif g_t[i,j] < 0 and g_t[i,j] > -1:y = g_t[i, j]first = eight_neiborhood[1, 2] + (eight_neiborhood[1, 2] - eight_neiborhood[2, 2]) * yy = -ysecond = eight_neiborhood[1, 0] + (eight_neiborhood[0, 0] - eight_neiborhood[1, 0]) * yif not (g_l[i, j] > first and g_l[i, j] > second):isLocalExtreme = Falseif isLocalExtreme:result[i,j] = g_l[i,j] #非极大值抑制result[result>=threshold[1]] = 255result[result<=threshold[0]] = 0# print(np.unique(result))# plt.imshow(result.astype(np.uint8), cmap='gray')# plt.axis('off')# plt.show()result = judgeConnect(result,threshold)result[result!=255] = 0return result

实现结果如下图所示。

4. 边缘的噪声点过多

问题

如下图红色区域部分，这些细碎的部分太多了，是不需要的噪声点，怎么去除？

解决方法

图像的边缘一般相对较长，是一段长距离的图像变化，如果短距离的我们可以认为是噪声点，想到的一种方式就是通过高斯模糊平滑掉这些短距离的图像变化点。

#coding:utf-8import numpy as npimport matplotlib.pyplot as pltimport osimport mathimport cv2def convolve(filter,mat,padding,strides):''':param filter:卷积核，必须为二维(2 x 1也算二维) 否则返回None:param mat:图片:param padding:对齐:param strides:移动步长:return:返回卷积后的图片。(灰度图，彩图都适用)@author:bilibili-会飞的吴克'''result = Nonefilter_size = filter.shapemat_size = mat.shapeif len(filter_size) == 2:if len(mat_size) == 3:channel = []for i in range(mat_size[-1]):pad_mat = np.pad(mat[:,:,i], ((padding[0], padding[1]), (padding[2], padding[3])), 'constant')temp = []for j in range(0,mat_size[0],strides[1]):temp.append([])for k in range(0,mat_size[1],strides[0]):val = (filter*pad_mat[j:j+filter_size[0],k:k+filter_size[1]]).sum()temp[-1].append(val)channel.append(np.array(temp))channel = tuple(channel)result = np.dstack(channel)elif len(mat_size) == 2:channel = []pad_mat = np.pad(mat, ((padding[0], padding[1]), (padding[2], padding[3])), 'constant')for j in range(0, mat_size[0], strides[1]):channel.append([])for k in range(0, mat_size[1], strides[0]):val = (filter * pad_mat[j:j + filter_size[0],k:k + filter_size[1]]).sum()channel[-1].append(val)result = np.array(channel)return resultdef linear_convolve(filter,mat,padding=None,strides=[1,1]):''':param filter:线性卷积核:param mat:图片:param padding:对齐:param strides:移动步长:return:返回卷积后的图片。(灰度图，彩图都适用) 若不是线性卷积核，返回None@author:bilibili-会飞的吴克'''result = Nonefilter_size = filter.shapeif len(filter_size) == 2 and 1 in filter_size:if padding == None or len(padding) < 2:if filter_size[1] == 1:padding = [filter_size[0]//2,filter_size[0]//2]elif filter_size[0] == 1:padding = [filter_size[1]//2,filter_size[1]//2]if filter_size[0] == 1:result = convolve(filter,mat,[0,0,padding[0],padding[1]],strides)elif filter_size[1] == 1:result = convolve(filter, mat, [padding[0],padding[1],0,0], strides)return resultdef _2_dim_divided_convolve(filter,mat):''':param filter: 线性卷积核,必须为二维(2 x 1也算二维) 否则返回None:param mat: 图片:return: 卷积后的图片,(灰度图，彩图都适用) 若不是线性卷积核，返回None'''result = Noneif 1 in filter.shape:result = linear_convolve(filter,mat)result = linear_convolve(filter.T,result)return resultdef judgeConnect(m2,threshold):e = 0.01s = []cood = []for i in range(m2.shape[0]):cood.append([])for j in range(m2.shape[1]):cood[-1].append([i,j])if abs(m2[i,j] - 255) < e:s.append([i,j])cood = np.array(cood)while not len(s) == 0:index = s.pop()jud = m2[max(0, index[0] - 1):min(index[0] + 2, m2.shape[1]), max(0, index[1] - 1):min(index[1] + 2, m2.shape[0])]jud_i = cood[max(0, index[0] - 1):min(index[0] + 2, cood.shape[1]), max(0, index[1] - 1):min(index[1] + 2, cood.shape[0])]jud = (jud > threshold[0])&(jud < threshold[1])jud_i = jud_i[jud]for i in range(jud_i.shape[0]):s.append(list(jud_i[i]))m2[jud_i[i][0],jud_i[i][1]] = 255return m2def DecideAndConnectEdge(g_l,g_t,threshold = None):if threshold == None:lower_boundary = g_l.mean()*0.5threshold = [lower_boundary,lower_boundary*3]result = np.zeros(g_l.shape)for i in range(g_l.shape[0]):for j in range(g_l.shape[1]):isLocalExtreme = Trueeight_neiborhood = g_l[max(0,i-1):min(i+2,g_l.shape[0]),max(0,j-1):min(j+2,g_l.shape[1])]if eight_neiborhood.shape == (3,3):if g_t[i,j] <= -1:x = 1/g_t[i,j]first = eight_neiborhood[0,1] + (eight_neiborhood[0,1] - eight_neiborhood[0,0])*xx = -xsecond = eight_neiborhood[2,1] + (eight_neiborhood[2,2] - eight_neiborhood[2,1])*xif not (g_l[i,j] > first and g_l[i,j] > second):isLocalExtreme = Falseelif g_t[i,j] >= 1:x = 1 / g_t[i, j]first = eight_neiborhood[0, 1] + (eight_neiborhood[0, 2] - eight_neiborhood[0, 1]) * xx = -xsecond = eight_neiborhood[2, 1] + (eight_neiborhood[2, 1] - eight_neiborhood[2, 0]) * xif not (g_l[i, j] > first and g_l[i, j] > second):isLocalExtreme = Falseelif g_t[i,j] >= 0 and g_t[i,j] < 1:y = g_t[i, j]first = eight_neiborhood[1, 2] + (eight_neiborhood[0, 2] - eight_neiborhood[1, 2]) * yy = -ysecond = eight_neiborhood[1, 0] + (eight_neiborhood[1, 0] - eight_neiborhood[2, 0]) * yif not (g_l[i, j] > first and g_l[i, j] > second):isLocalExtreme = Falseelif g_t[i,j] < 0 and g_t[i,j] > -1:y = g_t[i, j]first = eight_neiborhood[1, 2] + (eight_neiborhood[1, 2] - eight_neiborhood[2, 2]) * yy = -ysecond = eight_neiborhood[1, 0] + (eight_neiborhood[0, 0] - eight_neiborhood[1, 0]) * yif not (g_l[i, j] > first and g_l[i, j] > second):isLocalExtreme = Falseif isLocalExtreme:result[i,j] = g_l[i,j] #非极大值抑制result[result>=threshold[1]] = 255result[result<=threshold[0]] = 0# print(np.unique(result))# plt.imshow(result.astype(np.uint8), cmap='gray')# plt.axis('off')# plt.show()result = judgeConnect(result,threshold)result[result!=255] = 0return resultdef OneDimensionStandardNormalDistribution(x,sigma):E = -0.5/(sigma*sigma)return 1/(math.sqrt(2*math.pi)*sigma)*math.exp(x*x*E)if __name__ == '__main__':# Gaussian_filter_3 = 1.0/16*np.array([(1,2,1),(2,4,2),(1,2,1)]) #Gaussian smoothing kernel when sigma = 0.8, size: 3x3# Gaussian_filter_5 = 1.0/159*np.array([#[2,4,5,4,2],#[4,9,12,9,4],#[5,12,15,12,5],#[4,9,12,9,4],#[2,4,5,4,2]# ]) #Gaussian smoothing kernel when sigma = 1.4, size: 5x5sobel_kernel_y = np.array([[-1, -2, -1], [0, 0, 0], [1, 2, 1]])sobel_kernel_x = np.array([[-1, 0, 1], [-2, 0, 2], [-1, 0, 1]])pic_path = './img/'# pics = os.listdir(pic_path)pics = ["lena.png"]for i in pics:if i[-4:] == '.png'or i[-4:] == '.jpg' or i[-5:] == '.jpeg':filename = pic_path + iimg = plt.imread(filename)if i[-4:] == '.png':img = img*255img = img.mean(axis=-1) #this is a way to get a gray image.sigma = 1.52dim = int(np.round(6*sigma+1))if dim % 2 == 0:dim += 1linear_Gaussian_filter = [np.abs(t - (dim//2)) for t in range(dim)]linear_Gaussian_filter = np.array([[OneDimensionStandardNormalDistribution(t,sigma) for t in linear_Gaussian_filter]])linear_Gaussian_filter = linear_Gaussian_filter/linear_Gaussian_filter.sum()img2 = _2_dim_divided_convolve(linear_Gaussian_filter,img)# img2 = convolve(Gaussian_filter_5, img, [2, 2, 2, 2], [1, 1])plt.imshow(img2.astype(np.uint8), cmap='gray')plt.axis('off')plt.show()img3 = convolve(sobel_kernel_x,img2,[1,1,1,1],[1,1])img4 = convolve(sobel_kernel_y,img2,[1,1,1,1],[1,1])gradiant_length = (img3**2+img4**2)**(1.0/2)img3 = img3.astype(np.float64)img4 = img4.astype(np.float64)img3[img3==0]=0.00000001gradiant_tangent = img4/img3plt.imshow(gradiant_length.astype(np.uint8), cmap='gray')plt.axis('off')plt.show()#lower_boundary = 50final_img = DecideAndConnectEdge(gradiant_length,gradiant_tangent)cv2.imshow('edge',final_img.astype(np.uint8))cv2.waitKey(0)

实现结果为：

5. canny的完整流程