101210113-宋浪

1、2010级物理学专业 光电图像处理课程设计 光电图像课程设计报告书课题名称基于matlab的滤波器设计姓 名宋浪学 号101210113院 部通信与电子工程学院专 业物理学指导教师熊翠秀 2013年 06 月 07日一、设计任务及要求：设计任务：如何用matlab来实现滤波。要 求： 1. 选图像，并输入；2. 采用butter命令，设计一个巴特沃斯低通滤波器，以此实现信号的滤波；3. 在同一面上显示滤波前后的所有的图像和程序；4. 使用MATLAB分别实现IIR滤波器和FIR滤波器； 指导教师签名： 年 月 日 二、指导教师评语：指导教师签名： 年 月 日 三、成绩 验收盖章 年 月 日 图

2、像复原的MATLAB实现1课程设计目的（1）了解基于matlab的滤波器处理及其基本操作；（2）学习MATLAB在滤波器中的使用；（3）提高学习与解决问题的能力。2课程设计基本内容2.1滤波器的基本原理2.1.1设计数字滤波器的任务就是寻求一个因果稳定的线性时不变系统，并使系统函数H（z）具有指定的频率特性。 数字滤波器从实现的网络结构或者从单位冲激响应分类，可以分成无限单位冲激响应（IIR）数字滤波器和有限长单位冲激响应（FIR）数字滤波器。 数字滤波器频率响应的三个要素：（1） 幅度平方响应（2） 相位响应（3） 群时延响应 IIR数字滤波器： IIR数字滤波器的系统函数为有理分数，即II

3、R数字滤波器的逼近问题就是求解滤波器的系数和，使得在规定的物理意义上逼近所要求的特性的问题。如果是在s平面上逼近，就得到模拟滤波器，如果是在z平面上逼近，则得到数字滤波器。 FIR数字滤波器： 设FIR的单位脉冲响应h（n）为实数，长度为N，则其z变换和频率响应分别为按频域采样定理FIR数字滤波器的传输函数H(z)和单位脉冲响应h（z）可由它的N歌频域采值H(k)唯一确定。 MATLAB中提供了几个函数，分别用于实现IIR滤波器和FIR滤波器。 （1） 卷积函数conv，调用格式为，c=conv（a，b） 该格式可以计算两向量a和b的卷积，可以直接用于对有限长信号采用FIR滤波器的滤波。（2）

4、 函数filter的调用格式为，y=filter（b，a，x） 该格式采用数字滤波器对数据进行滤波，既可以用于IIR滤波器，也可以用于FIR滤波器。其中向量b和a分别表示系统函数的分子，分母多项式的系数，若a=1，此时表示FIR滤波器，否则就是IIR滤波器。该函数就是利用给出的向量b和a，对x中的数据进行滤波，结果放入向量y。（3） 函数fftfilt的调用格式为，y=fftfilt（b，x） 该格式是利用基于FFT的重叠相加法对数据进行滤波，这种频域滤波技术只对FIR滤波器有效。该函数是通过向量b描述的滤波器对x数据进行滤波。关于用butter函数求系统函数分子与分母的几种形式。 b，a=b

5、utter（N,wc，high）：设计N阶高通滤波器，wc为它的3dB边缘频率，以Hz为单位。 b，a=butter（N,wc）：当wc为具有两个元素的矢量wc=w1，w2时，它设计2N阶带通滤波器，3dB通带，以Hz为单位。b，a=butter（N,wc，stop）：若wc=w1，w2，则它设计2N阶带阻滤波器，3dB通带，以Hz为单位。 如果在这个函数输入变元的最后，加一个变元“s”，设计的是模拟滤波器。这里不作讨论。 为了设计任意的选项巴特沃斯滤波器，必须知道阶数N和3dB边缘频率矢量wc。这可以直接利用信号处理工具箱中的buttord函数来计算。如果已知滤波器指标，和，则调用格式为 N

6、,wc=buttord（wp，ws，Rp，As） 对于不同类型的滤波器，参数wp和ws有一些限制：对于低通滤波器，wp<ws;对于高通滤波器，wp>ws;对于带通滤波器，wp和ws分别为具有两个元素的矢量，wp=wp1,wp2和ws=ws1,ws2，并且ws1<wp1<wp2<ws2;对于带阻滤波器wp1<ws1<ws2<wp2。 2.2滤波器与MATLAB简介（1） 滤波器示例：在这里为了说明如何用matlab来实现滤波，特举出一个简单的函数信号滤波实例(对信号x(n)=sin(n/4)+5cos(n/2)进行滤波，信号长度为500点)，从中了

