Gabor滤波器

Gabor是一个用于边缘提取的线性滤波器，其频率和方向表达与人类视觉系统类似，能够提供良好的方向选择和尺度选择特性，而且对于光照变化不敏感，因此十分适合纹理分析。

图1是Gabor滤波器和脊椎动物视觉皮层感受野响应的比较

　　　　　　
　　　　　　　　　　　　　　　　图1. Gabor与脊椎动物视觉皮层感受野响应的比较

图1中第一行是脊椎动物的视觉响应，第二行是Gabor滤波器的响应，可以看到，二者相差极小。
基于以上特性，Gabor滤波器被广泛应用于人脸识别的预处理。

---

Gabor理论及公式

我们知道，傅里叶变换可以将信号从时域转换到频域，但无法获得频谱中不同频率之间的先后关系。
然而实际应用中我们更多的关心信号局部范围内的的特性，Gabor和小波变换突破了傅里叶变换的局限性。

Gabor变换是D.Gabor于1946年提出的，为了提取傅里叶变换的局部信息，引入了时间局部化的窗函数（把信号划分成许多小的时间间隔，用傅里叶变换分析每一个间隔）。因此Gabor变换又称为窗口傅里叶变换（短时傅里叶变换）。

二维Gabor滤波器

在空间域，一个二维的Gabor滤波器是一个正弦平面波和高斯核函数的乘积。前者是调谐函数，后者是窗口函数。

g(x,y;λ,θ,ψ,σ,γ)=e−x′2+γ2y′22σ2ei(2πx′λ+ψ)

可以分为实部与虚部的形式

⎧⎩⎨⎪⎪⎪⎪greal(x,y;λ,θ,ψ,σ,γ)=e−x′2+γ2y′22σ2cos(2πx′λ+ψ)gimag(x,y;λ,θ,ψ,σ,γ)=e−x′2+γ2y′22σ2sin(2πx′λ+ψ)

其中

