harris分块提取特征点（matlab）

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   Harris角点提取算法                                          %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% filename=‘G.tif‘;
% X是最后的结果影像 
%[posrposc]是兴趣点的坐标
%blocksize是我们要选择的分块大小
%filename是我们读取的原始影像%
function [Xposrposc]=Harris（blocksizefilename）
% clear;

 X= imread（filename）;     % 读取图像
Info=imfinfo（filename）;
filesize=length（X）/blocksize;
% if Info.BitDepth>8
%     f=rgb2gray（X）;
% f=X;
%end
cnt1 = 0;
for I=1:filesize
    for J=1:filesize
        f=X（blocksize*（I-1）+1:I*blocksizeblocksize*（J-1）+1:J*blocksize）;
%
% 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（fxori_im）;                % x方向滤波
% fy = [5 8 5;0 0 0;-5 -8 -5];          % 高斯函数一阶微分，y方向（用于改进的Harris角点提取算法）
fy = [-2;-1;0;1;2];                     % y方向梯度算子（用于Harris角点提取算法）
Iy = filter2（fyori_im）;                % y方向滤波
Ix2 = Ix.^2;
Iy2 = Iy.^2;
Ixy = Ix.*Iy;
clear Ix;
clear Iy;

h= fspecial（‘gaussian‘72）;        % 产生7*7的高斯窗函数，sigma=2

Ix2 = filter2（hIx2）;
Iy2 = filter2（hIy2）;
Ixy = filter2（hIxy）;

height = size（ori_im1）;
wi