最早接触霍夫变换(Hough Transform )[4,5,6,7,8]就是用在检测直线，大四的时候用来分割ROI用的，原理看起来很简单。

一．检测直线

直线y=k*x+b上一点(x0,y0)，满足方程y0=k*x0+b。那么换过来看，在坐标系(k,b)上，这就是一条直线，即b=- x0* k+ y0;再来直线上的一点(x1,y1)，同样满足y1=k*x1+b，那么在(k,b)坐标系上就是另外一条直线，即b=-x1*k+y1;这两条直线的交点就是(k,b)，即为(x,y)平面的斜率和截距。

(x,y)平面上的一条直线映射到(k,b)平面就变成一个点了。

为了解决k为无穷大的问题以及如果采用(k,b)矩阵，参数范围会比较庞大(k和b的范围不好取)，采用极坐标表示x*cos(theta) + y*sin(theta)=r，这样r就固定在x，y对角线范围内，角度可以设置为0到180°。程序关键在于步长怎么确定，太密了时间长，太稀疏了，可能找不到最大值。

1．  算法流程为:

1)       初始化累加器H为零

2)       遍历图像中每个边缘点

for theta = 0 to 180

r=x*cos(theta) + y*sin(theta)

H(theta,r) = H(theta,r)+1;

end

3)       H中最大的值对应的theta和r就是要求的直线(实际使用的时候可以设定最大值的阈值)

Code：

M = imread(‘sqr1.gif’);

if size(M,3) == 3

M1=rgb2gray(M);

else

M1=M;

end

BW = edge(M1,’canny’);

[sy,sx]=size(I);

[y,x]=find(I);

totalpix = length(x);

maxrho = round(sqrt(sx^2 + sy^2)); %Rou的最大坐标为对角线

HM = zeros(2*maxrho,180); %1°一个bin

for cnt = 1:totalpix;

cnt2 = 1;

for theta = -pi/2:pi/180:pi/2-pi/180

rho = round(x(cnt).*cos(theta) + y(cnt).*sin(theta));             HM(rho+maxrho,cnt2) = HM(rho+maxrho,cnt2) + 1;

cnt2 = cnt2 + 1;

end

end

theta = rad2deg(-pi/2:pi/180:pi/2-pi/180);

rho = -maxrho:maxrho-1;

imshow(uint8(HM),[],’xdata’,theta,’ydata’,rho);

xlabel(‘\theta’),ylabel(‘\rho’)

axis on, axis normal;

title(‘Hough Matrix’);

2．  改进算法

考虑到图像的边缘点是有梯度的，梯度的方向就是theta。因此在循环的时候只在该边缘点附近循环theta(边缘点的梯度方向可能会有误差，所以要加一个delta_theta)，同时给该点的梯度值增加一个阈值。即上式

for theta = -pi/2:pi/180:pi/2-pi/180

改为

for theta =theta(x(i)) – delta_theta:pi/180: theta(x(i)) + delta_theta；假如theta(x(i))为xi点的梯度方向

3．  噪声的影响

噪声等级越大，H矩阵中vote权值越来越小。均匀分布的噪声会产生一些假的直线，当这些噪声点变多的时候，也会加大vote的值。

解决办法：根据梯度值去掉一些小的边缘点。对相邻的bins进行加权(累加器)

4．  Example

直线

矩形和圆

二．检测圆[1]和椭圆

检测圆需要考虑三个参数，圆心和半径。如果半径已知，只需要构造一个H矩阵，索引为图上的点，H(a,b)，然后以已知的半径画圆，落在这个圆上的点就在累加器里面加1。如果半径不知，就变成三维问题，复杂度变多了。同理，考虑边缘点的梯度方向，边缘点的梯度方向都是指向圆心的，于是可以先检测圆心H(x0,y0)，然后再检测半径。或者干脆就循环半径r，H(x0,y0,r)，其中

X0=x+r *cos(gradient_theta)

Y0=y+r*sin(gradient_theta)

总体思路是减少不必要的循环。

椭圆就更恐怖了，五个参数，中心点，长轴，短轴，以及旋转角度。当然，如果用梯度的话，循环也会更少，只需要ab长短轴，再添加上角度的话也就是三个循环。

X0=x+a *cos(gradient_theta)

Y0=y+b*sin(gradient_theta)

Wikipedia上有改进算法

