zhuchuanwei
2008-04-03, 14:57
function [H,K]=compute_HKz(Z)
A=(1/70)*[2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2];
B=(1/70)*[2 2 2 2 2;
-1 -1 -1 -1 -1;
-2 -2 -2 -2 -2;
-1 -1 -1 -1 -1;
2 2 2 2 2];
C=(1/100)*[-4 -2 0 2 4;
-2 -1 0 1 2;
0 0 0 0 0;
2 1 0 -1 -2;
4 2 0 -2 -4];
D=(1/50)*[-2 -1 0 1 2;
-2 -1 0 1 2;
-2 -1 0 1 2;
-2 -1 0 1 2;
-2 -1 0 1 2];
E=(1/50)*[2 2 2 2 2;
1 1 1 1 1;
0 0 0 0 0;
-1 -1 -1 -1 -1;
-2 -2 -2 -2 -2];
F=(1/175)*[-13 2 7 2 -13;
2 17 22 17 2;
7 22 27 22 7;
2 17 22 17 2;
-13 2 7 2 -13];
A_coef=filter2(A,Z);
B_coef=filter2(B,Z);
C_coef=filter2(C,Z);
D_coef=filter2(D,Z);
E_coef=filter2(E,Z);
F_coef=filter2(F,Z);
[M,N]=size(A_coef);
z_xy=C_coef;
z_xx=2*A_coef;
z_yy=2*B_coef;
for i=1:M
for j=1:N
z_x(i,j)=2*A_coef(i,j)*i+C_coef(i,j)*j+D_coef(i,j);
z_y(i,j)=2*B_coef(i,j)*j+C_coef(i,j)*i+E_coef(i,j);
z(i,j)=A_coef(i,j)*i^(2)+ B_coef(i,j)*j^(2)+ C_coef(i,j)*i*j+ D_coef(i,j)*i+E_coef(i,j)*j+F_coef(i,j);
err(i,j)=uint8(Z(i,j))-uint8(z(i,j));%%%
H(i,j)=((1+z_x(i,j)^(2))*z_yy(i,j)+(1+z_y(i,j)^(2))*z_xx(i,j)-2*z_x(i,j)*z_y(i,j)*z_xy(i,j))/(2*(1+z_x(i,j)^(2)+z_y(i,j)^(2))^(3/2));
K(i,j)=(z_xx(i,j)*z_yy(i,j)-z_xy(i,j)^(2))/((1+z_x(i,j)^(2)+z_y(i,j)^(2))^(2));
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
err_sum=zeros(M,N);
temp1=zeros(5,5);
temp2=zeros(5,5);
for i=3:M-2
for j=3:N-2
err_sum(i,j)=sum(sum(err(i-2:i+2,j-2:j+2)));
end
end
err_threshold=100;%需要调整
for i=1:M
for j=1:N
if err_sum(i,j)>=err_threshold%err_threshold
temp1=err_sum(i-2:i+2,j-2:j+2);
[minValue,minIndex]=min(temp1(:));
[m,n]=ind2sub([5 5],minIndex);
u=3-m;
v=3-n;
temp2=Z(i-u-2:i-u+2,j-v-2:j-v+2);%出现问题[这里有错误]
z_x_H=(1/700)*[40*u-28*v-28,-20*u-14*v-14,-40*u,-20*u+14*v+14,40*u+28*v+28;
40*u-14*v-28,-20*u-7*v-14,-40*u,-20*u+7*v+14,40*u+14*v+28;
40*u-28,-20*u-14,-40*u,-20*u+14,40*u+28;
40*u+14*v-28,-20*u+7*v-14,-40*u,-20*u-7*v+14,40*u-14*v+28;
40*u+28*v-28,-20*u+14*v-14,-40*u,-20*u-14*v+14,40*u-28*v+28];
z_y_H=(1/700)*[-28*u+40*v+28,-14*u+40*v+28,40*v+28,14*u+40*v+28,28*u+40*v+28;
-14*u-20*v+14,-7*u-20*v+14,-20*v+14,7*u-20*v+14,14*u-20*v+14;
-40*v,-40*v,-40*v,-40*v,-40*v;
14*u-20*v-14,7*u-20*v-14,-20*v-14,-7*u-20*v-14,-14*u-20*v-14;
28*u+40*v-28,14*u+40*v-28,40*v-28,-14*u+40*v-28,-28*u+40*v-28];
z_xx_H=(1/35)*[2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2];
z_yy_H=(1/35)*[2 2 2 2 2;
-1 -1 -1 -1 -1;
-2 -2 -2 -2 -2;
-1 -1 -1 -1 -1;
2 2 2 2 2];
z_xy_H=(1/100)*[-4 -2 0 2 4;
-2 -1 0 1 2;
0 0 0 0 0;
2 1 0 -1 -2;
4 2 0 -2 -4];
z_x(i,j)=sum(sum(z_x_H.*temp2));
z_y(i,j)=sum(sum(z_y_H.*temp2));
z_xx(i,j)=sum(sum(z_xx_H.*temp2));
z_yy(i,j)=sum(sum(z_yy_H.*temp2));
z_xy(i,j)=sum(sum(z_xy_H.*temp2));
H(i,j)=((1+z_x(i,j)^(2))*z_yy(i,j)+(1+z_y(i,j)^(2))*z_xx(i,j)-2*z_x(i,j)*z_y(i,j)*z_xy(i,j))/(2*(1+z_x(i,j)^(2)+z_y(i,j)^(2))^(3/2));
K(i,j)=(z_xx(i,j)*z_yy(i,j)-z_xy(i,j)^(2))/((1+z_x(i,j)^(2)+z_y(i,j)^(2))^(2));
end
end
end
A=(1/70)*[2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2];
B=(1/70)*[2 2 2 2 2;
-1 -1 -1 -1 -1;
-2 -2 -2 -2 -2;
-1 -1 -1 -1 -1;
2 2 2 2 2];
C=(1/100)*[-4 -2 0 2 4;
-2 -1 0 1 2;
0 0 0 0 0;
2 1 0 -1 -2;
4 2 0 -2 -4];
D=(1/50)*[-2 -1 0 1 2;
-2 -1 0 1 2;
-2 -1 0 1 2;
-2 -1 0 1 2;
-2 -1 0 1 2];
E=(1/50)*[2 2 2 2 2;
1 1 1 1 1;
0 0 0 0 0;
-1 -1 -1 -1 -1;
-2 -2 -2 -2 -2];
F=(1/175)*[-13 2 7 2 -13;
2 17 22 17 2;
7 22 27 22 7;
2 17 22 17 2;
-13 2 7 2 -13];
A_coef=filter2(A,Z);
B_coef=filter2(B,Z);
C_coef=filter2(C,Z);
D_coef=filter2(D,Z);
E_coef=filter2(E,Z);
F_coef=filter2(F,Z);
[M,N]=size(A_coef);
z_xy=C_coef;
z_xx=2*A_coef;
z_yy=2*B_coef;
for i=1:M
for j=1:N
z_x(i,j)=2*A_coef(i,j)*i+C_coef(i,j)*j+D_coef(i,j);
z_y(i,j)=2*B_coef(i,j)*j+C_coef(i,j)*i+E_coef(i,j);
z(i,j)=A_coef(i,j)*i^(2)+ B_coef(i,j)*j^(2)+ C_coef(i,j)*i*j+ D_coef(i,j)*i+E_coef(i,j)*j+F_coef(i,j);
err(i,j)=uint8(Z(i,j))-uint8(z(i,j));%%%
H(i,j)=((1+z_x(i,j)^(2))*z_yy(i,j)+(1+z_y(i,j)^(2))*z_xx(i,j)-2*z_x(i,j)*z_y(i,j)*z_xy(i,j))/(2*(1+z_x(i,j)^(2)+z_y(i,j)^(2))^(3/2));
K(i,j)=(z_xx(i,j)*z_yy(i,j)-z_xy(i,j)^(2))/((1+z_x(i,j)^(2)+z_y(i,j)^(2))^(2));
end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
err_sum=zeros(M,N);
temp1=zeros(5,5);
temp2=zeros(5,5);
for i=3:M-2
for j=3:N-2
err_sum(i,j)=sum(sum(err(i-2:i+2,j-2:j+2)));
end
end
err_threshold=100;%需要调整
for i=1:M
for j=1:N
if err_sum(i,j)>=err_threshold%err_threshold
temp1=err_sum(i-2:i+2,j-2:j+2);
[minValue,minIndex]=min(temp1(:));
[m,n]=ind2sub([5 5],minIndex);
u=3-m;
v=3-n;
temp2=Z(i-u-2:i-u+2,j-v-2:j-v+2);%出现问题[这里有错误]
z_x_H=(1/700)*[40*u-28*v-28,-20*u-14*v-14,-40*u,-20*u+14*v+14,40*u+28*v+28;
40*u-14*v-28,-20*u-7*v-14,-40*u,-20*u+7*v+14,40*u+14*v+28;
40*u-28,-20*u-14,-40*u,-20*u+14,40*u+28;
40*u+14*v-28,-20*u+7*v-14,-40*u,-20*u-7*v+14,40*u-14*v+28;
40*u+28*v-28,-20*u+14*v-14,-40*u,-20*u-14*v+14,40*u-28*v+28];
z_y_H=(1/700)*[-28*u+40*v+28,-14*u+40*v+28,40*v+28,14*u+40*v+28,28*u+40*v+28;
-14*u-20*v+14,-7*u-20*v+14,-20*v+14,7*u-20*v+14,14*u-20*v+14;
-40*v,-40*v,-40*v,-40*v,-40*v;
14*u-20*v-14,7*u-20*v-14,-20*v-14,-7*u-20*v-14,-14*u-20*v-14;
28*u+40*v-28,14*u+40*v-28,40*v-28,-14*u+40*v-28,-28*u+40*v-28];
z_xx_H=(1/35)*[2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2;
2 -1 -2 -1 2];
z_yy_H=(1/35)*[2 2 2 2 2;
-1 -1 -1 -1 -1;
-2 -2 -2 -2 -2;
-1 -1 -1 -1 -1;
2 2 2 2 2];
z_xy_H=(1/100)*[-4 -2 0 2 4;
-2 -1 0 1 2;
0 0 0 0 0;
2 1 0 -1 -2;
4 2 0 -2 -4];
z_x(i,j)=sum(sum(z_x_H.*temp2));
z_y(i,j)=sum(sum(z_y_H.*temp2));
z_xx(i,j)=sum(sum(z_xx_H.*temp2));
z_yy(i,j)=sum(sum(z_yy_H.*temp2));
z_xy(i,j)=sum(sum(z_xy_H.*temp2));
H(i,j)=((1+z_x(i,j)^(2))*z_yy(i,j)+(1+z_y(i,j)^(2))*z_xx(i,j)-2*z_x(i,j)*z_y(i,j)*z_xy(i,j))/(2*(1+z_x(i,j)^(2)+z_y(i,j)^(2))^(3/2));
K(i,j)=(z_xx(i,j)*z_yy(i,j)-z_xy(i,j)^(2))/((1+z_x(i,j)^(2)+z_y(i,j)^(2))^(2));
end
end
end