Gabor filter 实现的几种方式

转载自：http://blog.csdn.net/watkinsong[+]

1. 方式一

 

Sx,Sy在公式里分别表示Guass函数沿着x,y轴的标准差，相当于其他的gabor函数中的 sigma. 同时也用Sx,Sy指定了gabor滤波器的大小。（滤波器矩阵的大小）

这里没有考虑到相位偏移.

 

[csharp]  view plain copy

1. %%%%%%%VERSION 2  
2. %%ANOTHER DESCRIBTION OF GABOR FILTER  
3.   
4. %The Gabor filter is basically a Gaussian (with variances sx and sy along x and y-axes respectively)  
5. %modulated by a complex sinusoid (with centre frequencies U and V along x and y-axes respectively)   
6. %described by the following equation  
7. %%  
8. %                            -1     x' ^     y'  ^               
9. %%% G(x,y,theta,f) =  exp ([----{(----) 2+(----) 2}])*cos(2*pi*f*x');  
10. %                             2    sx'       sy'  
11. %%% x' = x*cos(theta)+y*sin(theta);  
12. %%% y' = y*cos(theta)-x*sin(theta);  
13.   
14. %% Describtion :  
15.   
16. %% I : Input image  
17. %% Sx & Sy : Variances along x and y-axes respectively  
18. %% f : The frequency of the sinusoidal function  
19. %% theta : The orientation of Gabor filter  
20.   
21. %% G : The output filter as described above  
22. %% gabout : The output filtered image  
23.   
24.   
25.   
26. %%  Author : Ahmad poursaberi  e-mail : a.poursaberi@ece.ut.ac.ir  
27. %%          Faulty of Engineering, Electrical&Computer Department,Tehran  
28. %%          University,Iran,June 2004  
29.   
30. function [G,gabout] = gaborfilter1(I,Sx,Sy,f,theta)  
31.   
32. if isa(I,'double')~=1   
33.     I = double(I);  
34. end  
35.   
36. for x = -fix(Sx):fix(Sx)  
37.     for y = -fix(Sy):fix(Sy)  
38.         xPrime = x * cos(theta) + y * sin(theta);  
39.         yPrime = y * cos(theta) - x * sin(theta);  
40.         G(fix(Sx)+x+1,fix(Sy)+y+1) = exp(-.5*((xPrime/Sx)^2+(yPrime/Sy)^2))*cos(2*pi*f*xPrime);  
41.     end  
42. end  
43.   
44. Imgabout = conv2(I,double(imag(G)),'same');  
45. Regabout = conv2(I,double(real(G)),'same');  
46.   
47. gabout = sqrt(Imgabout.*Imgabout + Regabout.*Regabout);  

这个gaobor函数生成的gabor滤波器的图像如图：

不同的参数也会导致不同的结果。

 

2. 方式二

 

这个gabor滤波器的实现增加了尺度和方向变换，其他的参数以及意义都和上面的一样。

 

[csharp]  view plain copy

1. %%%%%%%VERSION 1  
2.   
3. %The Gabor filter is basically a Gaussian (with variances sx and sy along x and y-axes respectively)  
4. %modulated by a complex sinusoid (with centre frequencies U and V along x and y-axes respectively)   
5. %described by the following equation  
6. %%  
7. %               1                -1     x  ^    y  ^  
8. %%% G(x,y) = ---------- * exp ([----{(----) 2+(----) 2}+2*pi*i*(Ux+Vy)])  
9. %            2*pi*sx*sy           2    sx       sy  
10.   
11. %% Describtion :  
12.   
13. %% I : Input image  
14. %% Sx & Sy : Variances along x and y-axes respectively  
15. %% U & V : Centre frequencies  along x and y-axes respectively  
16.   
17. %% G : The output filter as described above  
18. %% gabout : The output filtered image  
19.   
20. %%  Author : Ahmad poursaberi  e-mail : a.poursaberi@ece.ut.ac.ir  
21. %%          Faulty of Engineering, Electrical&Computer Department,Tehran  
22. %%          University,Iran,June 2004  
23.   
24. function [G,gabout] = gaborfilter(I,Sx,Sy,U,V);  
25.   
26. if isa(I,'double')~=1   
27.     I = double(I);  
28. end  
29.   
30. for x = -fix(Sx):fix(Sx)  
31.     for y = -fix(Sy):fix(Sy)  
32.         G(fix(Sx)+x+1,fix(Sy)+y+1) = (1/(2*pi*Sx*Sy))*exp(-.5*((x/Sx)^2+(y/Sy)^2)+2*pi*i*(U*x+V*y));  
33.     end  
34. end  
35.   
36. Imgabout = conv2(I,double(imag(G)),'same');  
37. Regabout = conv2(I,double(real(G)),'same');  
38.   
39. gabout = uint8(sqrt(Imgabout.*Imgabout + Regabout.*Regabout));  

调用代码：

[csharp]  view plain copy

1. ori=imread('C:\Users\watkins\Pictures\cartoon.jpg');  
2. grayimg=rgb2gray(ori);  
3. gim=im2double(grayimg);   
4.   
5. Sx=32;  
6. Sy=32;  
7. f=sqrt(8);  
8. theta=pi/2;  
9. u=4;  
10. v=4;  
11. %[G,gabout] = gaborfilter1(gim,Sx,Sy,f,theta);  
12. [G,gabout] = gaborfilter1(gim,Sx,Sy,u,v);  
13.   
14. imshow(real(G));  
15. %imshow(real(gabout));  

 

滤波器图片：

 

3. 方式三

 

[csharp]  view plain copy

1. %%%%%%%VERSION 3  
2. %%ANOTHER DESCRIBTION OF GABOR FILTER  
3.   
4. %The Gabor filter is basically a Gaussian (with variances sx and sy along x and y-axes respectively)  
5. %modulated by a complex sinusoid (with centre frequencies U and V along x and y-axes respectively)   
6. %described by the following equation  
7. %%  
8. %               1                -1     x  ^    y  ^  
9. %%% Gi(x,y) = ---------- * exp ([----{(----) 2+(----) 2}])*Mi(x,y,f);   
10. %            2*pi*sx*sy           2    sx       sy  
11. %%% i =1,2  
12. %%% M1(x,y,f) = cos[2*pi*f*sqrt(x^2+y^2)];  
13. %%% M2(x,y,f) = cos[2*pi*f*(x*cos(theta) + y*sin(theta)];  
14.   
15. %% Describtion :  
16.   
17. %% I : Input image  
18. %% Sx & Sy : Variances along x and y-axes respectively  
19. %% f : The frequency of the sinusoidal function  
20. %% theta : The orientation of Gabor filter  
21.   
22. %% G1 & G2 : The output filters as described above  
23. %% gabout1 & gabout2 : The output filtered images  
24.   
25. %%  Author : Ahmad poursaberi  e-mail : a.poursaberi@ece.ut.ac.ir  
26. %%          Faulty of Engineering, Electrical&Computer Department,Tehran  
27. %%          University,Iran,June 2004  
28.   
29. function [G1,G2,gabout1,gabout2] = gaborfilter2(I,Sx,Sy,f,theta)  
30.   
31. if isa(I,'double')~=1   
32.     I = double(I);  
33. end  
34.   
35. for x = -fix(Sx):fix(Sx)  
36.     for y = -fix(Sy):fix(Sy)  
37.         M1 = cos(2*pi*f*sqrt(x^2+y^2));  
38.         M2 = cos(2*pi*f*(x*cos(theta)+y*sin(theta)));  
39.         G1(fix(Sx)+x+1,fix(Sy)+y+1) = (1/(2*pi*Sx*Sy)) * exp(-.5*((x/Sx)^2+(y/Sy)^2))*M1;  
40.         G2(fix(Sx)+x+1,fix(Sy)+y+1) = (1/(2*pi*Sx*Sy)) * exp(-.5*((x/Sx)^2+(y/Sy)^2))*M2;  
41.     end  
42. end  
43.   
44. Imgabout1 = conv2(I,double(imag(G1)),'same');  
45. Regabout1 = conv2(I,double(real(G1)),'same');  
46.   
47. Imgabout2 = conv2(I,double(imag(G2)),'same');  
48. Regabout2 = conv2(I,double(real(G2)),'same');  
49.   
50. gabout1 = sqrt(Imgabout1.*Imgabout1 + Regabout1.*Regabout1);  
51. gabout2 = sqrt(Imgabout2.*Imgabout2 + Regabout2.*Regabout2);  

 

调用代码：

[csharp]  view plain copy

1. ori=imread('C:\Users\watkins\Pictures\cartoon.jpg');  
2. grayimg=rgb2gray(ori);  
3. gim=im2double(grayimg);   
4.   
5. Sx=16;  
6. Sy=16;  
7. f=sqrt(2);  
8. theta=pi/2;  
9. u=8;  
10. v=0;  
11. %[G,gabout] = gaborfilter1(gim,Sx,Sy,f,theta);  
12. %[G,gabout] = gaborfilter1(gim,Sx,Sy,u,v);  
13. [G1,G2,gabout1,gabout2] = gaborfilter2(gim,Sx,Sy,f,theta);  
14.   
15. R=real(G2);  
16. k=127.5/max(max(abs(R)));  
17. imshow(uint8(k*R+127.5));  
18. %imshow(real(G2));  
19. %imshow(abs(real(gabout)));  

 

生成的滤波器图片：