数字图像处理 第4章-频率域滤波
4.1 一些基础知识
4.1.1 傅里叶级数
可以参考高等数学知识,或者考研数学一的内容。大致了解一下。
4.1.2 复数
复数的主要定义是:
C
=
R
+
j
I
C=R+jI
C=R+jI
其中R和I是实数,R表示复数的实部,I是复数的虚部。实数是I=0的复数的子集。
i在我们中学数学中有提到过,并且
i
=
−
1
i=\sqrt{-1}
i=−1
一个复数的C的共轭表示为C*,其定义为
C
∗
=
R
−
j
I
C^*=R-jI
C∗=R−jI
4.1.3 冲激及其取样特性
线性系统和傅里叶变换研究的核心是冲激及其取样特性。连续变量t在t=0处的单位冲激表示为δ(t),其定义是
δ
(
t
)
=
{
∞
t
=
0
0
t
≠
0
δ(t)=\begin{cases} ∞ & t=0 \\ 0 & t≠0 \end{cases}
δ(t)={∞0t=0t=0
物理上,如果我们把,解释为时间,那么一个冲激可看成是幅度无限、持续时间为0、具有单位面积的尖峰信号。一个冲激具有关于如下积分的所谓取样(sifting)特性:
∫
−
∞
∞
f
(
t
)
δ
(
t
)
d
t
=
f
(
0
)
∫^∞_{-∞}f(t)δ(t)dt=f(0)
∫−∞∞f(t)δ(t)dt=f(0)
假设f(t)在t=0处是连续的,这是实际中满足的一个典型条件。取样特性得到函数f(t)在冲激位置(也就是,在前面的公式中为原点t=0)的值。取样特性的一种更为一般的说明涉及位于任意点t0的冲激,表示为δ(t-t0)。在这种情况下,取样特性变为:
∫
−
∞
∞
f
(
t
)
δ
(
t
−
t
0
)
d
t
=
f
(
t
0
)
∫^∞_{-∞}f(t)δ(t-t_0)dt=f(t_0)
∫−∞∞f(t)δ(t−t0)dt=f(t0)
冲激串定义为无限多个分离的周期冲激单元ΔT之和:
s
Δ
T
(
t
)
=
∑
n
=
−
∞
∞
δ
(
t
−
n
Δ
T
)
s_{ΔT}(t)=∑^∞_{n=-∞}δ(t-nΔT)
sΔT(t)=n=−∞∑∞δ(t−nΔT)
冲激可以是连续的或者离散的。
4.1.4 卷积
在第三章介绍过这个概念。
4.2 取样和取样函数的傅里叶变换
4.2.1 取样
用计算机处理之前,连续函数必须转换为离散值序列。正如2.4节介绍的那样,这是用取样和量化来完成的。
模拟取样的一种方法是用一个ΔT单位间隔的冲激串作为取样函数去乘以f(t):
其中等式最左边(后续用f~(t)表示)表示取样后的函数。
4.2.2 取样函数的傅里叶变换
令F(μ)代表连续函数f(t)的傅里叶变换。
取样后的相应函数f~(t)是f(t)与一个冲激串的乘积。由卷积定理我们可知,空间域两个函数乘积的傅里叶变换是两个函数在频率域的卷积。
这样,取样后的函数f~(t)的傅里叶变换F(μ)是:
4.2.3 取样定理
对于以原点为中心的有限区间(带宽)[-μmax,μmax]之外的频率值,其傅里叶变换为零的函数f(t)称为带限函数。
如果以超过函数最高频率的两倍的取样率来获得样本,连续的带限函数可以完全地从它的样本集来恢复。这个结论就是众所周知的取样定理。
通过频率范围低端的频率,且会消除所有较高的频率的滤波器叫做低通滤波器。也被叫做理想低通滤波器,因为它在幅度上无限快速的过渡。
此外还有滤波器是从函数的取样恢复原始函数的手段叫做重建滤波器。
4.2.4 混淆
如果一个带限函数用低于其最高频率的两倍取样率取样将会发生什么情况?这相当于前一节中讨论的欠取样情况。
图4.9(a)说明了这一情况。以低于奈奎斯特取样率取样的最终效果是周期重叠,并且不管使用什么滤波器,都不可能分离出变换的一个单周期。例如,使用图4.9(b)中的理想低通滤波器,将会得到如图4.9(c)所示的一个变换,该变换已被来自邻近周期的频率破坏了。
然后,反变换会产生t的一个破坏了的函数。由函数欠取样导致的这种效果就是周知的频率混淆,或简称为混淆。
4.2.5 有取样后的数据重建(复原)函数
由一组样本集合来重建函数实际上可以减少样本间的内插。
即使简单地用显示介质显示一幅图像的行为都要求来自其样本的图像重建。
因此,理解取样数据重建的基础是很重要的。
卷积是这一理解的核心,因而再度显示了这一概念的重要性。
4.3 单变量的离散傅里叶变换(DFT)
离散傅里叶变换的表达式为:
F
m
=
∑
n
=
0
M
−
1
f
n
e
−
j
2
π
m
n
/
M
,
m
=
0
,
1
,
2
,
⋅
⋅
⋅
,
M
−
1
F_m=∑^{M-1}_{n=0}f_ne^{-j2πmn/M},m=0,1,2,···,M-1
Fm=n=0∑M−1fne−j2πmn/M,m=0,1,2,⋅⋅⋅,M−1
若给定一个由f(t)的M个样本组成的集合{fn},可以从上式得出一个与输入样本集合离散傅里叶变换相对应的M个复数离散值的样本集合{Fm}。反之,给定{Fm},我们可以用傅里叶反变换(IDFT)复原样本集{fn}:
f
n
=
1
M
∑
m
=
0
M
−
1
F
m
e
j
2
π
m
n
/
M
,
n
=
0
,
1
,
2
,
⋅
⋅
⋅
,
M
−
1
f_n={1\over M}∑^{M-1}_{m=0}F_me^{j2πmn/M},n=0,1,2,···,M-1
fn=M1m=0∑M−1Fmej2πmn/M,n=0,1,2,⋅⋅⋅,M−1
4.3.2 取样和频率间隔间的关系
如果f(x)由函数f(t)以ΔT为单位间隔取样后的M个样本组成,则包含集合{f(x)},x=0,1,2,···,M-1的记录的持续时间是
T
=
M
Δ
T
T=M ΔT
T=MΔT
离散频率域中的相应间隔Δu来自式:
Δ
u
=
1
M
Δ
T
=
1
T
Δu={1\over M ΔT}={1\over T}
Δu=MΔT1=T1
由DFT的M个分量跨越的整个频率范围是
Ω
=
M
Δ
u
=
1
Δ
T
Ω=M Δu={1\over ΔT}
Ω=MΔu=ΔT1
这样,可以看出,DFT的频率分辨率Δu取决于连续函数f(t)被取样的持续时间T,并日 DFT跨越的频率范围取决于取样间隔ΔT。很明显,两个表达式给出了T和ΔT的反转关系。
4.4 两个变量的函数的扩展
4.4.1 二维冲激及其取样特性
如在一维情况那样,对于坐标(x0,y0),二维冲激在积分下呈现了取样特性:
∫
−
∞
∞
∫
−
∞
∞
f
(
x
,
y
)
δ
(
x
−
x
0
,
y
−
y
0
)
d
x
d
y
=
f
(
x
0
,
y
0
)
∫^∞_{-∞}∫^∞_{-∞}f(x,y)δ(x-x_0,y-y_0)dxdy=f(x_0,y_0)
∫−∞∞∫−∞∞f(x,y)δ(x−x0,y−y0)dxdy=f(x0,y0)
4.4.2 二维连续傅里叶变换对
令f(t,z)是两个连续变量r和z的连续函数。则其二维连续傅里叶变换对由以下两个表达式给出:
F
(
μ
,
v
)
=
∫
−
∞
∞
∫
−
∞
∞
f
(
t
,
z
)
e
−
j
2
π
(
μ
t
+
v
z
)
d
t
d
z
F(μ,v)=∫^∞_{-∞}∫^∞_{-∞}f(t,z)e^{-j2π(μt+vz)}dtdz
F(μ,v)=∫−∞∞∫−∞∞f(t,z)e−j2π(μt+vz)dtdz
f
(
t
,
z
)
=
∫
−
∞
∞
∫
−
∞
∞
F
(
μ
,
v
)
e
j
2
π
(
μ
t
+
v
z
)
d
μ
d
v
f(t,z)=∫^∞_{-∞}∫^∞_{-∞}F(μ,v)e^{j2π(μt+vz)}dμdv
f(t,z)=∫−∞∞∫−∞∞F(μ,v)ej2π(μt+vz)dμdv
其中,μ和v是频率变量。当涉及图像时,t和z解释为连续空间变量。就像一维情况那样,变量μ和v的域定义了连续频率域。
4.4.3 二维取样和二维取样定理
类似于一维情况中的方式,二维取样可用取样函数(二维冲激串)建模:
s
Δ
T
Δ
Z
(
t
,
z
)
=
∑
m
=
−
∞
∞
∑
n
=
−
∞
∞
δ
(
t
−
m
Δ
T
,
z
−
n
Δ
Z
)
s_{ΔTΔZ}(t,z)=∑^{∞}_{m=-∞}∑^{∞}_{n=-∞}δ(t-mΔT,z-nΔZ)
sΔTΔZ(t,z)=m=−∞∑∞n=−∞∑∞δ(t−mΔT,z−nΔZ)
4.4.4 二维离散傅里叶变换及其反变换
类似于前面,二维傅里叶变换(DFT)如下:
F
(
u
,
v
)
=
∑
x
=
0
M
−
1
∑
y
=
0
N
−
1
f
(
x
,
y
)
e
−
j
2
π
(
u
x
/
M
+
v
y
/
N
)
F(u,v)=∑^{M-1}_{x=0}∑^{N-1}_{y=0}f(x,y)e^{-j2π(ux/M+vy/N)}
F(u,v)=x=0∑M−1y=0∑N−1f(x,y)e−j2π(ux/M+vy/N)
其中,f(x,y)是大小为MxN的数字图像。与一维中的情况一样,必须对离散变量u和v在u=0,1,2,…,M-1和v=0,1,2,…N-1范围内求值。
给出变换F(u,v),我们可以使用傅里叶反变换(IDFT)得到f(x,y):
f
(
x
,
y
)
=
1
M
N
∑
u
=
0
M
−
1
∑
v
=
0
N
−
1
F
(
u
,
v
)
e
j
2
π
(
u
x
/
M
+
v
y
/
N
)
,
n
=
0
,
1
,
2
,
⋅
⋅
⋅
,
M
−
1
f(x,y)={1\over MN}∑^{M-1}_{u=0}∑^{N-1}_{v=0}F(u,v)e^{j2π(ux/M+vy/N)},n=0,1,2,···,M-1
f(x,y)=MN1u=0∑M−1v=0∑N−1F(u,v)ej2π(ux/M+vy/N),n=0,1,2,⋅⋅⋅,M−1
其中,x=0,1,2…,M -1,y=0,1,2…,N-1。
两个式子构成了二维离散傅里叶变换对。本章的其余部分是以这两个公式的特性和它们在频率域图像滤波的应用为基础的。
4.5 二维离散傅里叶变换的一些性质
4.5.1 空间和频率间隔的关系
空间取样和相应频率域间隔的关系如前面中解释的那样。频率域样本间的间隔与空间样本间的间距和样本数成反比。
4.5.2 平移和旋转
用指数项乘以 f(x,y)将使 DFT 的原点移到点(u0,v0);反之,用负指数乘以F(u,v)将使f(x,y)的原点移到点(x0,y0)。
在极坐标下进行旋转,使用极坐标
x
=
r
c
o
s
θ
,
y
=
r
s
i
n
θ
,
u
=
ω
c
o
s
φ
,
v
=
φ
s
i
n
φ
x=rcosθ,y=rsinθ,u=ωcosφ,v=φsinφ
x=rcosθ,y=rsinθ,u=ωcosφ,v=φsinφ
可得到下列变换对:
f
(
r
,
θ
+
θ
0
)
↔
F
(
ω
,
φ
+
φ
0
)
f(r,θ+θ_0)↔F(ω,φ+φ_0)
f(r,θ+θ0)↔F(ω,φ+φ0)
它指出,若 f(x,y)旋转θ0角度,则 F(u,v)也旋转相同的角度。反之,若 F(u,v)旋转一个角度,则f(x,y)也旋转相同的角度。
4.5.3 周期性
如在一维情况中那样,二维傅里叶变换及其反变换在u方向和v方向是无限周期的.
4.5.4 对称性
如表4.1所示
4.5.5 傅里叶谱和相角
因为通常二维 DFT 一般是复函数,因此可使用极坐标形式来表示:
F
(
u
,
v
)
=
∣
F
(
u
,
v
)
∣
e
j
ϕ
(
u
,
v
)
F(u,v)=|F(u,v)|e^{jϕ(u,v)}
F(u,v)=∣F(u,v)∣ejϕ(u,v)
其中,幅度:
∣
F
(
u
,
v
)
∣
=
[
R
2
(
u
,
v
)
+
F
2
(
u
,
v
)
]
1
/
2
|F(u,v)|=[R^2(u,v)+F^2(u,v)]^{1/2}
∣F(u,v)∣=[R2(u,v)+F2(u,v)]1/2
称为傅里叶谱。
ϕ
(
u
,
v
)
=
a
r
c
t
a
n
[
I
(
u
,
v
)
R
(
u
,
v
)
]
ϕ(u,v)=arctan[{I(u,v)\over R(u,v)} ]
ϕ(u,v)=arctan[R(u,v)I(u,v)]
称为相角。
4.5.6 二维卷积定理
将一维的扩展至两个变量可得到下面的二维循环卷积的表达式:
f
(
x
,
y
)
★
h
(
x
,
y
)
=
∑
m
=
0
M
−
1
∑
n
=
0
N
−
1
f
(
m
,
n
)
h
(
x
−
m
,
y
−
n
)
f(x,y)★h(x,y)=∑^{M-1}_{m=0}∑^{N-1}_{n=0}f(m,n)h(x-m,y-n)
f(x,y)★h(x,y)=m=0∑M−1n=0∑N−1f(m,n)h(x−m,y−n)
其中,x=0,1,2…,M -1,y=0,1,2…,N-1。
4.6 频率域滤波基础
4.6.1 频率域滤波基础
频率域滤波由修改一幅图像的傅里叶变换然后计算其反变换得到处理后的结果组成,若给定一幅大小为M*N的数字图像f(x,y),则基本滤波公式为:
g
(
x
,
y
)
=
ℑ
−
1
[
H
(
u
,
v
)
F
(
u
,
v
)
]
,
g(x,y)=\Im ^{-1}[H(u,v)F(u,v)],
g(x,y)=ℑ−1[H(u,v)F(u,v)],
其中,
ℑ
−
1
\Im ^{-1}
ℑ−1为IDFT,
F
(
u
,
v
)
F(u,v)
F(u,v)为输入图像的DFT,
H
(
u
,
v
)
H(u,v)
H(u,v)是滤波函数。
4.7 使用频率域滤波器平滑图像
4.7.1 理想低通滤波器(ILPF)
在以原点为圆心、以D0为半径的圆内,无衰减地通过所有频率,而在该圆外“切断”所有频率的二维低通滤波器,称为理想低通滤波器(ILPF) 。
4.7.2 布特沃斯低通滤波器(BLPF)
截止频率位于距原点D0处的n阶布特沃斯低通滤波器(BLPF)的传递函数定义为:
H
(
u
,
v
)
=
1
1
+
[
D
(
u
,
v
)
/
D
0
]
2
n
H(u,v)={1\over {1+[D(u,v)/D_0]^{2n}}}
H(u,v)=1+[D(u,v)/D0]2n1
与ILPF不同,BLPF传递函数并没有在通过频率和滤除频率之间给出明显截止的尖锐的不连续性
4.7.3 高斯低通滤波器(GLPF)
表达式如下:
H
(
u
,
v
)
=
e
−
D
2
(
u
,
v
)
/
2
σ
2
H(u,v)=e^{-D^2(u,v)/2σ^2}
H(u,v)=e−D2(u,v)/2σ2
4.7.4 低通滤波的其他例子
可以用于机器感知领域,即字符识别的应用。
还可用于印刷和出版业。
还可以用到卫星图像和航空图像的处理。
4.8 使用频率域滤波器锐化图像
4.8.1 理想高通滤波器(IHPF)
IHPF与ILPF是相对的,在这种情况下,IHPF把以半径为D0的圆内的所有频率置零,而毫无衰减地通过圆外的所有频率。如同ILPF那样,IHPF也是物理上无法实现的。
4.8.2 布特沃斯高通滤波器(BHPF)
截止频率为D0的n阶布特沃斯高通滤波器(BHPF)定义为
H
(
u
,
v
)
=
1
1
+
[
D
0
/
D
(
u
,
v
)
]
2
n
H(u,v)={1\over {1+[D_0/D(u,v)]^{2n}}}
H(u,v)=1+[D0/D(u,v)]2n1
4.8.3 高斯高通滤波器(GHPF)
截止频率处在距频率矩形中心距离为D0的高斯高通滤波器(GHPF)的传递函数由下式给出:
H
(
u
,
v
)
=
1
−
e
−
D
2
(
u
,
v
)
/
2
D
0
2
H(u,v)=1-e^{-D^2(u,v)/2D_0^2}
H(u,v)=1−e−D2(u,v)/2D02
4.8.4 频率域的拉普拉斯算子
拉普拉斯算子可使用如下滤波器在频率城实现:
H
(
u
,
v
)
=
−
4
π
2
(
u
2
+
v
2
)
H(u,v)=-4π^2(u^2+v^2)
H(u,v)=−4π2(u2+v2)
然后,拉普拉斯图像由下式得到:
▽
2
f
(
x
,
y
)
=
ℑ
−
1
{
H
(
u
,
v
)
F
(
u
,
v
)
}
▽²f(x,y)=\Im ^{-1} \{H(u,v)F(u,v)\}
▽2f(x,y)=ℑ−1{H(u,v)F(u,v)}
F(u,v)是f(u,v)的傅里叶变换,增强可使用下式实现:
g
(
x
,
y
)
=
f
(
x
,
y
)
+
c
▽
2
f
(
x
,
y
)
g(x,y)=f(x,y)+c▽^2 f(x,y)
g(x,y)=f(x,y)+c▽2f(x,y)
整理一下,增强可以使用下式实现:
g
(
x
,
y
)
=
ℑ
−
1
{
F
(
u
,
v
)
−
H
(
u
,
v
)
F
(
u
,
v
)
}
=
ℑ
−
1
{
[
1
−
H
(
u
,
v
)
]
F
(
u
,
v
)
}
=
ℑ
−
1
{
[
1
+
4
π
2
D
2
(
u
,
v
)
]
F
(
u
,
v
)
}
g(x, y)= \Im ^{-1}\{F(u,v)-H(u,v)F(u,v)\}=\Im ^{-1}\{[1-H(u,v)]F(u,v)\}=\Im ^{-1}\{[1+ 4π²D^2(u,v)]F(u,v)\}
g(x,y)=ℑ−1{F(u,v)−H(u,v)F(u,v)}=ℑ−1{[1−H(u,v)]F(u,v)}=ℑ−1{[1+4π2D2(u,v)]F(u,v)}
4.8.4 滤波
以称之为同态系统的一类系统的特殊情况为基础,这种滤波称之为同态滤波,步骤如下图,方法的关键在于照射分量和反射分量的分离。
4.9 选择性滤波
选择性滤波指的是处理指定的频段或频率矩形的小区域。
第一类滤波器分别称为带阻滤波器或带通滤波器。第二类滤波器称为陷波滤波器。
4.10 实现
f(x,y)的二维DFT可通过计算f(x,y)的每一行的一维变换,然后,沿着计算结果的每一列计算一维变换来得到。
这是一个重要的简化,因为我们必须一次仅处理一个变量。类似的扩展可用一维IDFT来计算二维IDFT。