%Fisherface Code
clc;
clear all;
!dir *.mat
% eigenfeature.mat file contains the feature vectors
% you may skip the delay of loading many files by using this file directly
yload=input('If want to create a new eigenfeature.mat file, enter 1: ');
if isempty(yload), yload=0; end
if yload == 1,
dirname=[1:25]; % directory names are s1 to s25
nsubject=length(dirname);
filename=[1:10]; % 1.pgm to 10.pgm
nx=112; ny=92; % size of the individual image
veclen=nx*ny; % vector length after reshape reduced image to vector
imagearray = [];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Step 1. Preprocessing
% A. load a file
% B. low pass the image, and subsample the pixel at 3:1 ratio at each side
% Note that you should come up with your own preprocessing steps.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% construct a low pass Gaussian filter with 5 x 5 mask, zero mean and variance=1
h=fspecial('gaussian',5,1);
nxr=ceil(nx/3); nyr=ceil(ny/3); % size of reduced image
veclen=nxr*nyr; % vector length after reshape reduced image to vector
meanface=zeros(veclen,1); % averaged face image
xmi=zeros(veclen,nsubject); % averaged face image of i-th subject
data=[];
classavg = zeros(veclen, nsubject);
for i=1:length(dirname),
classa = zeros(veclen,1);
for j=1:length(filename),
fname=['s' int2str(dirname(i)) '/' int2str(filename(j)) '.pgm'];
tmp=imread(fname);
tmp1=filter2(h,tmp,'same'); % low pass anti-aliasing filter
% sub-sample image at 3:1 ratio at each dimension and convert to vector
vec=reshape(tmp1(1:3:112,1:3:92),veclen,1); %
imagearray = [imagearray,vec];
data = [data; vec' dirname(i)]; % raw data file, last column is subject id
xmi(:,i)=xmi(:,i)+vec; % accumulate each subject's averaged face
classa = classa + vec;
end
classarray = [];
for j=1:length(filename)
classarray = [classarray,classa];
end
%figure;plot(classarray);hold on;
classavg(:,(i-1)*j+1:i*j) = classarray / j;
end
clear tmp tmp1 vec;
xmi=xmi/length(filename); % averaged face of each subject
meanface=mean(xmi')'; % averaged face of all subjects
betweenarray = zeros(size(meanface,2), size(meanface,2));
for k = 1:j:size(classavg,2)%%25
betweenarray = betweenarray + (classavg(:,k) - meanface)*(classavg(:,k) - meanface)';
end
betweenarray = betweenarray/(length(dirname));
innervec = imagearray - classavg;
innerarray = (1/length(dirname))*(innervec*innervec');
figure;
for i=1:nsubject% nsubject=25,
eval(['subplot(5,5,' int2str(i) '),']);
imagesc(reshape(innerarray(:,i),38,31));
colormap('gray');
end
figure;
for i=1:nsubject% nsubject=25,
eval(['subplot(5,5,' int2str(i) '),']);
imagesc(reshape(betweenarray(:,i),38,31));
colormap('gray');
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Step 2. Feature Extraction
% here we uses Fisherfaces approach. You may also use
% other approaches
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
vec=data(:,1:1178)'; % 1178x 250
label=data(:,1179); % 250 x 1
nsample=length(label); % # of feature vectors = 250
vec2 = data(:,1:1178)';
vec2 = vec2 - classavg;
vec=vec-meanface*ones(1,nsample); % subtract mean value
[u,s] = eig(pinv(innerarray)*betweenarray, 'nobalance');
s = diag(s);
[unused, index] = sort(-real(s));
s = s(index);
u = u(:,index);
s=s/max(s); % normalized to 1
%[u,s,v]=svd(vec,0); clear v; % singular vavlue decomposition
% save some memory space making signular values a vector
save eigenfeature.mat
else
load eigenfeature.mat;
end
% now we need to choose # of basis (eigenfaces). We use a heuristic
% criterion so that the largest ns singular values will be chosen.
% the sum of singular values has a relation to the percentage energy
% of the signal that are retained.
% We let user to choose this number.
figure,clf, stairs(s),xlabel('n'),ylabel('% energy')
title('% energy of first n singular values');
ns=input('Choose eigenface dimension (< 250): ');
eface=u(:,1:ns); clear u % keep eigen faces
nfig=min(ns,25);
figure,clf,imagesc(reshape(meanface,nxr,nyr)),title(['average face']);
colormap('gray');
figure,clf
colormap('gray')
for i=1:nfig,
eval(['subplot(5,5,' int2str(i) '),']);
imagesc(reshape(eface(:,i),38,31));
end
feature=vec'*eface; % nsample x ns feature vectors
clear vec;
figure;
plot(feature);
featuregraph = zeros(2000,2000);
if size(feature,2) > 1
for k = 1:size(feature,1)
if feature(k,1) > -3000 & feature(k,2) > -3000 & feature(k,1) < 3000 & feature(k,2) < 3000
xpos = round(feature(k,1)/3)+1000;
ypos = round(feature(k,2)/3)+1000;
featuregraph(xpos-4:xpos+4,ypos-4:ypos+4) = ceil(k/10);
end
end
end
figure;
imagesc(featuregraph);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Step 3. Testing
% testing is done following testing protocols used in DARPA Ferret face
% recognition contest. See references:
% Phillips, P.J.; Hyeonjoon Moon; Rizvi, S.A.; Rauss, P.J.
% "The FERET evaluation methodology for face-recognition algorithms,"
% IEEE Trans. Pattern Analysis and Machine Intelligence, vol. 22, no. 10,
% pp. 1090- 1104, Oct. 2000.
% NOTE TESTIING PROTOCOL MUST BE FOLLOWED EXACTLY.
% You may choose your own distance measure. But the
% We use closed-set identification test. We separate the data into test
% set and query set. For each feature vector in the query set, a distance
% between it and each of the testing set feature vector will be measured.
% these distances will be ranked (smallest distance rank = 1). Then, we
% will compute percentage of correct identification at top n ranks.
% Since each subject has 10 poses, we rotate to use each post as test image
% the remaining nine as query images.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
idx=[1:nsample]; % original index of each feature vector
for j=1:10,
tidx{j}=idx(j:10:nsample); % test set indices 1 x 25
qidx{j}=setdiff(idx,tidx{j}); % query set indices 1 x 225
% calculate distance from query set to test set
dmat{j}=dist(feature(qidx{j},:),feature(tidx{j},:)); % 225 x 25
[tmp,tmpr]=sort(dmat{j}'); % tmpr is 25 x 225
[tmp1,tmp2]=sort(tmpr);
rank{j}=tmp2';% rank{j} is 225 x 25
end
% testing results
% compute cumulative face recognition test results
% rank{1}:rank{N}: M x nsubject matrix, M: # of query samples, N: # of partitions
% nsubject: # of test samples, each sample is a subject
% dirname: name of the subjects
% qidx{1}:qidx{N}: indices of query samples,
% tidx{1}:tidx{N}: indices of test samples
% label: labels of all samples including query and test samples
% each row of rank{j} gives ranks of distance from q to each of the 25 test images.
% 1 is closest and 25 is furtherest according to the distance measure.
% We need to check the rank of the correct label. If the label of q is sk, then
% we check the k-th entry of that row. For this, we first get the positions of these
% entries:
cnt=[];
for j=1:length(rank),
pos=[label(qidx{j})*ones(1,nsubject)==ones(length(qidx{j}),1)*label(tidx{j})'];
% consists of 1 at the column where label
% of q is the same as label of test image. O elsewhere. 225 x 25
rk=
- 1
- 2
前往页