频率域滤波的MATLAB设计与实现课程设计.docx
《频率域滤波的MATLAB设计与实现课程设计.docx》由会员分享，可在线阅读，更多相关《频率域滤波的MATLAB设计与实现课程设计.docx（31页珍藏版）》请在冰点文库上搜索。
频率域滤波的MATLAB设计与实现课程设计
综合课程设计
设计题目频率域滤波地MATLAB设计与实现
专业名称
班级学号
学生姓名
指导教师
设计时间
摘要
图像处理主要应用于对图像视觉效果地改善，如去噪处理、图像增强、几何校正等.而本次主要是进行频率域增强技术，它是增强技术地重要组成部分.本次设计重点用MATLAB对一幅图像作Fourier变换，然后对变换后得到地频谱图像进行分别进行理想低通滤波、二阶Butterworth低通滤波和高斯低通滤波，将原图和变换图放在一起进行比较，观察其图像，对其结果进行分谚辞調担鈧谄动禪泻類。
关键词：
Matlab，Fourier变换，Butterworth，高斯（Gaussian），低通滤波嘰觐詿缧铴嗫偽純铪锩。
1.数字图像处理
1.1发展简况：
数字图像处理最早出现于20世纪50年代，当时地电子计算机已经发展到一定水平，人们开始利用计算机来处理图形和图像信息.数字图像处理作为一门学科大约形成于20世纪60年代初期.早期地图像处理地目地是改善图像地质量，它以人为对象，以改善人地视觉效果为目地.图像处理中，输入地是质量低地图像，输出地是改善质量后地图像，常用地图像处理方法有图像增强、复原、编码、压缩等.首次获得实际成功应用地是美国喷气推进实验室（JPL）.他们对航天探测器徘徊者7号在1964年发回地几千张月球照片使用了图像处理技术，如几何校正、灰度变换、去除噪声等方法进行处理，并考虑了太阳位置和月球环境地影响，由计算机成功地绘制出月球表面地图，获得了巨大地成功.随后又对探测飞船发回地近十万张照片进行更为复杂地图像处理，以致获得了月球地地形图、彩色图及全景镶嵌图，获得了非凡地成果，为人类登月创举奠定了坚实地基础，也推动了数字图像处理这门学科地诞生.熒绐譏钲鏌觶鷹緇機库。
1.2关键技术：
1）图像变换：
由于图像阵列很大，直接在空间域中进行处理，涉及计算量很大.因此，往往采用各种图像变换地方法，如傅立叶变换、沃尔什变换、离散余弦变换等间接处理技术，将空间域地处理转换为变换域处理，不仅可减少计算量，而且可获得更有效地处理（如傅立叶变换可在频域中进行数字滤波处理）.目前新兴研究地小波变换在时域和频域中都具有良好地局部化特性，它在图像处理中也有着广泛而有效地应用.鶼渍螻偉阅劍鲰腎邏蘞。
2）图像编码压缩：
图像编码压缩技术可减少描述图像地数据量（即比特数），以便节省图像传输、处理时间和减少所占用地存储器容量.压缩可以在不失真地前提下获得，也可以在允许地失真条件下进行.编码是压缩技术中最重要地方法，它在图像处理技术中是发展最早且比较成熟地技术.纣忧蔣氳頑莶驅藥悯骛。
3）图像增强和复原：
按照特定地需要突出一幅图像中地某些信息或强化某些感兴趣地特征，将原来不清晰地图片变得清晰，使之改善图像质量和丰富信息量，提高图像地视觉效果和图像成分地清晰度，加强图像判读和识别效果地图像处理地方法.图像增强不考虑图像降质地原因，突出图像中所感兴趣地部分：
如强化图像高频分量，可使图像中物体轮廓清晰，细节明显；如强化低频分量可减少图像中噪声影响[1].颖刍莖蛺饽亿顿裊赔泷。
4）图像分割：
是将图像分为若干个特定地、具有独特性质地区域，其中每一个区域都是像素地一个连续集合.它是图像处理到图像分析地关键步骤.常用地分割方法主要分一下几类：
基于阀值地分割方法、基于区域地分割方法、基于边缘地分割方法.虽然近年来提出了很多新地分割方法，但并没有一种适用于所有图像地分割方法.在实际应用中，通常将多种分割算法有效地结合在一起使用以获得更好地图像分割效果[2].濫驂膽閉驟羥闈詔寢賻。
5）图像描述：
将图像分割为区域后，接下来通常要将分割区域加以表示与描述，以方便计算机处理.图像描述也是图像识别地必要前提.作为最简单地二值图像可采用其几何特性描述物体地特性，一般图像地描述方法采用二维形状描述，它有边界描述和区域描述两类方法[1].对于特殊地纹理图像可采用二维纹理特征描述.随着图像处理研究地深入发展，已经开始进行三维物体描述地研究，提出了体积描述、表面描述，广义圆柱体描述等方法[4].銚銻縵哜鳗鸿锓謎諏涼。
6）灰度变换：
其运算简单，可以达到增强对比度去除噪声地效果.但是该方法对于原图中所有灰度级地变换是相同地.然而在实际情况中，人们更关心图像中地目标，对背景不太关注，希望对目标地灰度级进行增强，对背景地灰度级进行压缩.这样，不但能提高目标地对比度，还可以更清楚地显示目标内部地细节变化，并且忽略了人们不关心地背景地部分细节.即使原灰度级地范围较大，该方法也可以得到满意地效果[3].挤貼綬电麥结鈺贖哓类。
7）空间滤波：
一种采用滤波处理地影响增强方法.其理论基础是空间卷积.目地是改善影像质量，包括出去高频噪声与干扰，及影像边缘增强、线性增强以及去模糊等.赔荊紳谘侖驟辽輩袜錈。
2．频率域滤波地产生背景及意义
2.1傅立叶级数和变换简史：
法国数学家傅立叶指出任何周期函数都可以表示为不同频率地正弦和或余弦之和地形式，每个正弦项和/或余弦乘以不同地系数（现在称该和为傅立叶级数）.无论函数多么复杂，只要它是周期地，并且满足某些适度地数学条件，都可以用这样地和来表示.我们现在认为这是理所当然地，但在当时，这个概念第一次出现之后，一个复杂函数可以表示为简单地正弦和余弦之和地概念一点也不直观，所以傅立叶思想遭到怀疑是不足为奇地.塤礙籟馐决穩賽釙冊庫。
甚至非周期函数也可用正弦和/或余弦乘以加权函数地积分来表示.在这种情况下地公式就是傅立叶变换，其作用在多数理论和应用学科中甚至远大于傅立叶级数.裊樣祕廬廂颤谚鍘羋蔺。
用傅立叶级数或变换表示地函数特征完全可以通过傅立叶反变换来重建，而不会丢失任何信息.这是这种表示方法地最重要特征之一，因为它可以使我们工作于“傅立叶域”，而且在返回到函数地原始域时不会丢失任何信息.总之，傅立叶级数和变换是解决实际问题地工具，它作为基础工具被广泛地学习和使用.仓嫗盤紲嘱珑詁鍬齊驁。
傅立叶概念地最初应用是在热扩散领域，在该领域，人们考虑用微分方程来表示热流动，并且使用这种方法第一次获得了结论.在过去一个世纪，特别是后50年，傅立叶地思想使整个工业和学术界都空前繁荣.早在20世纪60年代，数字计算地出现和快速傅立叶变换算法地“发现”在信号处理领域产生了巨大变革.这两种核心技术第一次允许人们对医学监视器和扫描仪到现代电子通信地异常重要地信号进行实际处理.绽萬璉轆娛閬蛏鬮绾瀧。
2.2频率域滤波地意义：
滤波法对于模糊图像地复原，边缘地强化和噪声地去除方面都有明显地效果.熟练地应用不同地滤波图像处理方法对试听资料地检验有很大地帮助，达到改善图像质量地目地.图像地平滑主要应用在去除图像地噪声上；图像地锐化则可以增强图像地细节和边缘；而图像地复原则可对运动模糊图像，离焦模糊图像等进行复原.骁顾燁鶚巯瀆蕪領鲡赙。
近年来，计算机技术地迅速发展，为数字图像处理提供了强大地软件和硬件支持，也促进了数字图像处理理论和方法地不断深入，使其迅速成为一门新兴技术，并被广泛应用于各个领域.滤波法新理论地不断提出，新算法地不断应用都使这门技术不断前进成为图像处理技术地重要一员.瑣钋濺暧惲锟缟馭篩凉。
3.频率域滤波地常用方法
频率域处理法[4]是在图像地某种变换域内，对变换域地系数进行运算，然后在反变换到原来地空域得到增强地图像，这是一种间接处理方法.比如，先对图像进行傅里叶变换，再对图像地频域进行滤波处理，最后将频域处理后地图像变换值反变换到空间域，从而得到增强后地图像[5].频率域数字图像滤波是通过不同地滤波器在频率域对图像进行处理地方法[6].它可分为频率域平滑，频率域锐化和同态滤波增强.对于一幅图像来说，高频部分大致对应图像中地边缘细节，低频部分大致对应着图像中过渡比较平缓地部分.每一类滤波法根据具体地算法不同都可以分为很多种，如中值滤波，现在就推广出许多种中值滤波器，如加权中值滤波器（weightedmedianfilter）,多级中值滤波器（multistagemedianfilerner）,组合滤波器（hybridfiler或称L1滤波器），堆滤波器（stackfilter）和置换滤波器（permutationfilter）等.鎦诗涇艳损楼紲鯗餳類。
3.1低通滤波
◆理想低通滤波
在以原点为圆心、以D0为半径地圆内，无衰减地通过所有频率，而在该圆外“切断”所有频率地二维低通滤泼器，称为理想低通滤波器（ILPF）；它由一下函数决定：
栉缏歐锄棗鈕种鵑瑶锬。
其中D0是一个正常数，D（u，v）是频率域中点（u，v）与频率矩形中心地距离：
D（u，v）=[（u-P/2）2+（v-Q/2）2]1/2
3.1.1理想低通滤波器地截面图
如下图3.1所示；
图3.1理想低通滤波器截面图
理想低通滤波器具有一个特性——振铃效果，效果图如下3.2所示
理想低通滤波后图
原图
图3.2理想低通滤波器地振铃效果
◆巴特沃斯低通滤波
物理上可实现（理想低通滤波器在数学上定义得很清楚，在计算机模拟中也可实现，但在截断频率处直上直下地理想低通滤波器是不能用实际地电子器件实现地）减少振铃效应，高低频率间地过渡比较光滑，辔烨棟剛殓攬瑤丽阄应。
n阶Butterworth低通滤波器地传递函数为：
其中，D0为截止频率.
=1时，
=0.5，它地特性是传递函数比较平滑，连续衰减，而不像理想滤波器那样陡峭变化，即明显地不连续.因此采用该滤波器滤波在抑制噪声地同时，图像边缘地模糊程度大大减小，没有振铃效应产生，滤波效果如图3.3所示.峴扬斕滾澗辐滠兴渙藺。
图3.3Butterworth低通滤波效果
巴特沃斯滤波器中阶数对振铃现象地影响：
阶数越高，越明显，如下图3.4所示：
图3.4巴特沃斯滤波器阶数对振铃现象地影响
◆高斯低通滤波
高斯（Gaussian）低通滤波器地传递函数为：
其中，s为标准偏差.令s＝D0,我们可以根据截止参数D0得到表达式：
当D（u,v）=D0时，滤波器H（u,v）由最大值1下降为0.607.GLPF没有振铃现象，但与阶数为2地BLPF相比，其通带要宽些，这样对应地空间滤波器地灰度级轮廓更窄些，因而平滑效果要差些.詩叁撻訥烬忧毀厉鋨骜。
对于巴特沃斯低通滤波器和高斯低通滤波器，振铃现象从严重到无，但平滑效果从好到差，BLPF可以看成ILPF和GLPF地过渡，阶为1时与GLPF差不多，阶数越高越接近BPLG.则鯤愜韋瘓賈晖园栋泷。
如下图3.5表示出了高斯低通滤波器对于不同D0值地滤波效果；
图3.5高斯（Gaussian）低通滤波器对于不同地D0值地滤波效果
3.2高通滤波
◆理想高通滤波
一个理想高通滤波器（IHPF）定义为：
其中，D0是截止频率，D（u，v）由下面公式给出：
D（u，v）=[（u-P/2）2+（v-Q/2）2]1/2
如同ILPF一样，IHPF在物理上也是无法实现地，但是IHPF可以用于解释空间域地振铃等现象.
下图3.6即为理想高通滤波器地滤波效果：
图3.6理想高通滤波器地滤波效果
◆巴特沃斯高通滤波
巴特沃斯n阶截止频率为DO地巴特沃斯高通滤波器地传递函数为：
其中：
D（u，v）=[（u-P/2）2+（v-Q/2）2]1/2
由频域滤波模型Q（U，v）=F（U，v）H（U，v）知，F（U，v）中地低频（小于D0）成分，因乘上一个远小于1地H（U，v）值而被衰减.而高频成分却被乘以一个接近于1地H（U，v）值而保留，这即是所谓地高通滤波地原理.当截止频率D0越大，滤掉地低频成分越多，同样损失地高频成份也越多.胀鏝彈奥秘孫戶孪钇賻。
如下图3.7即为巴特沃斯高通滤波器地滤波效果：
图3.7n阶巴特沃斯高通滤波器滤波效果
◆高斯高通滤波
其中，s为标准偏差.通过令s＝D0,我们可以根据截止参数D0得到s地值.
下图3.8即为高斯高通滤波器地滤波效果：
图3.8高斯高通滤波器滤波效果
3.3带阻滤波
带阻滤波器阻止一定频率范围内地信号通过而允许其它频率范围内地信号通过.
◆理想带阻滤波器
理想带阻滤波器地传递函数：
这里，W是频带地宽度，D0是频带地中心半径.
◆巴特沃斯带阻滤波器：
n阶地巴特沃思带阻滤波器地表达式为：
◆高斯带阻滤波器
高斯带阻滤波器地表达式为：
下图3.9是理想带阻滤波器、阶数为1地巴特沃斯带阻滤波器和高斯带阻滤波器地透视图：
图3.9理想滤波器、巴特沃思滤波器（阶数为1）和高斯带阻滤波器地透视图
3.4带通滤波
带通滤波器执行与带阻滤波器相反地操作，带通滤波器地传递函数据相应地带阻滤波器地传递函数并应用下式得到地：
鳃躋峽祷紉诵帮废掃減。
4．原理及实现
4.1频率域增强基本理论
不对Fourier变换（FT）和图像地频率域处理技术有所了解，就不可能完全理解图像增强这个最基本地图像处理任务.稟虛嬪赈维哜妝扩踴粜。
频域增强指在图像地频率域内，对图像地变换系数（频率成分）直接进行运算，然后通过Fourier逆变换以获得图像地增强效果.陽簍埡鲑罷規呜旧岿錟。
一般来说，图像地边缘和噪声对应Fourier变换中地高频部分，所以低通滤波能够平滑图像、去除噪声.
图像灰度发生聚变地部分与频谱地高频分量对应，所以采用高频滤波器衰减或抑制低频分量，能够对图像进行锐化处理.沩氣嘮戇苌鑿鑿槠谔應。
卷积理论是频域技术地基础，设函数f（x,y）与算子h（x,y）地卷积结果是g（x, y），即g（x,y）=h（x,y）*f（x,y），那么根据卷积定理在频域有：
钡嵐縣緱虜荣产涛團蔺。
其中G（u,v），H（u,v），F（u,v）分别是g（x,y），h（x,y），f（x,y）地傅立叶（或其它）变换，H（u,v）是转移函数.在具体增强应用中，f（x,y）是给定地（所以F（u,v）可利用变换得到），需要确定地是H（u,v），这样具有所需特性地g（x,y）就可算出G（u,v）而得到：
懨俠劑鈍触乐鹇烬觶騮。
4.2傅立叶变换
傅里叶变换是将时域信号分解为不同频率地正弦信号或余弦函数叠加之和.傅立叶变换是数字图像处理技术地基础，其通过在时空域和频率域来回切换图像，对图像地信息特征进行提取和分析，简化了计算工作量，被喻为描述图像信息地第二种语言，广泛应用于图像变换，图像编码与压缩，图像分割，图像重建等.因此，对涉及数字图像处理地工作者，深入研究和掌握傅立叶变换及其扩展形式地特性，是很有价值得.把傅立叶变换地理论通其物理解释相结合，将有助于解决大多数图像处理问题.傅里叶变换可分为连续傅里叶变换、离散傅里叶变换、快速傅里叶变换.謾饱兗争詣繚鮐癞别瀘。
快速傅里叶变换（FFT）是计算离散傅里叶变换（DFT）地快速算法.离散傅里叶变换运算量巨大，计算时间长，即运算时间很长.而快速傅里叶变换地提出将傅里叶变换地复杂度由降到了，很大程度上减少了计算量.呙铉們欤谦鸪饺竞荡赚。
令
，
，u=0，1，2，…，M-1
则
，
4.3频率域理想低通（ILPF）滤波器
一个二维地理想低通滤波器（ILPF）地转换函数满足（是一个分段函数）
其中：
D0为截止频率D（u,v）为距离函数D（u,v）=（u2+v2）1/2
4.3.1理想低通滤波器地截面图
理想低通滤波器截面图如图4.1所示：
图4.1理想低通滤波器截面图
振铃效果——理想低通滤波器地特性
理想低通滤波后图
原图
图4.2理想低通滤波器振铃效果
4.3频率域巴特沃兹（Butterworth）低通滤波器
物理上可实现（理想低通滤波器在数学上定义得很清楚，在计算机模拟中也可实现，但在截断频率处直上直下地理想低通滤波器是不能用实际地电子器件实现地）减少振铃效应，高低频率间地过渡比较光滑n阶Butterworth低通滤波器地传递函数为：
莹谐龌蕲賞组靄绉嚴减。
式中，D0为截止频率.
=1时，
=0.5，它地特性是传递函数比较平滑，连续衰减，而不像理想滤波器那样陡峭变化，即明显地不连续.因此采用该滤波器滤波在抑制噪声地同时，图像边缘地模糊程度大大减小，没有振铃效应产生，如图2.1所示.麸肃鹏镟轿騍镣缚縟糶。
图4.3Butterworth低通滤波效果
阶数对振铃现象地影响：
阶数越高，越明显，如图2.2所示
图4.4阶数对振铃现象地影响
4.4频率域高斯（Gaussian）低通滤波器
高斯（Gaussian）低通滤波器地传递函数为：
其中，
为标准偏差.通过令
＝D0,我们可以根据截止参数D0得到表达式：
当D（u,v）=D0时，滤波器H（u,v）由最大值1下降为0.607.
GLPF没有振铃现象，但与阶数为2地BLPF相比，其通带要宽些，这样对应地空间滤波器地灰度级轮廓更窄些，因而平滑效果要差些.納畴鳗吶鄖禎銣腻鰲锬。
以上两种滤波器，振铃现象从严重到无，但平滑效果从好到差，BLPF可以看成ILPF和GLPF地过渡，阶为1时与GLPF差不多，阶越高越接近BPLG..風撵鲔貓铁频钙蓟纠庙。
图4.5高斯（Gaussian）低通滤波器对于不同地D0值地滤波效果
5．程序设计
5.1算法设计（程序设计流程图）
算法流程图如图3.1所示：
图3.1程序设计流程图
5.2对灰度图像进行Fourier变换地程序
I=imread（'Testpicture.jpg'）。
f=double（I）。
%数据类型转换，MATLAB不支持图像地无符号整型地计算
g=fft2（f）。
figure,imshow（log（1+abs（fftshift（g）））,[]）%将频谱显示为一幅图像灭嗳骇諗鋅猎輛觏馊藹。
title（'原图Fourier变换后频谱图'）。
5.3频率域理想低通滤波器
I=imread（'Test1.jpg'）。
I=rgb2gray（I）。
Do=input（'请输入Do值:
'）。
subplot（221）,imshow（I）。
title（'原图像'）。
s=fftshift（fft2（I））。
subplot（223）,
imshow（abs（s）,[]）。
title（'图像傅里叶变换所得频谱'）。
subplot（224）,
imshow（log（abs（s））,[]）。
title（'图像傅里叶变换取对数所得频谱'）。
[a,b]=size（s）。
a0=round（a/2）。
b0=round（b/2）。
fori=1:
a
forj=1:
b
distance=sqrt（（i-a0）^2+（j-b0）^2）。
ifdistance<=Doh=1。
elseh=0。
end。
s（i,j）=h*s（i,j）。
end。
end。
s=uint8（real（ifft2（ifftshift（s））））。
subplot（222）,
imshow（s）。
title（sprintf（'%s%d','理想低通滤波结果Do=',Do））
5.4二阶巴特沃斯（Butterworth）低通滤波程序
此函数地功能为：
使用巴特沃斯低通滤波原理实现图像柔和平滑，为便比较，通过设置不同地d0值运行程序,如果对各个d0一起运行，运行速度奇慢，所以对各个d0单独运行.铹鸝饷飾镡閌赀诨癱骝。
I=imread（'a.jpg'）。
imshow（I）。
title（'原图'）。
s=fftshift（fft2（I））。
[M,N]=size（s）。
%分别返回s地行数到M中，列数到N中
n=2。
%对n赋初值
%BLPF滤波，d0=15，30，80
d0=input（'请输入d0值:
'）。
%初始化d0,d0依次取15，30，80
n1=floor（M/2）。
%对M/2进行取整
n2=floor（N/2）。
%对N/2进行取整
fori=1:
M
forj=1:
N
d=sqrt（（i-n1）^2+（j-n2）^2）。
%点（i,j）到傅立叶变换中心地距离
h=1/（1+（d/d0）^（2*n））。
%BLPF滤波函数
s（i,j）=h*s（i,j）。
%BLPF滤波后地频域表示
end
end
s=ifftshift（s）。
%对s进行反FFT移动
%对s进行二维反离散地Fourier变换后，取复数地实部转化为无符号8位整数
s=uint8（real（ifft2（s）））。
figure。
%创建图形图像对象
imshow（s）。
%显示BLPF滤波处理后地图像ttitle（sprintf（'%s%d%s','二阶Butterworth低通滤波结果（D0=',d0,'）'））攙閿频嵘陣澇諗谴隴泸。
5.5高斯（Gaussian）低通滤波程序
此函数地功能为：
使用高斯低通滤波原理实现图像越柔和平滑，为便比较，通过设置不同地d0值运行程序,如果对各个d0一起运行,运行速度奇慢，所以对各个d0单独运行.趕輾雏纨颗锊讨跃满賺。
%高斯低通滤波程序
ima=imread（'a.jpg'）。
ima=rgb2gray（ima）。
imshow（ima）。
title（'原图'）。
ima0=im2double（ima）。
Do=input（'请输入d0值:
'）。
%初始化d0,d0依次取15，30，80
[Row,Col]=size（ima）。
ima=im2double（ima）。
forix=1:
Row%滤波中心化
foriy=1:
Col
ima（ix,iy）=（-1）^（ix+iy）*ima（ix,iy）。
end
end
FourIma=fft2（ima）。
%中心变换后地傅立叶变换
[FRow,FCol]=size（FourIma）。
foru=1:
FRow%Do地高斯低通滤波器
forv=1:
FCol
D（u,v）=[（u-Row/2）^2+（v-Col/2）^2]^0.5。
factc=-D（u,v）^2/（2*Do^2）。
H（u,v）=exp（factc）。
G（u,v）=H（u,v）*FourIma（u,v）。
end
end
IFourIma=ifft2（G）。
%高斯低通滤波后地反傅立叶变换
forix=1:
FRow%反傅立叶变换后*（-1）^（x+y）并取实部
foriy=1:
FCol
LowPass（ix,iy）=（-1）^（ix+iy）*IFourIma（ix,iy）。
end
end
LowPass=real（LowPass）。
figure。
imshow（LowPass）。
%Do=15高斯低通滤波结果
title（sprintf（'%s%d%s','Do=',Do,'高斯低通滤波结果'））
6．结果与分析
我们还可以通过函数colorbar给一个坐标轴对象添加加一条色带，该色带可以指示出该图像中不同颜色地数据值，以显示出图像地灰度级.夹覡闾辁駁档驀迁锬減。
输入以下代码：
f=imread（'Test1.jpg'）。
imshow（f）。
colorbar。
可得到添加色带后地图像如图6.1所示，从图6.3可知该图像是数据类型为uint
升级会员 