基于OpenCV的Gabor变换及特征提取

一、Gabor变换概述

 　　Gabor变换是一种加窗短时Fourier变换（Window Fourier transform or Short Time Fourier Transform）。Fourier变换是整体上将信号分解为不同的频率分量（任何信号都可分解为复正弦信号之和），对确定性信号及平稳信号使用。其缺点为缺乏时间的局部性信息，并且对时变信号、非平稳信号的分析存在严重不足，(1)无法告知某些频率成分发生在哪些时间内；(2)无法告知某个时刻信号频谱的分布情况。

　　Gabor函数可以在频域不同尺度、不同方向上提取相关的特征。另外Gabor函数与人眼的生物作用相仿，所以经常用作纹理识别上，并取得了较好的效果，Gabor变换是短时Fourier变换中当窗函数取为高斯函数时的一种特殊情况。Gabor变换的本质实际上还是对二维图像求卷积。因此二维卷积运算的效率就直接决定了Gabor变换的效率。

　　下面简要说一下二维卷积运算的计算过程：

　　A可以理解成是待处理的笔迹纹理，B可以理解成Gabor变换的核函数，现在要求A与B的离散二维叠加卷积，我们首先对A的右边界和下边界填充0（zero padding)，然后将B进行水平翻转和垂直翻转，如下图：

　　

然后用B中的每个值依次乘以A中相对位置处的值并进行累加，结果填入相应位置处（注意红圈位置）。通常二维卷积的结果比A、B的尺寸要大，如下图：

 

二、Gabor函数表达式

　　通过频率参数和高斯函数参数的选取，Gabor变换可以选取很多纹理特征，但是Gabor是非正交的，不同特征分量之间有冗余。 Gabor是有高斯核函数与复正弦函数调制而成，如下图所示：

 

   图(a)为偏移x轴30°的正弦函数，图(b)为高斯核，图(c)为对应的Gabor filter。可以看出正弦函数是如何在空间上具有局部性的。

        二维Gabor函数的第一种表示形式：

 

       为正弦函数的波长，核函数方向，为相位偏移，为高斯标准差，为x、y两个方向的纵横比（指定了Gabor函数的椭圆率）。

       二维Gabor函数的第二种形式：

 

   v的取值决定了Gabor滤波的波长，u的取值表示Gabor核函数的方向，K表示总的方向数。参数决定了高斯窗口的大小，这里取

三、Gabor特征提取

　　提取Gabor特征之前，我们可以先对需要处理的图像I(x,y)进行实数形式的Gabor变换，得到处理后的图像。举个例子，我们需要处理的图像大小为128x128，可以先通过变换得到处理后的图像，大小也是128x128.但是我们不方便直接提取特征，因为这样提取出来的特征维数太高，不利于我们进行后续处理。这里我们对图像分块，分别水平和垂直方向取16等分，这样就可以将整个图像分成64个16x16大小的子图像块。如下图所示:

 

接着计算每一块对应的能量，第k块的能量定义如下：

 

这样计算之后就可以形成如下图所示的联合空间频率能量矩阵Energy。

 

然后，我们将能量矩阵降维成1x64的行向量，作为我们的原始图像在某一方向和尺度变换后的特征向量，如下：

 

接着就可以利用得到的特征向量进行下一步的处理，比如相似度判别或者聚类等等。当然一般情况下，需要提取多个方向和尺度特征，变换参数重复上述过程就可以了。

四、基于OpenCV 2.x的实现

　　下面的代码是根据上面Gabor变换的第二种形式来的，原始代码来源于这里，由于原始代码使用c语言实现，里面很多函数是基于OpenCV1.x的，比如使用了大量的IplImage类型，这种是比较老的图片处理类型，而且要手动进行内存管理，加之不方便进行矩阵运算，所以我在代码的基础上进行了修改，使用了Mat类型代替了原始的函数参数，并对原函数中的一些不必要的函数进行了精简。整个函数基于C++实现，已经调试通过。

