function result=center_rectangle(sdata)
%求取圖像的中心矩

size(sdata);
[M,N]=size(sdata);
%***************
%m00,miu00,yibu00分別與公式中數字相互對應，并取其音譯
%***************
m00=0.0;
m10=0.0;
m01=0.0;

for i=1:M
    for j=1:N
        m00=m00+double(sdata(i,j));
        m10=m10+i*double(sdata(i,j));
        m01=m01+j*double(sdata(i,j));
    end
end
ave_x=m10/m00;
ave_y=m01/m00;
m00;
miu20=0.0;
miu02=0.0;
miu11=0.0;
miu12=0.0;
miu21=0.0;
miu30=0.0;
miu03=0.0;
%求取7個中心矩
for i=1:M
    for j=1:N
        miu20=miu20+(i-ave_x).^2.*(j-ave_y).^0.*double(sdata(i,j));
        miu02=miu02+(i-ave_x).^0.*(j-ave_y).^2.*double(sdata(i,j));
        miu11=miu11+(i-ave_x).^1.*(j-ave_y).^1.*double(sdata(i,j));
        miu12=miu12+(i-ave_x).^1.*(j-ave_y).^2.*double(sdata(i,j));
        miu21=miu21+(i-ave_x).^2.*(j-ave_y).^1.*double(sdata(i,j));
        miu30=miu30+(i-ave_x).^3.*(j-ave_y).^0.*double(sdata(i,j));
        miu03=miu03+(i-ave_x).^0.*(j-ave_y).^3.*double(sdata(i,j));
    end
end
yibu20=miu20/(m00.^((2+0)/2+1));
yibu02=miu02/(m00.^((0+2)/2+1));
yibu11=miu11/(m00.^((1+1)/2+1));
yibu12=miu12/(m00.^((1+2)/2+1));
yibu21=miu21/(m00.^((2+1)/2+1));
yibu30=miu30/(m00.^((3+0)/2+1));
yibu03=miu03/(m00.^((0+3)/2+1));

fai1=yibu20+yibu02;
fai2=(yibu20-yibu02).^2+4*yibu11.*yibu11;
fai3=(yibu30-3*yibu12).^2+(3*yibu21-yibu03).^2;
fai4=(yibu30+yibu12).^2+(yibu21+yibu03).^2;
fai5=(yibu30-3*yibu12).*(yibu30+yibu12).*((yibu30+yibu12)*(yibu30+yibu12)-3*(yibu21+yibu03)*(yibu21+yibu03))+(3*yibu21-yibu03)*(yibu21+yibu03)*(3*(yibu30+yibu12).^2-(yibu21+yibu03)*(yibu21+yibu03));
fai6=(yibu20-yibu02).*((yibu30+yibu12).^2-(yibu21+yibu03).^2)+4*yibu11*(yibu30+yibu12)*(yibu21+yibu03);
fai7=(3*yibu21-yibu03)*(yibu30+yibu12)*((yibu30+yibu12)^2-3*(yibu21+yibu03)*(yibu21+yibu03))+(3*yibu12-yibu30)*(yibu21+yibu03)*(3*(yibu30+yibu12)^2-(yibu21+yibu03)^2);

% fai1
% fai2
% fai3
% fai4
% fai5
% fai6
% fai7
format long g

result=[fai1,fai2,fai3,fai4,fai5,fai6,fai7];