clear all;
close all;
clc;

img=double(imread('lena.jpg'));
[m n]=size(img);

w=4;
sigma=1.2;
[x y]=meshgrid(-w:w,-w:w);
%高斯函数对应的二阶偏导，如何求导数请自行脑补
Dxx = 1/(2*pi*sigma^4)*(x.^2/sigma^2-1)*exp(-(x.^2+x.^2)/(2*sigma^2));     
Dyy = 1/(2*pi*sigma^4)*(y.^2/sigma^2-1)*exp(-(x.^2+y.^2)/(2*sigma^2));
Dxy = 1/(2*pi*sigma^6)*(x.*y)*exp(-(x.^2+y.^2)/(2*sigma^2));

Ixx=imfilter(img,Dxx,'replicate');
Iyy=imfilter(img,Dyy,'replicate');
Ixy=imfilter(img,Dxy,'replicate');

H=cell(m,n);
for i=1:m
   for j=1:n 
    H{i,j}=[Ixx(i,j) Ixy(i,j);Ixy(i,j) Iyy(i,j)];    
   end
end