终于找到了源代码了：Harris角点特征提取MATLAB代码_夏安

http://blog.sina.com.cn/u/1789812875

首页博文目录关于我

个人资料

微博

加好友发纸条

写留言加关注

博客等级：
博客积分：

博客访问：
关注人气：
获赠金笔：0支
赠出金笔：0支
荣誉徽章：

正文字体大小：大中小

终于找到了源代码了：Harris角点特征提取MATLAB代码

(2012-10-03 14:26:21)

标签：

高斯

提取

算法

梯度

图像

杂谈

分类：技术贴

Harris角点特征提取MATLAB代码

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% Harris角点提取算法 %

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear;

%filename = 'Lena.jpg';

filename='object0024.view01.png';

X = imread(filename); % 读取图像

% imshow(X);

Info = imfinfo(filename); %获取图像相关信息

if (Info.BitDepth > 8)

f = rgb2gray(X);

end

%《基于特征点的图像配准与拼接技术研究》

%计算图像亮度f(x,y)在点(x,y)处的梯度-----------------------------------------------

% fx = [5 0 -5;8 0 -8;5 0 -5]; % 高斯函数一阶微分，x方向(用于改进的Harris角点提取算法)

ori_im = double(f) / 255; %unit8转化为64为双精度double64

fx = [-2 -1 0 1 2]; % x方向梯度算子(用于Harris角点提取算法)

Ix = filter2(fx, ori_im); % x方向滤波

% fy = [5 8 5;0 0 0;-5 -8 -5]; % 高斯函数一阶微分，y方向(用于改进的Harris角点提取算法)

fy = [-2; -1; 0; 1; 2]; % y方向梯度算子(用于Harris角点提取算法)

Iy = filter2(fy, ori_im); % y方向滤波

%构造自相关矩阵---------------------------------------------------------------

Ix2 = Ix .^ 2;

Iy2 = Iy .^ 2;

Ixy = Ix .* Iy;

clear Ix;

clear Iy;

h= fspecial('gaussian', [7 7], 2); % 产生7*7的高斯窗函数，sigma=2

Ix2 = filter2(h,Ix2);

Iy2 = filter2(h,Iy2);

Ixy = filter2(h,Ixy);

%提取特征点---------------------------------------------------------------

height = size(ori_im, 1);

width = size(ori_im, 2);

result = zeros(height, width); % 纪录角点位置，角点处值为1

R = zeros(height, width);

Rmax = 0; % 图像中最大的R值

k = 0.06; %k为常系数，经验取值范围为0.04~0.06

for i = 1 : height

for j = 1 : width

M = [Ix2(i, j) Ixy(i, j); Ixy(i, j) Iy2(i, j)]; % auto correlation matrix

R(i,j) = det(M) - k * (trace(M)) ^ 2; % 计算R

if R(i,j) > Rmax

Rmax = R(i, j);

end;

%T = 0.01 * Rmax;%固定阈值，当R(i, j) > T时，则被判定为候选角点

T = 0.1 * Rmax;%固定阈值，当R(i, j) > T时，则被判定为候选角点

%在计算完各点的值后，进行局部非极大值抑制-------------------------------------

cnt = 0;

for i = 2 : height-1

for j = 2 : width-1

% 进行非极大抑制，窗口大小3*3

if (R(i, j) > T && R(i, j) > R(i-1, j-1) && R(i, j) > R(i-1, j) && R(i, j) > R(i-1, j+1) && R(i, j) > R(i, j-1) && ...

R(i, j) > R(i, j+1) && R(i, j) > R(i+1, j-1) && R(i, j) > R(i+1, j) && R(i, j) > R(i+1, j+1))

result(i, j) = 1;

cnt = cnt+1;

end;

i = 1;

for j = 1 : height

for k = 1 : width

if result(j, k) == 1;

corners1(i, 1) = j;

corners1(i, 2) = k;

i = i + 1;

end;

[posc, posr] = find(result == 1);

figure,imshow(ori_im);

hold on;

plot(posr, posc, 'r+');

阅读┊ 收藏 ┊ 喜欢 ▼ ┊打印┊举报/Report

前一篇：只有把自己打点好，才有以后

后一篇：zz：爱拍照，来做操——以横沔为例谈谈我是怎样拍古镇的

新浪BLOG意见反馈留言板　欢迎批评指正