您的位置: 首页  >  文章  >  性别识别（中）

性别识别（中）

分类: 文章 • 2025-03-24 11:32:10

基于上一篇博文提取的人脸图像，本文运用Gabor小波对人脸图像提取特征，再运用PCA对其降维处理，最后用SVM训练一个性别分类器。这里不再深入介绍Gabor小波和PCA的原理。很多博文讲解的很详细。

可参考：http://blog.****.net/jinshengtao/article/details/17797641

           http://www.cnblogs.com/emouse/p/3611256.html

           http://blog.****.net/yutianzuijin/article/details/10823985

 

基于Gabor+PCA+SVM的性别识别(1):  http://www.cnblogs.com/xiaoming123abc/p/5078411.html

基于Gabor+PCA+SVM的性别识别(3):  http://www.cnblogs.com/xiaoming123abc/p/5079116.html 

这里只简单说一下自己的理解：

一维Gabor小波的实质是一个带通滤波器，具有频率选择性。

对于二维Gabor小波，它不仅具有频率选择性，还有频率方向的选择性；这时，可以理解Gabor是个空域的模板，匹配与模板频率大小和频率方向相同的成分。也就是说，图像中的频率大小和频率方向与该Gabor的频率大小和频率方向相同，则经过Gabor滤波后，响应会比较大，这样就把特征提取出来了。

二维Gabor小波是由二维Gabor滤波器函数 G(ω,θ)(Gabor滤波器有很多参数，为了方便理解，这里只写了两个)通过尺度伸缩和旋转生成的一组滤波器，其参数的选择通常在频率空间进行考虑。为了对一幅图像的整个频域进行采样，可以采用具有多个中心频率和方向的Gabor滤波器组来描述图像。参数 ω,θ的不同选择分别体现了二维Gabor小波在频率和方向空间的采样方式。

整个频率空间可以是0到无穷大的任意值，由于一幅图像实际的频率分布是有限的范围，所以，对于图像的局部特征来说，参数 ω只能在一个很小的范围内选取。由于图像的纹理是随机分布的， θ的取值范围为0到2*pi ，考虑到Gabor滤波器的对称性， 的实际取值范围为0到pi 。本文采用5个中心频率和8个方向组成40个Gabor滤波器。

 

性别识别（中）

40个Gabor滤波器，由于每一个Gabor滤波器都由实部和虚部组成。所以一共有80个Gabor模板。

Gabor小波提取特征的过程就是拿着这些Gabor模板与图像卷积，得出40个实部响应图和40个实部响应图。实部与虚部合并，形成40个总响应图(18X21)，最后把每一个响应图(18X21)拉直，既18X21的矩阵变成1X378的向量（这只是本人用的方法，其他方法也应该可以）。这样每幅人脸(18X21)的特征个数(40X378)，运用PCA对(40X378)降维。

 

性别识别（中）

 

 

PCA

PCA（Principal Component Analysis）是一种常用的数据分析方法。PCA通过线性变换将原始数据变换为一组各维度线性无关的表示，可用于提取数据的主要特征分量，常用于高维数据的降维。

PCA的算法步骤：

设有m条n维数据。

1)   将原始数据按列组成n行m列矩阵X；

2)   将X的每一行（代表一个属性字段）进行零均值化，即减去这一行的均值；

3)   求出协方差矩阵(样本协方差矩阵) ；求出协方差矩阵的特征值及对应的特征向量；

4)   将特征向量按对应特征值大小从上到下按行排列成矩阵，取前k行组成矩阵P；

5)   Y=PX即为降维到k维后的数据。

本程序实现了Gabor特征提取，PCA降维，SVM分类。由于每一幅图像要与80个Gabor模板卷积。man样本与woman样本一共600多个。可能需要运行一段时间，最后得到一个SVM分类器。

main.cpp

    

  1. 

    

    #include <iostream>
    

    2. 

  2. 

    

    #include <fstream>
    

    3. 

  3. 

    

    #include <opencv2/opencv.hpp>
    

    4. 

  4. 

    

    #include <opencv2/ml/ml.hpp>
    

    5. 

  5. 

    

    #include "GaborFR.h"
    

    6. 

  6. 

    

     
    

    7. 

  7. 

    

    
using namespace std;
    

    8. 

  8. 

    

    
using namespace cv;
    

    9. 

  9. 

    

     
    

    10. 

  10. 

    

    #define manNO 409     //man样本个数
    

    11. 

  11. 

    

    #define womanNO 287    //woman样本个数
    

    12. 

  12. 

    

    
