Gabor滤波器与特征提取

一、Gabor滤波器
Gabor滤波器，最主要使用优势体现在对物体纹理特征的提取上。

二维Gabor基函数能够很好地描述哺乳动物初级视觉系统中一对简单视觉神经元的感受野特性。随着小波变换和神经生理学的发展，Gabor变换逐渐演变成二维Gabor小波的形式。Gabor滤波器对于图像的亮度和对比度变化以及人脸姿态变化具有较强的健壮性，并且它表达的是对人脸识别最为有用的局部特征，故在计算机视觉及纹理分析中得到广泛的应用。

如果从Fourier变换的角度来看，Gobor变换就是窗函数取高斯窗时的短时Fourier变换。
如果从小波变换的角度来看，Gabor变换就是小波基函数取Gabor基的小波变换。

$f(x,y) =\frac{1}{2\pi*\sigma_x*\sigma_y\sqrt {1-\rho^2}}e^{-\frac{1}{2(1-\rho^2)}[({\frac{x-a}{\sigma_x}})^2+({\frac{y-b}{\sigma_y}})^2-({\frac{2\rho*(x-a)*(y-b)}{\sigma_x*\sigma_y}})]}$

$其中a,b,\sigma_1,\sigma_2,\rho都是常数，我们称(X,Y)服从参数为a,b,\sigma_1,\sigma_2,\rho 的**二维正态分布**，常把这个分布记作N(a,b,\sigma_1,\sigma_2,\rho)。$
$a,b,\sigma_1,\sigma_2,\rho的范围分为\sigma_1>0,\sigma_2>0,-\infty<a<+\infty,-\infty<b<+\infty,-1<\rho<1。$
$这个函数在三维空间中的图像好像是一个椭圆切面的钟倒扣在Oxy平面上，其中心在（a,b）点$

若$a=b=0,\rho=0$，
则$f(x,y) =\frac{1}{2\pi*\sigma_x*\sigma_y}e^{-\frac{1}{2}[({\frac{x}{\sigma_x}})^2+({\frac{y}{\sigma_y}})^2]}$。

Gabor小波实质上是由一个复指数信号调制一个二维高斯函数得到。

$g(x,y) =\frac{1}{2\pi*\sigma_x*\sigma_y}e^{-\frac{1}{2}[({\frac{x^{'}}{\sigma_x}})^2+({\frac{y^{'}}{\sigma_y}})^2]}e^{i2\pi *(\mu_0*x^{'}+\nu_0*y^{'})}$

$\left[ \begin{array}{cc} x^{'}\\ y^{'} \end{array} \right]=\left[ \begin{array}{cc} cos\theta&sin\theta\\ -sin\theta&cos\theta \end{array} \right] \left[ \begin{array}{cc} x\\ y \end{array} \right]$

$\theta：Gabor核函数的方向$
$\sigma_x,\sigma_y：高斯函数的标准差$

