harris角点检测改进算法

Harris角点检测算法优化

一、综述

用 Harris 算法进行检测，有三点不足：（1 ）该算法不具有尺度不变性；（2 ）该算法提取的角点是像素级的；（3 ）该算法检测时间不是很令人满意。

基于以上认识，我主要针对第（3 ）点对Harris 角点检测算法提出了改进。

二、改进 Harris 算法原理

在介绍我的方法之前，我先提出如下概念：图像区域像素的相似度。我们知道， Harris 角点检测是基于图像像素灰度值变化梯度的， 灰度值图像的角点附近，是其像素灰度值变化非常大的区域，其梯度也非常大。换句话说，在非角点位置邻域里，各点的像素值变化不大，甚至几乎相等，其梯度相对也比较小。从这个角度着眼，我提出了图像区域像素的相似度的概念，它是指检测窗口中心点灰度值与其周围n 邻域内其他像素点灰度值的相似程度，这种相似程度是用其灰度值之差来描述的。如果邻域内点的灰度值与中心点Image (i,j) 的灰度值之差 的绝对值 在一个阈值t 范围 内，那我就认为这个点与中心点是相似的。与此同时，属于该 Image (i,j) 点的相似点计数器nlike(i,j) 也随之加一。在 Image (i,j) 点的n 邻域 全部被遍历一边之后，就能得到 在这个邻域范围内 与中心点相似的 点个数的统计值nlike(i,j) 。根据nlike(i,j) 的大小，就可以判断这个中心点是否可能为角点。

由于我选择3*3 的检测窗口，所以， 对于中心像素点 , 在下面的讨论中只考虑其8 邻域内像素点的相似度。 计算该范围的像素点与中心像素点的灰度值之差的绝对值 ( 记为 Δ ) , 如果该值小于等于设定的阈值 ( 记为 t) , 则认为该像素点与目标像素点相

似 。

nlike(i,j)=sum(R(i+x,j+y))

( - 1 ≤ x ≤ 1 , - 1 ≤ y ≤ 1 , 且 x ≠ 0 , y ≠ 0) ,

其中 :

R(i+x, j+y)=1 , 当Δ( i + x , j + y) ≤ t

R(i+x, j+y)= 0 , 当Δ( i + x , j + y) > t

从定义中可以看出 : 0 ≤ nlike ( i , j ) ≤ 8 。 现在讨论 nlike( i , j) 值的含义 。

(1) nlike (i , j) = 8 , 表示当前中心像素点的 8 邻域范围内都是与之相似的像素点 , 所以该像素点邻域范围内的梯度不会很大 , 因此角点检测时 , 应该排除此类像素点，不将其作为角点的候选点 。

(2) nlike (i , j) = 0 , 表示当前中心像素点的 8 邻域范围内没有与之相似的像素点 , 所以该像素点为孤立像素点或者是噪声点 , 因此角点检测时 , 也应该排除此类像素点 。

(3) nlike (i , j) = 7 , 可以归结为以下的两者情况 , 其他情形都可以通过旋转来得到 ( 图中黑色区域仅表示与中心像素相似 , 而两个黑色区域像素可能是相似的 , 也可能不相似 ) 。 对于图 1 (a) 中 , 可能的角点应该是中心像素点的正上方的那个像素点 , 1(b) 图中可能的角点应该是中心像素点右上方的那个像素点 , 故这种情况下 , 中心像素点不应该作为角点的候选点 。

(4) nlike (i , j) = 1 , 可以归结为图 2 中的两种情况 ( 图中白色区域仅表示与中心像素不相似 , 而两个白色区域像素可能是相似的 , 也可能不相似 ) , 在这两种情况下 , 中心像素点也不可能为角点 。

(5) 2 ≤ nlike ( i , j) ≤ 6 , 情况比较复杂 , 无法确认像素点准确的性质。我采取的方法是先将其列入候选角点之列，对其进行计算CRF 等后续操作。

附：matlab源代码

%%%Prewitt Operator Corner Detection.m

%%%时间优化--相邻像素用取差的方法

