Harris角点检测

一、harris角点检测的原理和过程
参考课件：http://download.csdn.net/detail/zhangchen1003/9071467

---

二、代码实现

clear;clc;I=imread('10.tif');I=double(I);%计算差分Ix和Iysobelx=[-1 0 1;-2 0 2; -1 0 1];sobely=[-1 -2 -1;0 0 0;1 2 1];h=fspecial('gaussian',[3,3],2);Ix=imfilter(I,sobelx);Iy=imfilter(I,sobely);Ix2=Ix.^2;Iy2=Iy.^2;Ixy=Ix.*Iy;%乘以窗口Ix2=imfilter(Ix2,h);Iy2=imfilter(Iy2,h);Ixy=imfilter(Ixy,h);%[m,n]=size(I);R=zeros(m,n);Rmax=0;for i=1:m    for j=1:n        %构造矩阵        M=[Ix2(i,j),Ixy(i,j);Ixy(i,j),Iy2(i,j)];        R(i,j)=det(M)-0.04*(trace(M))^2;%0.04代表一个参数        if R(i,j)>Rmax            Rmax=R(i,j);        end    endend%局部非极大值抑制处理cnt=0;mark=zeros(m,n);for i=2:m-1    for j=2:n-1        if R(i,j)>0.01*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)               mark(i,j)=1;               cnt=cnt+1;        end    endend[x,y]=find(mark==1);figure,imshow(uint8(I));hold onplot(y,x,'w*');hold offcnt