1

直线检测问题

在纸上画一条直线，用手机拍下照片，把照片交给计算机识别。

计算机是如何知道这张照片中的这条直线的？

存在直线吗？

直线在哪里？

点、线、面是基本的几何元素。

欧几里得在《几何原本》中写道：

直线是点朝一个方向以及反方向的无限平铺。

一条直线在图像中，每个人一眼就能看出来。

我们不仅能够说出图像中存在一条直线，还可以说出直线经过哪个点，以及直线的方向。

如果有多条直线，我们还能说出直线之间的关系，相交还是(近似)平行，交点在哪里……

计算机怎么识别图像中的直线？

计算机的常规方法是提取边缘点，之后检测直线。

2

灰度矩阵与直线的灰度特征

计算机接收到一张照片，获得一个对应图像的灰度矩阵。

矩阵是一格一格的，一格称为一个像素(pixel)。

矩阵的横、纵坐标是离散化的。

矩阵中的元素是量化的值，从0到255，该值称为灰度。

灰度，也称亮度，是该点受到光照亮度的对应数字。

如果是彩色的照片，对于计算机来说就是3张灰度图，即3个灰度矩阵。

按顺序分别是红的、绿的和蓝的灰度图。

因为在摄像机成像元件(CCD或CMOS)前面嵌有拜尔滤镜(如下图所示)，使得三种颜色(即波长)的光进入到不同的灰度图中。

显示时，三种颜色同时显示，由于红、绿、蓝挨得很近，远处看就是彩色的。

1975年柯达美国分公司申请了拜尔滤镜发明专利(US3971065)，1976年获得授权。

下面是一张灰度图：

计算机的数据用Matlab软件显示出来：

白色直线——灰度呈山峰形貌(黑色直线——灰度呈裂缝形貌)。

下面是两块钢板拼在一起的灰度图，中间有条直线(接缝)，如何进行识别直线？

3

边缘检测原理

首先对灰度矩阵进行边缘提取。

边缘就是图像中灰度值变化大的点，凡是差值大于一定程度就当作边缘点。

在灰度矩阵中沿着水平方向从左到右移动当前位置，

从第二个像素开始，将该位置左右两个像素作差分，后者减去前者的灰度差存在另一个结果矩阵中对应当前位置处，直到倒数第二个像素为止，这样就获得了水平方向的边缘点差分结果矩阵。

沿着垂直方向也如此，获得垂直方向的边缘点差分结果矩阵。

有人觉得这样还不能精准找到边缘点，所以在上述思想上进行了修正，给出了索贝尔(Sobel)方法。

Sobel方法是将如下的左边小矩阵扣到灰度矩阵上，将上下对应项相乘再相加，得数就是这个中间位置处的水平差分结果，放入到水平差分结果矩阵中存储，然后移动该小矩阵向右，直到一行都完成这样的运算，再下移一行进行……直到全部算完，获得水平结果矩阵。

下方的右边小矩阵如法炮制，可以获得垂直结果矩阵。

这是在图像上进行计算的典型做法——邻域计算，邻域计算所用到的上述小矩阵被称为“算子”。

上面给出的是Sobel算子，分为x方向和y方向两种算子。

算子不同，得到的处理不同。

下面是Sobel算子与Prewitt算子对比：

在边缘检测领域，有roberts算子、canny算子、log算子、laplacian算子……

用不同的算子会得到不同的结果。

如果边缘点的方向是倾斜的(既不水平也不垂直)，在上述两个矩阵中就容易不明显，为了找到这些边缘点，还需要将上述两个矩阵进行一个整合。

整合办法是：

将上述两个矩阵的对应元素平方和、相加，再开根号，

这样就得到了全方向的边缘点差分结果矩阵。

整合公式：

包括了任意方向的灰度差分结果矩阵。

该差分结果矩阵将边缘点凸显出来。

接着进行二值化：

将上述边缘点差分结果矩阵进行二值化：

凡是边缘点差分大于某个事先预定的值(称为阈值)就当作边缘点，设其值为1，非边缘点设值为0。

从而获得一个元素值只有1或0的黑白图像。

下图为采用Sobel方法进行边缘提取得到的黑白图像。

阈值的取值需要慎重：

如果阈值取得过小，会将大量不是边缘的点当作边缘点。

如果阈值取得过大，会将大量是边缘的点忽略掉。

阈值取得过小的结果：

4

最小二乘法

接下来的任务是：

在上述边缘点黑白图像中识别出直线。