%%

clear;

for nfigure=1:6

switch 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(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,Ix2); % 用高斯窗口差分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;

%我并没有全部检测图像每个点，而是除去了边界上boundary个像素，

%因为我们感兴趣的角点并不出现在边界上

boundary=8;

for i=boundary:nrow-boundary+1

for j=boundary:ncol-boundary+1

nlike=0; %相似点个数

if Image(i-1,j-1)>Image(i,j)-t && Image(i-1,j-1)nlike=nlike+1;

end

if Image(i-1,j)>Image(i,j)-t && Image(i-1,j)nlike=nlike+1;

end

if Image(i-1,j+1)>Image(i,j)-t && Image(i-1,j+1)nlike=nlike+1;

end

if Image(i,j-1)>Image(i,j)-t && Image(i,j-1)nlike=nlike+1;

end

if Image(i,j+1)>Image(i,j)-t && Image(i,j+1)nlike=nlike+1;

end

if Image(i+1,j-1)>Image(i,j)-t && Image(i+1,j-1)nlike=nlike+1;

end

if Image(i+1,j)>Image(i,j)-t && Image(i+1,j)nlike=nlike+1;

end

if Image(i+1,j+1)>Image(i,j)-t && Image(i+1,j+1)nlike=nlike+1;

end

if nlike>=2 && nlike<=6

Corner(i,j)=1;%如果周围有0，1，7，8个相似与中心的（i,j）

%那(i,j)就不是角点，所以，直接忽略

end;

end;

end;

CRF = zeros(nrow,ncol); % CRF用来保存角点响应函数值,初值全零

CRFmax = 0; % 图像中角点响应函数的最大值，作阈值之用

t=0.05;

% 计算CRF

%工程上常用CRF(i,j) =det(M)/trace(M)计算CRF，那么此时应该将下面第103行的

%比例系数t设置大一些，t=0.1对采集的这几幅图像来说是一个比较合理的经验值

for i = boundary:nrow-boundary+1

for j = boundary:ncol-boundary+1

if Corner(i,j)==1 %只关注候选点

M = [A(i,j) C(i,j);

C(i,j) B(i,j)];

CRF(i,j) = det(M)-t*(trace(M))^2;

if CRF(i,j) > CRFmax

CRFmax = CRF(i,j);

end;

end

end;

end;

count = 0; % 用来记录角点的个数

t=0.01;

% 下面通过一个3*3的窗口来判断当前位置是否为角点

for i = boundary:nrow-boundary+1

for j = boundary:ncol-boundary+1

if Corner(i,j)==1 %只关注候选点的八邻域

if CRF(i,j) > t*CRFmax && CRF(i,j) >CRF(i-1,j-1) ......

&& CRF(i,j) > CRF(i-1,j) && CRF(i,j) > CRF(i-1,j+1) ......

&& CRF(i,j) > CRF(i,j-1) && CRF(i,j) > CRF(i,j+1) ......

&& CRF(i,j) > CRF(i+1,j-1) && CRF(i,j) > CRF(i+1,j)......

&& CRF(i,j) > CRF(i+1,j+1)

count=count+1;%这个是角点，count加1

else % 如果当前位置（i,j）不是角点，则在Corner(i,j)中删除对该候选角点的记录

Corner(i,j) = 0;

end;

end;

end;

end;

% disp(''角点个数'');

% disp(count)

figure,imshow(Image); % display Intensity Image

hold on;

% toc(t1)

for i=boundary:nrow-boundary+1

for j=boundary:ncol-boundary+1

column_ave=0;

row_ave=0;

k=0;

if Corner(i,j)==1

for x=i-3:i+3 %7*7邻域

for y=j-3:j+3

if Corner(x,y)==1

% 用算数平均数作为角点坐标，如果改用几何平均数求点的平均坐标，对角点的提取意义不大

row_ave=row_ave+x;

column_ave=column_ave+y;

k=k+1;

end

end

end

end

if k>0 %周围不止一个角点

plot( column_ave/k,row_ave/k ,''r.'');

end

end;

end;

end