tlda_1.m

人脸识别中TensorLDA算法的matlab程序

人脸识别中TesorLDAA算法Matlab程序
clc;clear;
cd E:\人脸识别\ORL人脸图库\样本：测试集\5_5\train
   fdir=dir('*.bmp');
   fnum1=size(fdir,1);
   nbcol=255;
   class=15;
%    M_Xtran=0;
   for     i=1:fnum1
              imgname=fdir(i,1).name;
              d=imread(imgname);
              a=double(d);
              [m,n]=size(a);            
              Xtran(:,:,i)=a;              
%              Xtran(:,i)=reshape(a,m*n,1);
%              M_Xtran=M_Xtran+a;
   end
     M_Xtran=sum(Xtran(:,:,:),3);
     M_X_mean=M_Xtran/fnum1;%求均值矩阵M（X）
     
 
     
  %计算类内均值矩阵M（X）

    for     i=1:class
              Xtran_class_mean(:,:,i)=sum(Xtran(:,:,5*(i-1)+1:5*(i-1)+5),3)/5;
    end
  
  %计算Sw,Sv
  Sw_U=0;
  Sb_U=0;
  Sw_V=0;
  Sb_V=0;
    for     i=1:class
             d=Xtran_class_mean(:,:,i);
        for   j=5*(i-1)+1:5*(i-1)+5
                 h=Xtran(:,:,j);
                 Sw_U=Sw_U+(h-d)'*(h-d);
                 Sw_V=Sw_V+(h-d)*(h-d)';
       end
    end
  C_k=fnum1/class;
    for    i=1:class
              d= Xtran_class_mean(:,:,i);
              Sb_U=Sb_U+C_k*(d-M_X_mean)'*(d-M_X_mean);
              Sb_V=Sb_V+C_k*(d-M_X_mean)*(d-M_X_mean)';
    end
    
    %计算投影矩阵U，V
    SS_V=inv(Sw_V)*(Sb_V);
    SS_U=inv(Sw_U)*(Sb_U);
    [U11,D11]=eig(SS_V);
    [V11,E11]=eig(SS_U);
    u1=U11(:,1);
    v1=V11(:,1);
    d1=D11(:,1);
    e1=E11(:,1);
    U=u1;
    V=v1;
    D=d1;
    E=e1;
    Au=eye(m);
    Av=eye(n);
    for    i=2:m
            Bv=U'*inv(Sw_V)*U;
            Mv=(Au-inv(Sw_V)*U*inv(Bv)*U')*SS_V;
            [UKK,DKK]=eig(Mv);
            uk=UKK(:,1);
            dk=DKK(:,1);            
            U=[U uk];
            D=[D dk];
    end
     for   i=2:n
              Bu=V'*inv(Sw_U)*V;
              Mu=(Av-inv(Sw_U)*V*inv(Bu)*V')*SS_U;
              [VKK,EKK]=eig(Mu);
              vk=VKK(:,1);
              ek=EKK(:,1);
              V=[V vk];
              E=[E ek];
     end
     
%降维     
   Tr_D=sum(D(1,:));
    for   i=1:m
            Tr_d=sum(D(1,1:i));
            Per_d_D=Tr_d/Tr_D;
          if Per_d_D>=0.90
              k=i;
              break,end
    end
      
    Tr_E=sum(E(1,:));
     for   i=1:n
             Tr_e=sum(E(1,1:i));
             Per_e_E=Tr_e/Tr_E;
           if  Per_e_E>=0.90
               l=i;
               break,end
     end
    u=U(:,1:k);
    v=V(:,1:l);    
     
% % %查看矩阵U，V的正交性
     S_U_1=U*U';
     S_V_1=V*V';

% %还原图像检测投影矩阵的可行性
%     figure(1)
%     t1=Xtran(:,:,1); 
%     imshow(uint8(t1));
%     grid on
%     figure(2)
%     t2=U*U'*t1*V*V';
%     imshow(uint8(t2));
%     grid on
%     figure(3)
%     t3=u*u'*t1*v*v';
%     imshow(uint8(t3));
%     grid on


% %叠代
% %计算Sw,Sv
  Sw_U=0;
  Sb_U=0;
  Sw_V=0;
  Sb_V=0;
    for     i=1:class
             d= Xtran_class_mean(:,:,i);
       for   j=5*(i-1)+1:5*(i-1)+5
                 h=Xtran(:,:,j); 
                 Sw_U=Sw_U+(h-d)'*u*u'*(h-d);
                 Sw_V=Sw_V+(h-d)*v*v'*(h-d)';
       end
    end
  C_k=fnum1/class;
    for    i=1:class
              d= Xtran_class_mean(:,:,i);
              Sb_U=Sb_U+C_k*(d-M_X_mean)'*u*u'*(d-M_X_mean);
              Sb_V=Sb_V+C_k*(d-M_X_mean)*v*v'*(d-M_X_mean)';
    end
 % %二次计算投影矩阵U，V
    SS_V=inv(Sw_V)*(Sb_V);
    SS_U=inv(Sw_U)*(Sb_U);
    [U11,D11]=eig(SS_V);
    [V11,E11]=eig(SS_U);
    u1=U11(:,1);
    v1=V11(:,1);
    d1=D11(:,1);
    e1=E11(:,1);
%     fprintf('d1=%d\n',d1);
%     fprintf('e1=%d\n',e1);
    U=u1;
    V=v1;
    D=d1;
    E=e1;
    Au=eye(m);
    Av=eye(n);
    for    i=2:m
            Bv=U'*inv(Sw_V)*U;
            Mv=(Au-inv(Sw_V)*U*inv(Bv)*U')*SS_V;
            [UKK,DKK]=eig(Mv);
            uk=UKK(:,1);
            dk=DKK(:,1);            
            U=[U uk];
            D=[D dk];
    end
     for   i=2:n
              Bu=V'*inv(Sw_U)*V;
              Mu=(Av-inv(Sw_U)*V*inv(Bu)*V')*SS_U;
              [VKK,EKK]=eig(Mu);
              vk=VKK(:,1);
              ek=EKK(:,1);
              V=[V vk];
              E=[E ek];
     end
     
     
%      Tr_d1_D=trace(d1)/trace(D);
%      Tr_e1_E=trace(e1)/trace(E);
%  %查看矩阵U，V的正交性
      S_U_2=U*U';
      S_V_2=V*V';
% % 
% % %还原图像检测投影矩阵的可行性   
%     t4=U*U'*t1*V*V';
%     figure(4)
%     imshow(uint8(t4));
%     grid on


% % 降维

% %选取占特征值综合的X%的特征值对应的特征向量组成投影矩阵,
 
   Tr_D=sum(D(1,:));
    for   i=1:m
            Tr_d=sum(D(1,1:i));
            Per_d_D=Tr_d/Tr_D;
          if Per_d_D>=0.90
              k=i;
              break,end
    end
      
    Tr_E=sum(E(1,:));
     for   i=1:n
             Tr_e=sum(E(1,1:i));
             Per_e_E=Tr_e/Tr_E;
           if  Per_e_E>=0.90
               l=i;
               break,end
     end
    u=U(:,1:k);
    v=V(:,1:l);
    
% % %查看矩阵u，v的正交性
%       S_u=u*u';
%       S_v=v*v';
% 
% % %还原图像检测投影矩阵的可行性   
%     t5=u*u'*t4*v*v';
%     figure(5)
%     imshow(uint8(t5));
%     grid on
%      figure(5)
%     z=u*u'*M_X_mean*v*v';
%     imshow(uint8(z));
%     grid on



%把图像矩阵影射到张量子空间中
   for  i=1:fnum1
            d=Xtran(:,:,i);
            e=u'*d*v;
            Ytran(:,:,i)=e;
   end
   
% % %查看投影到张量子空间后的图像形状
%    c=Ytran(:,1);
%    d=reshape(c,k,l);
%    figure(6)
%    imshow(uint8(d));
%    grid on

%计算投影到张量子空间中的图像类内均值
%   for  i=1:class
%              Ytran_class_mean(:,:,i)=sum(Ytran(:,:,10*(i-1)+1:10*(i-1)+10),3)/10;
%   end
   


%人脸识别
cd E:\人脸识别\ORL人脸图库\样本：测试集\5_5\test
     fdir=dir('*.bmp');
     fnum2=size(fdir,1);
   for     i=1:fnum2
              imgname=fdir(i,1).name;
              d=imread(imgname);
              a=double(d);
              [m,n]=size(a);
              Xtest(:,:,i)=a;
    end
    for    i=1:fnum2
             f=Xtest(:,:,i);
             g=u'*f*v;
             Ytest(:,:,i)=g;
    end
    
    
    correct=0;
    error=0;
   for     i=1:fnum2
             h=Ytest(:,:,i);
         for  j=1:fnum1
                  f=Ytran(:,:,j);                 
%                  dis(i,j) =sqrt(sum(sum((f-h)*(f-h)')));
%                   r=(f-h)*(f-h)';
%                   [P,Q]=eig(r);
%                   dis(i,j)=sqrt(trace(Q*Q'));
% %识别率为X%=95%,  Y=163,Z=81.5
                   dis(i,j) = sqrt(sum(sum((f-h).^2)));
  %识别率为 X%=95%,Y=168 ,Z=84
  %         X%=90%,Y=168,Z=84
         end
   end
    for i=1:fnum2
       [mind,ind] = min(dis');
       if ceil(ind(i)/5) == ceil(i/5)
            correct = correct + 1;
       else
          error = error + 1;
       end
    end

correct_detection=correct/(correct+error);
 
fprintf('correct=%d\n',correct);
fprintf('error=%d\n',error);