% The nearest neighbor classifier with the features extracted by PCA,
% 2D-PCA, LDA AND 2D-LDA.
% 2DPCA,then return a matrix contined all the data 

p=120;                 % the output number of principal componants in 1D-PCA
d=12;                 % the output number of principal componants in 2D-PCA
L=6;                  % the output number of principal componants in LDA
l=12;
nx=10;                 % the training sample number of semi-special-person
np=10;                 % the number of person in the Imageset
nsig=11;              % the number of sigle expression of one person


% load sequencedata32
% [Npix,Nperson]=size(Imageset);
% npix=sqrt(Npix);


% load cohndetseq32
load staticdata32
% load grame8data32     % in this dataset detimage=abs(I-II) and the
%                         interval between frames is 8
% load middlenonabs      % in this dataset detimage=II-I
Imageset=data.X;
Imagelab=data.y;
[Npix,dummy]=size(Imageset);
npix=sqrt(Npix);

%% get the unitary image
% Imagemax=max(Imageset);
% Imageset=Imageset./(ones(Npix,1)*Imagemax);


%  %%%% special-person recognizition
% Xspec=zeros(Npix,nx*np*7);
% Yspec=zeros(Npix,(nsig-nx)*np*7);
% standard=kron(1:7,ones((nsig-nx)*np,1));
% 
% for i=1:np*7
%     Xspec(:,(i-1)*nx+1:i*nx)=Imageset(:,(i-1)*nsig+1:(i-1)*nsig+nx);
%     Yspec(:,(i-1)*(nsig-nx)+1:i*(nsig-nx))=Imageset(:,(i-1)*nsig+nx+1:i*nsig);
% end
% % %%% 1D-PCA
% U=pcafeature(Xspec,p,1);
% trainset=U'*Xspec;                     % get the pca training set
% testset=U'*Yspec;                      % the test image after translate
% 
% % %% PCA+LDA
% % Upca=pcafeature(Xspec,p,1);
% % Z=Upca'*Xspec;
% % label=kron(1:7,ones(1,nx*np));
% % data=struct('X',Z,'y',label);
% % U=ldafeature(data,L);
% % trainset=U'*Z;
% % testset=U'*(Upca'*Yspec);
% 
% detvect=sum((kron(ones(1,(nsig-nx)*np*7),trainset)-kron(testset,ones(1,nx*np*7))).^2);
% detmat=reshape(detvect,nx*np*7,(nsig-nx)*np*7);
% [dummy,labspec]=min(detmat);
% labspec=ceil(labspec/(nx*np));
% specialmat=reshape(labspec,(nsig-nx)*np,7);
% recgspec=sum(specialmat==standard);
% recgspecratio=recgspec/((nsig-nx)*np)
% total=mean(recgspecratio)


%%% do oneself
% Xset=zeros(Npix,6*(18-1));
% Yset=zeros(Npix,6);
% labcohn=zeros(18,6);
% labY=zeros(18,6);

%% get the unitary image
% Imagemax=max(Imageset);
% Imageset=Imageset./(ones(Npix,1)*Imagemax);

% for i=1:18
%     Xset=[Imageset(:,1:(i-1)*6) Imageset(:,i*6+1:end)];
%     labX=[Imagelab(1:(i-1)*6) Imagelab(i*6+1:end)];
%     Yset=Imageset(:,(i-1)*6+1:i*6);
%     labY(i,:)=Imagelab((i-1)*6+1:i*6);
%     
%     %% 2D-LDA
%     data=struct('X',Xset,'y',labX);
%     intialR=rand(npix);
%     options=struct('R0',intialR,'low1',l,'low2',l,'ITE',5);
%     [U,B,V]=twodlda_modify(data,options,2,1000);
%     for j=1:6
%         Y=reshape(Yset(:,j),npix,npix);
%         testimage=reshape(U'*Y*V,l*l,1);
%         distance=sum((testimage*ones(1,6*(18-1))-B).^2);
%         [dummy,labindex]=min(distance);
%         labcohn(i,j)=labX(labindex);
%     end
%     
% end
% 
% recgcohn=sum(labcohn==labY);
% recgcohnratio=recgcohn/18


%% non-special-person recgnizition by leave-one method
Xnon=zeros(Npix,nsig*7*(np-1));
Ynon=zeros(Npix,nsig*7);
labX=zeros(1,nsig*7*(np-1));
recgnon=zeros(np,7);
recgnon2=zeros(np,7);
recgnon2dlda=zeros(np,7);
recgnon3=zeros(np,7);
standard=kron(1:7,ones(nsig,1));
err=zeros(1,np);
distribution=zeros(7,7,np);
for i=1:np       % select the ith person to test
    for j=1:7
        Xnon(:,(j-1)*(np-1)*nsig+1:j*(np-1)*nsig)=[Imageset(:,(j-1)*np*nsig+1:(j-1)*np*nsig+(i-1)*nsig) Imageset(:,(j-1)*np*nsig+i*nsig+1:j*np*nsig)];
        labX((j-1)*(np-1)*nsig+1:j*(np-1)*nsig)=[Imagelab((j-1)*np*nsig+1:(j-1)*np*nsig+(i-1)*nsig) Imagelab((j-1)*np*nsig+i*nsig+1:j*np*nsig)];
        Ynon(:,(j-1)*nsig+1:j*nsig)=Imageset(:,(j-1)*np*nsig+(i-1)*nsig+1:(j-1)*np*nsig+i*nsig);        
    end

%      %%% 1D-PCA+LDA
%      Upca=pcafeature(Xnon,p,1);
%      Z=Upca'*Xnon;
%      data=struct('X',Z,'y',labX);
%      U=ldafeature(data,L);
%      % train and test **********************************
%      trainset=U'*Z;                     % get the pca training set
%      testset=U'*(Upca'*Ynon);                      % the test image after translate
%      labnon=zeros(1,nsig*7);
%      for j=1:nsig*7
%          detvect=sum((trainset-testset(:,j)*ones(1,(np-1)*nsig*7)).^2);
%          detvect=reshape(detvect,(np-1)*nsig,7);
%          [dist,dummy]=min(detvect);
%          [dummy,templabel]=min(dist);
%          labnon(j)=templabel;
%      end
%      nonmat=reshape(labnon,nsig,7);
%      recgnon(i,:)=sum(nonmat==standard);
%      %% calculate the distribution of recognition
%      for j=1:7
%          distribution(j,:,i)=sum(nonmat==j*ones(nsig,7));
%      end
     
