pca.m

matlab获取PCA特征的小程序

% this function is used to obtain avereface and eigenface from the training set

function s = pca(sig)

affile = strcat('E:\FYC_FERET\datas\para\pcaaf_',sig,'.dat');
effile = strcat('E:\FYC_FERET\datas\para\pcaef_',sig,'.dat');
enfile = strcat('E:\FYC_FERET\datas\para\pcafg_',sig,'.dat');
afimg  = strcat('E:\FYC_FERET\datas\para\pcaaf_',sig,'.bmp');

mask  = 'E:\FYC_FERET\pimages\mask.bmp';
mimg  = imread(mask,'bmp');

imgpath = strcat('E:\FYC_FERET\pimages\para_',sig,'\');
inum = 456;

dim = 3648;
num = inum;
srcface=zeros(num,dim);
averface=zeros(1,dim);
differface=zeros(num,dim);
midmatrix=zeros(num,num);
mideigenvec=zeros(num,num);
mideigenval=zeros(num,num);
eigenface=zeros(dim,num);

%create the file path and name and read out source data
   for j=1:inum
      j
      filename=int2str(j);
      imnam = strcat(imgpath,filename,'.bmp')   
      immat = imread(imnam,'bmp');
      [imhei,imwid,cdim] = size(immat);
      immat = double(immat)/255;
      count = 0;
      for ii=1:imhei
         for jj=1:imwid
            if ~( immat(ii,jj,1)==0 & immat(ii,jj,2)==0 & immat(ii,jj,3)==1 )
               count = count+1;
               srcface(j,count) = immat(ii,jj,1);
            end
         end
      end
      if count~=dim
         fprintf(1,'mask error...');
         [ count dim ]
         pause;
      end
   end

%caculate averface and display and save the image
averface = sum(srcface);
averface = averface/num;
averfaceimg = double(mimg)/255;
count = 0;
for ii=1:imhei
   for jj=1:imwid
      if ~( mimg(ii,jj,1)==0 & mimg(ii,jj,2)==0 & mimg(ii,jj,3)==255 )
         count = count+1;
         averfaceimg(ii,jj,1) = averface(count);
         averfaceimg(ii,jj,2) = averface(count);
         averfaceimg(ii,jj,3) = averface(count);
      end
   end
end
imwrite(averfaceimg,afimg,'bmp');
figure;
imshow(averfaceimg);
fprintf(1,'Aver face is obtained...');
%pause;

% calculate dispersion matrix of samples
for i=1:num
   differface(i,:)=srcface(i,:)-averface;
end
midmatrix=differface*differface';

% obtain eigenface
[mideigenvec, LATENT, EXPLAINED] = PCACOV(midmatrix);
eigenface=differface'*mideigenvec;
    
%show eigenface here for check
figure;
title('First 40 eigenface');
eigenfaceimg = mimg;
neigen = min(min(eigenface));
xeigen = max(max(eigenface));
for i=1:40
   i
   eigenfaceimg=double(mimg)/255;
   count = 0;
   for ii=1:imhei
   	for jj=1:imwid
      	if ~( mimg(ii,jj,1)==0 & mimg(ii,jj,2)==0 & mimg(ii,jj,3)==255 )
         	count = count+1;
         	eigenfaceimg(ii,jj,1) = (eigenface(count,i)-neigen)/(xeigen-neigen);
	         eigenfaceimg(ii,jj,2) = (eigenface(count,i)-neigen)/(xeigen-neigen);
   	      eigenfaceimg(ii,jj,3) = (eigenface(count,i)-neigen)/(xeigen-neigen);
      	end
      end
   end
   subplot(5,8,i);
   imshow(eigenfaceimg);
end 
fprintf(1,'Eigen face is obtained...');
%pause;

%save averface
fid=fopen(affile,'w+t');
for i=1:dim
   i
   fprintf(fid,'%12.6f ',averface(i));
end
fclose(fid);
%save eigen faces
fid=fopen(effile,'w+t');
for i=1:dim
   for j=1:num
      [i j]
      fprintf(fid,'%12.6f ',eigenface(i,j));
   end
   fprintf(fid,'\n');
end
fclose(fid);
%save energy
fid=fopen(enfile,'w+t');
for j=1:num
    j
    fprintf(fid,'%12.6f %12.6f\n',LATENT(j), EXPLAINED(j));
end
fclose(fid);
fprintf(1,'Finishing store the result...');

s = 1;