zoukankan      html  css  js  c++  java

  • Gabor的OpenCV代码

    转自http://blog.csdn.net/guoming0000/article/details/7839917

    最近弄人脸识别，用到Gabor卷积核，但网上的代码似乎没有和我心意的，于是参考了自己写了下！参考了Zhou Mian以及matlab的Gabor实现代码的代码。虽然OpenCV的imporc下面有个gabor.cpp，但那个是一般形式的公式，不是用来做人脸识别的，可以参考文献A review on Gabor wavelets for face recognition，又说到。上代码和链接地址！下载地址~

       目前代码未经过更多的测试，不少功能为加入，但可以满足许多人的使用和参考了吧，很多人肯定非常非常需要，先开源下，欢迎指出错误之处。

    [cpp] view plaincopyprint?
    1. //GaborFR.h  
    [cpp] view plaincopyprint?
    1. #pragma once  
    2. #include "opencv2opencv.hpp"  
    3. #include <vector>  
    4. using namespace std;  
    5. using namespace cv;  
    6. class GaborFR  
    7. {  
    8. public:  
    9.     GaborFR();  
    10.     static Mat  getImagGaborKernel(Size ksize, double sigma, double theta,   
    11.                                     double nu,double gamma=1, int ktype= CV_32F);  
    12.     static Mat  getRealGaborKernel( Size ksize, double sigma, double theta,   
    13.                                     double nu,double gamma=1, int ktype= CV_32F);  
    14.     static Mat  getPhase(Mat &real,Mat &imag);  
    15.     static Mat  getMagnitude(Mat &real,Mat &imag);  
    16.     static void getFilterRealImagPart(Mat& src,Mat& real,Mat& imag,Mat &outReal,Mat &outImag);  
    17.     static Mat  getFilterRealPart(Mat& src,Mat& real);  
    18.     static Mat  getFilterImagPart(Mat& src,Mat& imag);  
    19.   
    20.     void        Init(Size ksize=Size(19,19), double sigma=2*CV_PI,  
    21.                     double gamma=1, int ktype=CV_32FC1);  
    22. private:  
    23.     vector<Mat> gaborRealKernels;  
    24.     vector<Mat> gaborImagKernels;  
    25.     bool isInited;  
    26. };  

    [cpp] view plaincopyprint?
    1. #include "stdafx.h"  
    2. #include "GaborFR.h"  
    3. GaborFR::GaborFR()  
    4. {  
    5.     isInited = false;  
    6. }  
    7. void GaborFR::Init(Size ksize, double sigma,double gamma, int ktype)  
    8. {  
    9.     gaborRealKernels.clear();  
    10.     gaborImagKernels.clear();  
    11.     double mu[8]={0,1,2,3,4,5,6,7};  
    12.     double nu[5]={0,1,2,3,4};  
    13.     int i,j;  
    14.     for(i=0;i<5;i++)  
    15.     {  
    16.         for(j=0;j<8;j++)  
    17.         {  
    18.             gaborRealKernels.push_back(getRealGaborKernel(ksize,sigma,mu[j]*CV_PI/8,nu[i],gamma,ktype));  
    19.             gaborImagKernels.push_back(getImagGaborKernel(ksize,sigma,mu[j]*CV_PI/8,nu[i],gamma,ktype));  
    20.         }  
    21.     }  
    22.     isInited = true;  
    23. }  
    24. Mat GaborFR::getImagGaborKernel(Size ksize, double sigma, double theta, double nu,double gamma, int ktype)  
    25. {  
    26.     double  sigma_x     = sigma;  
    27.     double  sigma_y     = sigma/gamma;  
    28.     int     nstds       = 3;  
    29.     double  kmax        = CV_PI/2;  
    30.     double  f           = cv::sqrt(2.0);  
    31.     int xmin, xmax, ymin, ymax;  
    32.     double c = cos(theta), s = sin(theta);  
    33.     if( ksize.width > 0 )  
    34.     {  
    35.         xmax = ksize.width/2;  
    36.     }  
    37.     else//这个和matlab中的结果一样，默认都是19 !  
    38.     {  
    39.         xmax = cvRound(std::max(fabs(nstds*sigma_x*c), fabs(nstds*sigma_y*s)));  
    40.     }  
    41.     if( ksize.height > 0 )  
    42.     {  
    43.         ymax = ksize.height/2;  
    44.     }  
    45.     else  
    46.     {  
    47.         ymax = cvRound(std::max(fabs(nstds*sigma_x*s), fabs(nstds*sigma_y*c)));  
    48.     }  
    49.     xmin = -xmax;  
    50.     ymin = -ymax;  
    51.     CV_Assert( ktype == CV_32F || ktype == CV_64F );  
    52.     float*  pFloat;  
    53.     double* pDouble;  
    54.     Mat kernel(ymax - ymin + 1, xmax - xmin + 1, ktype);  
    55.     double k        =   kmax/pow(f,nu);  
    56.     double scaleReal=   k*k/sigma_x/sigma_y;  
    57.     for( int y = ymin; y <= ymax; y++ )  
    58.     {  
    59.         if( ktype == CV_32F )  
    60.         {  
    61.             pFloat = kernel.ptr<float>(ymax-y);  
    62.         }  
    63.         else  
    64.         {  
    65.             pDouble = kernel.ptr<double>(ymax-y);  
    66.         }  
    67.         for( int x = xmin; x <= xmax; x++ )  
    68.         {  
    69.             double xr = x*c + y*s;  
    70.             double v = scaleReal*exp(-(x*x+y*y)*scaleReal/2);  
    71.             double temp=sin(k*xr);  
    72.             v   =  temp*v;  
    73.             if( ktype == CV_32F )  
    74.             {  
    75.                 pFloat[xmax - x]= (float)v;  
    76.             }  
    77.             else  
    78.             {  
    79.                 pDouble[xmax - x] = v;  
    80.             }  
    81.         }  
    82.     }  
    83.     return kernel;  
    84. }  
    85. //sigma一般为2*pi  
    86. Mat GaborFR::getRealGaborKernel( Size ksize, double sigma, double theta,   
    87.     double nu,double gamma, int ktype)  
    88. {  
    89.     double  sigma_x     = sigma;  
    90.     double  sigma_y     = sigma/gamma;  
    91.     int     nstds       = 3;  
    92.     double  kmax        = CV_PI/2;  
    93.     double  f           = cv::sqrt(2.0);  
    94.     int xmin, xmax, ymin, ymax;  
    95.     double c = cos(theta), s = sin(theta);  
    96.     if( ksize.width > 0 )  
    97.     {  
    98.         xmax = ksize.width/2;  
    99.     }  
    100.     else//这个和matlab中的结果一样，默认都是19 !  
    101.     {  
    102.         xmax = cvRound(std::max(fabs(nstds*sigma_x*c), fabs(nstds*sigma_y*s)));  
    103.     }  
    104.   
    105.     if( ksize.height > 0 )  
    106.         ymax = ksize.height/2;  
    107.     else  
    108.         ymax = cvRound(std::max(fabs(nstds*sigma_x*s), fabs(nstds*sigma_y*c)));  
    109.     xmin = -xmax;  
    110.     ymin = -ymax;  
    111.     CV_Assert( ktype == CV_32F || ktype == CV_64F );  
    112.     float*  pFloat;  
    113.     double* pDouble;  
    114.     Mat kernel(ymax - ymin + 1, xmax - xmin + 1, ktype);  
    115.     double k        =   kmax/pow(f,nu);  
    116.     double exy      =   sigma_x*sigma_y/2;  
    117.     double scaleReal=   k*k/sigma_x/sigma_y;  
    118.     int    x,y;  
    119.     for( y = ymin; y <= ymax; y++ )  
    120.     {  
    121.         if( ktype == CV_32F )  
    122.         {  
    123.             pFloat = kernel.ptr<float>(ymax-y);  
    124.         }  
    125.         else  
    126.         {  
    127.             pDouble = kernel.ptr<double>(ymax-y);  
    128.         }  
    129.         for( x = xmin; x <= xmax; x++ )  
    130.         {  
    131.             double xr = x*c + y*s;  
    132.             double v = scaleReal*exp(-(x*x+y*y)*scaleReal/2);  
    133.             double temp=cos(k*xr) - exp(-exy);  
    134.             v   =   temp*v;  
    135.             if( ktype == CV_32F )  
    136.             {  
    137.                 pFloat[xmax - x]= (float)v;  
    138.             }  
    139.             else  
    140.             {  
    141.                 pDouble[xmax - x] = v;  
    142.             }  
    143.         }  
    144.     }  
    145.     return kernel;  
    146. }  
    147. Mat GaborFR::getMagnitude(Mat &real,Mat &imag)  
    148. {  
    149.     CV_Assert(real.type()==imag.type());  
    150.     CV_Assert(real.size()==imag.size());  
    151.     int ktype=real.type();  
    152.     int row = real.rows,col = real.cols;  
    153.     int i,j;  
    154.     float*  pFloat,*pFloatR,*pFloatI;  
    155.     double* pDouble,*pDoubleR,*pDoubleI;  
    156.     Mat     kernel(row, col, real.type());  
    157.     for(i=0;i<row;i++)  
    158.     {  
    159.         if( ktype == CV_32FC1 )  
    160.         {  
    161.             pFloat = kernel.ptr<float>(i);  
    162.             pFloatR= real.ptr<float>(i);  
    163.             pFloatI= imag.ptr<float>(i);  
    164.         }  
    165.         else  
    166.         {  
    167.             pDouble = kernel.ptr<double>(i);  
    168.             pDoubleR= real.ptr<double>(i);  
    169.             pDoubleI= imag.ptr<double>(i);  
    170.         }  
    171.         for(j=0;j<col;j++)  
    172.         {  
    173.             if( ktype == CV_32FC1 )  
    174.             {  
    175.                 pFloat[j]= sqrt(pFloatI[j]*pFloatI[j]+pFloatR[j]*pFloatR[j]);  
    176.             }  
    177.             else  
    178.             {  
    179.                 pDouble[j] = sqrt(pDoubleI[j]*pDoubleI[j]+pDoubleR[j]*pDoubleR[j]);  
    180.             }  
    181.         }  
    182.     }  
    183.     return kernel;  
    184. }  
    185. Mat GaborFR::getPhase(Mat &real,Mat &imag)  
    186. {  
    187.     CV_Assert(real.type()==imag.type());  
    188.     CV_Assert(real.size()==imag.size());  
    189.     int ktype=real.type();  
    190.     int row = real.rows,col = real.cols;  
    191.     int i,j;  
    192.     float*  pFloat,*pFloatR,*pFloatI;  
    193.     double* pDouble,*pDoubleR,*pDoubleI;  
    194.     Mat     kernel(row, col, real.type());  
    195.     for(i=0;i<row;i++)  
    196.     {  
    197.         if( ktype == CV_32FC1 )  
    198.         {  
    199.             pFloat = kernel.ptr<float>(i);  
    200.             pFloatR= real.ptr<float>(i);  
    201.             pFloatI= imag.ptr<float>(i);  
    202.         }  
    203.         else  
    204.         {  
    205.             pDouble = kernel.ptr<double>(i);  
    206.             pDoubleR= real.ptr<double>(i);  
    207.             pDoubleI= imag.ptr<double>(i);  
    208.         }  
    209.         for(j=0;j<col;j++)  
    210.         {  
    211.             if( ktype == CV_32FC1 )  
    212.             {  
    213. //              if(pFloatI[j]/(pFloatR[j]+pFloatI[j]) > 0.99)  
    214. //              {  
    215. //                  pFloat[j]=CV_PI/2;  
    216. //              }  
    217. //              else  
    218. //              {  
    219. //                  pFloat[j] = atan(pFloatI[j]/pFloatR[j]);  
    220.                 pFloat[j] = asin(pFloatI[j]/sqrt(pFloatR[j]*pFloatR[j]+pFloatI[j]*pFloatI[j]));  
    221. /*              }*/  
    222. //              pFloat[j] = atan2(pFloatI[j],pFloatR[j]);  
    223.             }//CV_32F  
    224.             else  
    225.             {  
    226.                 if(pDoubleI[j]/(pDoubleR[j]+pDoubleI[j]) > 0.99)  
    227.                 {  
    228.                     pDouble[j]=CV_PI/2;  
    229.                 }  
    230.                 else  
    231.                 {  
    232.                     pDouble[j] = atan(pDoubleI[j]/pDoubleR[j]);  
    233.                 }  
    234. //              pDouble[j]=atan2(pDoubleI[j],pDoubleR[j]);  
    235.             }//CV_64F  
    236.         }  
    237.     }  
    238.     return kernel;  
    239. }  
    240. Mat GaborFR::getFilterRealPart(Mat& src,Mat& real)  
    241. {  
    242.     CV_Assert(real.type()==src.type());  
    243.     Mat dst;  
    244.     Mat kernel;  
    245.     flip(real,kernel,-1);//中心镜面  
    246. //  filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_CONSTANT);  
    247.     filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_REPLICATE);  
    248.     return dst;  
    249. }  
    250. Mat GaborFR::getFilterImagPart(Mat& src,Mat& imag)  
    251. {  
    252.     CV_Assert(imag.type()==src.type());  
    253.     Mat dst;  
    254.     Mat kernel;  
    255.     flip(imag,kernel,-1);//中心镜面  
    256. //  filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_CONSTANT);  
    257.     filter2D(src,dst,CV_32F,kernel,Point(-1,-1),0,BORDER_REPLICATE);  
    258.     return dst;  
    259. }  
    260. void GaborFR::getFilterRealImagPart(Mat& src,Mat& real,Mat& imag,Mat &outReal,Mat &outImag)  
    261. {  
    262.     outReal=getFilterRealPart(src,real);  
    263.     outImag=getFilterImagPart(src,imag);  
    264. }  

    main

    [cpp] view plaincopyprint?
    1. // Win32TestPure.cpp : 定义控制台应用程序的入口点。  
    2. #include "stdafx.h"  
    3. #include <vector>  
    4. #include <deque>  
    5. #include <iomanip>  
    6. #include <stdexcept>  
    7. #include <string>  
    8. #include <iostream>  
    9. #include <fstream>  
    10. #include <direct.h>//_mkdir()  
    11. #include "opencv2opencv.hpp"  
    12. #include "GaborFR.h"  
    13. using namespace std;  
    14. using namespace cv;  
    15. int main()  
    16. {   
    17.     //Mat M = getGaborKernel(Size(9,9),2*CV_PI,u*CV_PI/8, 2*CV_PI/pow(2,CV_PI*(v+2)/2),1,0);  
    18.     Mat saveM;  
    19.     //s8-4  
    20.     //s1-5  
    21.     //s1中年男人  
    22.     Mat I=imread("H:\pic\s1-5.bmp",-1);  
    23.     normalize(I,I,1,0,CV_MINMAX,CV_32F);  
    24.     Mat showM,showMM;Mat M,MatTemp1,MatTemp2;  
    25.     Mat line;  
    26.     int iSize=50;//如果数值比较大，比如50则接近论文中所述的情况了！估计大小和处理的源图像一样！  
    27.     for(int i=0;i<8;i++)  
    28.     {  
    29.         showM.release();  
    30.         for(int j=0;j<5;j++)  
    31.         {  
    32.             Mat M1= GaborFR::getRealGaborKernel(Size(iSize,iSize),2*CV_PI,i*CV_PI/8+CV_PI/2, j,1);  
    33.             Mat M2 = GaborFR::getImagGaborKernel(Size(iSize,iSize),2*CV_PI,i*CV_PI/8+CV_PI/2, j,1);  
    34.             //加了CV_PI/2才和大部分文献的图形一样，不知道为什么！  
    35.             Mat outR,outI;  
    36.             GaborFR::getFilterRealImagPart(I,M1,M2,outR,outI);  
    37. //          M=GaborFR::getPhase(M1,M2);  
    38. //          M=GaborFR::getMagnitude(M1,M2);  
    39. //          M=GaborFR::getPhase(outR,outI);  
    40. //          M=GaborFR::getMagnitude(outR,outI);  
    41.  //         M=GaborFR::getMagnitude(outR,outI);  
    42. //          MatTemp2=GaborFR::getPhase(outR,outI);  
    43. //          M=outR;  
    44.              M=M1;  
    45.             //      resize(M,M,Size(100,100));  
    46.             normalize(M,M,0,255,CV_MINMAX,CV_8U);  
    47.             showM.push_back(M);  
    48.             line=Mat::ones(4,M.cols,M.type())*255;  
    49.             showM.push_back(line);  
    50.         }  
    51.         showM=showM.t();  
    52.         line=Mat::ones(4,showM.cols,showM.type())*255;  
    53.         showMM.push_back(showM);  
    54.         showMM.push_back(line);  
    55.     }  
    56.     showMM=showMM.t();  
    57. //  bool flag=imwrite("H:\out.bmp",showMM);  
    58.     imshow("saveMM",showMM);waitKey(0);  
    59.     return 0;  
    60. }//endof   main()  

    一下图片可能和程序实际运行结果有点不同，图片只是示意图，代码暂时没问题。需要考虑的是iSize大小问题，首先iSize要用奇数，然后大部分文献iSize都比较大，好像是100左右，但没看到他们描述过卷积核的大小。




  • 相关阅读:
    11 MySQL视图
    10 MySQL索引选择与使用
    08 MySQL存储引擎
    09 MySQL字符集
    06 MySQL运算符
    07 MySQL常用内置函数
    05 MySQL数据类型的选择与使用
    04 MySQL数据类型
    js 当前日期后7天
    md5加密
  • 原文地址：https://www.cnblogs.com/walccott/p/4956946.html
Copyright © 2011-2022 走看看