int iSize=20;//  Gabor的scale

    

    13. 

  13. 

    

    int main()
    

    14. 

  14. 

    

    {
    

    15. 

  15. 

    

        	int DescriptorDim=200;  
    

    16. 

  16. 

    

    		
    

    17. 

  17. 

    

    		string ImgName;//图片名(绝对路径)

    

    18. 

  18. 

    

    		ifstream finMan("man.txt");//man样本图片的文件名列表

    

    19. 

  19. 

    

    		ifstream finWoman("woman.txt");//woman样本图片的文件名列表

    

    20. 

  20. 

    

    		 Mat ROI;
    

    21. 

  21. 

    

    		Mat Gabor_feature;
    

    22. 

  22. 

    

    		Mat PCA_feature;
    

    23. 

  23. 

    

    		Mat sampleFeatureMat;//所有训练样本的特征向量组成的矩阵，行数等于所有样本的个数，列数等于Gabor描述子维数	

    

    24. 

  24. 

    

    		Mat sampleLabelMat; //训练样本的类别向量，行数等于所有样本的个数，列数等于1；1表示man，-1表示woman

    

    25. 

  25. 

    

    		Vector<Mat> GaborReal;
    

    26. 

  26. 

    

    		Vector<Mat> GaborImag;
    

    27. 

  27. 

    

    	//生成Gabor 模板

    

    28. 

  28. 

    

    		for(int i=0;i<8;i++)
    

    29. 

  29. 

    

    	        {
    

    30. 

  30. 

    

    		      for(int j=0;j<5;j++)
    

    31. 

  31. 

    

    		     {
    

    32. 

  32. 

    

    			   Mat M1= GaborFR::getRealGaborKernel(Size(iSize,iSize),
    

    33. 

  33. 

    

    				                                 2*CV_PI,
    

    34. 

  34. 

    

    												 i*CV_PI/8+CV_PI/2,
    

    35. 

  35. 

    

    												 j,
    

    36. 

  36. 

    

    												 1);
    

    37. 

  37. 

    

    		       Mat M2 = GaborFR::getImagGaborKernel(Size(iSize,iSize),
    

    38. 

  38. 

    

    				                                 2*CV_PI,
    

    39. 

  39. 

    

    												 i*CV_PI/8+CV_PI/2,
    

    40. 

  40. 

    

    												 j,
    

    41. 

  41. 

    

    												 1);   
    

    42. 

  42. 

    

    			   GaborReal.push_back(M1);
    

    43. 

  43. 

    

    			   GaborImag.push_back(M2);
    

    44. 

  44. 

    

     
    

    45. 

  45. 

    

    			}
    

    46. 

  46. 

    

    		}
    

    47. 

  47. 

    

    	//依次读取man样本图片

    

    48. 

  48. 

    

    		for(int num=0; num<manNO && getline(finMan,ImgName); num++)
    

    49. 

  49. 

    

    		{
    

    50. 

  50. 

    

    	        cout<<"处理："<<ImgName<<endl;
    

    51. 

  51. 

    

    		    ImgName = "D:\\Mycode\\mantrain_cut\\" + ImgName+".png";
    

    52. 

  52. 

    

    		    ROI= imread(ImgName,0);//读取图片

    

    53. 

  53. 

    

                if(ROI.data ==0)
    

    54. 

  54. 

    

    	        {
    

    55. 

  55. 

    

                  printf("[error] 没有图片\n");
    

    56. 

  56. 

    

                  return -5;
    

    57. 

  57. 

    

                }
    

    58. 

  58. 

    

           //Gabor变换提取特征*************************

    

    59. 

  59. 

    

    			for(int i=0;i<8;i++)
    

    60. 

  60. 

    

    	        {
    

    61. 

  61. 

    

    		      for(int j=0;j<5;j++)
    

    62. 

  62. 

    

    		     {
    

    63. 

  63. 

    

    			  Mat outR,outI,M_Magnitude;
    

    64. 

  64. 

    

    		
    

    65. 

  65. 

    

    			 GaborFR::getFilterRealImagPart(ROI,GaborReal[i*5+j],GaborImag[i*5+j],outR,outI);
    

    66. 

  66. 

    

    	         M_Magnitude=GaborFR::getMagnitude(outR,outI);
    

    67. 

  67. 

    

    		     //cartToPolar( outR, outI, M_Magnitude, Angle, false );   //计算幅值和相角

    

    68. 

  68. 

    

    			  normalize(M_Magnitude,M_Magnitude,0,1,CV_MINMAX,CV_32F);
    

    69. 

  69. 

    

    			 Mat line= M_Magnitude.reshape(0,1);
    

    70. 

  70. 

    

    			
    

    71. 

  71. 

    

    		    Gabor_feature.push_back(line);
    

    72. 

  72. 

    

    		}
    

    73. 

  73. 

    

    	  }
    

    74. 

  74. 

    

           //PCA降维************************************************************

    

    75. 

  75. 

    

                PCA pca(Gabor_feature, cv::Mat(), CV_PCA_DATA_AS_ROW,5);
    

    76. 

  76. 

    

    	        Mat dst=pca.project(Gabor_feature) ;
    

    77. 

  77. 

    

    	        Mat line =dst.reshape(0,1);
    

    78. 

  78. 

    

     
    

    79. 

  79. 

    

                //处理第一个样本时初始化特征向量矩阵和类别矩阵，因为只有知道了特征向量的维数才能初始化特征向量矩阵

    

    80. 

  80. 

    

    			if( 0 == num )
    

    81. 

  81. 

    

    			{
    

    82. 

  82. 

    

    			//初始化所有训练样本的特征向量组成的矩阵，行数等于所有样本的个数，列数等于Gabor降维后描述子维数sampleFeatureMat

    

    83. 

  83. 

    

    			sampleFeatureMat = Mat::zeros(manNO+womanNO, DescriptorDim, CV_32FC1);
    

    84. 

  84. 

    

    			//初始化训练样本的类别向量，行数等于所有样本的个数，列数等于1；1表示有人，0表示无人

    

    85. 

  85. 

    

    			sampleLabelMat = Mat::zeros(manNO+womanNO, 1, CV_32FC1);
    

    86. 

  86. 

    

    			}
    

    87. 

  87. 

    

    			//将计算好的Gabor降维后描述子复制到样本特征矩阵sampleFeatureMat

    

    88. 

  88. 

    

    		    for(int i=0; i<DescriptorDim; i++) 
    

    89. 

  89. 

    

    			sampleFeatureMat.at<float>(num,i) = line.at<float>(0,i);//第num个样本的特征向量中的第i个元素

    

    90. 

  90. 

    

    			sampleLabelMat.at<float>(num,0) = 1;//man样本类别为1

    

    91. 

  91. 

    

    		}
    

    92. 

  92. 

    

     
    

    93. 

  93. 

    

           //依次读取woman样本图片，生成Gabor描述子并降维

    

    94. 

  94. 

    

    		for(int num=0; num<womanNO && getline(finWoman,ImgName); num++)
    

    95. 

  95. 

    

    		{
    

    96. 

  96. 

    

    			cout<<"处理："<<ImgName<<endl;
    

    97. 

  97. 

    

    		    ImgName = "D:\\Mycode\\womantrain_cut\\" + ImgName+".png";
    

    98. 

  98. 

    

    		    ROI= imread(ImgName,0);//读取图片

    

    99. 

  99. 

    

                if(ROI.data ==0)
    

    100. 

  100. 

    

    	        {
    

    101. 

  101. 

    

                  printf("[error] 没有图片\n");
    

    102. 

  102. 

    

                  return -5;
    

    103. 

  103. 

    

                }
    

    104. 

  104. 

    

    			//Gabor变换提取特征*************************

    

    105. 

  105. 

    

    			for(int i=0;i<8;i++)
    

    106. 

  106. 

    

    	        {
    

    107. 

  107. 

    

    		      for(int j=0;j<5;j++)
    

    108. 

  108. 

    

    		     {
    

    109. 

  109. 

    

    			   
    

    110. 

  110. 

    

    			
    

    111. 

  111. 

    

    			Mat outR,outI,M_Magnitude;
    

    112. 

  112. 

    

     
    

    113. 

  113. 

    

    			GaborFR::getFilterRealImagPart(ROI,GaborReal[i*5+j],GaborImag[i*5+j],outR,outI);
    

    114. 

  114. 

    

    	        M_Magnitude=GaborFR::getMagnitude(outR,outI);
    

    115. 

  115. 

    

    			//cartToPolar( outR, outI, M_Magnitude, Angle, false );   //计算幅值和相角

    

    116. 

  116. 

    

    			normalize(M_Magnitude,M_Magnitude,0,1,CV_MINMAX,CV_32F);
    

    117. 

  117. 

    

    			Mat line= M_Magnitude.reshape(0,1);	
    

    118. 

  118. 

    

    		    Gabor_feature.push_back(line);
    

    119. 

  119. 

    

     
    

    120. 

  120. 

    

    		   }
    

    121. 

  121. 

    

    	      }
    

    122. 

  122. 

    

     
    

    123. 

  123. 

    

    //PCA降维************************************************************
    

    124. 

  124. 

    

                PCA pca( Gabor_feature, cv::Mat(), CV_PCA_DATA_AS_ROW,5);
    

    125. 

  125. 

    

    	        Mat dst=pca.project(Gabor_feature) ;
    

    126. 

  126. 

    

    	        Mat line =dst.reshape(0,1);
    

    127. 

  127. 

    

     
    

    128. 

  128. 

    

    			for(int i=0; i<DescriptorDim; i++)
    

    129. 

  129. 

    

    			sampleFeatureMat.at<float>(num+manNO,i) = line.at<float>(0,i);//第PosSamNO+num个样本的特征向量中的第i个元素

    

    130. 

  130. 

    

    		
    

    131. 

  131. 

    

    			sampleLabelMat.at<float>(num+manNO,0) = -1;//woman样本类别为-1

    

    132. 

  132. 

    

    		}
    

    133. 

  133. 

    

     
    

    134. 

  134. 

    

    		//训练SVM分类器

    

    135. 

  135. 

    

    		//迭代终止条件，当迭代满1000次或误差小于FLT_EPSILON时停止迭代

    

    136. 

  136. 

    

    		CvSVM svm;          //SVM分类器

    

    137. 

  137. 

    

    		CvTermCriteria criteria = cvTermCriteria(CV_TERMCRIT_ITER+CV_TERMCRIT_EPS, 1000, FLT_EPSILON);
    

    138. 

  138. 

    

    		//SVM参数：SVM类型为C_SVC；线性核函数；松弛因子C=0.01

    

    139. 

  139. 

    

    		CvSVMParams param(CvSVM::C_SVC, CvSVM::LINEAR, 0, 1, 0, 0.01, 0, 0, 0, criteria);
    

    140. 

  140. 

    

    		cout<<"开始训练SVM分类器"<<endl;
    

    141. 

  141. 

    

    		svm.train(sampleFeatureMat, sampleLabelMat, Mat(), Mat(), param);//训练分类器

    

    142. 

  142. 

    

    		cout<<"训练完成"<<endl;
    

    143. 

  143. 

    

    		svm.save("SVM_PCA.xml");//将训练好的SVM模型保存为xml文件

    

    144. 

  144. 

    

    	     waitKey();//注意：imshow之后必须加waitKey，否则无法显示图像

    

    145. 

  145. 

    

    	    system("pause");
    

    146. 

  146. 

    

    		
    

    147. 

  147. 

    

    }
    

    148.

