1、Harris 角点检测算法优化 一、综述 用 Harris 算法进行检测，有三点不足：（1 ）该算法不具有尺度不变性；（2 ）该算法提取的角点是像素级的；（3 ）该算法检测时间不是很令人满意。 基于以上认识，我主要针对第（3 ）点对 Harris 角点检测算法提出了改进。 二、改进 Harris 算法原理 在介绍我的方法之前，我先提出如下概念：图像区域像素的相似度。我们知道， Harris 角点检测是基于图像像素灰度值变化梯度的， 灰度值图像的角点附近，是其像素灰度值变化非常大的区域，其梯度也非常大。换句话说，在非角点位置邻域里，各点的像素值变化不大，甚至几乎相等，其梯度相对也比较小。从这个角
2、度着眼，我提出了图像区域像素的相似度的概念，它是指检测窗口中心点灰度值与其周围 n 邻域内其他像素点灰度值的相似程度，这种相似程度是用其灰度值之差来描述的。如果邻域内点的灰度值与中心点 Image (i,j) 的灰度值之差 的绝对值 在一个阈值 t 范围 内，那我就认为这个点与中心点是相似的。与此同时，属于该 Image (i,j) 点的相似点计数器 nlike(i,j) 也随之加一。在 Image (i,j) 点的 n 邻域 全部被遍历一边之后，就能得到 在这个邻域范围内 与中心点相似的 点个数的统计值 nlike(i,j) 。根据 nlike(i,j) 的大小，就可以判断这个中心点是否可能
3、为角点。 由于我选择 3*3 的检测窗口，所以， 对于中心像素点 , 在下面的讨论中只考虑其 8 邻域内像素点的相似度。 计算该范围的像素点与中心像素点的灰度值之差的绝对值 ( 记为 ) , 如果该值小于等于设定的阈值 ( 记为 t) , 则认为该像素点与目标像素点相似 。 nlike(i,j)=sum(R(i+x,j+y)( - 1 x 1 , - 1 y 1 , 且 x 0 , y 0) , 其中 : 1 , ( i + x , j + y) t R(i+x, j+y)= 0 , ( i + x , j + y) t 从定义中可以看出 : 0 nlike ( i , j ) 8 。 现在讨
4、论 nlike( i , j) 值的含义 。 (1) nlike (i , j) = 8 , 表示当前中心像素点的 8 邻域范围内都是与之相似的像素点 , 所以该像素点邻域范围内的梯度不会很大 , 因此角点检测时 , 应该排除此类像素点，不将其作为角点的候选点 。 (2) nlike (i , j) = 0 , 表示当前中心像素点的 8 邻域范围内没有与之相似的像素点 , 所以该像素点为孤立像素点或者是噪声点 , 因此角点检测时 , 也应该排除此类像素点 。 (3) nlike (i , j) = 7 , 可以归结为以下的两者情况 , 其他情形都可以通过旋转来得到 ( 图中黑色区域仅表示与中心
5、像素相似 , 而两个黑色区域像素可能是相似的 , 也可能不相似 ) 。 对于图 1 (a) 中 , 可能的角点应该是中心像素点的正上方的那个像素点 , 1(b) 图中可能的角点应该是中心像素点右上方的那个像素点 , 故这种情况下 , 中心像素点不应该作为角点的候选点 。 图 1 （ a ） 图 1(b) (4) nlike (i , j) = 1 , 可以归结为图 2 中的两种情况 ( 图中白色区域仅表示与中心像素不相似 , 而两个白色区域像素可能是相似的 , 也可能不相似 ) , 在这两种情况下 , 中心像素点也不可能为角点 。 图 2 （ a ） 图 2 （ b ） (5) 2 nlike
6、 ( i , j) 6 , 情况比较复杂 , 无法确认像素点准确的性质。我采取的方法是先将其列入候选角点之列，对其进行计算 CRF 等后续操作。 三、改进 Harris 算法效率性能分析 作为对比，我用未加改进的 Harris 算法计算了图像中的每一个像素点( 离边界 boundary 个像素之内的点除外，boundary 为程序中定义的一个整数，名同其意) 的 CRF 值。Harris 算法采用 9*9 的矩形高斯窗口( 共 81 个像素) , 方向导数采用 Prewitt 算子 , 则在计算每个像素的自相关矩阵 M 过程中乘法次数为: （9*9+1 ）*3=246 次，加法次数为: （9*
7、9-1 ）*3=240 次 ，因此，计算点（i,j ）的 CRF 就要花费（1 ）乘法 246+2=248 次，（2 ）加法 240+1=241 次，（3 ）行列式计算 det(M) 一次，（4 ）对角和计算 trace(M) 一次。对于分辨率为 Height*Width 的图像, 不考虑边界的影响, 其算法时间复杂度为 248* (Height-boundary)*(Width-boundary) 次乘法 + 241* (Height-boundary)*(Width-boundary) 次加法+1* (Height-boundary)* (Width-boundary)+ det(M) 一
8、次+ trace(M) 一次，由于一次乘法的时间要远多于一次加法 , 因此在实时性要求较高的场合 , Harris 算法就显得速度较慢了。 改进 Harris 算法主要是在计算 CRF 值之前对像素点做了一个初始选择 , 这一部分主要是加法运算 , 运算量为 8 *(Height-boundary)*(Width-boundary) 次 。由于不涉及乘法运算 , 虽然在此增加了计算量 , 但为后续的 焦点提取剔除了大量的非角点像素点 , 而不必计算每个像素点的 CRF 值 , 大大减少了乘法的计算次数 , 计算效率在分辨率高的图像中可以明显地体现出来。 对于我们所采集的图像 , 预处理后 ,
9、剔除像素点比例占总像素点的 60 % 以上 ，乘法计算量一般都能降低一半。表 1 为改进前的传统 Harris 算法与改进 后的 Harris 算法在对一组 6 幅图像进行角点提取时，于核心部分（计算CRF 与确认角点）所耗费的平均计算时间的比较。从表中可以看出改进算法的时间约为原算法的 21.58 % , 在不使用其他技巧的情况下，大大缩短了计算时间 , 有利于在实时处理中的应用。 表 1 Harris 算法与改进 Harris 算法在核心部分 平均消耗消耗时间上的对比( 单位 : 秒) 图像名 Z1.jpg Z2.jpg Z3.jpg Z4.jpg Z5.jpg Z6.jpg 改进前 7.
10、466712 7.456286 7.655664 7.581580 7.501732 7.469117 改进后 1.270748 1.241037 1.698029 2.272551 2.025614 1.230464 四、其他 这样检测出来的角点，有可能会在真正的角点附近小的邻域里同时存在几个“ 伪角点” 。 我采用的解决方法是，对以该点为中心 7*7 的邻域内进行搜索，看是否存在这样的伪角点；如果存在，则用算数平均的方法求出全部这些点的平均坐标 row_ave 和 column_ave （见程序 137 和 138 行），并以此作为真正的角点坐标并将其用绿色的点标志出来。 但是，这样处理方
11、法的一个缺点是，对于图像边缘上（程序会检测出比较多的“角点”）杂点比较多的地方，会由于“角点”的聚集而导致更多的角点出现。 与此同时，削弱了噪声的影响（噪声点的 CRF 一般会比较大），会使图像 CRFmax 值变小，这也是导致用这种方法产生更多角点的一个原因。 表 2 给出了传统 Harris 算法与改进 Harris 算法在结果性能（检测出的角点个数）方面的对比。表 3 给出了传统 Harris 算法与改进 Harris 算法在计算图像 CRFmax 上的对比。 表 2 Harris 算法与改进 Harris 算法在结果性能方面的 对比（检测出的角点个数）（单位：个） 图像名 Z1.jpg
12、 Z2.jpg Z3.jpg Z4.jpg Z5.jpg Z6.jpg 改进前 77 71 76 71 70 76 改进后 79 72 99 86 98 78 表 3 改进前后 CRFmax 对比（单位 1 ） 图像名Z1.jpg Z2.jpg Z3.jpg Z4.jpg Z5.jpg Z6.jpg 改进前3.2947e+08 2.8773e+08 1.5712e+09 1.3835e+09 1.5645e+08 2.4432e+09 改进后3.2947e+08 2.8773e+08 1.5598e+09 1.3389e+09 1.5645e+08 2.4432e+09 附：matlab 源代
13、码 %Prewitt Operator Corner Detection.m%时间优化-相邻像素用取差的方法%clear;for nfigure=1:6switch nfigure %选择图片case 1 t=z1.jpg;case 2 t=z2.jpg;case 3 t=z3.jpg;case 4 t=z4.jpg;case 5 t=z5.jpg;case 6 t=z6.jpg;end% t1 = tic; %测算时间FileInfo = imfinfo(t); % 保存图像的所有信息Image = imread(t); % 读取图像% 转为灰度值图像(Intensity Image)if(
14、strcmp(truecolor,FileInfo.ColorType) = 1) %转为灰度值图像Image = im2uint8(rgb2gray(Image); end dx = -1 0 1;-1 0 1;-1 0 1; %dx：横向 Prewitt 差分模版Ix2 = filter2(dx,Image).2; Iy2 = filter2(dx,Image).2; Ixy = filter2(dx,Image).*filter2(dx,Image);%生成 9*9 高斯窗口。窗口越大，探测到的角点越少。h= fspecial(gaussian,9,2); A = filter2(h,I
15、x2); % 用高斯窗口差分 Ix2 得到 A B = filter2(h,Iy2); C = filter2(h,Ixy); nrow = size(Image,1); ncol = size(Image,2); Corner = zeros(nrow,ncol); %矩阵 Corner 用来保存候选角点位置,初值全零，值为 1 的点是角点%真正的角点在 137 和 138 行由(row_ave,column_ave)得到%参数 t:点(i,j)八邻域的“相似度”参数，只有中心点与邻域其他八个点的像素值之差在%（-t,+t）之间，才确认它们为相似点，相似点不在候选角点之列t=20;%我并没有
16、全部检测图像每个点，而是除去了边界上 boundary 个像素，%因为我们感兴趣的角点并不出现在边界上boundary=8;for i=boundary:nrow-boundary+1 for j=boundary:ncol-boundary+1nlike=0; %相似点个数if Image(i-1,j-1)Image(i,j)-t end; endend; end; count = 0; % 用来记录角点的个数t=0.01; % 下面通过一个 3*3 的窗口来判断当前位置是否为角点for i = boundary:nrow-boundary+1 for j = boundary:ncol-b
17、oundary+1if Corner(i,j)=1 %只关注候选点的八邻域if CRF(i,j) t*CRFmax %这个是角点，count 加 1else % 如果当前位置（i,j）不是角点，则在 Corner(i,j)中删除对该候选角点的记录Corner(i,j) = 0; end;end; end; end; % disp(角点个数);% disp(count)figure,imshow(Image); % display Intensity Imagehold on; % toc(t1)for i=boundary:nrow-boundary+1 for j=boundary:ncol-boundary+1column_ave=0;row_ave=0;k=0;if Corner(i,j)=1for x=i-3:i+3 %7*7 邻域for y=j-3:j+3if Corner(x,y)=1% 用算数平均数作为角点坐标，如果改用几何平均数求点的平均坐标，对角点的提取意义不大row_ave=row_ave+x;column_ave=column_ave+y;k=k+1;endendendendif k0 %周围不止一个角点 plot( column_ave/k,row_ave/k ,r.);endend; end; end