　　头文件(cvgabor.h)：

 1 #ifndef CVGABOR_H 2 #define CVGABOR_H 3  4 #include <iostream> 5  6  7 #include <opencv2/core/core.hpp> 8 #include <opencv2/highgui/highgui.hpp> 9 #include <opencv2/imgproc/imgproc.hpp>10 11 const double PI = 3.14159265;12 const int CV_GABOR_REAL = 1;13 const int CV_GABOR_IMAG = 2;14 const int CV_GABOR_MAG  = 3;15 const int CV_GABOR_PHASE = 4;16 17 using namespace cv;18 using namespace std;19 20 class CvGabor{21 public:22     CvGabor();23     ~CvGabor();24 25     CvGabor(int iMu, int iNu);26     CvGabor(int iMu, int iNu, double dSigma);27     CvGabor(int iMu, int iNu, double dSigma, double dF);28     CvGabor(double dPhi, int iNu);29     CvGabor(double dPhi, int iNu, double dSigma);30     CvGabor(double dPhi, int iNu, double dSigma, double dF);31 32     void Init(int iMu, int iNu, double dSigma, double dF);33     void Init(double dPhi, int iNu, double dSigma, double dF);34 35     bool IsInit();36     bool IsKernelCreate();37     int mask_width();38     int get_mask_width();39 40     void get_image(int Type, Mat& image);41     void get_matrix(int Type, Mat& matrix);42 43     void conv_img(Mat& src, Mat& dst, int Type);44 45 protected:46     double Sigma;47     double F;48     double Kmax;49     double K;50     double Phi;51     bool bInitialised;52     bool bKernel;53     int Width;54     Mat Imag;55     Mat Real;56   57 private:58     void creat_kernel();59     60 };61 62 #endif

 

