为什么要进行解析,因为自带的gabor函数有个小坑, 转opencv的时候,因为没有完全理解自带的gabor源码被小小的坑了一下, 所以做一下记录, 以方便后人。‘
版本是2016B
文章目录
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这个参数计算得到
σ = λ π l n 2 2 2 B W + 1 2 B W − 1 \sigma = \frac{\lambda}{\pi} \sqrt{\frac{ln2}{2}} \frac{2^{BW}+1}{2^{BW}-1} σ=πλ2ln22BW−12BW+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));%gaussian
hGaborEven = 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
注意: 这才是实际上用来和图像做卷积运算的gabor核。
\qquad \color{Red}{\text{注意: 这才是实际上用来和图像做卷积运算的gabor核。}}
注意: 这才是实际上用来和图像做卷积运算的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
总结
最后,如果要在opencv上取得相同的效果, 需要增加gabor的虚数部分, 如果要加快速度,需要采用A.K.J的实现方式~