%      %% 2D-PCA
%      labnon2=zeros(1,nsig*7);
%      U=pcafeature(Xnon,d,2);
%      trainset=zeros(npix,d,nsig*(np-1)*7);
%      Xnon2=reshape(Xnon,npix,npix,nsig*(np-1)*7);
%      Ynon2=reshape(Ynon,npix,npix,nsig*7);
%      for j=1:nsig*(np-1)*7
%          trainset(:,:,j)=Xnon2(:,:,j)*U;
%      end
%      for j=1:nsig*7
%          testimage=Ynon2(:,:,j)*U;
%          distmat=sum((reshape(trainset,npix,d*nsig*(np-1)*7)-kron(ones(1,nsig*(np-1)*7),testimage)).^2);
%          distvect=sum(reshape(distmat,d,nsig*(np-1)*7));
%          [dummy,templab]=min(distvect);
%          labnon2(j)=ceil(templab/(nsig*(np-1)));
%      end
%      nonmat2=reshape(labnon2,nsig,7);
%      recgnon2(i,:)=sum(nonmat2==standard);
%      %% calculate the distribution of recognition
%      for j=1:7
%          distribution(j,:,i)=sum(nonmat2==j*ones(nsig,7));
%      end
     
     %% 2D-LDA
     laby=kron(1:7,ones(1,nsig*(np-1)));
     intialR=rand(npix);
     data=struct('X',Xnon,'y',laby);
     options=struct('R0',intialR,'low1',l,'low2',l,'ITE',5);
     [U,B,V]=twodlda_modify(data,options,2,1000);
     labnon2dlda=zeros(1,nsig*7);
     for j=1:nsig*7
         Y=reshape(Ynon(:,j),npix,npix);
         testimage=reshape(U'*Y*V,l*l,1);
         distance=sum((testimage*ones(1,nsig*(np-1)*7)-B).^2);
         [dummy,templab]=min(distance);
         labnon2dlda(j)=ceil(templab/(nsig*(np-1)));
     end
     nonmat2d=reshape(labnon2dlda,nsig,7);
     recgnon2dlda(i,:)=sum(nonmat2d==standard);
     %% calculate the distribution of recognition
     for j=1:7
         distribution(j,:,i)=sum(nonmat2d==j*ones(nsig,7));
     end
     
% %  
% % 
%      % GLRAM
%      for j=1:nsig*(np-1)*7
%          A{j}=reshape(Xnon(:,j),npix,npix);
%      end
% %      L0=eye(npix);
%      L0=rand(npix,npix);
%      [U,trainset,V,err(i)]=GLRAM(A,L0,l,l,10);
%      labglram=zeros(1,nsig*7);
%      for j=1:nsig*7
%          testimage=U'*reshape(Ynon(:,j),npix,npix)*V;
%          distance=zeros(nsig*(np-1)*7,1);
%          for k=1:nsig*(np-1)*7
%              temp=(testimage-trainset{k}).^2;
%              distance(k)=sum(temp(:));
%          end
%          [dummy,templab]=min(distance);
%          labglram(j)=ceil(templab/(nsig*(np-1)));
%      end
%      nonmat3=reshape(labglram,nsig,7);
%      recgnon3(i,:)=sum(nonmat3==standard);
% %               
end
% 
% recgnonratio=mean(recgnon/nsig)
% recgnonratio2=mean(recgnon2/nsig)
distrib=sum(distribution,3)
% recgnonratioglram=mean(recgnon3/nsig)
recgnonratio2dlda=mean(recgnon2dlda/nsig)