　GarborFR.hpp　

    

  1. 

    

    #include "opencv2\opencv.hpp"
    

    2. 

  2. 

    

    #include <vector>
    

    3. 

  3. 

    

    
using namespace std;
    

    4. 

  4. 

    

    
using namespace cv;
    

    5. 

  5. 

    

    class GaborFR
    

    6. 

  6. 

    

    {
    

    7. 

  7. 

    

    
public:
    

    8. 

  8. 

    

    	GaborFR();
    

    9. 

  9. 

    

    	static Mat	getImagGaborKernel(Size ksize, double sigma, double theta, 

    

    10. 

  10. 

    

    									double nu,double gamma=1, int ktype= CV_32F);
    

    11. 

  11. 

    

    	static Mat	getRealGaborKernel( Size ksize, double sigma, double theta, 

    

    12. 

  12. 

    

    									double nu,double gamma=1, int ktype= CV_32F);
    

    13. 

  13. 

    

    	static Mat	getPhase(Mat &real,Mat &imag);
    

    14. 

  14. 

    

    	static Mat	getMagnitude(Mat &real,Mat &imag);
    

    15. 

  15. 

    

    	static void getFilterRealImagPart(Mat& src,Mat& real,Mat& imag,Mat &outReal,Mat &outImag);
    

    16. 

  16. 

    

    	static Mat	getFilterRealPart(Mat& src,Mat& real);
    

    17. 

  17. 

    

    	static Mat	getFilterImagPart(Mat& src,Mat& imag);
    

    18. 

  18. 

    

      /*

    

    19. 

  19. 

    

    	void		Init(Size ksize=Size(19,19), double sigma=2*CV_PI,
    

    20. 

  20. 

    

    					double gamma=1, int ktype=CV_32FC1);
    

    21. 

  21. 

    

      */
    

    22. 

  22. 

    

    
private:
    

    23. 

  23. 

    

    	//vector<Mat> gaborRealKernels;

    

    24. 

  24. 

    

    	//vector<Mat> gaborImagKernels;

    

    25. 

  25. 

    

    	bool isInited;
    

    26. 

  26. 

    

    };
    

    27.