　　函数实现(cvgabor.cpp)：

  1 #include "cvgabor.h"  2   3 CvGabor::CvGabor()  4 {  5 }  6   7   8 CvGabor::~CvGabor()  9 { 10 } 11  12  13 /*! 14  15 Parameters: 16 iMu        The orientation iMu*PI/8, 17 iNu         The scale, 18 dSigma         The sigma value of Gabor, 19 dPhi        The orientation in arc 20 dF        The spatial frequency 21  22 */ 23  24 CvGabor::CvGabor(int iMu, int iNu) 25 { 26     double dSigma = 2*PI; 27     F = sqrt(2.0); 28     Init(iMu, iNu, dSigma, F); 29 } 30  31 CvGabor::CvGabor(int iMu, int iNu, double dSigma) 32 {  33     F = sqrt(2.0); 34     Init(iMu, iNu, dSigma, F); 35 } 36  37 CvGabor::CvGabor(int iMu, int iNu, double dSigma, double dF) 38 { 39     Init(iMu, iNu, dSigma, dF); 40 } 41  42 CvGabor::CvGabor(double dPhi, int iNu) 43 { 44     Sigma = 2*PI; 45     F = sqrt(2.0); 46     Init(dPhi, iNu, Sigma, F); 47 } 48  49 CvGabor::CvGabor(double dPhi, int iNu, double dSigma) 50 { 51     F = sqrt(2.0); 52     Init(dPhi, iNu, dSigma, F); 53 } 54  55 CvGabor::CvGabor(double dPhi, int iNu, double dSigma, double dF) 56 { 57     Init(dPhi, iNu, dSigma,dF); 58 } 59  60 /*! 61 Parameters: 62 iMu     The orientations which is iMu*PI.8 63 iNu     The scale can be from -5 to infinit 64 dSigma     The Sigma value of gabor, Normally set to 2*PI 65 dF     The spatial frequence , normally is sqrt(2) 66  67 Initilize the.gabor with the orientation iMu, the scale iNu, the sigma dSigma, the frequency dF, it will call the function creat_kernel(); So a gabor is created. 68 */ 69 void CvGabor::Init(int iMu, int iNu, double dSigma, double dF) 70 { 71     //Initilise the parameters 72     bInitialised = false; 73     bKernel = false; 74  75     Sigma = dSigma; 76     F = dF; 77  78     Kmax = PI/2; 79  80     //Absolute value of K 81     K = Kmax / pow(F, (double)iNu); 82     Phi = PI*iMu/8; 83     bInitialised = true; 84  85     Width = mask_width(); 86     creat_kernel(); 87 } 88  89 /*! 90 Parameters: 91 dPhi     The orientations 92 iNu     The scale can be from -5 to infinit 93 dSigma     The Sigma value of gabor, Normally set to 2*PI 94 dF     The spatial frequence , normally is sqrt(2) 95  96 Initilize the.gabor with the orientation dPhi, the scale iNu, the sigma dSigma, the frequency dF, it will call the function creat_kernel(); So a gabor is created.filename     The name of the image file 97 file_format     The format of the file, 98 */ 99 void CvGabor::Init(double dPhi, int iNu, double dSigma, double dF)100 {101 102     bInitialised = false;103     bKernel = false;104     Sigma = dSigma;105     F = dF;106 107     Kmax = PI/2;108 109     // Absolute value of K110     K = Kmax / pow(F, (double)iNu);111     Phi = dPhi;112     bInitialised = true;113 114     Width = mask_width();115     creat_kernel();116 }117 118 /*!119 Returns:120 a boolean value, TRUE is created or FALSE is non-created.121 122 Determine whether a gabor kernel is created.123 */124 125 bool CvGabor::IsInit()126 {127     return bInitialised;128 }129 130 bool CvGabor::IsKernelCreate()131 {132     return bKernel;133 }134 135 /*!136 Returns:137 The long type show the width.138 139 Return the width of mask (should be NxN) by the value of Sigma and iNu.140 */141 int CvGabor::mask_width()142 {143     int lWidth;144     if (IsInit() == false)145     {146         cerr << "Error: The Object has not been initilised in mask_width()!\n" << endl;147         return 0;148     }149     else150     {151         //determine the width of Mask152         double dModSigma = Sigma/K;153         int dWidth = cvRound(dModSigma*6 + 1);154 155         //test whether dWidth is an odd.156         if((dWidth % 2) == 0)157         {158             lWidth = dWidth + 1;159         }160         else161         {162             lWidth = dWidth;163         }164         return lWidth;165     }166 }167 168 /*!169 170 Returns:171 Pointer to long type width of mask.172 173 */174 int CvGabor::get_mask_width()175 {176     return Width;177 }178 179 /*!180 \fn CvGabor::creat_kernel()181 Create gabor kernel182 183 Create 2 gabor kernels - REAL and IMAG, with an orientation and a scale 184 */185 void CvGabor::creat_kernel()186 {187 188     if (IsInit() == false)189     {190         cerr << "Error: The Object has not been initilised in creat_kernel()!" << endl;191     }192     else193     {194         Mat mReal(Width, Width, CV_32FC1);195         Mat mImag(Width, Width, CV_32FC1);196 197         /**************************** Gabor Function ****************************/ 198         int x, y;199         double dReal;200         double dImag;201         double dTemp1, dTemp2, dTemp3;202 203         for (int i = 0; i < Width; i++)204         {205             for (int j = 0; j < Width; j++)206             {207                 x = i-(Width-1)/2;208                 y = j-(Width-1)/2;209                 dTemp1 = (pow(K,2)/pow(Sigma,2))*exp(-(pow((double)x,2)+pow((double)y,2))*pow(K,2)/(2*pow(Sigma,2)));210                 dTemp2 = cos(K*cos(Phi)*x + K*sin(Phi)*y) - exp(-(pow(Sigma,2)/2));211                 dTemp3 = sin(K*cos(Phi)*x + K*sin(Phi)*y);212                 dReal = dTemp1*dTemp2;213                 dImag = dTemp1*dTemp3;214 215                 mReal.row(i).col(j) = dReal;216                 mImag.row(i).col(j) = dImag;217             }218         }219         /**************************** Gabor Function ****************************/220         bKernel = true;221 222         mReal.copyTo(Real);223         mImag.copyTo(Imag);224     }225 }226 227 228 /*!229 \fn CvGabor::get_image(int Type)230 Get the speific type of image of Gabor231 232 Parameters:233 Type        The Type of gabor kernel, e.g. REAL, IMAG, MAG, PHASE   234 235 Returns:236 Pointer to image structure, or NULL on failure    237 238 Return an Image (gandalf image class) with a specific Type   "REAL"    "IMAG" "MAG" "PHASE"  239 */240 void CvGabor::get_image(int Type, Mat& image)241 {242     if(IsKernelCreate() == false)243     { 244         cerr << "Error: the Gabor kernel has not been created in get_image()!" << endl;245         return;246     }247     else248     {  249         Mat re(Width, Width, CV_32FC1);250         Mat im(Width, Width, CV_32FC1);251         Mat temp;252 253         switch(Type)254         {255         case 1:  //Real256             temp = Real.clone();257             normalize(temp, temp, 255.0, 0.0, NORM_MINMAX);258             break;259         case 2:  //Imag260             temp = Imag.clone();261             break; 262         case 3:  //Magnitude263             re = Real.clone();264             im = Imag.clone();265 266             pow(re, 2, re);267             pow(im, 2, im);268             add(im, re, temp);269             pow(temp, 0.5, temp);270             break;271         case 4:  //Phase272             ///@todo273             break;274         }275 276         convertScaleAbs(temp, image, 1, 0);277     }278 }279 280 /*!281 \fn CvGabor::get_matrix(int Type)282 Get a matrix by the type of kernel283 284 Parameters:285 Type        The type of kernel, e.g. REAL, IMAG, MAG, PHASE286 287 Returns:288 Pointer to matrix structure, or NULL on failure.289 290 Return the gabor kernel.291 */292 void CvGabor::get_matrix(int Type, Mat& matrix)293 {294     if (!IsKernelCreate())295     {296         cerr << "Error: the gabor kernel has not been created!" << endl;297         return;298     }299     switch (Type)300     {301     case CV_GABOR_REAL:302         matrix = Real.clone();303         break;304     case CV_GABOR_IMAG:305         matrix = Imag.clone();306         break;307     case CV_GABOR_MAG:308         break;309     case CV_GABOR_PHASE:310         break;311     }312 }313 314 /*!315 \fn CvGabor::conv_img_a(IplImage *src, IplImage *dst, int Type)316 */317 void CvGabor::conv_img(Mat &src, Mat &dst, int Type)318 {319     Mat mat = src.clone();320 321     Mat rmat(src.rows, src.cols, CV_32FC1);322     Mat imat(src.rows, src.cols, CV_32FC1);323 324     switch (Type)325     {326     case CV_GABOR_REAL:327         filter2D(mat, mat, 1, Real, Point((Width-1)/2, (Width-1)/2));328         break;329 330     case CV_GABOR_IMAG:331         filter2D(mat, mat, 1, Imag, Point( (Width-1)/2, (Width-1)/2));332         break;333 334     case CV_GABOR_MAG:335         /* Real Response */336         filter2D(mat, rmat, 1, Real, Point((Width-1)/2, (Width-1)/2));337 338         /* Imag Response */339         filter2D(mat, imat, 1, Imag, Point( (Width-1)/2, (Width-1)/2));340 341         /* Magnitude response is the square root of the sum of the square of real response and imaginary response */342         pow(rmat, 2, rmat);343         pow(imat, 2, imat);344         add(rmat, imat, mat);345         pow(mat, 0.5, mat);346         break;347 348     case CV_GABOR_PHASE:349         break;350     }351 352 //    cvNormalize(mat, mat, 0, 255, CV_MINMAX, NULL);353     mat.copyTo(dst);354 }

 

 函数的使用方法如下：