二、Gabor函数结合代码分析：
$g(x,y;\lambda，\theta，\phi，\sigma,\gamma) =e^{-\frac{1}{2}[({\frac{x^{'}}{\sigma}})^2+({\frac{\gamma y^{'}}{\sigma}})^2]}e^{i（ (\frac{2\pi*x^{'}}{\lambda}+\phi)）}$
实数：
$g(x,y;\lambda，\theta，\phi，\sigma,\gamma) =e^{-\frac{1}{2}[({\frac{x^{'}}{\sigma}})^2+({\frac{\gamma y^{'}}{\sigma}})^2]}cos{（ (\frac{2\pi*x^{'}}{\lambda}+\phi)）}$
虚数：
$g(x,y;\lambda，\theta，\phi，\sigma,\gamma) =e^{-\frac{1}{2}[({\frac{x^{'}}{\sigma}})^2+({\frac{\gamma y^{'}}{\sigma}})^2]}sin{（ (\frac{2\pi*x^{'}}{\lambda}+\phi)）}$

波长$（\lambda）$：它的值以像素为单位指定，通常大于等于2.但不能大于输入图像尺寸的五分之一。
方向$（\theta）$：这个参数指定了Gabor函数并行条纹的方向，它的取值为0到$2\pi$。
相位偏移$（\phi）$：它的取值范围为$-\pi到\pi$。其中，0、$\pi$度分别对应中心对称的center-on函数和center-off函数，而-$\pi/2$和$\pi/2$对应反对称函数。
长宽比$（\gamma）$：空间纵横比，决定了Gabor函数形状的椭圆率$（ellipticity）$。
当$γ= 1$时，形状是圆的。当$γ< 1$时，形状随着平行条纹方向而拉长。通常该值为0.5。
带宽$（b）$：$Gabor$滤波器的半响应空间频率带宽b和σ/ λ的比率有关，其中$σ$表示$Gabor$函数的高斯因子的标准差，如下：
$g(x,y;\lambda，\theta，\phi，\sigma,\gamma) =e^{-\frac{1}{2}[({\frac{x^{'}}{\sigma}})^2+({\frac{\gamma y^{'}}{\sigma}})^2]}sin{（ (\frac{2\pi*x^{'}}{\lambda}+\phi)）}$

$b =log_2\frac{\frac{\sigma*\pi}{\lambda}+\sqrt{\frac{ln2}{2}}}{\frac{\sigma*\pi}{\lambda}-\sqrt{\frac{ln2}{2}}}$

$\frac{\sigma}{\lambda}=\frac{1}{\pi}\sqrt{\frac{ln2}{2}}*\frac{2^b+1}{2^b-1}$

$σ$的值不能直接设置，它仅随着带宽$b$变化。带宽值必须是正实数，通常为1，此时，标准差和波长的关系为：$\sigma= 0.56 \lambda$。
带宽越小，标准差越大，$Gabor$形状越大，可见平行和抑制区条纹数量越多。

    //构造Gabor函数
    for (int k = 0;k < N;k++)                   //定义N方向
    {
        Theta = PI*((double)k/N);      //0*PI/4,1*PI/4,2*PI/4,3*PI/4
        for (int i = 0;i < gabor_height;i++)    //定义模版大小
        {
            for (int j = 0;j < gabor_width;j++)
            {
                x = j - gabor_width/2; //-5,-4,-3,-2,-1,0,+1,+2,+3,+4
                y = i - gabor_height/2;//-5,-4,-3,-2,-1,0,+1,+2,+3,+4

                xtmp = (double)x*cos(Theta) + (double)y*sin(Theta);
                ytmp = (double)y*cos(Theta) - (double)x*sin(Theta);

                tmp1 = /*(1/(2*PI*sigma_x*sigma_y))**/exp(-(pow(xtmp,2)/pow(sigma_x,2)+pow(ytmp,2)/pow(sigma_y,2))/2);
                tmp2 = cos(2*PI*xtmp/Lambda);
                tmp3 = sin(2*PI*xtmp/Lambda);

                re = tmp1*tmp2;
                im = tmp1*tmp3;

                rel[i*gabor_width+j] = re;
                img[i*gabor_width+j] = im;
                //printf("%f,%f\n",re,im);
            }
        }
        //用不同方向的Gabor函数对图像滤波并保存图片 
        g_filer(src_data,src->height,src->width,rel,img,10,10,k);
    }

参考资料：
基于测地线活动区域模型的非监督式纹理分割
何 源, 罗予频, 胡东成

三、参考书籍资料研究：
《人脸识别与人体动作识别技术与应用》 曹林著

四、参考网上代码：

// gabor.cpp : 定义控制台应用程序的入口点。

//#include "stdafx.h"
#include<iostream>
#include <opencv2/core/core.hpp>  
#include <opencv2/highgui/highgui.hpp>  
#include <opencv2/imgproc/imgproc.hpp>
#include "math.h"
#define  PI 3.1415926
#define  N 8
using namespace std;
using namespace cv;

void m_filer(double *src,int height,int width,double *mask_rel,double *mask_img,int mW,int mH,int k)
{
    IplImage *tmp;
    double a,b,c;
    char res[20];       //保存的图像名称

    tmp = cvCreateImage(cvSize(width,height),IPL_DEPTH_8U,1);

    for (int i = 0;i < height;i++)
    {
        for (int j = 0;j < width;j++)
        {
            a = 0.0;
            b = 0.0;
            c = 0.0;
            //去掉靠近边界的行列，无法滤波，超出范围 卷积
            if (i > int(mH/2) && i < height - int(mH/2) && j > int(mW) && j < width - int(mW/2))
            {
                for (int m = 0;m < mH;m++)
                {
                    for(int n = 0;n < mW;n++)
                    {
                        //printf("%f\n",src[(i+m-int(mH/2))*width+(j+n-int(mW))]);
                        a += src[(i+m-int(mH/2))*width+(j+n-int(mW))]*mask_rel[m*mW+n];
                        b += src[(i+m-int(mH/2))*width+(j+n-int(mW))]*mask_img[m*mW+n];
                        //printf("%f,%f\n",a,b);
                    }
                }
            }
            c = sqrt(a*a+b*b);
            c /= mW*mH;
            tmp->imageData[i*width+j] = (unsigned char)c;
        }
    }
    sprintf(res,"result%d.jpg",k);
    cvSaveImage(res,tmp);
    cvReleaseImage(&tmp);
}

int main(int argc, char* argv[])
{
    IplImage *src;
    double *rel,*img,*src_data,xtmp,ytmp,tmp1,tmp2,tmp3,re,im;
    double Theta,sigma,Gamma,Lambda,Phi;        //公式中的5个参数
    int gabor_height,gabor_width,x,y;

    src = cvLoadImage("1.jpg",CV_LOAD_IMAGE_GRAYSCALE);
    gabor_height = 10;
    gabor_width = 10;
    Gamma = 1.0;
    Lambda = 10.0;
    sigma = 100;
    Phi = 0;

    rel = (double *)malloc(sizeof(double)*gabor_width*gabor_height);//实数部分
    img = (double *)malloc(sizeof(double)*gabor_width*gabor_height);//虚数部分
    src_data = (double *)malloc(sizeof(double)*src->widthStep*src->height); //图像数据 

    for (int i=0;i<src->height;i++)
    {
        for (int j=0;j<src->widthStep;j++)
        {
            src_data[i*src->widthStep+j]=(unsigned char)src->imageData[i*src->widthStep+j];
            //printf("%f\n",src_data[i*src->widthStep+j]);
        }
    }

    //构造gabor函数
    for (int k = 0;k < N;k++)                   //定义N方向
    {
        Theta = PI*((double)k/N);      //0*PI/8,1*PI/8,2*PI/8,3*PI/8,4*PI/8,5*PI/8,6*PI/8,7*PI/8
        for (int i = 0;i < gabor_height;i++)    //定义模版大小
        {
            for (int j = 0;j < gabor_width;j++)
            {
                x = j - gabor_width/2; //-5,-4,-3,-2,-1,0,+1,+2,+3,+4
                y = i - gabor_height/2;//-5,-4,-3,-2,-1,0,+1,+2,+3,+4

                xtmp = (double)x*cos(Theta) + (double)y*sin(Theta);
                ytmp = (double)y*cos(Theta) - (double)x*sin(Theta);

                tmp1 = exp(-(pow(xtmp,2)+pow(ytmp*Gamma,2))/(2*pow(sigma,2)));
                tmp2 = cos(2*PI*xtmp/Lambda + Phi);
                tmp3 = sin(2*PI*xtmp/Lambda + Phi);

                re = tmp1*tmp2;
                im = tmp1*tmp3;

                rel[i*gabor_width+j] = re;
                img[i*gabor_width+j] = im;
                //printf("%f,%f\n",re,im);
            }
        }
        //用不同方向的GABOR函数对图像滤波并保存图片
        m_filer(src_data,src->height,src->width,rel,img,10,10,k);
    }

    free(rel);
    free(img);
    free(src_data);
    return 0;
}

Matlab研究分析Gabor滤波器
http://blog.csdn.net/forest_world/article/details/51018493