频域滤波(1)--基础与低通滤波器

                          
               

       之前的博文主要介绍了空间域内的滤波器,本文主要从频域的角度进行分析。主要使用傅里叶变换,将空间域的图像转换到频域内,在频域内进行数字图像处理。这部分的内容及其重要,频域内的处理可以解决空间域内无法完成的图像增强。本文首先从数学角度,对图像的频域内的性质进行分析,然后在着重介绍滤波器在频域内的性质。

        1.傅里叶变换与频域

        在之前的文中,我们已经进行过一些基本的图像处理。比如,使用低通滤波可以将图像模糊,也有些许降噪的作用。这些都是在空间域内进行的滤波处理,这个处理主要是依靠卷积来进行计算的。首先,从连续的一维卷积入手,如下所示。


       将上式进行傅里叶变换,可以得到如下结果。


        从这个式子,我们可以得到一个重要的结论。也就是,函数卷积的傅里叶变换所得到的结果,是函数的傅里叶变换的乘积。再将其总结得简单易懂一些,有如下结论。


        在将其扩展到二维的形况下,假设尺寸为MxN的图像,如下关系是成立的。


       其实到这,基本的原理就明了的。我们所看到的图像,均为空间域内的表现形式,我们无法辨识出频域内的图像。要进行频域内的滤波器处理,首先就需要进行傅里叶变换,然后直接进行滤波处理,最后再用反傅里叶变换倒回到空间域内。

       到此,已经可以开始空间域内的滤波处理了。但是,还有一点需要注意的地方。使用某个一维信号来举例子,一维信号的傅里叶变换是以2π为周期的函数。所以,我们常常使用的范围[-π,π]来表示这个信号的傅里叶变换,如下所示。


        这样做的好处是,靠近0的成分就是低频,靠近-π与π的成分就表示高频。而对于图像而言,在Matlab中,我们使用fft2()这个函数来求取图像的傅里叶变换。

g = fft2(f);   
        上面这个代码求取的,其实是 范围[0,π]内的傅里叶变换。为了方便理解,下图画出了本行代码所求取的图像的傅里叶变换的范围(右)和与其等效的一维傅里叶变换的范围(左)。


       很显然,这并不是希望的范围,下面这个代码可以求取[0,2π]内的傅里叶变换。

P = 2*M;
Q = 2*N;
F = fft2(f,P,Q);
       下图画出了本行代码所求取的图像的傅里叶变换的范围(右)和与其等效的一维傅里叶变换的范围 (左) 。所得到的图像F(u,v)的尺寸为PxQ。


         我们需要对其移动一下,如下图所示,我们需要的是粉色范围的区域。


下面,从数学上分析一下,如何获得这个部分的频谱。对于傅里叶变换,有如下性质。


这个特性称为平移特性,粉色部分的频谱,将带入上式,我们可以得到如下式子。


