Matlab gabor函数解析

涉及函数：gabor ，imgaborfilt。注意，这两个函数是一套，相辅相成，其中gabor是一个类，这加大了源码解读难度。这两个函数都是在matlab 2015b中引入的。

1 gabor基本公式

这里首先贴上gabor基本公式(wikipedia)。
值得注意的是

matlab自带gabor函数利用的是公式1，也就是含有实虚部的gabor公式；
opencv 3.1自带的gabor使用的是公式2，仅使用了gabor滤波器的实部，目的是减少计算量并带来相似的效果。

2 matlab的gabor类参数解析

根据原始gabor公式，拥有以下参数：

λ（lambda, 波长） θ（theta, 角度），σ（sigma）, γ（gamma）， φ（phi，相移）。具体的每个参数的作用以及设置技巧，网上太多资料，自己收集，这里就不多做说明。

根据matlab gabor函数，拥有以下输入参数：
wavelength ， orientation， SpatialFrequencyBandwidth， SpatialAspectRation
其中，

wavelength 就是原始公式中的λ(lambda)，
orientation就是原始公式中的θ(theta，这里是角度制)
SpatialAspectRatio 就是原始公式中的 , γ（gamma）
SpatialFrequencyBandwidth(简称BW)则比较特殊，并不对应于其中的任何一个参数，但是它和σ(sigma)与λ(lambda)有以下关系。 也就是参数σ(sigma)可通过BW这个参数计算得到
$\sigma = \frac{\lambda}{\pi} \sqrt{\frac{ln2}{2}} \frac{2^{BW}+1}{2^{BW}-1}$

另外： **参数 φ(phi)**在gabor函数中默认定义为0

总结

参数	λ（lambda, 波长）	θ（theta, 角度）	γ（gamma）	σ（sigma）	φ（phi，相移）
Matlab里的gabor参数	wavelength	orientation（角度制）	SpatialAspectRatio	由参数SpatialFrequencyBandwidth和λ（lambda）转化得到	默认为0

3 源码解读

详细完整的gabor和imgaborfilt源码，自行在matlab中查看，这里只截取各自最为重要的部分；

3.1 gabor源码

3.1.1 空域gabor

$\qquad$ $\color{Red}{\text{注意，大坑！！}}$ ， get.SpatialKernel()返回的gabor核是你所能看到的gabor核，接下来imgaborfilt函数并不用此gabor核对你想要滤波的图进行卷积！
$\qquad$ 简单的说，更改此函数里面的内容，可视化gabor核会发现不同，但是会发现imgaborfilt函数所输出的gabor响应图并没有什么不同。

function h = get.SpatialKernel(self)% Parameterization of spatial kernel frequency includes Phi as% an independent variable. We use a constant of 0.phi = 0 ;%φ默认为0% 这里把σ分为 σ_x, σ_y.与公式对应的关系是： %σ_x = σ；  σ_y = σ/γsigmax = self.SigmaX;sigmay = self.SigmaY;%获取gabor核的大小， 也就是gabor滤波器的大小（正方形）r = max(self.Rx, self.Ry);[X,Y] = meshgrid(-r:r,-r:r);%可以和公式对应上的Xprime = X .*cosd(self.Orientation) - Y .*sind(self.Orientation);Yprime = X .*sind(self.Orientation) + Y .*cosd(self.Orientation);%可以和公式对应上hGaussian = exp( -1/2*( Xprime.^2 ./ sigmax^2 + Yprime.^2 ./ sigmay^2));%gaussianhGaborEven = hGaussian.*cos(2*pi.*Xprime ./ self.Wavelength+phi);%实部hGaborOdd  = hGaussian.*sin(2*pi.*Xprime ./ self.Wavelength+phi);%虚部h = complex(hGaborEven, hGaborOdd);%复数形式end

3.1.2 频域gabor

$\qquad \color{Red}{\text{注意：这才是实际上用来和图像做卷积运算的gabor核。}}$ 或许不该叫gabor核，并不是我们所事先预料的那种生成gabor核然后和图像卷积的方式。而是直接依据gabor函数的输入参数在频域构建出与原图大小相同的gabor核的频率图！！！
$\qquad$ 而且这个函数(或者按照c++叫法，类成员函数)是在imgaborfilt被调用的。这种方式会显著加快速度，算法依据是A，K，J的论文(大牛恐怖如斯啊~)；

function H = makeFrequencyDomainTransferFunction(self, imageSize, classA)% Directly construct frequency domain transfer function of% Gabor filter. (Jain, Farrokhnia, "Unsupervised Texture% Segmentation Using Gabor Filters", 1999)% Directly construct frequency domain transfer function of% Gabor filter. (Jain, Farrokhnia, "Unsupervised Texture% Segmentation Using Gabor Filters", 1999)M = imageSize(1);N = imageSize(2);u = cast(images.internal.createNormalizedFrequencyVector(N),classA);v = cast(images.internal.createNormalizedFrequencyVector(M),classA);[U,V] = meshgrid(u,v);Uprime = U .*cosd(self.Orientation) - V .*sind(self.Orientation);Vprime = U .*sind(self.Orientation) + V .*cosd(self.Orientation);sigmau = 1/(2*pi*self.SigmaX);sigmav = 1/(2*pi*self.SigmaY);freq = 1/self.Wavelength;A = 2*pi*self.SigmaX*self.SigmaY;H = A.*exp(-0.5*( ((Uprime-freq).^2)./sigmau^2 + Vprime.^2 ./ sigmav^2) );end

3.2 imgaborfilt 源码

得到响应图的过程调用了gabor类的成员函数，以加快速度

%这段代码来源于matlab
for p = 1:length(GaborBank)%GaborBank是生成的gabor kernel的集合%sizeAPadded是对原始图像A（待卷积图）边界处理以后的尺寸；%class(A)是原始图像A的数据类型%函数makeFrequencyDomainTransferFunction是在类gabor的成员函数H = makeFrequencyDomainTransferFunction(GaborBank(p),sizeAPadded,class(A));outPadded = ifft2(A .* ifftshift(H));out(:,:,p) = unpadSlice(outPadded,padSize,outSize);
end

上述代码等价于下面代码现，下面的代码是我们通常所理解的方式

%这段代码与上面代码相同效果，使我们通常所理解的实现方式
for p = 1:length(GaborBank)outPadded = conv2(a, GaborBank(p).SpatialKernel, 'same');out(:,:,p) = unpadSlice(outPadded,padSize,outSize);end