用人脸图像测试PCA

利用上一个帖子的程序和一个很小的人脸来测试PCA的效果。

#include <opencv/cv.h>#include <opencv/highgui.h>#include <stdio.h>#include <stdlib.h>using namespace cv;using namespace std;#define WIDTH100#define HEIGHT100#define SAMPLE_NUM10#define DIMENTIONS(WIDTH*HEIGHT)#define NUM_EIGENFACES10#define PCA_MEAN"mean"#define PCA_EIGEN_VECTOR"eigen_vector"const char pSamplePath[SAMPLE_NUM][256] = {"YALE_s1.bmp","YALE_s17.bmp","YALE_s25.bmp","YALE_s38.bmp","YALE_s49.bmp","YALE_s65.bmp","YALE_s71.bmp","YALE_s79.bmp","YALE_s126.bmp","YALE_s159.bmp"};const char pTestPath[256] = "faces\\YALE_s1.bmp";int main(){//load imagesMat SampleSet(SAMPLE_NUM, DIMENTIONS, CV_32FC1);for (int s=0; s<SAMPLE_NUM; ++s){IplImage *img = cvLoadImage(pSamplePath[s]);IplImage *img_g = cvCreateImage(cvSize(img->width,img->height),IPL_DEPTH_8U, 1);cvCvtColor(img, img_g, CV_BGR2GRAY);//Change it to gray levelMat frame(img_g);for (int row=0; row<img->height; ++row){for (int col=0; col<img->width; ++col){float f = (float)(frame.at<uchar>(row, col));SampleSet.at<float>(s, (row*(img->width) + col) ) = f;}}cvReleaseImage(&img);cvReleaseImage(&img_g);}//TrainingPCA *pca = new PCA(SampleSet, Mat(), CV_PCA_DATA_AS_ROW);/////////////////Change format into eigenFaceMat EigenFace[NUM_EIGENFACES];for (int e=0; e<NUM_EIGENFACES; ++e){EigenFace[e].create(HEIGHT, WIDTH, CV_32FC1);float max_value=-1., min_value=1000.;for (int row=0; row<HEIGHT; ++row){for (int col=0; col<WIDTH; ++col){float f = (float)(pca->eigenvectors.at<float>(e, (row*WIDTH + col)));EigenFace[e].at<float>(row,col) = f;if (f > max_value)max_value = f;if(f < min_value)min_value = f;}}for (int row=0; row<HEIGHT; ++row){for (int col=0; col<WIDTH; ++col){float f = (float)(pca->eigenvectors.at<float>(e, (row*WIDTH + col)));EigenFace[e].at<float>(row,col) = (f-min_value)*256/(max_value-min_value);}}IplImage imgEigenFace0 = EigenFace[e];char pEigenFacePath[256];sprintf(pEigenFacePath,"EigenFace_%d.bmp",e);cvSaveImage(pEigenFacePath, &imgEigenFace0);}//calculate the decreased dimensions  int index;  float sum=0, sum0=0, ratio;  for (int d=0; d<pca->eigenvalues.rows; ++d)  {  sum += pca->eigenvalues.at<float>(d,0);  }  for (int d=0; d<pca->eigenvalues.rows; ++d)  {  sum0 += pca->eigenvalues.at<float>(d,0);  ratio = sum0/sum;  if(ratio > 0.9){  //0.9 is the thresholdindex = d;  break;  }  }  Mat eigenvetors_d;  eigenvetors_d.create((index+1), DIMENTIONS, CV_32FC1);//eigen values of decreased dimension  for (int i=0; i<(index+1); ++i)  {  pca->eigenvectors.row(i).copyTo(eigenvetors_d.row(i));  }  //cout << "eigenvectors" <<endl << eigenvetors_d << endl;  FileStorage fs_w("config.xml", FileStorage::WRITE);//write mean and eigenvalues into xml file  fs_w << PCA_MEAN << pca->mean;  fs_w << PCA_EIGEN_VECTOR << eigenvetors_d;  fs_w.release();  //EncodePCA *pca_encoding = new PCA();  //Read configFileStorage fs_r("config.xml", FileStorage::READ);  fs_r[PCA_MEAN] >> pca_encoding->mean;  fs_r[PCA_EIGEN_VECTOR] >> pca_encoding->eigenvectors;  fs_r.release();  //Mat encoded(SAMPLE_NUM, NUM_EIGENFACES, CV_32FC1);Mat encoded(SAMPLE_NUM, pca_encoding->eigenvectors.rows, CV_32FC1);for (int s=0; s<SAMPLE_NUM; ++s){Mat in = SampleSet.row(s);Mat out = encoded.row(s);//pca->project(in, out);pca_encoding->project(in, out);cout << "Encode" << s << encoded.row(s) << endl;}//DecodePCA *pca_decoding = new PCA();  //Read configFileStorage fs_d("config.xml", FileStorage::READ);  fs_d[PCA_MEAN] >> pca_decoding->mean;  fs_d[PCA_EIGEN_VECTOR] >> pca_decoding->eigenvectors;  fs_d.release();  for (int s=0; s<SAMPLE_NUM; ++s){Mat decode(1, DIMENTIONS,CV_32FC1);Mat in = encoded.row(s);//pca->backProject(in, decode);pca_decoding->backProject(in, decode);Mat mat_decode(HEIGHT, WIDTH, CV_32FC1);for (int row=0; row<HEIGHT; ++row){for (int col=0; col<WIDTH; ++col){float f = (float)(decode.at<float>(0, (row*WIDTH + col)));mat_decode.at<float>(row,col) = f;}}IplImage DecodedImgEigenFace = mat_decode;char pDecodedEigenFacePath[256];sprintf(pDecodedEigenFacePath,"DecodedEigenFace_%d.bmp",s);cvSaveImage(pDecodedEigenFacePath, &DecodedImgEigenFace);}//Test new inputIplImage *img_test = cvLoadImage(pTestPath);IplImage *img_test_g = cvCreateImage(cvSize(img_test->width,img_test->height),IPL_DEPTH_8U, 1);cvCvtColor(img_test, img_test_g, CV_BGR2GRAY);//Change it to gray levelMat mat_img(img_test_g);//inputMat mat_test(1, DIMENTIONS, CV_32FC1);for (int row=0; row<HEIGHT; ++row){for (int col=0; col<WIDTH; ++col){float f = (float)(mat_img.at<uchar>(row, col));mat_test.at<float>(0, (row*(WIDTH) + col) ) = f;}}Mat encoded_test(1, pca_encoding->eigenvectors.rows, CV_32FC1);pca_encoding->project(mat_test, encoded_test);cout << "test img encoded:" << encoded_test << endl;cvReleaseImage(&img_test);cvReleaseImage(&img_test_g);//Recogfloat min_sum = CV_MAX_ALLOC_SIZE;int min_index;for (int s=0; s<SAMPLE_NUM; ++s){float sum=0;for (int e=0; e<pca_encoding->eigenvectors.rows; ++e){float fs = encoded.at<float>(s,e);float fi = encoded_test.at<float>(0,e);sum += ((fs-fi)*(fs-fi));}if(sum < min_sum){min_sum = sum;min_index = s;}cout << s << "" << sum << endl;}cout << "The most similar one is :" << pSamplePath[min_index];delete pca;delete pca_encoding;delete pca_decoding;return 0;}

有个帖子【1】比较适合我这种初学者。回头有时间回来再总结下或者翻译下这篇帖子。

【1】 http://www.cognotics.com/opencv/servo_2007_series/part_4/page_3.html