Canny 边缘检测与内边界追踪

对目标图片进行边缘检测与边缘连接

原理

Canny 边缘检测算法边缘连接/ Inner boundary tracing 的算法过程如下.

高斯滤波

使用高斯滤波矩阵作用于图像, 用以平滑图像来消除噪点.

梯度值与方向

使用 Sobel 算子 (或 Roberts, Prewitt 等), 分别得到横向与纵向两个方向的梯度. 根据横向梯度 \(g_x\) 和纵向梯度 \(g_y\), 计算梯度值与梯度方向. \[ G = \sqrt{g_x^2 + g_y^2}\\ \theta = arctan \frac{g_y}{g_x} \]

非最大值抑制

在穿过某个点的四个直线上, 分别与相邻的两个梯度值进行比较, 只有它比相邻梯度值都大时将其保留.

这样能够将比较多的梯度值简化为较少的边界.

双阈值过滤

使用双阈值对图片的梯度值进行过滤, 若大于较大阈值则认为是边界; 若小于较小阈值则认为不是边界; 若位于双阈值之间, 且有超过较大阈值的梯度值相邻, 则也认为是边界, 否则认为不是边界. 伪代码:

1
2
3
4
5
6
if (image(row, col) > higher_threshold or
(image(row, col) > lower_threshold and
neighbor_greater_than_higher_threshold(row, col))):
image(row, col) = 1;
else
image(row, col) = 0;

内边界追踪

使用 Sobel 算子得到的图片有八连通边界, 以此为条件, inner boundary tracing 过程如下

  • 设置初始 dir 为 7, 不同 dir 对应的方向如图 5.13(b);

  • 逆时针搜索点的八近邻, 直到发现下一个边界点设为 dir, 且:

    • dir = (dir + 7) % 8, if dir is even (图 5.13(d))

    • dir = (dir + 6) % 8, if dir is odd (图 5.13(e))

  • 搜索直到无法找到下一个边界点或回到起点.

细节

Canny 边缘检测算法在维基百科中的介绍已经非常全面了, 但一些实现细节仍需探讨.

频域上的变化

考虑到傅立叶变换后对高频空间进行过滤也能够有效地平滑图像, 因此考虑进行快速傅立叶变换之后, 将高频部分置为 0.

当然, 将低频部分过滤也可能起到非常好的作用, 但难以使用同样的指标对不同图片都进行低频过滤, 因此放弃过滤低频.

双阈值的取值

维基百科并没有对双阈值的取值给出准确的推荐值, 根据 Quora: How to set the threshold in canny edge detection, Stackoverflow: Canny edge detector threshold, 认为有两种对双阈值取值的方法.

  • 根据 Lab Book Pages: Otsu Thresholding 计算 OTSU 阈值作为高阈值, 并取其 1/2 作为低阈值.

  • 原灰度图平均值的 1.33 倍作为高阈值, 0.66 倍作为低阈值.

经过实践, 认为第二种方法的效果会更好并更具有可控性. 两种算法在代码中都有体现.

效果

Pic A

Pic B

实现

完整代码: edge_detection.m edge_linking.m

高斯滤波

高斯滤波在 MatLab 中可以通过 imfilter 实现.

1
2
3
4
5
6
7
8
function gauss_ret = gauss_blur(gauss_img)
filter_matrix = 1/159 * [
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];
gauss_ret = imfilter(gauss_img, filter_matrix);

梯度值与方向

使用 Sobel 算子也可以使用 imfilter 函数. 在计算出 \(\theta\) 后, 为了避免 nan 引入的问题, 将这些地方设为 2.

1
2
3
4
5
6
7
8
9
function [sobel_g, sobel_theta] = sobel(sobel_img)
sobel_x = [-1 0 1; -2 0 2; -1 0 1];
sobel_y = [1 2 1; 0 0 0; -1 -2 -1];
x_ret = imfilter(sobel_img, sobel_x);
y_ret = imfilter(sobel_img, sobel_y);
sobel_g = sqrt(x_ret.^2 + y_ret.^2);
sobel_theta = atan(y_ret./x_ret) * 180 / pi;
sobel_theta = round(sobel_theta / 45) * 45;
sobel_theta(isnan(sobel_theta)) = 2;

非极大值抑制

非极大值抑制只需要在四个方向上进行比较并拷贝. 为什么需要拷贝而不是在原图上操作? 避免在原图上的新结果影响最终的结果值.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
edge = zeros(height, width);
for i = 2 : height - 1
for j = 2 : width - 1
if (...
(abs(theta(i, j)) == 90 &&...
g(i, j) >= g(i+1, j) &&...
g(i, j) >= g(i-1, j)) ||...
(theta(i, j) == 0 &&...
g(i, j) >= g(i, j + 1) &&...
g(i, j) >= g(i, j - 1)) ||...
(theta(i, j) == 45 &&...
g(i, j) >= g(i + 1, j + 1) &&...
g(i, j) >= g(i - 1, j - 1)) ||...
(theta(i, j) == -45 &&...
g(i, j) >= g(i + 1, j - 1) &&...
g(i, j) >= g(i - 1, j + 1))...
)
edge(i, j) = g(i, j);
else
edge(i, j) = 0;
end
end
end

双阈值过滤

双阈值的取值设为平均值的 1.33 倍与 0.66 倍.

1
2
3
4
5
6
upper_thresh = 1.33 * mean_value;
lower_thresh = 0.66 * mean_value;
edge((lower_thresh < edge) & (edge < upper_thresh)) = 0.5;
edge(edge < lower_thresh) = 0;
edge(edge > upper_thresh) = 1.0;
edge = filter_low_thresh(edge);

在 filter_low_thresh 函数中, 遍历所有点的八近邻, 如果当前点小于较大阈值大于较小阈值, 八近邻中有点大于较大阈值, 则视为边界.

内边界追踪

由于要将数字 0-7 转化为 8 个方向, 代码比较冗长, 因此这里省略部分代码.

在这部分, 需要明确:

  • 什么时候终止继续寻找边界, 回到起点或找不到下一个边界点

  • MatLab index 值为 1-8, 需要进行转换

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
output = [row col;];
image(row, col) = 0;
dir = 7;
while (1)
[row, col, next_found, dir] = find_next_dir(row, col, dir, image);
if (~next_found)
break;
end
if (output(1, 1) == row && output(1, 2) == col)
break;
end
image(row, col) = 0;
output = [output; row col;];
end

function [frh, frw, frfound, frp] = find_next_dir(fh, fw, fdir, fimage)
[height, width] = size(fimage);
signs = zeros(1, 8);
# ...
# set signs
# ...
frfound = 0;
frp = mod(fdir + 1, 8) + 1;
while (1)
if (signs(frp))
frfound = 1;
break;
end
frp = mod(frp, 8) + 1;
if (frp == mod(fdir + 1, 8) + 1)
break;
end
end
# ...
# set return frh & frw
# ...
frp = frp - 1;
if (mod(frp, 2) == 1)
frp = mod(frp + 6, 8);
else
frp = mod(frp + 7, 8);
end

参考资料

Lab Book Pages: Otsu Thresholding

Quora: How to set the upper and lower threshold in canny edge detection

Stackoverflow: Canny edge detector threshold

Wiki: Canny edge detector

Wiki: Sobel detector

Border Tracing

0%