candy算⼦python_Python实现Canny边缘检测算法
Canny 边缘检测算法由计算机科学家 John F. Canny 于 1986 年提出的。其不仅提供了算法，还带来了⼀套边缘检测的理论，分阶段的解释如何实现边缘检测。Canny 检测算法包含下⾯⼏个阶段：
灰度化
⾼斯模糊
计算图⽚梯度幅值
⾮极⼤值抑制
双阈值选取
灰度化
灰度化实际上是⼀种降维的操作，可以减少计算。如果算法不进⾏⾊彩相关的识别的话，不灰度化，也可以直接进⾏后⾯的阶段。
# 灰度化
def gray(lf, img_path):男生贾里全传
"""
计算公式：
Gray(i,j) = [R(i,j) + G(i,j) + B(i,j)] / 3
or :
Gray(i,j) = 0.299 * R(i,j) + 0.587 * G(i,j) + 0.114 * B(i,j)
"""
天猫新风尚# 读取图⽚
img = plt.imread(img_path)
# BGR 转换成 RGB 格式
img_rgb = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)
# 灰度化
img_gray = np.dot(img_rgb[...,:3], [0.299, 0.587, 0.114])
return img_gray
⾼斯模糊
在实际的图⽚中，都会包含噪声。但有时候，图⽚中的噪声会导致图⽚中边缘信息的消失。对此的解决⽅案就是使⽤⾼斯平滑来减少噪声，即进⾏⾼斯模糊操作。该操作是⼀种滤波操作，与⾼斯分布有关，下⾯是⼀个⼆维的⾼斯函数，其中 (x, y) 为坐标，σ 为标准差：
进⾏⾼斯滤波之前，需要先得到⼀个⾼斯滤波器(kernel)。如何得到⼀个⾼斯滤波器？其实就是将⾼斯函数离散化，将滤波器中对应的横纵坐标索引代⼊⾼斯函数，即可得到对应的值。不同尺⼨的滤波器，得到的值也不同，下⾯是 (2k+1)x(2k+1) 滤波器的计算公式 :
常⽤尺⼨为 5x5，σ=1.4 的⾼斯滤波器。下⾯是 5x5 ⾼斯滤波器的实现代码：
# 去除噪⾳ - 使⽤ 5x5 的⾼斯滤波器
def smooth(lf, img_gray):
苏子的功效与作用# ⽣成⾼斯滤波器
"""
要⽣成⼀个 (2k+1)x(2k+1) 的⾼斯滤波器，滤波器的各个元素计算公式如下：
H[i, j] = (1/(2*pi*sigma**2))*exp(-1/2*sigma**2((i-k-1)**2 + (j-k-1)**2))
"""
sigma1 = sigma2 = 1.4
gau_sum = 0
gaussian = np.zeros([5, 5])
for i in range(5):
for j in range(5):
gaussian[i, j] = p((-1/(2*sigma1*sigma2))*(np.square(i-3)
+ np.square(j-3)))/(2*math.pi*sigma1*sigma2)
gau_sum = gau_sum + gaussian[i, j]
很想
# 归⼀化处理
gaussian = gaussian / gau_sum
# ⾼斯滤波
W, H = img_gray.shape
new_gray = np.zeros([W-5, H-5])
for i in range(W-5):
for j in range(H-5):
new_gray[i, j] = np.sum(img_gray[i:i+5, j:j+5] * gaussian)
return new_gray
图⽚梯度幅值
边缘是图像强度快速变化的地⽅，可以通过图像梯度幅值，即计算图像强度的⼀阶导数来识别这些地⽅。由于图⽚是离散的，可以⽤有限导数来近似图⽚的梯度：
图⽚梯度幅值为：
梯度⽅向为：
实现代码如下：
# 计算梯度幅值
def gradients(lf, new_gray):
"""
:type: image which after smooth
:rtype:
dx: gradient in the x direction
dy: gradient in the y direction
M: gradient magnitude
theta: gradient direction
"""
W, H = new_gray.shape
dx = np.zeros([W-1, H-1])
dy = np.zeros([W-1, H-1])
M = np.zeros([W-1, H-1])
theta = np.zeros([W-1, H-1])
for i in range(W-1):
for j in range(H-1):
dx[i, j] = new_gray[i+1, j] - new_gray[i, j]
如何叠爱心dy[i, j] = new_gray[i, j+1] - new_gray[i, j]
# 图像梯度幅值作为图像强度值吃的笔顺笔画
M[i, j] = np.sqrt(np.square(dx[i, j]) + np.square(dy[i, j]))
# 计算 θ - artan(dx/dy)
theta[i, j] = math.atan(dx[i, j] / (dy[i, j] + 0.000000001))
return dx, dy, M, theta
⾮极⼤值抑制(NMS)
理想情况下，最终得到的边缘应该是很细的。因此，需要执⾏⾮极⼤值抑制以使边缘变细。原理很简单：遍历梯度矩阵上的所有点，并保留边缘⽅向上具有极⼤值的像素。
梯度⽅向与边缘⽅向相互垂直
下⾯说说 NMS 的细节内容。NMS 在 4 个⽅向上进⾏，分别是 0，90，45，135，没有⾓度包含两个领域，因此，⼀共⽤⼋个领域：上，下，左，右，左上，左下，右上，右下，如下图所⽰，C 周围的 8 个点就是其附近的⼋个领域。
这样做的好处是简单， 但是这种简化的⽅法⽆法达到最好的效果， 因为，⾃然图像中的边缘梯度⽅向不⼀定是沿着这四个⽅向的。因此，就有很⼤的必要进⾏插值，找出在⼀个像素点上最能吻合其所在梯度⽅向的两侧的像素值。
NMS 是要找出局部最⼤值，因此，需要将当前的像素的梯度，与其他⽅向进⾏⽐较。如下图所⽰，g1，g2，g3，g4 分别是 C ⼋个领域中的 4 个点，蓝线是 C 的梯度⽅向。如果 C 是局部最⼤值的话，C 点的梯度幅值就要⼤于梯度⽅向直线与 g1g2，g4g3 两个交点的梯度幅值，即⼤于点 dTemp1 和 dTemp2 的梯度幅值。上⾯提到这种⽅法⽆法达到最好的效果，因为 dTemp1 和 dTemp2 不是整像素，⽽是亚像素。亚像素的意思就是在两个物理像素之间还有像素。
那么，亚像素的梯度幅值怎么求？可以使⽤线性插值的⽅法，计算 dTemp1 在 g1，g2 之间的权重，就可以得到其梯度幅值。计算公式如下：
weight = |gx| / |gy| or |gy| / |gx|
dTemp1 = weight*g1 + (1-weight)*g2
dTemp2 = weight*g3 + (1-weight)*g4
下⾯两幅图是 y ⽅向梯度值⽐较⼤的情况，即梯度⽅向靠近 y 轴。所以，g2 和 g4 在 C 的上下位置，此时 weight = |gy| / |gx| 。左边的图是 x，y ⽅向梯度符号相同的情况，右边是 x，y ⽅向梯度符号相反的情况。
对于左边的图来说，以 C 点为当前位置 - d[i, j] ，那么 g2 在 C 的前⼀⾏，g4 在 C 的后⼀⾏，所以位置坐标是：
g2 = d[i-1, j]；g4 = d[i+1, j]。根据左图的位置关系可以得到：g1 = d[i-1, j-1]；g3 = d[i+1, j+1]。
同理，根据右图的位置关系可以得到：g1 = d[i-1, j+1]；g3 = d[i+1, j-1]。
大理丽江旅游下⾯两幅图是 x ⽅向梯度值⽐较⼤的情况，即梯度⽅向靠近 x 轴。所以，g2 和 g4 在 C 的左右位置，此时 weight = |gy| / |gx| 。左边的图是 x，y ⽅向梯度符号相同的情况，右边是 x，y ⽅向梯度符号相反的情况。
由上⾯可知，可以得到如下信息：g2 = d[i, j-1]；g4 = d[i, j+1]；
身陷囹圄
左图：g1 = d[i+1, j-1]；g3 = d[i-1, j+1]；
右图：g1 = d[i-1, j-1]；g3 = d[i+1, j+1]。
下⾯的这两幅图，可能会带来理解帮助：
然后，根据以上信息，代码实现如下：
def NMS(lf, M, dx, dy):
d = np.copy(M)
W, H = M.shape
NMS = np.copy(d)
NMS[0, :] = NMS[W-1, :] = NMS[:, 0] = NMS[:, H-1] = 0
for i in range(1, W-1):
for j in range(1, H-1):
# 如果当前梯度为0，该点就不是边缘点
if M[i, j] == 0:
NMS[i, j] = 0
el:
gradX = dx[i, j] # 当前点 x ⽅向导数
gradY = dy[i, j] # 当前点 y ⽅向导数
gradTemp = d[i, j] # 当前梯度点
# 如果 y ⽅向梯度值⽐较⼤，说明导数⽅向趋向于 y 分量
if np.abs(gradY) > np.abs(gradX):
weight = np.abs(gradX) / np.abs(gradY) # 权重
grad2 = d[i-1, j]
grad4 = d[i+1, j]
# 如果 x, y ⽅向导数符号⼀致
# 像素点位置关系
# g1 g2
# c
# g4 g3
if gradX * gradY > 0:
grad1 = d[i-1, j-1]
grad3 = d[i+1, j+1]
# 如果 x，y ⽅向导数符号相反
# 像素点位置关系
# g2 g1
# c
# g3 g4
el:
grad1 = d[i-1, j+1]
grad3 = d[i+1, j-1]
# 如果 x ⽅向梯度值⽐较⼤
el:
weight = np.abs(gradY) / np.abs(gradX) grad2 = d[i, j-1]
grad4 = d[i, j+1]
# 如果 x, y ⽅向导数符号⼀致
# 像素点位置关系
# g3
# g2 c g4
# g1
if gradX * gradY > 0:
grad1 = d[i+1, j-1]
grad3 = d[i-1, j+1]
# 如果 x，y ⽅向导数符号相反
# 像素点位置关系
# g1
# g2 c g4
# g3
el:
grad1 = d[i-1, j-1]
grad3 = d[i+1, j+1]
