数字摄影测量之——Matlab实现Harris算子
程序员文章站
2024-01-19 19:15:16
...
Harris算子理论上来说是Moravec算子的改进,不过实际在求解过程中与Moravec算子的算法有很大不同。首先是使用高斯滤波,对整幅图像进行处处理,以达到降噪的目的。由于Matlab中可以很容易地实现滤波操作,因此在对一阶梯度进行求解时,可以使用[-1 0 1]和[-1 0 1]T的滤波核对整幅图像进行滤波处理,得到的即是一阶梯度的数据。之后根据得到的一阶梯度,求解图像的自相关矩阵M,通过自相关矩阵M的迹和行列式的值,计算响应因子R,从而根据响应因子R来判定其是否为角点。
具体步骤:
1、对整幅图像做高斯滤波进行降噪
2、使用-1 0 1]和[-1 0 1]T作为滤波核,分别求出图像的水平和竖直方向的一阶梯度。
3、基于上一步的梯度结果,计算自相关矩阵M,进而计算响应因子R
4、设定阈值,此处选择阈值T=40000
5、进行抑制局部非最大操作,将符合条件的点显示在图像上。
close all;
clear all;
clc;
img=(imread('lena512color.tif'));
img =double( rgb2gray(img));
[X,Y]=size(img);
img2 = imread ('lena512color.tif');
img2 = rgb2gray(img2);
figure;
subplot (1,2,1);
imshow(img2);
imgn=zeros(X,Y);
% 定义高斯滤波
sigma = 1.6;
gausFilter = fspecial('gaussian',5,sigma);
% blur = imfilter(img,gausFilter,'conv');
% 对图像求一阶梯度
fx=[-1,0,1];
Ix=filter2(fx,img);
fy=[-1;0;1];
Iy=filter2(fy,img);
Ix2=Ix.^2;
Iy2=Iy.^2;
Ixy=Ix.*Iy;
gx=filter2(gausFilter,Ix2);
gy=filter2(gausFilter,Iy2);
gxy=filter2(gausFilter,Ixy);
% 遍历所有像元,计算响应因子R
count1 = 0 ;
for x=1:X
for y=1:Y
% 构建自相关矩阵M
M=[gx(x,y),gxy(x,y);gxy(x,y),gy(x,y)];
% 计算响应算子
k = 0.04;
R = det( M )-k*trace(M)^2;
imgn(x,y)=R ;
if(R>0)
count1 = count1 + 1;
end
end
end
% %设阈值,小于均置零
% % % % % % 多种阈值的显示%%%%%%
sortedarray= sort(imgn(:),'descend');%降序排序
% T = mean(sortarray(1:(count1*0.9)));
% T = sum(sortarray(1:(count1*0.2)))/count1;%求前百分之20的均值
% T=0.01*max(max(imgn)); %最大值的0.01为阈值
T=40000;%手动经验阈值
% T = mean(imgn(:)); %整幅图像均值
% T = mean(sortedarray(1:count1));
% T = 0;
ind=find(imgn<T);
imgn(ind)=0;
subplot(1,2,2);
imshow(img2);hold on;
m=3;
count = 0;
for x=1+m:X-m %选局部最大(窗口为11*11)且非零值作为特征点
for y=1+m:Y-m
window2=imgn(x-m:x+m,y-m:y+m);
if max(window2(:))==imgn(x,y) && imgn(x,y)>0
plot(y,x,'+','color','red');
count=count+1;
end
end
end
程序结果如下:
上一篇: C++ STL中 Stack讲解
下一篇: Scrapy原理讲解