%最小區域法求平面度誤差
function flatness=flat(X,Y,Z) %X,Y,Z是反算出來的世界坐標系的坐標值

z_abs=[ 0   2.66054447411648   9.47025280660865   9.02753734079175  14.32252332887071    18.68882664048553  18.68882664048553  20.62367902623812  23.13614541603834  26.009511403437125]

i_fenxi=6;

path_name1=strcat('E:\vision3d_2008\vision_20080822_02\result2\');
path_name1=strcat('E:\vision3d_2008\vision_3d_data\result\');
 
filename1=strcat(path_name1,int2str(i_fenxi),'\x3d.txt');
x3d=load(filename1);

filename1=strcat(path_name1,int2str(i_fenxi),'\y3d.txt');
y3d=load(filename1);

filename1=strcat(path_name1,int2str(i_fenxi),'\z3d.txt');
z3d=load(filename1);

X=x3d(:);
Y=y3d(:);
Z=z3d(:);
 
S1=ones(841,1);
S=[X,Y,S1];
A=inv(S'*S)*S'*Z;
W_matrix_x=reshape(X,29,29);
W_matrix_y=reshape(Y,29,29);
W_matrix_z=reshape(Z,29,29);
for i=1:29
    for j=1:29
        distance(i,j)=(W_matrix_z(i,j)-A(1)*W_matrix_x(i,j)-A(2)*W_matrix_y(i,j)-A(3))/sqrt(A(1)^2+A(2)^2+1);
    end
end

if 1
    % Z=Ax+By+C
    l1=length(X);
    M1=zeros(l1,1)+1;
    XY1=[ X Y M1]
    ABC=XY1\Z;

end

clear i j;
[a,b]=max(distance);%For matrices,MAX(X) is a row vector containing the maximum element from each column
[c,d]=max(a);
maxd=c;
clear a b c d;
[a,b]=min(distance);%For matrices,MIN(X) is a row vector containing the minimum element from each column
[c,d]=min(a);
mind=c;
flatness=maxd-mind;

lx=29;ly=29
stpx=1:1:lx;
stpy=1:1:ly;
[X_axi,Y_axi]=meshgrid(stpy,stpx);
surface(X_axi,Y_axi,distance)
distance;
end