%矢量X的一些特征值計算
figure(1);
sig=imread('lena1.bmp');
sig=double(sig)/255;
[m_sig,n_sig]=size(sig);
siz_word=16;
siz_book=1024;
imshow(sig);

load codebook16_1024;
%k=256;
%N=1024;

%重新排列圖像
num=m_sig/siz_word;
ss=siz_word*siz_word;
nn=num*num;
re_sig=[];
for i=1:m_sig
    for j=1:m_sig
        f1=floor(i./siz_word);
        m1=mod(i,siz_word);
        if m1==0
            m1=siz_word;
            f1=f1-1;
        end
        f2=floor(j./siz_word);
        m2=mod(j,siz_word);
        if m2==0
            m2=siz_word;
            f2=f2-1;
        end
        re_sig(num*f1+f2+1,siz_word*(m1-1)+m2)=sig(i,j);
    end 
end
x=re_sig';
HX=hadamard(k);
X=HX*x;

%搜索初始匹配碼字
tic
Sb=[];
D=[];
P=[];
for i=1:nn
    Sa=abs(X(1,i)-Y(1,1));
    p=1;
    for j=2:N
        arg=abs(X(1,i)-Y(1,j));
        if arg<Sa
            Sa=arg;
            p=j;
        end
    end
    Sb=[Sb,Sa];
    Dmin=0;
    for m=1:k
        Dmin=Dmin+(X(m,i)-Y(m,p))^2;
    end
    D=[D,Dmin];
    P=[P,p];
end

%編碼
for i=1:nn
    codenumber(i)=1;     
    T=P;
    for j=1:N
       if P(i)+j<=N
           sum6=0;
           for m=1:k
                sum6=sum6+(X(m,i)-Y(m,P(i)+j))^2;
                if D(i)<sum6
                    codenumber(i)=T(i);
                    break;                         
                end
           end
           if D(i)>sum6
                D(i)=sum6; 
                codenumber(i)=P(i)+j;
                T(i)=codenumber(i);
           else
                break;
           end
       else
           break;
       end
    end
       
   for j=1:N
       if P(i)-j>=1
           sum6=0;
           for m=2:k
                sum6=sum6+(X(m,i)-Y(m,P(i)-j))^2;
                if D(i)<sum6
                    codenumber(i)=T(i);
                    break;                         
                end
           end
           if D(i)>sum6
                D(i)=sum6; 
                codenumber(i)=P(i)-j;
                T(i)=codenumber(i);
           else
                break;
           end
       else
           break;
       end
    end
end
toc
t=toc

%figure(2);
%imshow(sig);