Yonghong Xie; Qiang Ji (2002). A new efficient ellipse detection method. pp. 957

考虑到梯度信息对噪声比较敏感，这篇paper没有使用梯度信息，而是利用椭圆长轴上的两个点，然后利用数学公式来进行短轴的估计运算。也有三个循环，不过复杂度应该降低了，循环只是在边缘点上。

三．广义霍夫变换Generalized Hough Transform[3]

其实这个是关键，主要用于图像匹配，识别等。GHT不考虑具体的曲线的数学表达式，而是将任意形状转化成点集。具体来说就是预先把参考图形边缘点对中心的径向量保存起来。线性时间！

假如有个参考图形，选取一个参考点(内部点，例如重心等XR)，连接边缘点和参考点，两点与x轴的夹角为alpha，向量为r，同时记住该点的梯度方向phi。遍历所有边缘点，形成R-table。这个表以梯度方向为索引(常用delta_phi的倍数表示)，表里面存的是二元组(r,alpha)，这样可能一个梯度phi对应多个二元组。(由参考点+边缘点-à二元组)

匹配的时候，根据待匹配的图像的边缘点的梯度方向phi，找到对应的一个或多个(r,alpha)，求得对应的参考点(二元组à得到参考点)。累加器为A[XR,YR]

公式：

XR=x+r *cos(alpha)

YR=y+r*sin(alpha)

考虑到图形的尺度变化和旋转变化，加上两个参数尺度s和选择tau，得到A[XR,YR,s,tau]

XR=x+r *s*cos(alpha+tau)

YR=y+r*s*sin(alpha+tau)

回到初始的霍夫变换：

半径已知的圆就是累加器去掉尺度s，因为圆旋转不变，tau也可以不要。如果圆的半径不知道，那么就得要加上s参数了。

对应椭圆的话参数就少很多了，保存在表中的数据依然是(r,alpha)，累加的时候也只需要两个参数。这里就可以看出来广义霍夫变换其实是不考虑曲线的数学表达式的，它只是一种转换，换个空间表达同样的内容。

四．图像识别

上述累加器其实就是计算参考点的权重，不仅仅能对参考点加权，也可以对一种转化加权。同样的idea可以用在feature上[2]。

Training：

1． 在兴趣点周围抽取图像patch，用聚类建立码本

2． 将兴趣点周围的图像patch映射到最近的码本[梯度方向]

3． 对映射后的码本，存储其相对于目标中心的所有的位置[r,alpha]

Testing：

1． 给一幅测试图像，抽取patchs，匹配其对应的码本(相当于梯度方向)

2． 对目标中心可能的位置进行加权[XR,YR]

3． 在加权空间寻找最大值

4． 还原目标

说明：

1． 监督学习，需要参考点以及分割mask

2． 加权空间是连续的，不是离散的，聚类算法需要找到最值

3． 对于尺度变化 a)搜索一系列的尺度，类似于不知道圆的半径一样b)使用尺度相关的兴趣点

4． Verification很重要，像素与像素之间的比较

五．总结

Pros：

1． 能够处理non-locality和遮挡问题

2． 能够处理一个模型的多个实例

3． 对噪声有点健壮性：噪声点不太可能对任何单一的bin都会有作用的

4． 能够发现图像中潜在的结构(尤其是边缘很杂乱)

Cons：

1． 随着模型参数的增加，搜索时间复杂度指数增加

2． 非目标形状也能够形成一些峰值(发现多个形状)

3． 网格大小不好确定

六．后记

霍夫变换和最小二乘找直线的区别：霍夫变换是寻找点最多的直线，而最小二乘是找误差最小的直线。有的时候可以将二者结合起来

首先,用霍夫变换剔除数据点集中的干扰点或噪声,并将分布在不同直线附近的点分离出来;然后,用最小二乘法拟合各直线.该方法既解决了直接使用最小二乘法拟合时,拟合直线易受干扰点或噪声的影响和数据点分布在多条直线附近而无法拟合的两个问题;同时也解决了直接使用霍夫变换时,拟合直线精度不高和直线段有效区间不容易控制的问题

Reference

[2]B. Leibe, A. Leonardis, and B. Schiele, Combined Object Categorization and Segmentation with an Implicit Shape Model, ECCV Workshop on Statistical Learning in Computer Vision 2004

[3]

-33.883000

151.217000

赞过：

赞 正在加载……

相关