1、第7章 频域处理第第7 7章章 频域处理频域处理7.1 频域与频域变换7.2 傅立叶变换7.3 频域变换的一般表达式7.4 离散余弦变换(DCT)7.5 离散沃尔什-哈达玛变换(WHT)7.6 频域中图像处理的实现7.7 用MatrixLIBC+库实现图像变换的Visual C+编程7.8 小波变换简介第7章 频域处理7.1 频域与频域变换频域变换的理论基础就是“任意波形都可以用单纯的正弦波的加权和来表示”。如图7-1(a)所示的任意波形，可分解为图7-1(b)、(c)、(d)所示的不同幅值、不同频率的正弦波的加权和。为便于理解，将图7-1(b)所示的正弦波取出来，如图7-2所示。如果将虚线表
2、示的振幅为1且初相位为0的正弦波作为基本正弦波，则实线表示的波形可由其振幅A和初相位确定。第7章 频域处理图7-1 任意波形可分解为正弦波的加权和 第7章 频域处理图7-2 正弦波的振幅A和相位 第7章 频域处理由此，图7-1(b)、(c)、(d)三个不同的正弦波形可以描述为图7-3所示的两幅图。其中图7-3(a)表示振幅与频率之间的关系，称为幅频特性;而图7-3(b)表示初相位与频率之间的关系，称为相频特性。这样便将图7-1(a)所示的时域波形f(x)变换到图7-3所示的频域F(f)。显然，不管波形多么复杂，均可将其变换到频域。第7章 频域处理图7-3 图7-1(a)波形的频域表示 第7章
3、频域处理时域和频域之间的变换可用数学公式表示如下：)(),()(ffAxf正变换逆变换(7-1)式中：A(f)、(f)分别为幅值和相位与频率f之间的关系。为能同时表示信号的振幅和相位，通常采用复数表示法。式(7-1)可用复数表示法表示为)()(fFxf正变换逆变换(7-2)式中：F(f)用复数表示幅值、相位与频率f之间的关系。第7章 频域处理7.2 傅立叶变换7.2.1 连续函数的傅立叶变换若一个一维信号f(x)满足狄里赫莱条件，即(1)具有有限个间断点；(2)具有有限个极值点；(3)绝对可积。则其傅立叶变换对(傅立叶变换和逆变换)一定存在。在实际应用中，这些条件一般总是可以满足的。第7章 频
4、域处理一维傅立叶变换对定义为 dxexfuFxfxuj2)()()(dueuFxfuFxuj21)()()(7-3)(7-4)式中：；x为时域变量；u为频域变量。1j第7章 频域处理以上一维傅立叶变换可以很容易推广到二维。如果二维函数f(x，y)满足狄里赫莱条件，则它的二维傅立叶变换对为 dxdyeyxfvuFyxfvyuxj)(2),(),(),(dudvevuFyxfvuFvyuxj)(21),(),(),(7-5)(7-6)式中：x，y为时域变量；u，v为频域变量。第7章 频域处理7.2.2 离散傅立叶变换在数字图像处理中应用傅立叶变换，还需要解决两个问题：一是在数学中进行傅立叶变换的f
5、(x)为连续(模拟)信号，而计算机处理的是数字信号(图像数据)；二是数学上采用无穷大概念，而计算机只能进行有限次计算。通常，将受此限制的傅立叶变换称为离散傅立叶变换(Discrete Fourier Transform,DFT)。第7章 频域处理定义：设f(x)|f(0)，f(1)，f(2)，f(N-1)为一维信号f(x)的N个抽样，其离散傅立叶变换对为102)()()(NxNxujexfuFxf1021)(1)()(NxNxujeuFNxfuF(7-7)(7-8)式中：x，u=0，1，2，N1。第7章 频域处理注意：式(7-8)中的系数1/N也可以放在式(7-7)中，有时也可在傅立叶正变换和
6、逆变换前分别乘上，这是无关紧要的，只要正变换和逆变换前系数乘积等于1/N即可。由欧拉公式可知：N1sincosjej(7-9)将式(7-9)代入式(7-7)，并利用cos()=cos()，可得：)2sin2cos)()(10NxujNxuxfuFNx(7-10)第7章 频域处理可见，离散序列的傅立叶变换仍是一个离散的序列，对每一个u对应的傅立叶变换结果是所有输入序列f(x)的加权和(每一个f(x)都乘以不同频率的正弦和余弦值)，u决定了每个傅立叶变换结果的频率。通常傅立叶变换为复数形式，即)()()(ujIuRuF(7-11)式中：R(u)和I(u)分别是F(u)的实部和虚部。式(7-11)也
7、可表示成指数形式：)()()(ujeuFuF(7-12)第7章 频域处理式中：)()()(22uIuRuF)()(arctan)(uRuIu(7-13)(7-14)通常，称|F(u)|为f(x)的频谱或傅立叶幅度谱，(u)为f(x)的相位谱。频谱的平方称为能量谱或功率谱，它表示为)()()()(222uIuRuFuE(7-15)第7章 频域处理 考虑到两个变量，就很容易将一维离散傅立叶变换推广到二维。二维离散傅立叶变换对定义为1010)(2),(),(),(MxNyNvyMuxjeyxfvuFyxf1010)(21),(1),(),(MuNvNvyMuxjevuFMNyxfvuF(7-16)(
8、7-17)式中：u，x=0，1，2，M-1；v，y=0，1，N-1；x，y为时域变量；u，v为频域变量。第7章 频域处理像一维离散傅立叶变换一样，系数1/(MN)可以在正变换或逆变换中，也可以分别在正变换和逆变换前分别乘上系数，只要两系数的乘积等于1/(MN)即可。二维离散函数的傅立叶频谱、相位谱和能量谱分别为MN1),(),(),(22vuIvuRvuF),(),(arctan),(vuRvuIvu),(),(),(22vuIvuRvuE(7-18)(7-19)(7-20)式中：R(u，v)和I(u，v)分别是F(u，v)的实部和虚部。第7章 频域处理7.2.3 离散傅立叶变换的性质二维离散
9、傅立叶变换的性质对图像的分析具有十分重要的作用，因此，有必要理解和掌握二维DFT的性质。二维离散傅立叶变换的主要性质如表7-1所示。第7章 频域处理表7-1 二维DFT的性质 第7章 频域处理 续表第7章 频域处理以上具有重要意义的两个性质的含义如下。1.可分离性由可分离性可知，一个二维傅立叶变换可分解为两步进行，其中每一步都是一个一维傅立叶变换。可先对f(x，y)按行进行傅立叶变换得到F(x，v)，再对F(x，v)按列进行傅立叶变换，便可得到f(x，y)的傅立叶变换结果F(u，v)，如图7-4所示。显然先按列进行傅立叶变换，再按行进行傅立叶变换也是可行的。第7章 频域处理图7-4 用两次一维
10、DFT计算二维DFT 第7章 频域处理同理，傅立叶变换的逆变换也具有可分离性。利用傅立叶变换的可分离性，可以简化傅立叶变换的软、硬件设计，用一维傅立叶变换软件或硬件便可实现二维傅立叶变换。第7章 频域处理2.平移性质平移性质表明只要将f(x，y)乘上因子(1)x+y再进行离散傅立叶变换，即可将图像的频谱原点(0，0)移动到图像中心(M/2，N/2)处。图7-5(a)是一简单方块图像，图(b)是其无平移的傅立叶频谱，图(c)是平移后的傅立叶频谱。可见，利用傅立叶变换的平移性质将图像频谱原点移动到图像中心，更便于分析和处理，特别是设计滤波器时更加方便。第7章 频域处理图7-5 傅立叶频谱平移示意图
11、 第7章 频域处理由表7-1中性质9可知，图像的频谱原点(0，0)代表的是图像灰度的平均值，是图像信号中的直流分量。因此，平移后的频谱中，图像能量的低频成分将集中到频谱中心，图像上的边缘、线条细节信息等高频成分将分散在图像频谱的边缘。第7章 频域处理7.2.4 快速离散傅立叶变换基本离散傅立叶变换计算量非常大，运算时间长。可以证明其运算次数正比于N2，特别是当N较大时，其运算时间将迅速增长，以至于无法容忍。为此，需要研究离散傅立叶变换的快速算法(Fast Fourier Transform，FFT)。1965年Cooley和Tukey首先提出了一种称为逐次加倍法的快速傅立叶变换算法(FFT)。
12、采用该FFT算法，其运算次数正比于N lbN，在N很大时计算量可以大大减少。例如，FFT的运算次数和DFT的运算次数之比，当N=1024时为1/102.4；当N=4096时可达1/341.3。第7章 频域处理Cooley-Tukey FFT算法的基本思想是将f(x)序列按x的奇偶进行分组计算，并充分利用傅立叶变换的周期性和对称性进行计算，采用迭代法，大大简化了程序设计的复杂度，提高了计算速度。Sande-Tukey FFT算法与Cooley-Tukey FFT算法类似，只不过它是将f(x)序列按中心位置点进行分组计算的。限于篇幅，FFT算法的详细内容不再冗述。当然对于计算机专业的学生而言，每个
13、人都应尝试编写快速傅立叶变换的程序。有关傅立叶变换的算法还有很多，网上的FFT算法源代码也非常多，但不建议大家拿来就用。当你得到类似的代码后，一定要认真分析其实现过程和思路，只有这样才能不断地提高编程水平。第7章 频域处理用傅立叶变换处理和分析信号，就像用三棱镜分解光线一样。让一束白光通过三棱镜，可将白光分解成七色的彩虹，若将分解开的七色光再次通过三棱镜，又可以得到白光。从形式上看这是由简单变换出了繁复，实则是将混合的东西分解成了基本的元素，通过对其基本元素的分析与处理，进而完成对信号的处理和分析。因此，傅立叶变换又有“数字棱镜”的美誉。第7章 频域处理7.3 频域变换的一般表达式7.3.1
14、可分离变换二维傅立叶变换可用通用关系式来表示：1010),(),(),(MxNyvuyxgyxfvuF1010),(),(),(MuNvvuyxhvuFyxf(7-21)(7-22)式中：x，u取0，1，2，M1；y，v取0，1，2，N1；g(x，y，u，v)和h(x，y，u，v)分别称为正向变换核和反向变换核。第7章 频域处理如果),(),(),(21vyguxgvuyxg),(),(),(21vyhuxhvuyxh(7-23)(7-24)则称正反变换核是可分离的。进一步，如果g1和g2，h1和h2在函数形式上一样，则称该变换核是对称的。第7章 频域处理二维傅立叶变换对是式(7-21)和式(
15、7-22)的一个特殊情况，它们的变换核为NvyjMuxjNvyMuxjeeevuyxg22)(2),(NvyjMuxjNvyMuxjeNeMeMNvuyxh22)(2111),(7-25)(7-26)可见，它们均为可分离的和对称的。第7章 频域处理如前所述，二维傅立叶变换可以利用变换核的可分离性，用两次一维变换来实现，即可先对f(x，y)的每一行进行一维变换得到F(x，v)，再沿F(x，v)每一列取一维变换得到变换结果F(u，v)。对其它的图像变换，只要其变换核是可分离的，同样可用两次一维变换来实现二维变换。若先对f(x，y)的每一列进行一维变换得到F(y，u)，再沿F(y，u)每一行取一维变
16、换得到F(u，v)，其最终结果相同。该结论对逆变换也适用。第7章 频域处理7.3.2 图像变换的矩阵表示数字图像都是实数矩阵，设f(x，y)为MN的图像灰度矩阵，通常，为了分析、推导方便，将可分离变换写成矩阵的形式：QfPF 11FQPf(7-27)(7-28)式中：F、f是二维MN的矩阵；P是MM矩阵；Q是NN矩阵。第7章 频域处理图像变换的矩阵表达式和代数表达式其本质相同，将式(7-27)写成代数表达式如下：1010),(),(),(),(MxNyvxQyxfuxPvuF(7-29)式中：u取0，1，2，M1；v取0，1，2，N1。对二维离散傅立叶变换，则有：MxujeuxguxP21),
17、(),(NvyjevygvyQ22),(),(7-30)(7-31)图像处理实践中，除了DFT变换之外，还可采用其它正交变换，例如离散余弦变换、沃尔什-哈达玛变换、K-L变换等。下面对常用的变换作简要介绍。第7章 频域处理7.4 离散余弦变换(DCT)离散余弦变换(Discrete Cosine Transform，DCT)的变换核为余弦函数，因其变换核为实数，所以，DCT计算速度比变换核为复数的DFT要快得多。DCT除了具有一般的正交变换性质外，它的变换阵的基向量能很好地描述人类语音信号、图像信号的相关特征。因此，在对语音信号、图像信号的变换中，DCT变换被认为是一种准最佳变换。近年颁布的一
18、系列视频压缩编码的国际标准中，均把DCT作为其中的一个基本处理模块。此外，DCT也是一种可分离的变换。第7章 频域处理7.4.1 一维离散余弦变换一维DCT的变换核定义为NuxNuCuxg2)12(cos2)(),(7-32)式中：x，u取0，1，2，N1；且其它1021)(uuC(7-33)第7章 频域处理设f(x)|x=0，1，N-1为离散的信号列，则一维DCT定义如下：102)12(cos)(2)()(NxNuxxfNuCuF(7-34)式中：u，x取0，1，2，N1。将变换式展开整理后，可以写成矩阵形式：GfF(7-36)第7章 频域处理其中：2)12)(1(cos)23)1(cos(
19、)2)1(cos(2 )22)12(cos()26cos()22cos(2 )2)12(cos()23cos()2cos(2 1111NNNNNNNNNNNNNNNNNNNG(7-36)第7章 频域处理一维DCT的逆变换IDCT(Inverse DCT)定义为 102)12(cos)()(2)(NuNuxuFuCNxf(7-37)式中：x，u取0，1，2，N1。可见一维DCT的逆变换核与正变换核是相同的。第7章 频域处理7.4.2 二维离散余弦变换考虑到两个变量，很容易将一维DCT的定义推广到二维DCT。其正变换核为NvyMuxvCuCMNvuyxg2)12(cos2)12(cos)()(2)
20、,(7-38)式中：C(u)、C(v)定义同式(7-33)；x，u取0，1，2，M1；y，v取0，1，2，N1。第7章 频域处理设f(x，y)为MN的二维离散信号，其二维DCT定义如下：10102)12(cos2)12(cos)()(),(2),(MxNyNvyMuxvCuCyxfMNvuF(7-39)式中：x，u取0，1，2，M1；y，v取0，1，2，N1。第7章 频域处理二维DCT逆变换定义如下：10102)12(cos2)12(cos),()()(2),(MuNvNvyMuxvuFvCuCMNyxf(7-40)式中，x，u取0，1，2，M1；y，v取0，1，2，N1。类似一维矩阵形式的D
21、CT，可以写出二维DCT的矩阵形式如下：TGfGF(7-41)第7章 频域处理同时，由式(7-40)和式(7-39)可知二维DCT的逆变换核与正变换核相同，且是可分离的，即),(),(),(21vyguxgvuyxgNvyvCNMuxuCM2)12(cos)(22)12(cos)(2(7-42)式中：C(u)、C(v)定义同式(7-33)；x，u取0，1，2，M1；y，v取0，1，2，N1。第7章 频域处理根据DCT的可分离性，二维DCT可用两次一维DCT来完成，其算法流程与DFT类似：),(),(),(),(),(),(),(vuFvuFvxFFvxFvxFyxfFyxfTTT转置列转置行(
22、7-43)第7章 频域处理离散余弦变换的计算量也相当大，在实用中非常不方便，也需要研究相应的快速算法。目前已有多种快速DCT(FCT)，其一是由FFT的思路发展起来的，它利用FFT和IFFT便可实现快速DCT(FCT)和快速IDCT算法(IFCT)。不过，由于FFT及IFFT中涉及到复数运算，所以，该FCT及IFCT算法并不是最佳的。限于篇幅，不再介绍其它的快速DCT及IDCT算法，请参考相关资料。最后，需要注意的是二维DCT的频谱分布与DFT相差一倍，如图7-6所示。由图可以看出，对于DCT而言，(0，0)点对应于频谱的低频成分，(N-1，N-1)点对应于高频成分；而同阶的DFT中，(N/2
23、，N/2)点对应于高频成分(注：此频谱图中未作频谱中心平移)。第7章 频域处理图7-6 DFT和DCT的频谱分布 第7章 频域处理 7.5 离散沃尔什-哈达玛变换 7.5.1 一维离散沃尔什-哈达玛变换1.沃尔什函数沃尔什函数是1923年由美国数学家沃尔什(Walsh)提出的。沃尔什函数系是一个完备正交函数系，其值只能取1和1。从排列次序上可将沃尔什函数分为三种定义方法：按照沃尔什排列来定义(按列率排序)；按照佩利排列来定义(按自然排序)；按照哈达玛排列来定义。由于哈达玛排序的沃尔什函数是由2n阶哈达玛矩阵(Hadamard Matrix)得到的，哈达玛矩阵的最大优点在于它具有简单的递推关系，
24、即高阶矩阵可用两个低阶矩阵的克罗内克积求得，因此，在此只介绍这种定义的沃尔什变换。第7章 频域处理N=2n阶哈达玛矩阵每一行的符号变化规律对应于某一个沃尔什函数的符号变化规律，即N=2n阶哈达玛矩阵的每一行对应于一个离散沃尔什函数，哈达玛矩阵与沃尔什函数系不同之处仅仅是行的次序不同。2n阶哈达玛矩阵有如下形式：11H11112H111111111111111122224HHHHH(7-44)(7-45)(7-46)第7章 频域处理 哈达玛矩阵的阶数是按N2n规律排列的(n0，1，2，)，高阶数哈达玛矩阵可以利用矩阵的克罗内克积运算，由低阶哈达玛矩阵递推得到：2222222222211111NN
25、NNNHHHHHHHHHHHHnnnnnn(7-47)第7章 频域处理矩阵的克罗内克积(Kronecker Product)运算用符号记作A B，设：ijiijjmnmmnnbbbbbbbbbBaaaaaaaaaA212222111211212222111211,第7章 频域处理则其运算规律如下：BaBaBaBaBaBaBaBaBaBAmnmmnn212222111211AbAbAbAbAbAbAbAbAbABijiijj212222111211(7-48)(7-49)利用式(7-48)和式(7-49)便可产生任意阶的哈达玛矩阵。第7章 频域处理2.离散沃尔什-哈达玛变换一维离散沃尔什变换定义
26、为10),()(1)(NxxuWalshxfNuW(7-50)一维离散沃尔什逆变换定义为10),()()(NuxuWalshuWxf(7-51)式中：Walsh(u，x)为沃尔什函数。第7章 频域处理若将Walsh(u，x)用哈达玛矩阵表示，并将变换表达式写成矩阵形式，则式(7-50)和式(7-51)为)1()1()0(1)1()1()0(NfffHNNWWWN)1()1()0()1()1()0(NWWWHNfffN(7-52)(7-53)式中：HN为N阶哈达玛矩阵。第7章 频域处理由哈达玛矩阵的特点可知，沃尔什-哈达玛变换的本质是将离散序列f(x)的各项值的符号按一定规律改变后，进行加减运算
27、，因此，它比采用复数运算的DFT和采用余弦运算的DCT要简单得多。第7章 频域处理7.5.2 二维离散沃尔什变换很容易将一维WHT的定义推广到二维WHT，二维WHT的正变换核和逆变换核分别为1010),(),(),(1),(MxNyyvWalshxuWalshyxfMNvuW1010),(),(),(),(MuNvyvWalshxuWalshvuWyxf(7-54)(7-55)式中：x,u取0,1,2,M1；y,v取0,1,2,N1。第7章 频域处理二维离散沃尔什变换的矩阵形式表达式为NMHfHMNW1NMWHHf(7-56)(7-57)式中：HM为M阶哈达玛矩阵；HN为N阶哈达玛矩阵。第7章
28、 频域处理例如，有两个二维数字图像信号矩阵分别为 13311331133113311f11111111111111112f和求这两个信号的二维WHT。第7章 频域处理根据题意，式(7-56)中M=N=4，其二维WHT变换核为11111111111111114H第7章 频域处理所以：00000000000010021111111111111111133113311331133111111111111111114121W00000000000000011111111111111111111111111111111111111111111111114122W第7章 频域处理再如，图7-7是一幅数字图
29、像及对其进行二维WHT变换的结果。图7-7 二维WHT结果 第7章 频域处理注意：图7-7中的结果是哈达玛变换后再用沃尔什排序后的结果。从以上的例子可看出，二维WHT具有能量集中的特性，而且原始数据中数字越是均匀分布，经变换后的数据越集中于矩阵的边角上。因此，二维WHT可用于压缩图像信息。第7章 频域处理7.5.3 快速沃尔什变换(FWHT)类似于FFT计算，WHT也有快速算法FWHT，也可将输入序列f(x)按奇偶进行分组，分别进行WHT。FWHT的基本关系为)()(21)2()()(21)(uWuWNuWuWuWuW奇偶奇偶(7-58)第7章 频域处理WHT的变换核是可分离和对称的，因此，二
30、维WHT也可分为两个一维的WHT，分别用FWHT进行变换而得到最终结果，由此便可实现二维的FWHT。综上所述，WHT将一个函数变换成取值为1或1的基本函数构成的级数，用它来逼近数字脉冲信号要比FFT有利。同时，WHT只需要进行实数运算，存储量比FFT要少得多，运算速度也快得多。因此，它在图像传输、通信技术和数据压缩中被广泛使用。第7章 频域处理 7.6 频域中图像处理的实现7.6.1 理解数字图像的频谱图在7.2.3节指出，数字图像平移后的傅立叶频谱中，图像的能量将集中到频谱中心(低频成分)，图像上的边缘、线条细节信息(高频成分)将分散在图像频谱的边缘。也就是说，频谱中低频成分代表了图像的概貌
31、，高频成分代表了图像中的细节。例如，一幅室内图像，墙和地板的灰度变化平缓，它们对应的是频谱中靠近中心的分量，当进一步远离频谱中心点时，较高的频率分量开始对应图像中变化急剧的灰度级，如墙和地板的交界、噪声等图像成分。第7章 频域处理图7-8是一幅图像及其傅立叶频谱。图7-8(a)是一幅放大了近2500倍的集成电路电子扫描显微镜图像，图中沿大约45方向存在强边缘，并且有两个因热感应不足而产生的白色氧化物突起。图7-8(a)的傅立叶频谱如图7-8(b)所示，在图(b)中沿着45的亮带对应了图像中的强边缘；沿着纵轴偏左的部分存在两个亮带对应氧化物突起，在此应注意亮带偏离纵轴的角度与白色氧化物突起的对应
32、关系。第7章 频域处理图7-8 图像及其傅立叶频谱 第7章 频域处理7.6.2 频域图像处理步骤在频域中进行图像处理的一般步骤如下：(1)用(-1)x+y乘以输入图像各像素值，以将图像频谱原点移动到频谱图中心；(2)计算图像的DFT，得到F(u，v)；(3)用滤波函数H(u，v)乘以F(u，v)，得到处理结果G(u，v)；(4)计算滤波后的IDFT；(5)取IDFT变换结果中的实部；第7章 频域处理(6)用(-1)x+y乘以IDFT变换结果的实部，得到处理后的图像。H(u，v)称为滤波器，它具有允许某些频率成分通过，而阻止其它频率成分通过的特性。该处理过程可表示为(7-59),(),(),(v
33、uHvuFvuGH和G的相乘是点乘运算，即H的第一个元素乘以F的第一个元素，H的第二个元素乘以F的第二个元素，依此类推。滤波后的图像可以由IDFT得到：(7-60),(),(1vuGyxg第7章 频域处理注意：最终的图像是取IDFT的实部，并乘以(-1)x+y得到的。图7-9给出了频域中图像处理的基本步骤。第7章 频域处理图7-9 频域处理的基本步骤 第7章 频域处理7.6.3 频域滤波如7.6.1所述，频谱中低频成分代表了图像的概貌，即灰度变化平缓的部分；高频成分代表了图像中的细节，即图像中的边缘或噪声等。因此，合理地选择滤波器，通过在频域中对相应频率成分进行抑制或增强，便可完成图像的增强处
34、理。频域滤波中，存在各种各样的滤波器。其基本类型有：低通滤波器、高通滤波器、带通滤波器和带阻滤波器。若将这些基本滤波器有机地组合在一起，便可形成各种各样的频域滤波器。第7章 频域处理1.低通滤波器顾名思义，低通滤波器允许低频成分通过，而抑制高频成分。因此，它能够去除图像中的噪声，实现图像平滑作用。当然，这必然会引起图像模糊。理想低通滤波器的滤波函数为(7-61)00),(0),(1),(DvuDDvuDvuH第7章 频域处理2.高通滤波器与低通滤波器相反，高通滤波器允许高频成分通过，而抑制低频成分。因此，它能够强化图像中目标的边缘，起锐化作用。但它同时也强化了图像中的噪声。理想高通滤波器的滤波
35、函数为(7-62)00),(1),(0),(DvuDDvuDvuH第7章 频域处理3.带通滤波器带通滤波器允许指定范围的频率成分通过，而抑制其它频率成分。理想带通滤波器的滤波函数为(7-63)其它0),(1),(21DvuDDvuH第7章 频域处理4.带阻滤波器带阻滤波器抑制指定范围的频率成分，而允许其它频率成分通过。理想带阻滤波器的滤波函数为(7-64)其它1),(0),(21DvuDDvuH式(7-61)式(7-64)中，D0、D1、D2是指定的非负值，称为截止频率；D(u，v)是(u，v)点到原点的距离。四种基本类型滤波器的频率响应特性如图7-10所示。第7章 频域处理图7-10 基本滤
36、波器的频率响应 第7章 频域处理图7-11是分别采用D010、D030、D060、D0160(以像素为单位)进行理想低通滤波的结果。在图7-11(c)中，存在着严重的模糊现象，这表明图像中多数的细节信息包含在被滤除掉的频率成分之中。随着滤波半径的增加，滤除的能量越来越少，图7-11(d)到图7-11(f)中的模糊现象也就越来越轻。当被滤除的高频成分减少时，图像的质量会逐渐变好，但其平滑作用也将逐渐减弱。一个值得注意的问题是在图7-11(c)到图7-11(e)中存在有明显的振铃现象，“振铃”现象产生的原因在此不再讨论，请参阅有关信号处理方面的资料。第7章 频域处理图7-11 理想低通滤波器处理结
37、果 第7章 频域处理由于理想滤波器存在明显的“振铃”现象，且其垂直的频率响应特性仅能用软件方法实现，无法用电路实现。因此，研究实用的滤波器是非常有价值的。一种常用的频域滤波器是巴特沃斯(Butterworth)滤波器。低通巴特沃斯滤波器的滤波函数为(7-65)nDvuDvuH20),(11),(第7章 频域处理高通巴特沃斯滤波器的滤波函数为(7-66)nvuDDvuH20),(11),(式中：；n为巴特沃斯滤波器乘阶数。巴特沃斯滤波器的频率响应特性如图7-12所示。22),(vuvuD第7章 频域处理图7-12 巴特沃斯滤波器的频率响应 第7章 频域处理图7-13是对图7-11(a)采用D0=
38、60，n=1的低通巴特沃斯滤波器的滤波结果。可以看出巴特沃斯滤波器可以有效地抑制“振铃”现象。第7章 频域处理图7-13 巴特沃斯滤波器及处理效果 第7章 频域处理注意：所有示例只是针对低通滤波而言的，对于高通滤波，其处理过程与低通滤波类似，在此不再赘述。在频域中，还可完成对指定频率成分的处理。图7-14是另一频域滤波实例。图7-14(a)是有条纹干扰的原图像，在图7-14(b)频谱图中，可以明显看到图像中存在高频噪声点。在频域中去除这些高频噪声点后如图7-14(c)所示，再经过IDFT变换，便可获得图7-14(d)所示的去除条纹干扰的图像。第7章 频域处理图7-14 频域滤波实例 第7章 频
39、域处理7.7 用MatrixLIBC+库实现图像变换 的Visual C+编程7.7.1 MatrixLIB简介及其与Visual C+工程的集成1.MatixLIB简介MatrixLibC+数学库是MathTools公司利用MatCom技术开发的一个面向专业从事工程技术和科学计算人员的矩阵运算动态链接库。该C+库提供了绝大多数的关于矩阵类、矩阵操作、数值计算等函数的定义，它提供了一个双精度Matrix类型Mm，该矩阵可以是复数矩阵、实数矩阵、稀疏矩阵甚至n维矩阵；同时还提供了近千个与矩阵运算相关的函数，并对如、*、/、()等操作符进行了重载。库中所提供的函数涉及线性代数、多项式数学、信号处理
40、、文件输入输出、绘图等，充分利用这些库函数，便可完成数字图像变换的各种操作。第7章 频域处理2.在C+工程中集成MatrixLIBC+数学库为了将MatrixLIBC+库集成到C+的工程文件中，需做如下操作：(1)将库文件v4501v.lib添加到C+的工程中；(2)包含库声明头文件matlib.h；(3)初始化库；(4)创建矩阵变量；(5)访问矩阵单元；(6)调用MatrixLIBC+库函数完成矩阵操作。第7章 频域处理为了使读者掌握如何在VC+中使用MatrixLIBC+数学库的方法，下面详细介绍从创建工程到最后调用MatrixLIBC+库函数的全过程。(1)创建一个新的VC+工程。在VC
41、中，采用默认选项创建一个MDI MFC工程。当然，MatrixLIBC+库可以用在任何其它类型的工程中，如Console、OWL、ActiveX等。(2)将库文件添加到工程中。为使C+编译器能够将MatrixLIBC+库链接到最后的执行文件中去，必须将MatrixLIBC+库加入到C+的工程文件中。用菜单命令Project/Add to Project/Files将C：matcom45libv4501v.lib(如果v4501v.lib存在不同的文件夹中，请选择正确的路径)加入到C+工程中。第7章 频域处理(3)包含matlib.h头文件。在要使用MatrixLIBC+库的.h或.cpp文件的
42、头部，将matlib.h头文件包含到工程中，即在使用MatrixLIBC+库的.h或.cpp文件头部添加如下代码：#include matlib.h同时用菜单命令Tools/Options/Directories设置matlib.h所在的路径，以便编译器能找到相应的头文件。(4)初始化MatrixLIBC+库。在调用MatrixLIBC+库函数之前，要用“initM(MATCOM_VERSION)”函数来初始化类库调用，并用“exitM()”函数来结束类库调用，故还应在.cpp 中加入如下代码：第7章 频域处理 initM(MATCOM_VERSION)；/调用MatrixLIBC+库函数以完
43、成矩阵操作 exitM()；当然可以在一个类的构造函数中添加：initM(MATCOM_VERSION)；以完成类库的初始化工作。可在该类的析构函数中添加：exitM()；以完成结束类库调用的工作。MATCOM_VERSION参数是一个在matlib.h头文件中定义的值为4501的常量(不同版本其值有所差异)，它保证了动态链接库与matlib.h的匹配。第7章 频域处理(5)创建矩阵变量。添加了MatrixLIBC+库后，可用关键字Mm来定义矩阵类型变量，例如：Mm a，b，x；该VC+语句在当前工程文件中创建了a、b、x三个矩阵变量，但a、b和x只是空的矩阵，它们的矩阵单元不包含任何值。(6
44、)访问矩阵单元。矩阵单元的访问包括读和写操作，主要有：通过.r()来访问矩阵的实部，通过.i()来访问矩阵的虚部；通过BR()函数把一个数据列表赋给矩阵各个单元；通过.addr()或.addi()返回矩阵变量的实部或虚部的内存指针，来完成对矩阵单元的访问。第7章 频域处理3.安装与MatrixLIBC+库对应的动态链接库当采用MatrixLIBC+数学库完成矩阵的运算时，要求系统中有相关的ago4501.dll、v4501v.dll、opengl32.dll和glu32.dll四个动态链接库文件。如果系统中安装有Matcom4.5，这些动态链接库文件将会自动安装在Windowsystem或Wi
45、ndowsystem32下；如果系统没有安装Matcom4.5，只需将这四个动态链接库文件拷贝到相应目录下即可。第7章 频域处理7.7.2 创建图像数据矩阵为利用矩阵运算完成图像变换，首先应将其图像数据赋给矩阵变量。设指针pBits指向读入内存中的图像数据首地址，可用如下代码将图像数据赋给矩阵变量(为简单其见，设图像为8位灰度图像)。(1)在类定义头文件.h中定义矩阵变量：Mm m_matBits；(2)在类实现文件.cpp中完成对矩阵变量的赋值操作：第7章 频域处理/图像数据总字节数，其中，nWidthBytes为每行图像/数据占用的字节数DWORD SizeImage=nWidthByte
46、s*nHeight；/通过创建一个1SizeImage的一维0矩阵，为矩阵/变量分配内存m_matBits=zeros(1，SizeImage)；/将矩阵由double型转换为unsigned char型，以便将/图像数据赋给矩阵m_matBits=muint8(m_matBits)；第7章 频域处理/首先利用m_matBits.addr()获得矩阵实部的内存指针，/然后利用memcpy()将图像数据一次性赋给矩阵变量memcpy(m_matBits.addr()，pBits，SizeImage)；/由于Mm型的矩阵在内存中是按先列后行的顺序存储/的，先用MatrixLIBC+的库函数resh
47、ape()将一维/矩阵m_matBits变维为nHeightnWidthBytes的二维矩/阵，再用rot90()函数将其逆时针旋转90度，将矩阵变/为nWidthBytesnHeight的二维矩阵m_matBits=rot90(reshape(m_matBits，nWidthBytes，nHeight)；第7章 频域处理/若图像宽度与其字节宽度不同，则将由系统补齐的/每行字节数为4的整数倍的各列0删除，以减小矩阵/大小加快处理速度if(nWidthBytes!=nWidth)m_matBits=m_matBits(c_p，colon(1，1，nWidth)；/将矩阵由unsigned cha
48、r型转换为double型，以便进行/运算m_matBits=mdouble(m_matBits)；第7章 频域处理7.7.3 将矩阵数据赋给图像数据区为使数据符合图像文件格式，还应将处理后的矩阵数据赋给图像数据区，其代码如下：/检查m_matBits矩阵是否为空if(isempty(m_matBits)=1)return FALSE；第7章 频域处理/将矩阵数据范围限定于0255之间m_matBits=mabs(m_matBits)；Mm L=m_matBits 255.0；Mm NotL=!L；m_matBits=times(m_matBits，NotL)+L*255.0；/将矩阵由doub
49、le型转换为unsigned char型，以便将/图像数据赋给矩阵m_matBits=muint8(m_matBits)；第7章 频域处理/计算nWidth/4的余数int nTmp=(int)rem(nWidth，4)；int nPadding；/为矩阵各列补0，以使其为4的整数倍if(nTmp 0)nPadding=4-nTmp；m_matBits=cat(2，m_matBits，第7章 频域处理 repmat(muint8(0)，(BR(size(m_matBits，1)，nPadding)；/将矩阵顺时针旋转90度m_matBits=rot90(m_matBits，-1)；/将图像数据
50、赋给矩阵memcpy(pBits，m_matBits.addr()，(nWidthBytes*nHeight)*sizeof(unsigned char)；第7章 频域处理7.7.4 利用矩阵运算进行图像变换建立图像数据矩阵后，便可利用MatrixLIBC+的库函数，通过各种矩阵运算完成图像变换处理。1.离散傅立叶变换MatrixLIBC+的库提供了dft()、fft()和fft2()函数，分别用于完成一维DFT、一维FFT和二维FFT运算，函数ifft()和ifft2()用于完成一维FFT和二维FFT的逆傅立叶变换。在图像处理中常用的是fft2()和ifft2()，利用这两个函数便可完成图像
