# Gabor filters for digital image processing

For the application of Gabor filters in feature-based (e.g. lines and edges in images) digital image processing, it is desired that the DC component is 0. Why is this so? What effect does it have mathematically if the DC component of the filter is not 0? Thanks in advance!

In short: If a filter have DC component it means that for a constant region of the image, the filter response will depend on the grey value of that region. Constant region with higher values ==> higher response.

One usually wants higher responses only in case of edges/corners and not in the case of bright but constant regions.

As a picture worth 1000 words, here is an example. I took a simple image consisting of two horisontal lines, each one in different background. I took the Gabor and LogGabor filter. LogGabor is a variant that have 0 DC response.

I filtered the image with a Gabor filter (tuned to respond to horizontal signals) and then, the same image, with LogGabor filter (tuned with exactly the same parameters as Gabor filter).

Original image: Gabor convoluted: LogGabor convoluted: Another example of Gabor filter, with a bandwidth of 0.7 octaves. Now, let's ignore the edges for a moment and focus on the constant regions. In the second image, you see different responses for the constant regions! Darker region have a lower response than the lighter one. For the LogGabor case (third image) you see that the response is the same regardless of the original intensity.

Note: The response values were scaled in 0-255 interval, so, neutral grey means 0 response, darker means negative values and lighter positive ones. The bandwidth is set to 2 octaves (except for the 4'th image)

EDIT

According to a paper by Boukerroui et al., for bandwidths around 0.7 the DC response is close to the computational errors ($10^{-6}$). The 4'th figure demonstrate this.

My personal advice to you is to use another bandpass filter. I strongly recommend Gauss Derivatives or Difference of Gaussians (in this order) if it is necessary for you to stay in spatial domain. If you don't have this restriction, I personally prefer LogGabor.

From my experience there is no net advantage or disadvantage in terms of end result between logGabor, Gaussian Derivatives, Difference of Gaussians or more exotic Cauchy filters. BUT this might vary, depending on your practical application.

As for the code, here is a snippet of my Gabor/LogGabor implementation. The code is in Matlab, the convolutions are done in Fourier space. The complexity is $O(n log n)$ and it depends only on the image dimensions, not on the filter dimensions. The code follows the guidelines in Boukerroui et al.. It is developed by me and released under MIT licence (not shown here) so feel free to use it as you see fit.

Img must be a MxNx1 of doubles, (preferably between 0 and 1 but it should work with 0-255 too); P0 is the tuning period (wavelength), that is, $P_0=1/f$, where f is tuning frequency; Orient is the filter orientation, expressed in radians; FBW frequency bandwidth in octaves; ABW the angle bandwidth.

function [ ReConv,ImConv ] = GaborConv( Img, P0,Orient, FBW, ABW,varargin )
ReConv=Img; ImConv=Img;
ImgFFt=fft2(Img);
lg2pi=sqrt(log(2)/pi);
A=P0*lg2pi*(2^FBW+1)/(2^FBW-1); B=P0*lg2pi*(1/tan(0.5*ABW));
K=1/(A*B); KF=2*pi/P0;
[N,M]=size(Img); u0=1/P0*cos(Orient); v0=1/P0*sin(Orient);
FFilt=zeros([M,N]);
stx=1/(M-1); sty=1/(N-1);
[U V]=meshgrid(-0.5:stx:0.5,-0.5:sty:0.5);
Ur=(U-u0)*cos(Orient)+(V-v0)*sin(Orient);
Vr=-(U-u0)*sin(Orient)+(V-v0)*cos(Orient);
FFilt=exp(-pi*(Ur.*Ur.*A.*A+Vr.*Vr.*B.*B)) ;
FFilt=ifftshift(FFilt);
normFact=sum(sum((abs(ifft2(FFilt)))));
RezFreq=FFilt.*ImgFFt/normFact; RezSpace=ifft2(RezFreq);
ReConv=real(RezSpace); ImConv=imag(RezSpace);
end


Log Gabor method:

function [  ReConv,ImConv ] = LogGaborConv( Img, P0,Orient, FBW, ABW,useCos ).
ReConv=Img;
ImConv=Img;
ImgFFt=fft2(Img);
%Compute the k_beta according to the desired frequency bandwidth
kBeta=exp(-1/4*sqrt(2*log(2))*FBW);
gaussSD=sqrt(ABW^2/(2*log(2))); %compute Gauss SD to ensure that Gauss(ABW,0,SD)==0.5
lnkB2=2*log(kBeta)^2;
nc=exp(-1/8*lnkB2/2)*sqrt(-2*sqrt(pi)/(1/P0*log(kBeta)));
[N M]=size(Img);
stx=1/(M-1);
sty=1/(N-1);
[Ug Vg]=meshgrid(-0.5:stx:0.5,-0.5:sty:0.5);
%prepare a low pass filter, as proposed by Kovesi
%Compute the bandpass filter
%Compute the angle component
Theta=atan2(Vg,Ug);
sintheta = sin(Theta);
costheta = cos(Theta);
ds = sintheta * cos(Orient) - costheta * sin(Orient);
dc = costheta * cos(Orient) + sintheta * sin(Orient);
dtheta = abs(atan2(ds,dc));
if ~exist('useCos','var')
useCos=0;
end
if useCos==0
%Use a gauss angle
AngFilter=exp(-dtheta.^2/(2*gaussSD^2));
else
%Use a cosine angle
AngFilter=cos(dtheta).^2.*heaviside(dtheta./pi+0.5).*heaviside(0.5-dtheta./pi);
end
logGaborFilt=ifftshift(logGaborFilt);
logGaborFilt(1,1)=0;
normFact=sum(sum((abs(ifft2(logGaborFilt)))));
RezFreq=logGaborFilt.*ImgFFt/normFact;
RezSpace=ifft2(RezFreq);
ReConv=real(RezSpace);
ImConv=imag(RezSpace);
end


As for how the DC response varies with various parameters and why, a comprehensive discussion is in Boukerroui et al., equations (19) and (20) at pages 75, 76. Take a peek at Appendix C (pag 101). I would like to redirect you to the paper rather than just copy/paste the relevant paragraphs.

Hope it helps!

• Hi @visoft, thanks for your answer! anyway is it possible to make the gabor filter 0 DC to improve its performance? – meta_warrior Feb 4 '14 at 13:21
• In theory, if you make the frequency bandwidth small enough (~0.7) the DC response will be low (close to the computational errors). I don't use it because lower bandwidth means larger spatial support (larger filter) I'll edit the post later with an example and some references. – visoft Feb 4 '14 at 13:28
• Can you show mathematically why does the DC response varies with the values of the image if it is nonzero? Also can you show me your code for Gabor and log Gabor respectively?:) – meta_warrior Feb 4 '14 at 15:07
• @freak_warrior I added more details in the answer. – visoft Feb 4 '14 at 16:31
• p.s. A good piece of code is in your previous question: dsp.stackexchange.com/questions/13907/… – visoft Feb 4 '14 at 16:34

