本节书摘来自异步社区《精通Matlab数字图像处理与识别》一书中的第6章,第6.5节,作者 张铮 , 倪红霞 , 苑春苗 , 杨立红,更多章节内容可以访问云栖社区“异步社区”公众号查看
6.5 频域低通滤波器
精通Matlab数字图像处理与识别
在频谱中,低频主要对应图像在平滑区域的总体灰度级分布,而高频对应图像的细节部分,如边缘和噪声。因此,图像平滑可以通过衰减图像频谱中的高频部分来实现,这就建立了空间域图像平滑与频域低通滤波之间对应关系。
6.5.1 理想低通滤波器及其实现
1.理论基础
最容易想到的衰减高频成份的方法就是在一个称为“截止频率”的位置“截断”所有的高频成份,将图像频谱中所有高于这一截止频率的频谱成份设为0,低于截止频率的成份设为保持不变。能够达到这种效果的滤波器如图6.16所示,我们称之为理想低通滤波器。如果图像的宽度为M,高度为N,那么理想低通频域滤波器可形式化地描述为
其中D 0表示理想低通滤波器的截止频率,滤波器的频率域原点在频谱图像的中心处,在以截止频率为半径的圆形区域之内的滤镜元素值全部为1,而该圆之外的滤镜元素值全部为0。理想低通滤波器的频率特性在截止频率处十分陡峭,无法用硬件实现,这也是我们称之为理想的原因,但其软件编程的模拟实现较为简单。
理想低通滤波器可在一定程度上去除图像噪声,但由此带来的图像边缘和细节的模糊效应也较为明显,其滤波之后的处理效果比较类似于5.3.1中的平均平滑。实际上,理想低通滤波器是一个与频谱图像同样尺寸的二维矩阵,通过将矩阵中对应较高频率的部分设为0,较低频率的部分(靠近中心)设为1,可在与频谱图像相乘后有效去除频谱的高频成份(由于是矩阵对应元素相乘,频谱高频成份与滤波器中的0相乘)。其中0与1的交界处即对应滤波器的截止频率。
2.Matlab实现
利用我们编写的imidealflpf函数可以得到截止频率为freq的理想低通滤波器。
imidealflpf函数的完整实现如下,该函数被封装在金羽图书论坛(http://bbs.book95.com)的“金羽图书与答疑”板块与本书同名的主题帖子附件的“chapter6/ code”目录下的imidealflpf.m文件中。
function out = imidealflpf(I, freq)
% imidealflpf函数 构造理想的频域低通滤波器
% I参数 输入的灰度图像
% freq参数 低通滤波器的截止频率
% 返回值:out – 指定的理想低通滤波器
[M,N] = size(I);
out = ones(M,N);
for i=1:M
for j=1:N
if (sqrt(((i-M/2)^2+(j-N/2)^2))>freq)
out(i,j)=0;
end
end
end
下面我们仍以第5章中图像平滑中曾使用过的婴儿老照片baby_noise.bmp为例,使用频域理想低通滤波器进行处理,相应的Matlab代码如例6.3所示。
[例6.3]理想低通滤波实现。
I = imread('../baby_noise.bmp'); %读入原图像
% 生成滤镜
ff = imidealflpf(I, 20);
% 应用滤镜
out = imfreqfilt(I, ff);
figure (1);
subplot(2,2,1);
imshow(I);
title('Source');
% 计算FFT并显示
temp = fft2(double(I));
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,1);
imshow(temp, []);
title('Source');
figure (1);
subplot(2,2,2);
imshow(out);
title('Ideal LPF, freq=20');
% 计算FFT并显示
temp = fft2(out);
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,2);
imshow(temp, []);
title(' Ideal LPF, freq=20');
% 生成滤镜
ff = imidealflpf(I, 40);
% 应用滤镜
out = imfreqfilt(I, ff);
figure (1);
subplot(2,2,3);
imshow(out);
title('Ideal LPF, freq=40');
% 计算FFT并显示
temp = fft2(out);
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,3);
imshow(temp, []);
title(' Ideal LPF, freq=40');
% 生成滤镜
ff = imidealflpf(I, 60);
% 应用滤镜
out = imfreqfilt(I, ff);
figure (1);
subplot(2,2,4);
imshow(out);
title('Ideal LPF, freq=60');
% 计算FFT并显示
temp = fft2(out);
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,4);
imshow(temp, []);
title(' Ideal LPF, freq=60');
上述程序的运行效果如图6.17和图6.18所示。
从图6.17和图6.18可见,当截止频率非常低时,只有非常靠近原点的低频成份能够通过,图像模糊严重;截止频率越高,通过的频率成份就越多,图像模糊的程度越小,所获得的图像也就越接近原图像。但可以看出,理想低通滤波器并不能很好地兼顾噪声滤除与细节保留两个方面,这和空域中采用平均模板时的情形比较类似。下面,我们将介绍频域的高斯低通滤波器,并比较它与理想低通滤波器的处理效果。
6.5.2 高斯低通滤波器及其实现
1.理论基础
高斯低通滤波器的频率域二维形式由下式给出。
高斯函数具有相对简单的形式,而且它的傅立叶变换和傅立叶反变换都是实高斯函数。为了简单起见,下面仅给出一个一维高斯函数的傅立叶变换和傅立叶反变换作为例子,式(6-52)告诉我们一维高斯函数的傅立叶反变换和正变换仍为高斯函数,该式的证明留给有兴趣的读者自己完成(提示:可以利用高斯分布的概率密度函数在定义域上积分为1的性质)。
式中:σ是高斯曲线的标准差。
频域和与之对应的空域一维高斯函数的图形如图6.19所示。
从图6.19可以发现,当σ增大时,H(u)的图像倾向于变宽,而h(x)的图像倾向于变窄和变高。这也体现了频率域和空间域的对应关系。频率域滤波器越窄,滤除的高频成份越多,图像就越平滑(模糊);而在空间域,对应的滤波器就越宽,相应的卷积模板越平坦,平滑(模糊)效果就越明显。
我们在图6.20中分别给出σ取值为3和1时的频域二维高斯滤波器的三维曲面表示,可以看出频域下的二维高斯滤波器同样具有上述一维情况时的特点。
2.Matlab实现
根据上面二维高斯低通滤波器的定义,可以编写高斯低通滤波器的生成函数如下,该函数被封装在金羽图书论坛(http://bbs.book95.com)的“金羽图书与答疑”板块与本书同名的主题帖子附件的“chapter6/code”目录下的imgaussflpf.m文件中。
f`
javascript
unction out = imgaussflpf(I, sigma)
% imgaussflpf函数 构造频域高斯低通滤波器
% I参数 输入的灰度图像
% sigma参数 高斯函数的Sigma参数
[M,N] = size(I);
out = ones(M,N);
for i=1:M
for j=1:N
out(i,j) = exp(-((i-M/2)^2+(j-N/2)^2)/2/sigma^2);
end
end
下面的例6.4给出了针对图像baby_noise.bmp,sigma不同取值时高斯低通滤波的Matlab程序。
[例6.4]高斯低通滤波实现。
I = imread('../baby_noise.bmp'); %读入原图像
% 生成滤镜
ff = imgaussflpf (I, 20);
% 应用滤镜
out = imfreqfilt(I, ff);
figure (1);
subplot(2,2,1);
imshow(I);
title('Source');
% 计算FFT并显示
temp = fft2(double(I));
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,1);
imshow(temp, []);
title('Source');
figure (1);
subplot(2,2,2);
imshow(out);
title('Gauss LPF, sigma=20');
% 计算FFT并显示
temp = fft2(out);
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,2);
imshow(temp, []);
title(' Gauss LPF, sigma=20');
% 生成滤镜
ff = imgaussflpf (I, 40);
% 应用滤镜
out = imfreqfilt(I, ff);
figure (1);
subplot(2,2,3);
imshow(out);
title('Gauss LPF, sigma =40');
% 计算FFT并显示
temp = fft2(out);
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,3);
imshow(temp, []);
title(' Gauss LPF, sigma =40');
% 生成滤镜
ff = imgaussflpf (I, 60);
% 应用滤镜
out = imfreqfilt(I, ff);
figure (1);
subplot(2,2,4);
imshow(out);
title('Gauss LPF, sigma =60');
% 计算FFT并显示
temp = fft2(out);
temp = fftshift(temp);
temp = log(1 + abs(temp));
figure (2);
subplot(2,2,4);
imshow(temp, []);
title(' Gauss LPF, sigma =60');
上述程序的运行后得到的滤波效果如图6.21所示。
![image](https://yqfile.alicdn.com/a9854718f6fa46b35dfba9024094d38f68761026.png)
图6.21中各幅图像所对应的频域图像如图6.22所示。显然,高斯滤波器的截止频率处不是陡峭的。
![image](https://yqfile.alicdn.com/fff89c1086a9a87c7ced978576b2dbcc9ad0fe75.png)
高斯低通滤波器在Sigma参数取40的时候可以较好地处理被高斯噪声污染的图像,而相比于理想低通滤波器而言,处理效果上的改进是显而易见的。高斯低通滤波器在有效抑制噪声的同时,图像的模糊程度更低,对边缘带来的混叠程度更小,从而使高斯低通滤波器在通常情况下获得了比理想低通滤波器更为广泛的应用。