{x′=xcosθ+ysinθy′=−xsinθ+ycosθ

　　　　　　　　　　　　表1. 二维Gabor滤波器参数解释

参数 物理意义 描述 λ 波长 直接影响滤波器的滤波尺度，通常大于等于2 θ 方向 滤波器的方向 ψ 相位偏移 调谐函数的相位偏移，取值-180到180 γ 空间纵横比 决定滤波器的形状，取1时为圆形，通常取0.5 σ 带宽 高斯滤波器的方差，通常取2π

---

实验

本节参考自【图像处理】Gabor滤波器、Gabor滤波器学习、Gabor的OpenCV代码

　　　　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　图2. 取不同参数的Gabor滤波器

图2对比了不同参数的Gabor滤波器，从第1到第5行λ取3,6,9,12,15，从第1到第8列θ取0,π8,π4,3π8,π2,5π8,3π4,7π8，其他参数保持不变：ψ=0,σ=2π,γ=0.5。

图4使用图2的滤波器对图3中两幅人脸做卷积
　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　图3. FERET中同一人在不同光照不同角度的人脸

　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　图4.1. 使用图2滤波器对图3左卷积结果

　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　图4.2. 使用图2滤波器对图3右卷积结果

可以看到，不同方向不同尺度的Gabor滤波器可以提取人脸中不同特征
此外还可以看到，Gabor滤波结果对不同光照也能保持较稳定结果

---

代码

// gaborfilter.cpp#include "gaborfilter.h"int CGabor::getFilterSize(float theta, float sigma, float gamma){    float sigma_x = sigma*sigma;    float sigma_y = sigma_x/(gamma*gamma);    float sqrt_sigma_y = sqrt(sigma_y);    float c_theta = cos(theta);    float s_theta = sin(theta);    // calculate filter size (3sigma)    int nstds = 3;    int xmax = max(abs(nstds*sigma*c_theta), abs(nstds*sqrt_sigma_y*s_theta));    int ymax = max(abs(nstds*sigma*s_theta), abs(nstds*sqrt_sigma_y*c_theta));    int half_filter_size = xmax>ymax?xmax:ymax;    int filter_size = 2*half_filter_size+1;    return filter_size;}void CGabor::getGaborFilter(float lambda, float theta, float fhi, float sigma, float gamma,     Mat& realGabor, Mat& imagGabor, Mat& mag){    if(abs(lambda-0.0f) < 1e-6)        lambda = 1.0f;    float sigma_x = sigma*sigma;    float sigma_y = sigma_x/(gamma*gamma);    float sqrt_sigma_y = sqrt(sigma_y);    float c_theta = cos(theta);    float s_theta = sin(theta);    int filter_size = getFilterSize(theta, sigma, gamma);    int half_filter_size = (filter_size-1)/2;    realGabor = Mat(filter_size, filter_size, CV_32F);    imagGabor = Mat(filter_size, filter_size, CV_32F);    mag       = Mat(filter_size, filter_size, CV_32F);    for(int i=0; i<filter_size; i++)    {        float* p_real = realGabor.ptr<float>(i);        float* p_imag = imagGabor.ptr<float>(i);        float* p_mag  = mag.ptr<float>(i);        int y = i - half_filter_size;        for(int j=0; j<filter_size; j++)        {            int x = j - half_filter_size;            float x_theta = x*c_theta + y*s_theta;            float y_theta = y*c_theta - x*s_theta;            float value_ = exp(-0.5*(x_theta*x_theta/sigma_x + y_theta*y_theta/sigma_y));            p_real[j] = value_ * cos(2*CV_PI*x_theta/lambda + fhi);            p_imag[j] = value_ * sin(2*CV_PI*x_theta/lambda + fhi);            p_mag[j]  = sqrt(p_real[j]*p_real[j] + p_imag[j]*p_imag[j]);        }    }}Mat CGabor::crossGaborFilter(Mat gaborFilter, Mat image){    CV_Assert(!gaborFilter.empty() && !image.empty());    int half_filter_size = (max(gaborFilter.rows, gaborFilter.cols)-1)/2;    Mat filtered_img(image.rows, image.cols, CV_32F);    for(int i=0; i<image.rows; i++)    {        float* p_fil = filtered_img.ptr<float>(i);        for(int j=0; j<image.cols; j++)        {            float sum_value = 0.0f;            for(int m=0; m<gaborFilter.rows; m++)            {                float* p_gab = gaborFilter.ptr<float>(m);                int img_i = i + m - half_filter_size;                img_i = img_i<0 ? 0 : img_i;                img_i = img_i>=image.rows ? (image.rows-1) : img_i;                uchar* p_img = image.ptr<uchar>(img_i);                for(int n=0; n<gaborFilter.cols; n++)                {                    int img_j = j + n - half_filter_size;                    img_j = img_j<0 ? 0 : img_j;                    img_j = img_j>=image.cols ? (image.cols-1) : img_j;                    sum_value += ((float)p_img[img_j] * p_gab[n]);                }            }            p_fil[j] = sum_value;        }    }    return filtered_img;}Mat CGabor::normalizer(Mat src){    CV_Assert(!src.empty());    float min_ = FLT_MAX;    float max_ = FLT_MIN;    for(int i=0; i<src.rows; i++)    {        float* p_src = src.ptr<float>(i);        for(int j=0; j<src.cols; j++)        {            if(p_src[j] > max_)                max_ = p_src[j];            if(p_src[j] < min_)                min_ = p_src[j];        }    }    Mat src_show(src.size(), CV_8UC1);    float scale = max_ - min_;    for(int i=0; i<src.rows; i++)    {        float* p_src = src.ptr<float>(i);        uchar* p_src_show = src_show.ptr<uchar>(i);        for(int j=0; j<src.cols; j++)        {            if(scale > 0.01)                p_src_show[j] = (uchar)((p_src[j]-min_)*255/scale);            else                p_src_show[j] = 255;        }    }    return src_show;}

\\ main.cpp#include "stdafx.h"#include "gaborfilter.h"void main(){    Mat img[7];    string path = "C:\\Users\\sq\\Desktop\\humanface\\FaceDataSet\\FERET_80_80\\FERET-001\\";    for(int i=0; i<7; i++)    {        char num[7];        sprintf(num, "0%d.tif", i+1);        string wholepath = num;        wholepath = path + wholepath;        Mat tmp = imread(wholepath, 0);        resize(tmp, img[i], Size(200, 200)); // 为了确保滤波器尺寸小于图片大小    }    // calc gabor    CGabor m_gabor;    float lambda = 5;    float theta = 0;    float fhi = 0;    float sigma = 2*CV_PI;    float gamma = 0.5;    Mat realGabor[40], imagGabor, magGabor;    int size_ = INT_MAX;    for(int i=0; i<8; i++)    {        theta = i*CV_PI/8;        int _size_ = m_gabor.getFilterSize(theta, sigma, gamma);        if(_size_ < size_)            size_ = _size_;    }    Mat showit((size_+3)*5-3, (size_+3)*8-3, CV_8UC1, Scalar::all(255));    for(int i=0; i<8; i++)    {        theta = i*CV_PI/8;        for(int j=0; j<5; j++)        {            lambda = (j+1)*3;            m_gabor.getGaborFilter(lambda, theta, fhi, sigma, gamma, realGabor[j*8+i], imagGabor, magGabor);            Mat saveimg = m_gabor.normalizer(realGabor[j*8+i]);            Mat saveimg_roi(saveimg, Rect(saveimg.rows/2-size_/2, saveimg.rows/2-size_/2, size_, size_));            int start_row = j*(size_+3);            int start_col = i*(size_+3);            for(int m=0; m<size_; m++)            {                uchar* p_img = saveimg_roi.ptr<uchar>(m);                uchar* p_show = showit.ptr<uchar>(start_row + m);                for(int n=0; n<size_; n++)                   p_show[start_col + n] = p_img[n];            }        }    }    imshow("gabor", showit);    imwrite("gabor.jpg", showit);    waitKey(1);    // filter face    for(int i=0; i<7; i++)    {        Mat face((img[0].rows+3)*5-3, (img[0].cols+3)*8-3, CV_8UC1, Scalar::all(255));        for(int m=0; m<5; m++)        {            int start_row = m*(img[0].rows+3);            for(int n=0; n<8; n++)            {                Mat filted = m_gabor.crossGaborFilter(realGabor[m*8+n], img[i]);                Mat filted_norm = m_gabor.normalizer(filted);                int start_col = n*(img[0].cols+3);                for(int ro=0; ro<img[0].rows; ro++)                {                    uchar* p_filted = filted_norm.ptr<uchar>(ro);                    uchar* p_face = face.ptr<uchar>(start_row + ro);                    for(int co=0; co<img[0].cols; co++)                        p_face[start_col + co] = p_filted[co];                }            }        }        char num[7];        sprintf(num, "_%d.jpg", i);        string filename = num;        filename = "face"+filename;        imwrite(filename, face);        imshow(filename, face);        waitKey(1);    }    waitKey(0);}

---