7、解滤波的实现过程。程序如下：Wn=0.2*pi;N=5;b,a=butter(N,Wn/pi);N=0:1:499;x=sin(pi*n/4)+5*cos(pi*n/2);X=fft(x,4096);Subplot(2,2,1);plot(x);title(滤波前信号的波形);Subplot(2,2,2);plot(X);title(滤波前信号的频谱);y=filter(b,a,x);Y=fft(y,4096);Subplot(2,2,3);plot(y);title(滤波后信号的波形);Subplot(2,2,4);plot(Y);title(滤波后信号的频谱);结果如图： 在这里，是采用了

8、butter命令，设计出一个巴特沃斯低通滤波器，从频谱图中可以很明显的看出来。下面，也就是本课题的主要内容，也都是运用到了butter函数，以便容易的得到系统函数的分子与分母系数，最终以此来实现信号的滤波。（2） N阶高通滤波器的设计（在这里，以5阶为例，其中wc为其3dB边缘频率，以Hz为但位），程序设计如下：x=wavread(ding.wav);Sound(x);N=5;wc=0.3;b,a=butter(N,wc,high);X=fft(x);Subplot(3,2,1);plot(x);title(滤波前信号的波形);Subplot(3,2,2);plot(X);title(滤波前信

9、号的频谱);y=filter(b,a,x);Y=fft(y);Subplot(3,2,3);plot(y);title(IIR滤波后信号的波形);Subplot(3,2,4);plot(Y);title(IIR滤波后信号的频谱);z=fftfilt(b,x);Z=fft(z);Subplot(3,2,5);plot(z);title(FIR滤波后信号的波形);Subplot(3,2,6);plot(Z);title(FIR滤波后信号的频谱);得到结果如图：（3）2N阶带通滤波器的设计（在这里，以10阶为例，其中wc为其3dB边缘频率，以Hz为单位，在wc=w1,w2,w1 wc w2）,程序设

10、计如下： x=wavread(ding.wav);Sound(x);N=5;wc=0.3,0.6;b,a=butter(N,wc);X=fft(x);Subplot(3,2,1);plot(x);title(滤波前信号的波形);Subplot(3,2,2);plot(X);title(滤波前信号的频谱);y=filter(b,a,x);Y=fft(y);Subplot(3,2,2);plot(X);title(滤波前信号的频谱);y=filter(b,a,x);Y=fft(y);Subplot(3,2,3);plot(y);title(IIR滤波后信号的波形);Subplot(3,2,4);p

11、lot(Y);title(IIR滤波后信号的频谱);z=fftfilt(b,x);Z=fft(z);Subplot(3,2,5);plot(z);title(FIR滤波后信号的波形);Subplot(3,2,6);plot(Z);title(FIR滤波后信号的频谱);得到的结果如图：（4）关于MATLAB MATLAB对于技术计算来说是一种高性能的语言。它以易于应用的环境集成了计算、可视化和编程，在该环境下，问题及其解以我们熟悉的数学表示法来表示。典型的应用包括如下方面：数学和计算/算法开发/数据获取/建模、模拟和原型设计/数据分析、研究和可视化/科学和工程图形/应用开发，包括图像用户界面构建

12、。 图2约束最小二乘滤波恢复对比图 如图2所示，利用振铃抑制恢复图像是3几种中恢复效果最好的，其他几种方法也可以恢复但是比较模糊，效果不是很明显。2.1.3使用Lucy-Richardson算法的迭代非线性复原 Lucy-Richardson(LR)算法假设图像服从Poission 分布，采用最大似然法进行估计，是一种基于贝叶斯分析的迭代算法。Lucy-Richardson算法是目前世界上应用最广泛的函数恢复技术之一，它是一种迭代非线性复原算法，能够按照泊松噪声统计标准求出与给定的点扩散函数PSF (Point Spread Function)卷积后最有可能成为输入模糊图像的图像。当PSF已知