为次,我们已经得到了粉色范围的频谱。越靠近傅里叶频谱图像中间的成分,代表了低频成分。其Matlab代码如下所示。


   
   
  1. [M,N] = size(f);
  2. P = 2*M;
  3. Q = 2*N;
  4. fc = zeros(M,N);
  5. for x = 1 : 1 :M
  6.     for y = 1 : 1 :N
  7.         fc(x,y) = f(x,y) * (- 1)^(x+y);
  8.     end
  9. end
  10. F = fft2(fc,P,Q);

        代码所得到的结果,如下图所示。


        接下来,我们总结一下频域滤波的步骤:

        ①:先将图像做频域内的水平移动,然后求原图像f(x,y)的DFT,得到其图像的傅里叶谱F(u,v)。


        ②:与频域滤波器做乘积,

        ③:求取G(u,v)的IDFT,然后再将图像做频域内的水平移动(移动回去),其结果可能存在寄生的虚数,此时忽略即可。


        ④:这里使用ifft2函数进行IDFT变换,得到的图像的尺寸为PxQ。切取左上角的MxN的图像,就能得到结果了。

        2.低通滤波器

        2.1理想的低通滤波器


       其中,D0表示通带的半径。D(u,v)的计算方式也就是两点间的距离,很简单就能得到。

       使用低通滤波器所得到的结果如下所示。低通滤波器滤除了高频成分,所以使得图像模糊。由于理想低通滤波器的过度特性过于急峻,所以会产生了振铃现象。


         

        2.2巴特沃斯低通滤波器

       同样的,D0表示通带的半径,n表示的是巴特沃斯滤波器的次数。随着次数的增加,振铃现象会越来越明显。

   


       2.3高斯低通滤波器


       D0表示通带的半径。高斯滤波器的过度特性非常平坦,因此是不会产生振铃现象的。

       3.实现代码


    
    
  1. close all;
  2. clear all;
  3. %% ---------Ideal Lowpass Filters (Fre. Domain)------------
  4. f = imread( 'characters_test_pattern.tif');
  5. f = mat2gray(f,[ 0 255]);
  6. [M,N] = size(f);
  7. P = 2*M;
  8. Q = 2*N;
  9. fc = zeros(M,N);
  10. for x = 1: 1:M
  11.     for y = 1: 1:N
  12.         fc(x,y) = f(x,y) * (- 1)^(x+y);
  13.     end
  14. end
  15. F = fft2(fc,P,Q);
  16. H_1 = zeros(P,Q);
  17. H_2 = zeros(P,Q);
  18. for x = (-P/ 2): 1:(P/ 2)- 1
  19.      for y = (-Q/ 2): 1:(Q/ 2)- 1
  20.         D = (x^ 2 + y^ 2)^( 0.5);
  21.         if(D <= 60)  H_1(x+(P/ 2)+ 1,y+(Q/ 2)+ 1) = 1; end   
  22.         if(D <= 160)  H_2(x+(P/ 2)+ 1,y+(Q/ 2)+ 1) = 1; end
  23.      end
  24. end
  25. G_1 = H_1 .* F;
  26. G_2 = H_2 .* F;
  27. g_1 = real(ifft2(G_1));
  28. g_1 = g_1( 1: 1:M, 1: 1:N);
  29. g_2 = real(ifft2(G_2));
  30. g_2 = g_2( 1: 1:M, 1: 1:N);        
  31. for x = 1: 1:M
  32.     for y = 1: 1:N
  33.         g_1(x,y) = g_1(x,y) * (- 1)^(x+y);
  34.         g_2(x,y) = g_2(x,y) * (- 1)^(x+y);
  35.     end
  36. end
  37. %% -----show-------
  38. figure();
  39. subplot( 1, 2, 1);
  40. imshow(f,[ 0 1]);
  41. xlabel( 'a).Original Image');
  42. subplot( 1, 2, 2);
  43. imshow(log( 1 + abs(F)),[ ]);
  44. xlabel( 'b).Fourier spectrum of a');
  45. figure();
  46. subplot( 1, 2, 1);
  47. imshow(H_1,[ 0 1]);
  48. xlabel( 'c).Ideal Lowpass filter(D=60)');
  49. subplot( 1, 2, 2);
  50. h = mesh( 1: 20:P, 1: 20:Q,H_1( 1: 20:P, 1: 20:Q));
  51. set(h, 'EdgeColor', 'k');
  52. axis([ 0 P 0 Q 0 1]);
  53. xlabel( 'u');ylabel( 'v');
  54. zlabel( '|H(u,v)|');
  55. figure();
  56. subplot( 1, 2, 1);
  57. imshow(log( 1 + abs(G_1)),[ ]);
  58. xlabel( 'd).Result of filtering using c');
  59. subplot( 1, 2, 2);
  60. imshow(g_1,[ 0 1]);
  61. xlabel( 'e).Result image');
  62. figure();
  63. subplot( 1, 2, 1);
  64. imshow(H_2,[ 0 1]);
  65. xlabel( 'f).Ideal Lowpass filter(D=160)');
  66. subplot( 1, 2, 2);
  67. h = mesh( 1: 20:P, 1: 20:Q,H_2( 1: 20:P, 1: 20:Q));
  68. set(h, 'EdgeColor', 'k');
  69. axis([ 0 P 0 Q 0 1]);
  70. xlabel( 'u');ylabel( 'v');
  71. zlabel( '|H(u,v)|');
  72. figure();
  73. subplot( 1, 2, 1);
  74. imshow(log( 1 + abs(G_2)),[ ]);
  75. xlabel( 'g).Result of filtering using e');
  76. subplot( 1, 2, 2);
  77. imshow(g_2,[ 0 1]);
  78. xlabel( 'h).Result image');
  79. close all;
  80. clear all;
  81. %% ---------Butterworth Lowpass Filters (Fre. Domain)------------
  82. f = imread( 'characters_test_pattern.tif');
  83. f = mat2gray(f,[ 0 255]);
  84. [M,N] = size(f);
  85. P = 2*M;
  86. Q = 2*N;
  87. fc = zeros(M,N);
  88. for x = 1: 1:M
  89.     for y = 1: 1:N
  90.         fc(x,y) = f(x,y) * (- 1)^(x+y);
  91.     end
  92. end
  93. F = fft2(fc,P,Q);
  94. H_1 = zeros(P,Q);
  95. H_2 = zeros(P,Q);
  96. for x = (-P/ 2): 1:(P/ 2)- 1
  97.      for y = (-Q/ 2): 1:(Q/ 2)- 1
  98.         D = (x^ 2 + y^ 2)^( 0.5);
  99.         D_0 = 100;
  100.         H_1(x+(P/ 2)+ 1,y+(Q/ 2)+ 1) = 1/( 1+(D/D_0)^ 2);  
  101.         H_2(x+(P/ 2)+ 1,y+(Q/ 2)+ 1) = 1/( 1+(D/D_0)^ 6);
  102.      end
  103. end
  104. G_1 = H_1 .* F;
  105. G_2 = H_2 .* F;
  106. g_1 = real(ifft2(G_1));
  107. g_1 = g_1( 1: 1:M, 1: 1:N);
  108. g_2 = real(ifft2(G_2));
  109. g_2 = g_2( 1: 1:M, 1: 1:N);        
  110. for x = 1: 1:M
  111.     for y = 1: 1:N
  112.         g_1(x,y) = g_1(x,y) * (- 1)^(x+y);
  113.         g_2(x,y) = g_2(x,y) * (- 1)^(x+y);
  114.     end
  115. end
  116. %% -----show-------
  117. figure();
  118. subplot( 1, 2, 1);
  119. imshow(f,[ 0 1]);
  120. xlabel( 'a).Original Image');
  121. subplot( 1, 2, 2);
  122. imshow(log( 1 + abs(F)),[ ]);
  123. xlabel( 'b).Fourier spectrum of a');
  124. figure();
  125. subplot( 1, 2, 1);
  126. imshow(H_1,[ 0 1]);
  127. xlabel( 'c)Butterworth Lowpass (D_{0}=100,n=1)');
  128. subplot( 1, 2, 2);
  129. h = mesh( 1: 20:P, 1: 20:Q,H_1( 1: 20:P, 1: 20:Q));
  130. set(h, 'EdgeColor', 'k');
  131. axis([ 0 P 0 Q 0 1]);
  132. xlabel( 'u');ylabel( 'v');
  133. zlabel( '|H(u,v)|');
  134. figure();
  135. subplot( 1, 2, 1);
  136. imshow(log( 1 + abs(G_1)),[ ]);
  137. xlabel( 'd).Result of filtering using c');
  138. subplot( 1, 2, 2);
  139. imshow(g_1,[ 0 1]);
  140. xlabel( 'e).Result image');
  141. figure();
  142. subplot( 1, 2, 1);
  143. imshow(H_2,[ 0 1]);
  144. xlabel( 'f).Butterworth Lowpass (D_{0}=100,n=3)');
  145. subplot( 1, 2, 2);
  146. h = mesh( 1: 20:P, 1: 20:Q,H_2( 1: 20:P, 1: 20:Q));
  147. set(h, 'EdgeColor', 'k');
  148. axis([ 0 P 0 Q 0 1]);
  149. xlabel( 'u');ylabel( 'v');
  150. zlabel( '|H(u,v)|');
  151. figure();
  152. subplot( 1, 2, 1);
  153. imshow(log( 1 + abs(G_2)),[ ]);
  154. xlabel( 'g).Result of filtering using e');
  155. subplot( 1, 2, 2);
  156. imshow(g_2,[ 0 1]);
  157. xlabel( 'h).Result image');
  158. close all;
  159. clear all;
  160. clc;
  161. %% ---------Gaussian Lowpass Filters (Fre. Domain)------------
  162. f = imread( 'characters_test_pattern.tif');
  163. f = mat2gray(f,[ 0 255]);
  164. [M,N] = size(f);
  165. P = 2*M;
  166. Q = 2*N;
  167. fc = zeros(M,N);
  168. for x = 1: 1:M
  169.     for y = 1: 1:N
  170.         fc(x,y) = f(x,y) * (- 1)^(x+y);
  171.     end
  172. end
  173. F = fft2(fc,P,Q);
  174. H_1 = zeros(P,Q);
  175. H_2 = zeros(P,Q);
  176. for x = (-P/ 2): 1:(P/ 2)- 1
  177.      for y = (-Q/ 2): 1:(Q/ 2)- 1
  178.         D = (x^ 2 + y^ 2)^( 0.5);
  179.         D_0 = 60;
  180.         H_1(x+(P/ 2)+ 1,y+(Q/ 2)+ 1) = exp(-(D*D)/( 2*D_0*D_0));  
  181.         D_0 = 160;
  182.         H_2(x+(P/ 2)+ 1,y+(Q/ 2)+ 1) = exp(-(D*D)/( 2*D_0*D_0));
  183.      end
  184. end
  185. G_1 = H_1 .* F;
  186. G_2 = H_2 .* F;
  187. g_1 = real(ifft2(G_1));
  188. g_1 = g_1( 1: 1:M, 1: 1:N);
  189. g_2 = real(ifft2(G_2));
  190. g_2 = g_2( 1: 1:M, 1: 1:N);        
  191. for x = 1: 1:M
  192.     for y = 1: 1:N
  193.         g_1(x,y) = g_1(x,y) * (- 1)^(x+y);
  194.         g_2(x,y) = g_2(x,y) * (- 1)^(x+y);
  195.     end
  196. end
  197. %% -----show-------
  198. close all;
  199. figure();
  200. subplot( 1, 2, 1);
  201. imshow(f,[ 0 1]);
  202. xlabel( 'a).Original Image');
  203. subplot( 1, 2, 2);
  204. imshow(log( 1 + abs(F)),[ ]);
  205. xlabel( 'b).Fourier spectrum of a');
  206. figure();
  207. subplot( 1, 2, 1);
  208. imshow(H_1,[ 0 1]);
  209. xlabel( 'c)Gaussian Lowpass (D_{0}=60)');
  210. subplot( 1, 2, 2);
  211. h = mesh( 1: 20:P, 1: 20:Q,H_1( 1: 20:P, 1: 20:Q));
  212. set(h, 'EdgeColor', 'k');
  213. axis([ 0 P 0 Q 0 1]);
  214. xlabel( 'u');ylabel( 'v');
  215. zlabel( '|H(u,v)|');
  216. figure();
  217. subplot( 1, 2, 1);
  218. imshow(log( 1 + abs(G_1)),[ ]);
  219. xlabel( 'd).Result of filtering using c');
  220. subplot( 1, 2, 2);
  221. imshow(g_1,[ 0 1]);
  222. xlabel( 'e).Result image');
  223. figure();
  224. subplot( 1, 2, 1);
  225. imshow(H_2,[ 0 1]);
  226. xlabel( 'f).Gaussian Lowpass (D_{0}=160)');
  227. subplot( 1, 2, 2);
  228. h = mesh( 1: 20:P, 1: 20:Q,H_2( 1: 20:P, 1: 20:Q));
  229. set(h, 'EdgeColor', 'k');
  230. axis([ 0 P 0 Q 0 1]);
  231. xlabel( 'u');ylabel( 'v');
  232. zlabel( '|H(u,v)|');
  233. figure();
  234. subplot( 1, 2, 1);
  235. imshow(log( 1 + abs(G_2)),[ ]);
  236. xlabel( 'g).Result of filtering using e');
  237. subplot( 1, 2, 2);
  238. imshow(g_2,[ 0 1]);
  239. xlabel( 'h).Result image');

      原文发于博客:http://blog.csdn.net/thnh169/ 

           
               
  • 0
    点赞
  • 3
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值