　　Gabor.cpp

    

  1. 

    

    //#include "StdAfx.h"
    

    2. 

  2. 

    

    #include "GaborFR.h"
    

    3. 

  3. 

    

    GaborFR::GaborFR()
    

    4. 

  4. 

    

    {
    

    5. 

  5. 

    

    	isInited = false;
    

    6. 

  6. 

    

    }
    

    7. 

  7. 

    

    /*
    

    8. 

  8. 

    

    void GaborFR::Init(Size ksize, double sigma,double gamma, int ktype)
    

    9. 

  9. 

    

    {
    

    10. 

  10. 

    

    	gaborRealKernels.clear();
    

    11. 

  11. 

    

    	gaborImagKernels.clear();
    

    12. 

  12. 

    

    	double mu[8]={0,1,2,3,4,5,6,7};
    

    13. 

  13. 

    

    	double nu[5]={0,1,2,3,4};
    

    14. 

  14. 

    

    	int i,j;
    

    15. 

  15. 

    

    	for(i=0;i<5;i++)
    

    16. 

  16. 

    

    	{
    

    17. 

  17. 

    

    		for(j=0;j<8;j++)
    

    18. 

  18. 

    

    		{
    

    19. 

  19. 

    

    			gaborRealKernels.push_back(getRealGaborKernel(ksize,sigma,mu[j]*CV_PI/8,nu[i],gamma,ktype));
    

    20. 

  20. 

    

    			gaborImagKernels.push_back(getImagGaborKernel(ksize,sigma,mu[j]*CV_PI/8,nu[i],gamma,ktype));
    

    21. 

  21. 

    

    		}
    

    22. 

  22. 

    

    	}
    

    23. 

  23. 

    

    	isInited = true;
    

    24. 

  24. 

    

    }
    

    25. 

  25. 

    

    

    26. 

  26. 

    

    */
    

    27. 

  27. 

    

    Mat GaborFR::getImagGaborKernel(Size ksize, 
    

    28. 

  28. 

    

    	                           double sigma, 
    

    29. 

  29. 

    

    							   double theta, 
    

    30. 

  30. 

    

    							   double nu,
    

    31. 

  31. 

    

    							   double gamma, 
    

    32. 

  32. 

    

    							   int ktype)   
    

    33. 

  33. 

    

    {
    

    34. 

  34. 

    

    	double	sigma_x		= sigma;
    

    35. 

  35. 

    

    	double	sigma_y		= sigma/gamma;
    

    36. 

  36. 

    

    	int		nstds		= 3;
    

    37. 

  37. 

    

    	double	kmax		= CV_PI/2;
    

    38. 

  38. 

    

    	double	f			= cv::sqrt(2.0);
    

    39. 

  39. 

    

    	int xmin, xmax, ymin, ymax;
    

    40. 

  40. 

    

    	double c = cos(theta), s = sin(theta);
    

    41. 

  41. 

    

    	if( ksize.width > 0 )
    

    42. 

  42. 

    

    	{
    

    43. 

  43. 

    

    		xmax = ksize.width/2;
    

    44. 

  44. 

    

    	}
    

    45. 

  45. 

    

    	else//这个和matlab中的结果一样，默认都是19 !

    

    46. 

  46. 

    

    	{
    

    47. 

  47. 

    

    		xmax = cvRound(std::max(fabs(nstds*sigma_x*c), fabs(nstds*sigma_y*s)));
    

    48. 

  48. 

    

    	}
    

    49. 

  49. 

    

    	if( ksize.height > 0 )
    

    50. 

  50. 

    

    	{
    

    51. 

  51. 

    

    		ymax = ksize.height/2;
    

    52. 

  52. 

    

    	}
    

    53. 

  53. 

    

    	else

    

    54. 

  54. 

    

    	{
    

    55. 

  55. 

    

    		ymax = cvRound(std::max(fabs(nstds*sigma_x*s), fabs(nstds*sigma_y*c)));
    

    56. 

  56. 

    

    	}
    

    57. 

  57. 

    

    	xmin = -xmax;
    

    58. 

  58. 

    

    	ymin = -ymax;
    

    59. 

  59. 

    

    	CV_Assert( ktype == CV_32F || ktype == CV_64F );
    

    60. 

  60. 

    

    	float*	pFloat;
    

    61. 

  61. 

    

    	double*	pDouble;
    

    62. 

  62. 

    

    	Mat kernel(ymax - ymin + 1, xmax - xmin + 1, ktype); //初始化gabor的尺寸

    

    63. 

  63. 

    

    	double k		=	kmax/pow(f,nu);
    

    64. 

  64. 

    

    	double scaleReal=	k*k/sigma_x/sigma_y;
    

    65. 

  65. 

    

    	for( int y = ymin; y <= ymax; y++ )
    

    66. 

  66. 

    

    	{
    

    67. 

  67. 

    

    		if( ktype == CV_32F )
    

    68. 

  68. 

    

    		{
    

    69. 

  69. 

    

    			pFloat = kernel.ptr<float>(ymax-y);
    

    70. 

  70. 

    

    		}
    

    71. 

  71. 

    

    		else

    

    72. 

  72. 

    

    		{
    

    73. 

  73. 

    

    			pDouble = kernel.ptr<double>(ymax-y);
    

    74. 

  74. 

    

    		}
    

    75. 

  75. 

    

    		for( int x = xmin; x <= xmax; x++ )
    

    76. 

  76. 

    

    		{
    

    77. 

  77. 

    

    			double xr = x*c + y*s;
    

    78. 

  78. 

    

    			double v = scaleReal*exp(-(x*x+y*y)*scaleReal/2);
    

    79. 

  79. 

    

    			double temp=sin(k*xr);
    

    80. 

  80. 

    

    			v	=  temp*v;
    

    81. 

  81. 

    

    			if( ktype == CV_32F )
    

    82. 

  82. 

    

    			{
    

    83. 

  83. 

    

    				pFloat[xmax - x]= (float)v;
    

    84. 

  84. 

    

    			}
    

    85. 

  85. 

    

    			else

    

    86. 

  86. 

    

    			{
    

    87. 

  87. 

    

    				pDouble[xmax - x] = v;
    

    88. 

  88. 

    

    			}
    

    89. 

  89. 

    

    		}
    

    90. 

  90. 

    

    	}
    

    91. 

  91. 

    

    	return kernel;
    

    92. 

  92. 

    

    }
    

    93. 

  93. 

    

    //sigma一般为2*pi
    

    94. 

  94. 

    

    Mat GaborFR::getRealGaborKernel( Size ksize,
    

    95. 

  95. 

    

    	                            double sigma,
    

    96. 

  96. 

    

    								double theta, 
    

    97. 

  97. 

    

    	                            double nu,
    

    98. 

  98. 

    

    								double gamma, 
    

    99. 

  99. 

    

    								int ktype)
    

    100. 

  100. 

    

    {
    

    101. 

  101. 

    

    	
    

    102. 

  102. 

    

    	double	sigma_x		= sigma;
    

    103. 

  103. 

    

    	double	sigma_y		= sigma/gamma;
    

    104. 

  104. 

    

    	int		nstds		= 3;
    

    105. 

  105. 

    

    	double	kmax		= CV_PI/2;
    

    106. 

  106. 

    

    	double	f			= cv::sqrt(2.0);
    

    107. 

  107. 

    

    	int xmin, xmax, ymin, ymax;
    

    108. 

  108. 

    

    	double c = cos(theta), s = sin(theta);
    

    109. 

  109. 

    

    	if( ksize.width > 0 )
    

    110. 

  110. 

    

    	{
    

    111. 

  111. 

    

    		xmax = ksize.width/2;
    

    112. 

  112. 

    

    	}
    

    113. 

  113. 

    

    	else//这个和matlab中的结果一样，默认都是19 !

    

    114. 

  114. 

    

    	{
    

    115. 

  115. 

    

    		xmax = cvRound(std::max(fabs(nstds*sigma_x*c), fabs(nstds*sigma_y*s)));
    

    116. 

  116. 

    

    	}
    

    117. 

  117. 

    

     
    

    118. 

  118. 

    

    	if( ksize.height > 0 )
    

    119. 

  119. 

    

    		ymax = ksize.height/2;
    

    120. 

  120. 

    

    	else

    

    121. 

  121. 

    

    		ymax = cvRound(std::max(fabs(nstds*sigma_x*s), fabs(nstds*sigma_y*c)));
    

    122. 

  122. 

    

    	xmin = -xmax;
    

    123. 

  123. 

    

    	ymin = -ymax;
    

    124. 

  124. 

    

    	CV_Assert( ktype == CV_32F || ktype == CV_64F );
    

    125. 

  125. 

    

    	float*	pFloat;
    

    126. 

  126. 

    

    	double*	pDouble;
    

    127. 

  127. 

    

    	Mat kernel(ymax - ymin + 1, xmax - xmin + 1, ktype);
    

    128. 

  128. 

    

    	double k		=	kmax/pow(f,nu);
    

    129. 

  129. 

    

    	double exy		=	sigma_x*sigma_y/2;
    

    130. 

  130. 

    

    	double scaleReal=	k*k/sigma_x/sigma_y;
    

    131. 

  131. 

    

    	int	   x,y;
    

    132. 

  132. 

    

    	for( y = ymin; y <= ymax; y++ )
    

    133. 

  133. 

    

    	{
    

    134. 

  134. 

    

    		if( ktype == CV_32F )
    

    135. 

  135. 

    

    		{
    

    136. 

  136. 

    

    			pFloat = kernel.ptr<float>(ymax-y);
    

    137. 

  137. 

    

    		}
    

    138. 

  138. 

    

    		else

    

    139. 

  139. 

    

    		{
    

    140. 

  140. 

    

    			pDouble = kernel.ptr<double>(ymax-y);
    

    141. 

  141. 

    

    		}
    

    142. 

  142. 

    

    		for( x = xmin; x <= xmax; x++ )
    

    143. 

  143. 

    

    		{
    

    144. 

  144. 

    

    			double xr = x*c + y*s;
    

    145. 

  145. 

    

    			double v = scaleReal*exp(-(x*x+y*y)*scaleReal/2);
    

    146. 

  146. 

    

    			double temp=cos(k*xr) - exp(-exy);
    

    147. 

  147. 

    

    			v	=	temp*v;
    

    148. 

  148. 

    

    			if( ktype == CV_32F )
    

    149. 

  149. 

    

    			{
    

    150. 

  150. 

    

    				pFloat[xmax - x]= (float)v;
    

    151. 

  151. 

    

    			}
    

    152. 

  152. 

    

    			else

    

    153. 

  153. 

    

    			{
    

    154. 

  154. 

    

    				pDouble[xmax - x] = v;
    

    155. 

  155. 

    

    			}
    

    156. 

  156. 

    

    		}
    

    157. 

  157. 

    

    	}
    

    158. 

  158. 

    

    	return kernel;
    

    159. 

  159. 

    

    }
    

    160. 

  160. 

    

     
    

    161. 

  161. 

    

     
    

    162. 

  162. 

    

    Mat GaborFR::getMagnitude(Mat &real,Mat &imag)
    

    163. 

  163. 

    

    {
    

    164. 

  164. 

    

    	CV_Assert(real.type()==imag.type());
    

    165. 

  165. 

    

    	CV_Assert(real.size()==imag.size());
    

    166. 

  166. 

    

    	int ktype=real.type();
    

    167. 

  167. 

    

    	int row = real.rows,col = real.cols;
    

    168. 

  168. 

    

    	int i,j;
    

    169. 

  169. 

    

    	float*	pFloat,*pFloatR,*pFloatI;
    

    170. 

  170. 

    

    	double*	pDouble,*pDoubleR,*pDoubleI;
    

    171. 

  171. 

    

    	Mat		kernel(row, col, real.type());
    

    172. 

  172. 

    

    	for(i=0;i<row;i++)
    

    173. 

  173. 

    

    	{
    

    174. 

  174. 

    

    		if( ktype == CV_32FC1 )
    

    175. 

  175. 

    

    		{
    

    176. 

  176. 

    

    			pFloat = kernel.ptr<float>(i);
    

    177. 

  177. 

    

    			pFloatR= real.ptr<float>(i);
    

    178. 

  178. 

    

    			pFloatI= imag.ptr<float>(i);
    

    179. 

  179. 

    

    		}
    

    180. 

  180. 

    

    		else

    

    181. 

  181. 

    

    		{
    

    182. 

  182. 

    

    			pDouble = kernel.ptr<double>(i);
    

    183. 

  183. 

    

    			pDoubleR= real.ptr<double>(i);
    

    184. 

  184. 

    

    			pDoubleI= imag.ptr<double>(i);
    

    185. 

  185. 

    

    		}
    

    186. 

  186. 

    

    		for(j=0;j<col;j++)
    

    187. 

  187. 

    

    		{
    

    188. 

  188. 

    

    			if( ktype == CV_32FC1 )
    

    189. 

  189. 

    

    			{
    

    190. 

  190. 

    

    				pFloat[j]= sqrt(pFloatI[j]*pFloatI[j]+pFloatR[j]*pFloatR[j]);
    

    191. 

  191. 

    

    			}
    

    192. 

  192. 

    

    			else

    

    193. 

  193. 

    

    			{
    

    194. 

  194. 

    

    				pDouble[j] = sqrt(pDoubleI[j]*pDoubleI[j]+pDoubleR[j]*pDoubleR[j]);
    

    195. 

  195. 

    

    			}
    

    196. 

  196. 

    

    		}
    

    197. 

  197. 

    

    	}
    

    198. 

  198. 

    

    	return kernel;
    

    199. 

  199. 

    

    }
    

    200. 

  200. 

    

    Mat GaborFR::getPhase(Mat &real,Mat &imag)
    

    201. 

  201. 

    

    {
    

    202. 

  202. 

    

    	CV_Assert(real.type()==imag.type());
    

    203. 

  203. 

    

    	CV_Assert(real.size()==imag.size());
    

    204. 

  204. 

    

    	int ktype=real.type();
    

    205. 

  205. 

    

    	int row = real.rows,col = real.cols;
    

    206. 

  206. 

    

    	int i,j;
    

    207. 

  207. 

    

    	float*	pFloat,*pFloatR,*pFloatI;
    

    208. 

  208. 

    

    	double*	pDouble,*pDoubleR,*pDoubleI;
    

    209. 

  209. 

    

    	Mat		kernel(row, col, real.type());
    

    210. 

  210. 

    

    	for(i=0;i<row;i++)
    

    211. 

  211. 

    

    	{
    

    212. 

  212. 

    

    		if( ktype == CV_32FC1 )
    

    213. 

  213. 

    

    		{
    

    214. 

  214. 

    

    			pFloat = kernel.ptr<float>(i);
    

    215. 

  215. 

    

    			pFloatR= real.ptr<float>(i);
    

    216. 

  216. 

    

    			pFloatI= imag.ptr<float>(i);
    

    217. 

  217. 

    

    		}
    

    218. 

  218. 

    

    		else

    

    219. 

  219. 

    

    		{
    

    220. 

  220. 

    

    			pDouble = kernel.ptr<double>(i);
    

    221. 

  221. 

    

    			pDoubleR= real.ptr<double>(i);
    

    222. 

  222. 

    

    			pDoubleI= imag.ptr<double>(i);
    

    223. 

  223. 

    

    		}
    

    224. 

  224. 

    

    		for(j=0;j<col;j++)
    

    225. 

  225. 

    

    		{
    

    226. 

  226. 

    

    			if( ktype == CV_32FC1 )
    

    227. 

  227. 

    

    			{
    

    228. 

  228. 

    

    // 				if(pFloatI[j]/(pFloatR[j]+pFloatI[j]) > 0.99)
    

    229. 

  229. 

    

    // 				{
    

    230. 

  230. 

    

    // 					pFloat[j]=CV_PI/2;
    

    231. 

  231. 

    

    // 				}
    

    232. 

  232. 

    

    // 				else
    

    233. 

  233. 

    

    // 				{
    

    234. 

  234. 

    

    //					pFloat[j] = atan(pFloatI[j]/pFloatR[j]);
    

    235. 

  235. 

    

    				pFloat[j] = asin(pFloatI[j]/sqrt(pFloatR[j]*pFloatR[j]+pFloatI[j]*pFloatI[j]));
    

    236. 

  236. 

    

    /*				}*/
    

    237. 

  237. 

    

    //				pFloat[j] = atan2(pFloatI[j],pFloatR[j]);
    

    238. 

  238. 

    

    			}//CV_32F

    

    239. 

  239. 

    

    			else

    

    240. 

  240. 

    

    			{
    

    241. 

  241. 

    

    				if(pDoubleI[j]/(pDoubleR[j]+pDoubleI[j]) > 0.99)
    

    242. 

  242. 

    

    				{
    

    243. 

  243. 

    

    					pDouble[j]=CV_PI/2;
    

    244. 

  244. 

    

    				}
    

    245. 

  245. 

    

    				else

    

    246. 

  246. 

    

    				{
    

    247. 

  247. 

    

    					pDouble[j] = atan(pDoubleI[j]/pDoubleR[j]);
    

    248. 

  248. 

    

    				}
    

    249. 

  249. 

    

    //				pDouble[j]=atan2(pDoubleI[j],pDoubleR[j]);
    

    250. 

  250. 

    

    			}//CV_64F

    

    251. 

  251. 

    

    		}
    

    252. 

  252. 

    

    	}
    

    253. 

  253. 

    

    	return kernel;
    

    254. 

  254. 

    

    }
    

    255. 

  255. 

    

    Mat GaborFR::getFilterRealPart(Mat& src,Mat& real)
    

    256. 

  256. 

    

    {
    

    257. 

  257. 

    

    	//CV_Assert(real.type()==src.type());

    

    258. 

  258. 

    

    	Mat dst;
    

    259. 

  259. 

    

    	Mat kernel;
    

    260. 

  260. 

    

    	flip(real,kernel,-1);//中心镜面

    

    261. 

  261. 

    

    //	filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_CONSTANT);
    

    262. 

  262. 

    

    	filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_REPLICATE);
    

    263. 

  263. 

    

    	return dst;
    

    264. 

  264. 

    

    }
    

    265. 

  265. 

    

    Mat GaborFR::getFilterImagPart(Mat& src,Mat& imag)
    

    266. 

  266. 

    

    {
    

    267. 

  267. 

    

    	//CV_Assert(imag.type()==src.type());

    

    268. 

  268. 

    

    	Mat dst;
    

    269. 

  269. 

    

    	Mat kernel;
    

    270. 

  270. 

    

    	flip(imag,kernel,-1);//中心镜面

    

    271. 

  271. 

    

    //	filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_CONSTANT);
    

    272. 

  272. 

    

    	filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_REPLICATE);
    

    273. 

  273. 

    

    	return dst;
    

    274. 

  274. 

    

    }
    

    275. 

  275. 

    

     
    

    276. 

  276. 

    

     
    

    277. 

  277. 

    

     
    

    278. 

  278. 

    

    
void GaborFR::getFilterRealImagPart(Mat& src,Mat& real,Mat& imag,Mat &outReal,Mat &outImag)
    

    279. 

  279. 

    

    {
    

    280. 

  280. 

    

    	outReal=getFilterRealPart(src,real);
    

    281. 

  281. 

    

    	outImag=getFilterImagPart(src,imag);
    

    282. 

  282. 

    

    }
    

    283.

　　本人能力有限，错误在所难免。敬请赐教

转自：https://blog.****.net/u012507022/article/details/50979003

相关推荐