首先显示核函数图像：

 1     //创建一个方向是PI/4而尺度是3的gabor 2     double Sigma = 2*PI; 3     double F = sqrt(2.0); 4     CvGabor gabor(PI/4, 3, Sigma, F); 5  6     //获得实部并显示它 7     Mat kernel(gabor.get_mask_width(), gabor.get_mask_width(), CV_8UC1); 8     gabor.get_image(CV_GABOR_REAL, kernel); 9     imshow("Kernel", kernel);10     cout << kernel.rows << endl;11     cout << kernel.cols << endl;12     cvWaitKey(0);

 显示效果如下：

接着，对输入的图像进行处理：

 1     //载入一个图像并显示 2     Mat img = imread( "test.jpg", CV_LOAD_IMAGE_GRAYSCALE ); 3     imshow("Original Image", img); 4     cvWaitKey(0); 5  6     //获取载入图像的gabor滤波响应的实部并且显示 7     Mat reimg(img.rows,img.cols, CV_32FC1); 8     gabor.conv_img(img, reimg, CV_GABOR_REAL); 9     imshow("After Image", reimg);10     cvWaitKey(0);

 

下面是显示效果:

接着就可以对得到的图像，按照上面体征提取部分的步骤进行特征提取。

 

参考链接：http://blog.csdn.net/renjinr/article/details/13768655

　　　　　http://blog.sina.com.cn/s/blog_75e063c10101455s.html

　　　　　http://blog.163.com/huai_jing@126/blog/static/171861983201172091718341/

参考文献：基于Gabor变换的特征提取及其应用

 

本文为原创内容，转载请注明出处！http://www.cnblogs.com/Jack-Lee/p/3649114.html