13、但图像噪声信息未知时，也可以使用这个函数进行有效的恢复。其最优估计以最大似然准则作为标准，即要使概率密度函数最大，推导出的迭代式为： (1) 其中，和分别为卷积运算和相关运算；为迭代次数，可以令进行迭代，可以证明，当噪声可以忽略当不断增大时会依概率收敛于，从而恢复出原始图像。 当噪声不可忽略时，把式(1.2)代入式(7.1)可得到： (2) 从上式可看出，若噪声不可忽略，则以上过程的收敛性将难以保证，即 LR 存在放大噪声的缺陷。因此，处理噪声项是 LR 算法应用于低信噪比图像复原的关键。 是未退化图像的估计，这个算法的迭代本质是显而易见的。它的非线性本质是载方程的右边用来除产生的。就像大多数

14、非线性方法一样，关于什么时候停止L-R算法通常是很难回答。通常，接下来的处理途径是对于给定的应用，在我们获得满意的结果时，观察输出并终止算法。在IPT中，L-R算法是由名为deconvlucy的函数完成的，该函数的语法为 fr=deconvlucy(g,PSF,NUMIT,DAMPAR,WEIGHT)其中，fr代表复原的图像，g代表退化的图像，PSF是点扩散函数，NUMIT为迭代次数（默认为10次），DAMPAR和WEIGHT定义如下。DAMPAR是一个标量，它指定了结果图像与原图像g之间的偏离阈值。当像素偏离原值的范围在DAMPAR之内时，就不用再迭代。这既抑制了这些像素上的噪声，又保存了必

15、要的图像细节。默认值为0（无衰减）。WEIGHT是一个与g同样大小的数组，它为每一个像素分配一个权重来反映其重量。当用一个指定的PSF来模拟模糊时，WEIGHT可以从计算像素中剔除那些来自图像边界的像素点，因此，PSF造成的模糊是不同的。若PSF的大小为，则在WEIGHT中用到的零边界的宽度是ceil（n/2）。默认值是同输入图像g同等大小的一个单位数组。若复原图像呈现出由算法中所用的离散傅里叶变换所引入的振铃，则在调用函数deconvlucy之前，要利用函数edgetaper。下面通过程序事例说明 Lucy-Richardson算法恢复图像：I=imread('G:1.jpg'

16、;);PSF=fspecial('gaussian',5,5) ;Blurred=imfilter(I,PSF,'symmetric','conv');V=.003;BN=imnoise(Blurred,'gaussian',0,V);luc=deconvlucy(BN,PSF,5);figuresubplot(2,2,1);imshow(I);title('原始图像');subplot(2,2,2);imshow (Blurred);title('模糊后的图像');subplot(2,2,3);i

17、mshow (BN);title('加噪后的图像');subplot(2,2,4);imshow (luc);title('恢复后的图像'); 图3 Lucy-Richardson算法恢复图像对比图 用Lucy-Richardson算法可以较好的恢复图像9。但由于迭代产生的噪声痕迹是最大化可能性数据逼近法的常见问题，在低信噪比条件下，恢复图像可能会出现一些斑点，这些斑点并不代表图像的真实结构，只不过是恢复图像过于逼近噪声所产生的结果。另外本方法存在一些较严重的缺陷，问题一是噪声放大问题，这也是这类方法（利用最大似然性求解）的通病。问题二是对于恢复图像中的不同部分

18、，分别执行多少迭代才合适的问题。因为图像中信噪比高的部分可能需要数百次迭代才能获得满意的结果；而另一些光滑的对象可能只需很少次数即可达到满意的结果，所以适当选择迭代次数对图像恢复也很重要。这两个问题若得不到解决，将会对最终结果产生不利影响从图3我们可以看到，经多次迭代，尤其是在低信噪比情况下，重建图像可能会出现一些斑点，这些斑点并不代表图像的真实结构，是输出图像过于逼近噪声所产生的结果。用Lucy-Richardson算法可以较好的恢复图像。但由于迭代产生的噪声痕迹是最大化可能性数据逼近法的常见问题，在低信噪比条件下，恢复图像可能会出现一些斑点，这些斑点并不代表图像的真实结构，只不过是恢复图像过于逼近噪声所产生的结总 结本文介绍了图像退化的原因并且简要介绍了当前主流的图像复原方法，并通过对各种复原方法的仿真，了解了各种方法的优劣性，为我们在实际生活提供依据。但是无论是哪一种方法都有所局限性，我们应该努力致力于研究新型的优秀的图像复原方法，来获得更好的图像复原效果。同时，我们知道总而言之，图像恢复的算法很多，无论何种算