有一种方法是用于找到多个点的拟合直线的，叫做最小二乘法。

最小二乘法就是找到一条直线，该直线到各点的距离之和最短。

但是这种方法不适用于边缘点黑白图中识别直线，

原因是边缘点有不少存在错误，即存在噪点，

一个偏离直线很远的点可能被错认为是边缘点，

而且在直线两侧的点的数量也不一样多，

都会导致最小二乘法找不出正确的直线。

噪点存在的原因是多方面的，很难没有一个噪点：

1)光线总是或多或少的照进摄像机；

2)镜头畸变，弧面精度有限，玻璃材质均匀程度偏差导致折射率有变化；

3)摄像机电路干扰、电磁影响；

4)边缘点提取时阈值选择不当会带来大量噪点……

因此最小二乘法无法在存在大量错误点的情况下找到直线。

5

哈夫变换

大量采用的直线检测的方法是哈夫变换。

哈夫变换(Hough Transform, 简称HT，又译作霍夫变换)可以找到经过点最多的直线。

在计算机视觉和图像处理领域，哈夫变换作为一种形状分析技术被广泛应用。

哈夫变换是1959年由Paul Hough申请的美国发明专利技术，1962年获得授权。

Hough发明专利的名称为：

Method and Means for Recognizing Complex Patterns(用于识别复杂图案的方法和手段)

拿穿羊肉串打个比方，就是要找到羊肉穿得最多的竹签。

经过边缘点最多的直线就是要找的直线。

这个方法可以有效克服噪点的影响，

因为噪点往往不在直线上，噪点的数量一般是少于经过直线上的点的。

如果噪点的数量多于经过直线上的点，是什么情况？

就是人眼也分辨不出该直线了，

这种情况出现在前面的边缘点提取环节，灰度跳变的阈值选小了。

哈夫变换方法怎么做的呢？

哈夫变换方法是：

利用直线参数方程：

建立一个对应的哈夫空间。

哈夫空间的横坐标是直线相对于坐标横轴的角度，

而纵坐标是坐标原点到直线的垂足长度。

这样一来，每个直角坐标空间的点都对应哈夫空间中的一条正弦曲线，

反之，每个哈夫空间中的点都对应直角坐标空间中的一条直线。

这个变换称为哈夫变换。

哈夫变换就是直角坐标空间(也称为笛卡尔空间)与极坐标空间的变换。

哈夫变换将每个边缘点换算成为哈夫空间的正弦曲线，

多个边缘点就有多条正弦曲线，

哈夫空间中的这些正弦曲线的交点，

对应于经过了最多点的一条直线。

下图为哈夫变换得到的哈夫空间结果：

存在一个高耸的塔尖，塔尖这一点的高度就是投票数，说明有最多的点在对应该点的直线上，塔尖这一点在哈夫空间中的角度值、垂足距离值可以首先得到，用参数方程就得到了对应的直线。

在原来边缘点图像中找直线的任务，就转换成在对应的哈夫空间中找最多正弦曲线的交点的任务。

计算机善于找最大值或最小值，

模模糊糊、大约多少的事情计算机不擅长。

如何找最多曲线的交点？

计算机采用投票的方式，

哈夫空间中，正弦曲线经过的每个小格子里都会投上一票，

当所有的边缘点对应的正弦曲线都画完，则投票结束，

数数得票数，就找到了得票最多的那个点。

有了这个点，其对应的直线参数也就找到了，

直线也就识别出来了。

直线识别的结果(白线)叠加在原图上：

之所以划分很多格子，

是因为计算机只能做量化计算。

每个格子的边长称为步长(或周期)。

一个是直线角度的步长，

另一个是垂足距离的步长。

步长越小，计算得越精确，

但是耗时间就越长。

不仅直线可以通过哈夫变换来检测，圆也可以通过哈夫变换来识别。

凡是能够有确定的解析式的曲线，也都可以使用哈夫变换的方法来检测。

哈夫变换还可以识别多条直线、多个圆弧。

哈夫变换是图像处理中从图像中识别几何形状的基本方法之一。

哈夫变换的基本原理在于利用点与线的对偶性，将原始图像空间的给定的曲线通过曲线表达形式变为参数空间的一个点，从而把原始图像中给定曲线的检测问题转化为寻找参数空间中的峰值问题。可以适用于直线、椭圆、圆等。

6

哈夫变换教程视频

来看看哈夫变换的教程视频：

7

哈夫变换的不足