Log-gabor filter is a very good one alternative to gabor filter has no DC component. Yet in ordinary gabor filter it is also possible to remove the DC component. See the Matlab code below:

function [GABOUT1,GABOUT2]=gaborfilter(I,S,F,W,P)

if isa(I,'double')~=1
I=double(I);
end
size=fix(1.5/S); % exp(-1.5^2*pi) < 0.1%

k=1;
for x=-size:size
for y=-size:size
G1(size+x+1,size+y+1)=k*exp(-pi*S^2*(x*x+y*y))*...
(exp(sqrt(-1)*(2*pi*F*(x*cos(W)+y*sin(W))+P)));
G2(size+x+1,size+y+1)=k*exp(-pi*S^2*(x*x+y*y))*...
(exp(sqrt(-1)*(2*pi*F*(x*cos(W)+y*sin(W))+P))-exp(-pi*(F/S)^2+sqrt(-1)*P));
end
end



In the code exp(-pi*(F/S)^2+sqrt(-1)*P)) is the DC component due to the overlapping. And GABOUT1 is the filtered image with DC component, while GABOUT2 is the filtered image with the compensation of DC component.

The issue with DC component is that one cannot achieve both broad spectral information and maximal spatial localization. As a result, you may over-estimate the low-frequency components but under-estimate the high-frequency components in a image. And this is especially reflected if you observe the magnitude of the filtered image. But there is no evident difference if you are observing the imag part of the image.

Check one example below:

I=imread('cat.jpg');
I=rgb2gray(I);
imshow(I),figure,
[G1,G2]=gaborfilter(I,0.5,0.1,0,0);
I1=abs(G1);
I1=I1/max(I1(:))*256;
imshow(uint8(I1))

I2=abs(G2);
I2=I2/max(I2(:))*256;
figure,imshow(uint8(I2))   Compared with the original image (left), the filtered one with DC component (middle) deteriorates at the high frequency, while the one without DC component (right) enhances the high-frequency details (tiny edges at various orientations) very much.

You can observe the imag part of the filtered image by changing the code to I1=imag(G1); and I2=imag(G2); and you may not see that much difference, i.e., the influence of DC component is not that much at the imag part of the filtered image.

• Hi lennon310: thanks for your answer. BTW u posted an answer for the prediction of fine-level wavelet coefficients. Could you explain more for the derivation of the formula for that question? – meta_warrior Feb 5 '14 at 1:23
• Thank you freak_warrior. I couldn't find that post, seems like being deleted for some reason. Actually that formula is from a publication, and I noticed later that they made some assumptions on that, so probably that's not the coefficient for general case. Sorry about that. – lennon310 Feb 5 '14 at 1:30
• @lennon310 The DC removal solution looks like Ronse et al. proposal. Am I right? – visoft Feb 7 '14 at 15:41
• @freak_warrior I see that nobody answered your remark for past 2 days, so here is a quick and dirty attempt: In the Appendix C of the reference paper I put in my answer, there are two methods that deal with DC response. Both are fairly complex and don't not completely solves the problem. IMHO whatever you do in the Fourier space or to the spatial filter itself will alter its properties in an unpredictable way and you will no longer have a Gabor filter. – visoft Feb 7 '14 at 15:43
• @visoft Movellan, J. R. - Tutorial on Gabor Filters. Tech. rep., 2002. – lennon310 Feb 7 '14 at 15:50