去除伪角点matlab,请教一下怎么把图像周围的角点给判断删掉

本帖最后由 没有圆周的钟 于 2018-10-27 10:10 编辑

clear;

clc;

I1= imread('D:\matlab练习\a.png');%读入灰度图像

I1=rgb2gray(I1);

[I2,rect]=imcrop(I1);%裁剪角点区域

ImageData=I2;

orimg=ImageData;

fx = [5 0 -5;8 0 -8;5 0 -5]; % 高斯函数一阶微分，x方向

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

fy = [5 8 5;0 0 0;-5 -8 -5];  % 高斯函数一阶微分，y方向

Iy = filter2(fy,orimg);  % 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(orimg,1);

width = size(orimg,2);

result = zeros(height,width);% 纪录角点位置

R = zeros(height,width);

Rmax = 0;

for i = 1:height

for j = 1:width

M = [Ix2(i,j) Ixy(i,j);Ixy(i,j) Iy2(i,j)];% 对M特征值进行分析

R(i,j) = det(M)-0.05*(trace(M))^2;% 角点响应函数

if R(i,j) > Rmax  %R取正整数为角点

Rmax = R(i,j);

end;

end;

end;

cnt = 0;

for i = 2:height-1

for j = 2:width-1

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

if R(i,j) > 0.005*Rmax && 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;

end;

end;

[y, x] = find(result == 1);

Rmax;

cnt; % 角点个数

imshow(orimg);

hold on;

plot(x,y,'r+');

%计算角点的亚像素坐标并输出结果

x_subpixel=zeros(1,cnt);

y_subpixel=zeros(1,cnt);

r=R;

for ii=1:cnt

pixelmatrix=[1,y(ii)-1,x(ii)-1,(y(ii)-1)^2,(y(ii)-1)*(x(ii)-1),(x(ii)-1)^2;

1,y(ii)-1, x(ii) ,(y(ii)-1)^2,  (y(ii)-1)*x(ii)  , x(ii)^2   ;

1,y(ii)-1,x(ii)+1,(y(ii)-1)^2,(y(ii)-1)*(x(ii)+1),(x(ii)+1)^2;

1,y(ii),x(ii)-1,y(ii)^2,y(ii)*(x(ii)-1),(x(ii)-1)^2;

1,y(ii), x(ii) ,y(ii)^2,  y(ii)*x(ii)  , x(ii)^2   ;

1,y(ii),x(ii)+1,y(ii)^2,y(ii)*(x(ii)+1),(x(ii)+1)^2;

1,y(ii)+1,x(ii)-1,(y(ii)+1)^2,(y(ii)+1)*(x(ii)-1),(x(ii)-1)^2;

1,y(ii)+1, x(ii) ,(y(ii)+1)^2,  (y(ii)+1)*x(ii)  , x(ii)^2   ;

1,y(ii)+1,x(ii)+1,(y(ii)+1)^2,(y(ii)+1)*(x(ii)+1),(x(ii)+1)^2;];

rmatrix = [r(y(ii)-1,x(ii)-1),r(y(ii)-1,x(ii)),r(y(ii)-1,x(ii)+1),r(y(ii),x(ii)-1),r(y(ii),x(ii)),r(y(ii),x(ii)+1),r(y(ii)+1,x(ii)-1),r(y(ii)+1,x(ii)),r(y(ii)+1,x(ii)+1)]';

amatrix = pixelmatrix\rmatrix; %矩阵除法求解

a0=amatrix(1,1);

a1=amatrix(2,1);

a2=amatrix(3,1);

a3=amatrix(4,1);

a4=amatrix(5,1);

a5=amatrix(6,1);

subx=(2*a2*a3-a1*a4)/(a4^2-4*a3*a5);

suby=(2*a1*a5-a2*a4)/(a4^2-4*a3*a5);  %计算每个像素点的亚像素坐标

x_subpixel(ii)=subx;

y_subpixel(ii)=suby;

disp([x(ii),y(ii)]); %像素级坐标

disp([x_subpixel(ii),y_subpixel(ii)]); %亚像素坐标

end

P=[x_subpixel(:),y_subpixel(:)];