哈夫变换方法存在一个问题：效率低。

因为哈夫变换方法需要将每个点对应的正弦曲线画在哈夫空间中，

一个点就要画出179根线(如果是步长为1度的话)，

每个点都这样做一遍，

计算量比较大，程序运行耗时长。

因此哈夫变换方法也被人取了一个绰号：咖啡方法。

一杯咖啡慢悠悠地喝完了，

哈夫变换的计算结果才姗姗来迟。

8

哈夫变换改进方法

人们针对哈夫变换方法提出了多种改进措施，

于是产生了随机哈夫变换、快速哈夫变换、随机快速哈夫变换等。

在产生了这些改进的哈夫变换后，

原来的哈夫变换就被定义为：

标准哈夫变换。

随机哈夫变换方法是随机找边缘点来进行哈夫变换，

在不断运行的过程中，就会找到交于一点的直线数量逐渐上升，

当直线数量达到预设的一个数值时，就终止程序，

认为此时已经找到了该直线。

这个方法有一定的风险性，

只要交于一点的直线数量足够的大就好，

风险虽然仍然存在，但是已经很小。

快速哈夫变换方法是将两个点一同拿来做哈夫变换。

因为两点有且仅有一条直线，

因此用该直线在哈夫空间去绘制一个点，

这样计算量就大大降低了。

原来标准哈夫变换方法需要画出整条正弦曲线，

两个点就是两条正弦曲线，

现在只用在哈夫空间画出两个点所在直线的一个点即可。

8

哈夫变换Matlab程序

事先拍摄一张照片，存为灰度图像，名为a1.bmp，将其放在下面的程序同一个文件夹中。打开Matlab，将下面的程序新建一个脚本文件，运行该文件。

%Hough detection

clear;          %清除上次运行的所有变量

close all;     %关闭上次运行的子窗口

%=====下面读取原图并显示

I1 = imread('a1.bmp');     %事先有张灰度图，名为a1.bmp

figure(1),imshow(I1)         %将图像显示出来

I2 = im2double(I1);          %将图片的数据类型从unit8变为double型，只有这样回叙才能进行处理

%=====下面进行Sobel方法的边缘提取并显示

I_edge = edge(I2,'sobel');  %边缘提取，采用Sobel方法

figure(2),imshow(I_edge)   %将结果图像显示出来

%=====输入源图像大小

image_y_max=570;

image_x_max = 760;

%=====输入搜索初始角度

alpha_initial = 100;

%=====输入搜索步长角度

thelta_step = 1;

%=====输入搜索总角度范围

total_angle_for_search = 180;

m = round(total_angle_for_search/thelta_step+1);

%=====初始化hough投票矩阵，未防止计算的lamda距离超出矩阵，设置一个大的矩阵

n=2000;

M_hough=zeros(n,m);

%=====防止距离为负

lamda_initial =round(n/2);

%=====对每个像素操作

for x = 1 : image_x_max - 1

for y = 1 : image_y_max - 1

%=====转换为矩阵下标的坐标系

u = -y + image_y_max;

v=x;

%=====如果有点，就投票

if I_edge(u,v) == 1

%======循环m次，直线旋转投票

for k = 1 : m-1

%=====计算角度

alpha1 = alpha_initial + k*thelta_step;

alpha1 = round(alpha1);

%=====计算距离

lamda1 = x*cos(alpha1*pi/180)+y*sin(alpha1*pi/180);

lamda1 = round(lamda1);

M_hough(lamda_initial + lamda1,k) = M_hough(lamda_initial + lamda1,k) + 1;

end

end

end

end

%=====显示M_hough矩阵

figure(5),mesh(M_hough);

%======取出最大值所在的行、列

[rm,km]=find(M_hough == max(max(M_hough)))

%=====计算该点原距离

rm = rm - lamda_initial;

%draw the line

%=====计算该点原角度

alpha_line = alpha_initial + km*thelta_step

alpha_line_rad = alpha_line*pi/180;

lamda_line = rm/sin(alpha_line_rad);

%=====在原图像上重叠画出该直线

for x=1:image_x_max-1

%=====直线方程式

y = -x*(1/tan(alpha_line_rad))+lamda_line;

%=====只画在图像中的直线

if y <= image_y_max - 1

if y >= 1

v = x;

u = -y + image_y_max;

u = round(u);

I2(u,v) = 255;

end

end

end

%=====显示检测结果

figure(6),imshow(I2);

关于图像上的直线检测，你是否还有更好的办法呢？