边缘检测——人体描边大师
微信公众号：命运探知之魔眼
如有问题或建议，请后台留言，我会尽力解决你的问题。
前言
我是人工智障，一名程序猿。做过嵌入式、爬虫，目前在自学计算机视觉 。注册 「命运探知之魔眼」（名字取自影视作品「命运石之门」）这个公号已有些日子，真正有心将它运营起来是看到朋友狗哥运营它的公号「一个优秀的废人」之后。注册这个号的初衷是分享我的 计算机视觉 学习笔记，希望更多的人加入我的行业。
边缘是什么
边缘可以定义为颜色剧变的区域。想象下面这样一幅图。
在中间的黑白交界处就是颜色剧变的地方，假如取其中一个像素观察，会发现它左边的像素值（255）比右边的像素值（0）大，左边差值的绝对值是255，可以理解为变化的大小或强度，变化的方向是沿着 x 方向的。这样有大小有方向的东西跟大学时学的梯度很像，因此最基本的边缘检测可以理解为求图像上每一点的梯度。
梯度的计算
我们把梯度分解为 x 方向的梯度 G_{x} 和 y 方向的梯度 G_{y} 。合并的梯度大小和方向就可以分别用以下的公式表示：
梯度大小：
|G| = \sqrt{ G_{x}^{2} + G_{y}^{2} }
梯度方向：
\theta = arctan( G_{y} / G_{x} )
某个像素点的 x 方向的梯度的计算可以通过这个像素点左右两边的像素值的差值的绝对值计算出来，而 y 方向的梯度可以通过该像素点上下两边的像素值的差值的绝对值计算。
根据计算方法的不同衍生了不同的边缘检测算法，其中一种就是用索贝尔（Sobel）算子与图像进行卷积。Sobel 算子是两个 3x3 的矩阵，分别负责计算 x 方向和 y 方向的梯度。
其中计算 x 方向梯度的算子是：
\left[ \begin{matrix} -1 & 0 & 1 \\ -2 & 0 & 2 \\ -1 & 0 & 1 \end{matrix} \right]
计算 y 方向梯度的算子是：
\left[ \begin{matrix} -1 & -2 & -1 \\ 0 & 0 & 0 \\ 1 & 2 & 1 \end{matrix} \right]
用这两个算子分别和图像做卷积就可以得到 x 方向和 y 方向的梯度。
图像卷积
假设我们用下面这样的卷积核：
\left[ \begin{matrix} -1 & 0 & 1 \\ -2 & 0 & 2 \\ -1 & 0 & 1 \end{matrix} \right]
对下面这样的图像进行卷积：
在结果矩阵中，卷积核中间的位置对应到结果矩阵上就等于 100x(-1) + 50x0 + 255x1 + 10x(-2) + 20x0 + 120x2 + 25x(-1) + 30x0 + 150x1 = 700 。
然后我们把卷积核往右移一格，再重复这个过程，直到遍历完整幅图。
这时候算出来的值是：285
得到的结果矩阵：
但是最左边一列和最上边一行的值我们算不出来，因为卷积核是 3x3 的，而图像外面没有值，我通常的处理方式是：如果是在算一幅比较大的图，就直接把边界值填充为 0 ，而如果是在做数值计算或者图像很小，去掉会很影响精度的情况下，就把图像往外扩充一个像素，并且把值填充为 0 ，这在卷积神经网络中也叫 padding 。
代码实现
下面是卷积的函数实现。
def filter2D(img, kernel):
img_h, img_w = img.shape[:2]
k_size = kernel.shape[0]
pad = (k_size - 1) // 2
result = np.zeros((img_h, img_w), dtype=np.float32)
for i in range(pad, img_h-pad):
for j in range(pad, img_w-pad):
result[i,j] = np.sum(img[i-pad:i+pad+1, j-pad:j+pad+1] * kernel)
return result
这里对边界的处理是直接填充为 0，注意 img[i-pad:i+pad+1, j-pad:j+pad+1] * kernel 这里，乘号两边是相同大小的矩阵，运算的结果是两个等大的矩阵对应位置的元素相乘。
在卷积完以后，要表示梯度大小，还要对卷积结果取绝对值，如果要用于显示，还要把所有值限制在 0 - 255 之间。
# x方向梯度
dx = np.abs(filter2D(img, dx_kernel))
dx[dx > 255] = 255
dx[dx < 0] = 0
# y方向梯度
dy = np.abs(filter2D(img, dy_kernel))
dy[dy > 255] = 255
dy[dy < 0] = 0
我们要处理的图片是这一张，原图有露点所以手动打了码。
lena.jpg
显示的效果如下：
x方向梯度
y方向梯度
合并的梯度按公式做是这样的：
# 合并梯度
grad = np.sqrt(dx**2 + dy**2)
grad[grad > 255] = 255
但在OpenCV实现的 Sobel 中其实是用两个方向梯度的平均，这样可以节省耗时的操作，而且不用再对超出范围的值做处理。
# 合并梯度
grad = dx / 2 + dy / 2
显示效果如下：
合并的梯度
完整的代码贴在这里：
# coding: utf-8
import cv2
import numpy as np
def filter2D(img, kernel):
img_h, img_w = img.shape[:2]
k_size = kernel.shape[0]
pad = (k_size - 1) // 2
result = np.zeros((img_h, img_w), dtype=np.float32)
for i in range(pad, img_h-pad):
for j in range(pad, img_w-pad):
result[i,j] = np.sum(img[i-pad:i+pad+1, j-pad:j+pad+1] * kernel)
return result
def main():
# 读取图像
img = cv2.imread('lena.jpg', 0)
# 索贝尔算子
# x方向算子
dx_kernel = np.array([[-1,0,1],[-2,0,2],[-1,0,1]])
# y方向算子
dy_kernel = np.array([[-1,-2,-1],[0,0,0],[1,2,1]])
# x方向梯度
dx = np.abs(filter2D(img, dx_kernel))
dx[dx > 255] = 255
dx[dx < 0] = 0
# y方向梯度
dy = np.abs(filter2D(img, dy_kernel))
dy[dy > 255] = 255
dy[dy < 0] = 0
# 合并梯度
grad = dx / 2 + dy / 2
# grad[grad > 255] = 255
# 显示
cv2.imshow('dx', dx.astype(np.uint8))
cv2.imshow('dy', dy.astype(np.uint8))
cv2.imshow('grad', grad.astype(np.uint8))
cv2.imwrite('dx.jpg', dx)
cv2.imwrite('dy.jpg', dy)
cv2.imwrite('grad.jpg', grad)
cv2.waitKey(0)
if __name__ == '__main__':
main()
上面的代码展示的是我们怎么手工实现 Sobel，但 OpenCV 有自带的实现方式。
def opencv_sobel(img):
x = cv2.Sobel(img, cv2.CV_16S, 1, 0)
y = cv2.Sobel(img, cv2.CV_16S, 0, 1)
absX = cv2.convertScaleAbs(x) # 转回uint8
absY = cv2.convertScaleAbs(y)
dst = cv2.addWeighted(absX,0.5,absY,0.5,0)
return dst
为什么边缘重要
边缘其实是很多其他更高层的特征如 HOG、SIFT 等图像特征的 building block 。我们要表达一个物体的轮廓时也是首先描绘边缘，就像画画时要先画线稿。一个更有力的科学证据就是 Hubel & Wiesel 在 1958 年对猫的视觉皮层做的实验（为什么科学家都喜欢虐猫？），他们把两个电极插进猫脑子里，探测初级视觉皮层的信号，并且发现猫只有在看到移动的边缘的时候初级视觉皮层才有信号，并由此逐步推导出视觉皮层对图像的理解是逐级往后传递的，越深的层级表达的信息越复杂，如初级只能理解边缘、到了后面就能理解基本形状、直到我们的人脸、汽车等。
后语
我不是大神，不是什么牛人，于 计算机视觉 领域来说，我是菜鸡，但谁刚开始接触一个领域的时候不是菜鸡呢。 写这个号的目的是为了记录我自学 计算机视觉 的笔记。
如果本文对你哪怕有一丁点帮助请右下角点赞，否则忽略就好。平时工作也是在做计算机视觉，希望大家多多指教。
我一直认为学习不能有所见即所得的想法，一看就会，一做就错。千万不要偷懒，所谓大神都是一个一个坑踩过来的。
最后，如果对 计算机视觉 感兴趣请长按二维码关注一波，我会努力带给你们价值，赞赏就不必了，能力没到，受之有愧。
