一、简介
我们知道,在一幅图像中,其低频成分对应者图像变化缓慢的部分,对应着图像大致的相貌和轮廓。而其高频成分则对应着图像变化剧烈的部分,对应着图像的细节(图像的噪声也属于高频成分)。
低频滤波器,顾名思义,就是过滤掉或者大幅度衰减图像的高频成分,让图像的低频成分通过。低频滤波器可以平滑图像,虑去图像的噪声。而与此相反的高频滤波器,则是过滤低频成分,通过高频成分,可以达到锐化图像的目的。
理想低通滤波器的滤波非常尖锐,而高斯低通滤波器的滤波则非常平滑。Butterworth低通滤波器则介于两者之间,当Butterworth低通滤波器的阶数较高时,接近于理想低通滤波器,阶数较低时,则接近于高斯低通滤波器。
下面我们所说的都是对应于二维图像处理的情况。
二、这三种滤波器的相关介绍
在这三种低通滤波器的表达式中,我们都用D0来表示其截止频率。D(u, v)表示距离频率矩形中心的距离。
我并没有说明振铃现象产生的原因,只是说明了什么情况下振铃现象较为明显。
1. 理想低通滤波器(ILPF)
理想低通滤波器在以原点为圆心、D0为半径的园内,通过所有的频率,而在圆外截断所有的频率。(圆心的频率最低,为变换的直流(dc)分量)。函数如下:
可以看出,理想低通滤波器的过渡非常急剧,会产生振铃现象。
2. Butterworth低通滤波器
函数表达式如下(在有些书中,Butterworth的函数的平方才等于右边的表达式,在这里我们按照课本的写法,计算较为方便),其中n称为Butterworth低通滤波器的阶数:
从滤波器的函数图中,我们可以看出其过渡没有理想低通滤波器那么剧烈,从图(c)中可以看出,阶数越高,滤波器的过度越剧烈,振铃现象将越明显。
3. 高斯低通滤波器(GLPF)
函数表达式如下:
高斯滤波器的过度特性非常平坦,因此不会产生振铃现象。
三、频域滤波步骤
1. 基本步骤
1. 给定一幅大小为m*n的图像f(x,y)。选择适当的填充参数P和Q,一般令P = 2m,Q = 2n。
2. 对图像f(x, y)填充0,填充后得到图像大小为P*Q的图像fp(x, y)。
3. 用(-1)^(x+y)乘以fp(x,y)将其移到变换中心(中心化)。
4. 计算fp(x, y)的DFT,得到F(u,v)。
5. 生成一个实的,对称的滤波函数H(u, v),大小为P*Q,中心在(P/2, Q/2)处。然后相乘(矩阵点乘)得到G(u,v) = H(u,v)F(u,v)。
6. 对G(u, v)反傅里叶变换,然后取实部,再乘以(-1)^(x+y)进行反中心变换最后得到gp(x,y)。
7. 提取gp(x,y)左上角的m*n区域,对提取的部分进行标准化处理,得到最终的结果图像g(x,y)。
2. 相关步骤说明
1. 对图像进行0填充,得到大小为P*Q的图像,主要是为了避免在循环卷积中出现的缠绕错误。当两个矩阵大小相同时,P≥2m-1,Q≥2n-1时可以避免环绕错误。由于对于偶数尺寸的矩阵计算其傅里叶变换较快,因此P取2m,Q取2n。
2. 图像乘以(-1)^(x+y),再对图像进行傅里叶变换可以得到将原点移到中心的傅里叶变换。这样,对于F(u,v)来说,中心的频率最低,四周的频率较高。每一点的值表示该频率对于的幅度。在matlab中,也可以不乘以(-1)^(x+y),直接对填充后的图像进行傅里叶变换,之后使用fftshift函数对其傅里叶变换进行中心化。得到的结果是一样的。
3. 由于图像是一个实函数,所以其傅里叶变换是一个旋转对称的傅里叶变换。
4. 关于对称中心。课本中指出,对于一个长度为M的一维序列,当M为偶数时,位置0和M/2呈现零的特性,当M为奇数时,只有位置0呈现零的特性。因此,在matlab中,由于矩阵下标是从1开始的,我们的P = 2m,Q = 2n。因此中心点的位置为(P/2+1, Q/2+1),即(m+1, n+1)。
5. 标准化处理。对最终得到的图像每一点的值进行处理。使其范围变为[0, 255]的uint8。或者[0, 1.0]的double值。
四、matlab代码实现Butterworth低通滤波器
在使用matlab代码的实现过程中,对于这三种低通滤波器,只是在实验低通滤波器函数H(u,v)的代码中有部分不同,其他部分一致。因此,在下面中,只给出实现Butterworth低通滤波器的代码,不给出其他两种滤波器的代码。
1. Butterworth滤波器的代码如下:
该函数为Bfilter,输入为需要进行Butterworth滤波的灰度图像,Butterworth滤波器的截止频率D0以及Butterworth滤波器的阶数n。输出为进行滤波之后的图像(图像的值已经归一化到[[0, 255])。
-
function [image_out] = Bfilter(image_in, D0, N)
-
%
Butterworth
滤波器,在频率域进行滤波
-
% 输入为需要进行滤波的灰度图像,
Butterworth
滤波器的截止频率
D0
,阶数
N
-
% 输出为滤波之后的灰度图像
-
-
[
m
,
n
] =
size
(image_in)
;
-
P =
2 * m;
-
Q =
2 * n;
-
-
fp = zeros(P, Q);
-
%对图像填充
0,并且乘以(-
1)^(x+y) 以移到变换中心
-
for i =
1 : m
-
for j =
1 : n
-
fp(i, j) = double(image_in(i, j)) * (-
1)^(i+j);
-
end
-
end
-
% 对填充后的图像进行傅里叶变换
-
F1 = fft2(fp);
-
-
% 生成Butterworth滤波函数,中心在(m+
1,n+
1)
-
Bw = zeros(P, Q);
-
a = D0^(
2 * N);
-
for u =
1 : P
-
for v =
1 : Q
-
temp = (u-(m+
1.0))^
2 + (v-(n+
1.0))^
2;
-
Bw(u, v) =
1 / (
1 + (temp^N) / a);
-
end
-
end
-
-
%进行滤波
-
G = F1 .* Bw;
-
-
% 反傅里叶变换
-
gp = ifft2(G);
-
-
% 处理得到的图像
-
image_out = zeros(m, n,
’uint8’);
-
gp = real(gp);
-
g = zeros(m, n);
-
for i =
1 : m
-
for j =
1 : n
-
g(i, j) = gp(i, j) * (-
1)^(i+j);
-
-
end
-
end
-
mmax = max(g(:));
-
mmin = min(g(:));
-
range = mmax-mmin;
-
for i =
1 : m
-
for j =
1 : n
-
image_out(i,j) = uint8(
255 * (g(i, j)-mmin) / range);
-
end
-
end
-
-
end
2. 测试代码如下
测试时Butterworth滤波器的阶数为2。截止频率分别为10,30,60,160,460。
-
clear all;
-
close all;
-
clc;
-
-
image1 = imread(
'3.bmp');
-
-
image2 = Bfilter(image1,
10,
2);
-
image3 = Bfilter(image1,
30,
2);
-
image4 = Bfilter(image1,
60,
2);
-
image5 = Bfilter(image1,
160,
2);
-
image6 = Bfilter(image1,
460,
2);
-
-
% 显示图像
-
subplot(
2,
3,
1), imshow(image1), title(
'原图像');
-
subplot(
2,
3,
2), imshow(image2), title(
'D0 = 10, n = 2');
-
subplot(
2,
3,
3), imshow(image3), title(
'D0 = 30, n = 2');
-
subplot(
2,
3,
4), imshow(image4), title(
'D0 = 60, n = 2');
-
subplot(
2,
3,
5), imshow(image5), title(
'D0 = 160, n = 2');
-
subplot(
2,
3,
6), imshow(image6), title(
'D0 = 460, n = 2');