Matlab和双线性变换的滤波器设计

1、题 目：MATLAB课程设计一基于Matlab和双线性变换的滤波器设计初始条件：MATLA仿真软件数字信号处理与图像处理基础知识双线性变换法的原理和算法巴特沃斯高通滤波器的性能指标要求完成的主要任务：利用MATLAB仿真软件系统结合双线性变换法设计一个数字巴特沃斯高通、低通和带 通IIR滤波器。说明：自己确定滤波器的不同参数得到不同的结果，最后自己分析原因。课程设计题目由指导教师提供，每人隶属一组完成任务，每组完成的内容不能雷同（按学号分组）时间安排：序 号阶段内容所需时间1方案设计1天2软件设计2天3系统调试1天4答辩1天合计5天指导教师签名：系主任（或责任教师）签名:目录摘要数字滤波器是数

2、字信号处理的基础，用来对信号进行过滤、检测与参数估计等处理， 在通信、图像、语音、雷达等许多领域都有着十分广泛的应用。尤其在图像处理、数据压 缩等方面取得了令人瞩目的进展和成就。鉴于此，数字滤波器的设计就显得尤为重要。此报告重点介绍了用双线性不变法设计 IIR 数字滤波器的基本流程，比较了各种设计 方法的优缺点，总结了模拟滤波器的性能特征。最后以双线性不变法设计了一个高通巴特 沃斯 FIR 数字滤波器，介绍了设计步骤，然后在 Matlab 环境下进行了仿真与调试，实现 了设计目标。关键字：数字滤波器 巴特沃斯 MATLAB 双线性变法AbstractMatlab for research an

3、d application of digital filter provides an intuitive, efficient and convenient weap on. He matrix operati ons based on the calculati on, visualizati on, program desig n in tegrated in to an in teractive worki ng en vir onment. In particular, Matlab toolbox to en able researchersin various fields ca

4、n be easy and intuitive application of scientific research and engin eeri ng.This report in troduced with emphasis of the basic flow of desig ning the IIR digit filter by the bili near political reform, compared with each kind of desig n method's good and bad poin ts, summarized an alog filter&#

5、39;s performa nce characteristic. Fin ally desig n one by the bili near political reform to pass Butterworth high IIR digit filter, introduced the design procedure, then has carried on the simulatio n and the debuggi ng un der the Matlab en vir onment, has achieved the project objective.Keyword: dig

6、ital filter Butterworth MATLAB1 数字滤波器1.1 数字滤波器的基本概念数字滤波器是数字信号处理的重要基础，是对信号都是过滤检测与参数估计等处理过 程中，它是使用最为广泛的一种线性系统。数字滤波器处理的对象是经由采样期间将模拟 信号转换而得到的数字信号。数字滤波器是指完成信号滤波处理功能的，用有限精度算法 实现的离散时间线性非时变系统。数字滤波器的输入是一组数字量。它本身既可以是用数 字硬件装配而成的一台用于完成给定运算的专用数字计算机，也可以是将所需的运算编写 的程序通过计算机来执行。数字滤波器具有稳定性高、精度高、灵活性大等突出优点。随 着数字技术的发展，用数

7、字技术实现滤波器的功能愈来愈受到人们的重视，并得到了广泛 的应用。数字滤波器的数学运算通常有两种实现方式。一种是频域法，即利用 FFT 快速运算方 法对输入信号进行离散傅里叶变换， 分析其频谱， 然后再根据所希望的频率特性进行滤波， 再利用傅里叶反变换得到时域信号。这种方法具有较好的频域选择特性和灵活性，并且由 于信号频率与所希望的频率特性是简单的相乘关系，所以它比计算等价的时域卷积要快得 多。另一种方法是时域法，这种方法是通过离散的抽样数据做差分数学运算来达到滤波目 的的。1.2 数字滤波器的分类数字滤波器按照不同的分类方法，可分为许多种，但总体来讲可以分成两大类。 一类称为经典滤波器，即一

8、般滤波器，特点是输入信号中有用的频率成分和希望滤除 的频率成分各占有不同的频带，通过一个合适的选频滤波器达到滤波的目的。但对于 一般滤波器，如果信号和干扰的频带互不重叠，则不能完成对干扰的有效滤除，这是 需要采用另一类所谓的现代滤波器，例如维纳滤波器、卡尔曼滤波器、自适应滤波器 等最佳滤波器，这些滤波器可按照随机信号内部的一些统计分布规律，从干扰中最佳 地提取信号。从功能上，一般数字滤波器可以分为低通、高通、带通、带阻和全通等，此种分 类方法是和模拟滤波器一样的1.3数字滤波器的MATLA实现MATLAE工具箱未滤波器的设计应用提供了丰富而简便的方法，如函数方法和图形工具 方法等，使原来非常繁

9、琐复杂的程序设计变成了简单的函数调用，为滤波器的设计和实现 开辟了广阔的天地。数字滤波器的一般设计过程为：按照实际需要，确定滤波器的性能要求； 用一个因果的、稳定的离散线性时不变系统，去逼近这一性能指标；用有限精度的运 算实现所设计的系统；通过模拟，验证所设计的系统是否符合给定性能要求。 Dd2双线性变换法2.1双线性变换法知识简介脉冲响应不变法的主要缺点是产生频率响应的混叠失真。这是因为从S平面到Z平面是多值的映射关系所造成的。为了克服这一缺点，可以采用非线性频率压缩方法，将整个 频率轴上的频率范围压缩到-n/Tn/T之间，再用z=e转换到Z平面上。也就是说，第 一步先将整个S平面压缩映射到

10、Si平面的-n/Tn/T 条横带里；第二步再通过标准变 换关系z=es1T将此横带变换到整个Z平面上去。这样就使S平面与Z平面建立了一一对应 的单值关系，消除了多值变换性，也就消除了频谱混叠现象。由图2-1看出，在零频率附近，模拟角频率Q与数字频率3之间的变换关系接近于线 性关系；但当Q进一步增加时，3增长得越来越慢，最后当QX时，3终止在折叠频率 3 =冗处，因而双线性变换就不会出现由于高频部分超过折叠频率而混淆到低频部分去的 现象，从而消除了频率混叠现象。2.2双线性变换法设计数字滤波器原理将S平面j 轴压缩变换到si平面j轴上的-/T到/T 一段，可以采用以下变换关系：tan iT2这样

11、变换到i -，0变换到i 0，可将上式写成j丄 j卫2 2. e ej.1Tj j-2 2 e e令j s，j 1 s,，解析延拓到整个s平面和s1平面，可得s,TqTr2c 2eesstst22ee再将$平面通过以下变换关系映射到Z平面，即S1Tz e从而得到s平面和z平面的单值映射关系为般来说，为了使模拟滤波器的某一频率与数字滤波器的任一频率有对应关系，可引入待定常数c，1c-1eT将z e5代入到上式，可得1c-1i z1 z在MATLAB，双线性Z变换可以通过bilinear 函数实现，其调用格式为:Bz，Az = bilinear（B ，A，Fs）;其中B，A为模拟滤波器传递函数 q

12、s）的分子分母多项式的系数向量，而 Bz，Az为数字滤 波器的传递函数Hz）的分子分母多项式的系数向量。3设计任务及方案选择3.1设计任务与要求利用MATLAB仿真软件系统结合双线性变换法设计一个数字巴特沃斯高通、低通 和带通FIR滤波器。说明：自己确定滤波器的不同参数得到不同的结果，最后自己分 析原因。假设一个信号：x(t) sin 2 f it0.5cos 2 f 2t其中f仁100Hz, f2=550Hz。试将该信号分别通过低通、高通、带通滤波器，对输出信号进行分析比较。Step1:Step2:Step3:Step4:3.2方案设计论证与选择本次滤波器设计的主要方法是先设计低通模拟滤波器

13、，然后转换为高通、带通或带阻 数字滤波器。对于其他如高通，带通，则通过频率变换转换为设计相应的高通，带通等。 在设计的全过程的各个步骤，matlab都提供相应的工具箱函数，使得IIR数字滤波器设计 变得非常简单。总的来说，我的设计思路主要有以下两种：方案一：从归一化模拟低通原型出发，先在模拟域内经频率变换成为所需类型的模拟 滤波器；然后进行双线性变换，由S域变换到Z域，而得到所需类型的数字滤波器。图3-1先频率变换再离散方案二：先进行双线性变换，将模拟低通原型滤波器变换成数字低通滤波器；然后在Z域内经数字频率变换为所需类型的数字滤波器图3-2先离散再频率变换以上两种方案都可以，我最后选择了第一

14、种方案进行设计，即先在模拟域内经频率变换成为所需类型的模拟滤波器；然后进行双线性变换，由S域变换到Z域，而得到所需类型的数字滤波器。因为高通，带通滤波器的设计方法是先将要设计的滤波器的技术特性指标 通过频率转换成模拟低通滤波器的技术指标，再根据这些性能指标设计出低通滤波器的传 递函数，为了设计程序的简洁，故选择了方案一。4程序设计与调试4.1设计要求若实验设计要求得巴特沃斯高通滤波器的性能指标为：通带边界频率fp: fp=400Hz;通带波纹Rp：Rp=1dB;阻带边界频率fs： fs=300Hz;阻带衰减Rs：Rs=30dB；采样频率Fs：Fs=2000Hz;要求巴特沃斯低通滤波器通带边界频

15、率与阻带边界频率与高通的正好相反模拟带通滤波器中心频率为 f0=200Hz，带宽为Bw=400Hz，4.2程序设计与说明模拟低通滤波首先应用双线性变换作为数字低通滤波器到模拟高通滤波器的频率预畸变处理:2 tanT 22p Twscot 222 fT利用上式确定的技术指标和已知的参数Rp、Rs，就可以进行模拟原型低通滤波器设计了。该部分的程序设计与注释如下：Fs=2000; Ts=1/Fs; wp=400*2*pi;ws=300*2*pi;Rp=1;Rs=30;wpl=Fs*2*ta n(ws*Ts/2);wsl=Fs*2*ta n(wp*Ts/2);N,Wn=buttord(wpl,wsl,

16、Rp,Rs,'s');%米样频率和米样周期%低通、高通数字滤波器截止频率%通带波纹和阻带衰减%按频率转换公式进行转换z,p,k=buttap(N);%计算模拟滤波器的最小阶数%设计模拟原型低通滤波器b,a=zp2tf(z,p,k);%变零极点增益形式为系统传递函数形式H,w=freqs(b,a); figure;%求模拟低通滤波器的频率响应plot(w,abs(H);grid; title('模拟低通原型');%绘制低通频率特性xlabel('w');ylabel(' 幅度');4.2.2 模拟到数字的转换(1)从低通向低通的转换

17、 b,a=lp2lp(Bap,Aap,Wn); 功能是把模拟滤波器原型转换成截至频率为 Wn的低通滤波器。再由得到的模拟低通滤波器通过函数 Bz, Az = bilinear(B ,A，Fs)得到数字低通滤波器。该部分的程序设计如下：e1,f1=lp2lp(b,a,300*2*pi) ;%原型低通转换为低通模拟滤波器h1,w1=freqs(e1,f1);figure;%求模拟低通滤波器的频率响应plot(w1/(2*pi),abs(h1); grid;title('模拟低通滤波器');%绘制高通频率特性xlabel('频率 Hz');ylabel('幅度

18、');r1,q1=bilinear(e1,f1,1000);%模拟低通滤波器转换为数字滤波器s1,w2=freqz(r1,q1); figure;%求数字低通滤波器的频率响应plot(w2/pi, 20*log10(abs(s1);grid;响应%绘制数字低通滤波器的频率title('数字低通滤波器');xlabel('w/pi');ylabel('幅度(dB)');(2)从低通向高通的转换 b,a=lp2hp(Bap,Aap,Wn); 功能是把模拟滤波器原型转换 成截至频率为 Wn 的模拟高通滤波器，再进一步转换成数字高通滤波器。程序设

19、计如下：e2,f2=lp2hp(b,a,300*2*pi) ;%原型低通模拟滤波器转换为高通模拟滤波器h2,w3=freqs(e2,f2);%求模拟高通滤波器的频率响应figure;plot(w3/(2*pi),abs(h2);grid;title(' 模拟高通滤波器 ');xlabel('频率 Hz');ylabel('幅度');r2,q2=bilinear(e2,f2,1000);s2,w4=freqz(r2,q2);figure;plot(w4/pi, 20*log10(abs(s2);grid;%绘制高通频率特性%模拟高通滤波器转换为数字

20、滤波器%求数字高通滤波器的频率响应%绘制数字高通滤波器的频率响应title('数字高通滤波器');xlabel('w/pi');ylabel('幅度(dB)');(3) 从低通向带通的转换 b,a=lp2bp(Bap,Aap,W0 ，Bw); 功能是把模拟滤波器原型 转换成中心频率为 W0带宽为Bw的模拟带通滤波器，在进一步转换成数字带通滤波器。 程序设计如下：W0=200*2*pi;Bw=200; e3,f3=lp2bp(b,a,W0,Bw) ;h3,w5=freqs(e3,f3); figure;plot(w5/(2*pi),abs(h3);

21、grid;title('模拟带通滤波器'); xlabel('频率 Hz');ylabel('幅度');r3,q3=bilinear(e3,f3,1000); s3,w6=freqz(r3,q3);figure;plot(w6/pi, 20*log10(abs(s3);grid; title('数字带通滤波器'); xlabel('w/pi');ylabel('幅度(dB)');%原型低通转换为带通模拟滤波器%求模拟带通滤波器的频率响应%绘制带通频率特性%模拟滤波器转换为数字滤波器%求数字带通滤波器

22、的频率响应%绘制数字带通滤波器的频率响应%时间序列%输入信号频率 %输入信号%绘制输入信号y=filter(bz,az,x)%对输入信号进行低通滤波%绘制输出信号filter 中参量改为 r2、q2 和 r3、q3 即可高通、带通滤波器，调用函数4.1.3 对输入信号进行滤波读取输入信号程序：n=0:180;dt=1/Fs;t=n*dt;f1=100;f2=550;x=sin(2*pi*f1*t)+cos(2*pi*f2*t)+sin(2*pi*1*t); figure;plot(t,x);grid;title('输入信号'); xlabel('时间 t');y

23、label('幅度');将输入信号分别通过得到的低通、 对输入信号进行滤波。y=filter(r1,q1,x);figure;plot(t,y);grid;title('低通滤波输出信号'); xlabel('时间 t');ylabel('幅度') 再进行高通、带通滤波只需将调用函数4.3 仿真结果与分析将设计好的程序在 matlab 软件下进行仿真，得到仿真波形仿真结果记录图4-1模拟低通原型幅频响应Figure gSH S恰目Fil s Edat View InsertDsskt盯f Window HslpFigure 2模拟

24、低通滤液器1 Si!：File Idi t Vi e* InEart TooLe Desktop If 1 ndo* Help 竺 Q 略殴Q紗®贬 0 - 数字低逋滤漩器100aw(豊 grson3-(a)低通模拟滤波器的频率特性曲线(b)巴特沃斯低通数字滤波器的频率特性曲线图4-2(a)高通模拟滤波器的频率特性曲线(b)巴特沃斯高通数字滤波器的频率特性曲线图4-3(b)巴特沃斯带通数字滤波器的频率特性曲线(a)带通模拟滤波器的频率特性曲线图4-4图4-5输入信号波形图4-6经过巴特沃斯低通滤波器滤波后的输出波形图4-7经过巴特沃斯高通滤波器滤波后的输出波形图4-8经过巴特沃斯高通

25、滤波器滤波后的输出波形432仿真结果分析由仿真结果可知，设计出的数字低通，高通，带通的幅频响应均符合各自的特性，从模拟滤波器特性曲线可以看出，所设计的高通滤波器在大于400Hz为通带，其衰减均小于1dB;小于300Hz为阻带，其衰减大于30dB，而设计的低通滤波器情况正好相反，符合本 次设计的要求。当滤波器输入100Hz和550Hz两种信号后，低通滤波器可以滤除频率为 550Hz的 高频信号，而低频信号100Hz的信号可以通过滤波器，验证了该低通滤波器的功能。 而高通滤波器可以滤除频率为100Hz的信号，而高频信号550Hz的信号可以通过滤波 器。程序调试中出现的问题带通滤波器的仿真结果虽然符

26、合其带通特性，但要求得到的性能指标存在较大误差， 由于时间紧迫，将在课后进一步寻找原因及解决方案。一个线性是不变系统的典型表达式有状态空间型、传递函数型、零极增益型、极点留 数型共四种，在开始的设计过程中我发现在把零极点增益型向传递函数型的转换过程中多 用了一步，即先把零极点增益型转换为状态空间型， 进行相应的频率变换和双线性变换后, 再由状态空间型转换为传递函数型，进而绘制出其频率响应图，但是经仿真分析后发现， 性能并不理想，所以并没有采用。5 心得体会本次课程设计是基于数字滤波器的 MATLAB 软件应用，在本学期中， 我在数字信 号处理课程中学习了关于本次课程设计的一些内容， 对于 MA

27、TLAB 软件在数字滤波器 的设计中的一些函数也有所认识，通过着一周的时间更是让自己在这个方面的学习更 进一步，使自己在这方面的学习深度有所加深。本次用双线性变换法对巴特沃斯低通，高通及带通滤波器的设计，要求我能够用 所学理论知识对巴特沃斯高通滤波器进行设计，并能够对设计的滤波器的频率特性进 行分析。在设计过程中，首先，我通过查找资料等途径，对关于所要设计的题目进行 更加全面深入的学习，认真分析其中的原理，弄明白设计的思路；其次，根据自己选 择的性能指标要求，结合自己通过学习掌握的知识和资料进行程序设计，根据自己提 出的方案利用软件实现；最后，运行程序，分析所得结果是否符合自己的设计要求， 对

28、结果进行分析，实现由理论到实践，再由实践结果验证理论的过程。总的来说，这次课程设计让我对 MATLAB有了更深刻的了解，对数字滤波器的设计流 程有了大致的了解，掌握了一些设计滤波器的基本方法，提高了理论用于实践的能力，掌 握了更多专业相关的使用知识与技能。同时，也暴露了我很多的不足，在以后的学习中， 将进一步发扬有点，克服缺点。参考文献1 袁润章，自蔓延高温合成技术研究进展 . 武汉: 武汉工业大学出版社 ,19942 Schacht E.Industrial polysaccharides. Amsterdam:Elsevier Science,19873 郭景坤.陶瓷材料的强化与增韧新途径

29、的探索 .无机材料学报 ,1998,13910;23-264 Spriggs G E.A history of fine grained hardmetal. Int J of Refractory Metal and Hard Material,1995,13:241-2555 王为民.二硼化钛陶瓷的自蔓延高温合成与加工 : 博士论文. 武汉工业大学 ,19886 吴中庆 , 张基田 .SiC 颗粒与液体铝的界石层 . 广州 :92 秋季中国材料科学研讨 会,C-MRS,1992.603-609.附录本次课设总程序如下：Fs=2000; Ts=1/Fs; wp=400*2*pi;ws=300

30、*2*pi;Rp=1;Rs=30; wpl=Fs*2*tan(ws*Ts/2); wsl=Fs*2*tan(wp*Ts/2);N,Wn=buttord(wpl,wsl,Rp,Rs,'s');z,p,k=buttap(N);b,a=zp2tf(z,p,k);H,w=freqs(b,a);figure;plot(w,abs(H);grid; title(' 模拟低通原型 '); xlabel('w');ylabel(' 幅度');e1,f1=lp2lp(b,a,300*2*pi) ;h1,w1=freqs(e1,f1);figure;

31、plot(w1/(2*pi),abs(h1);grid;title(' 模拟低通滤波器 '); xlabel('频率 Hz');ylabel('幅度');r1,q1=bilinear(e1,f1,1000); s1,w2=freqz(r1,q1);figure;plot(w2/pi, 20*log10(abs(s1);grid; 响应%采样频率和采样周期%低通、高通数字滤波器截止频率%通带波纹和阻带衰减%按频率转换公式进行转换%计算模拟滤波器的最小阶数%设计模拟原型低通滤波器%变零极点增益形式为系统传递函数形式%求模拟低通滤波器的频率响应%绘制低

32、通频率特性%原型低通转换为低通模拟滤波器%求模拟低通滤波器的频率响应%绘制高通频率特性%模拟低通滤波器转换为数字滤波器%求数字低通滤波器的频率响应%绘制数字低通滤波器的频率title(' 数字低通滤波器 '); xlabel('w/pi');ylabel(' 幅度 (dB)');e2,f2=lp2hp(b,a,300*2*pi) ;%原型低通模拟滤波器转换为高通模拟滤波器h2,w3=freqs(e2,f2);figure;%求模拟高通滤波器的频率响应plot(w3/(2*pi),abs(h2);grid;%绘制高通频率特性title('

33、模拟高通滤波器 ');xlabel('频率 Hz');ylabel('幅度');r2,q2=bilinear(e2,f2,1000);%模拟高通滤波器转换为数字滤波器s2,w4=freqz(r2,q2);figure;%求数字高通滤波器的频率响应plot(w4/pi, 20*log10(abs(s2);grid;title('数字高通滤波器'); xlabel('w/pi');ylabel('幅度(dB)');%绘制数字高通滤波器的频率响应W0=200*2*pi;Bw=200;e3,f3=lp2bp(b,a

34、,W0,Bw) ;%原型低通转换为带通模拟滤波器h3,w5=freqs(e3,f3);figure;%求模拟带通滤波器的频率响应plot(w5/(2*pi),abs(h3);grid;title('模拟带通滤波器');%绘制带通频率特性xlabel('频率 Hz');ylabel('幅度');r3,q3=bilinear(e3,f3,1000);%模拟滤波器转换为数字滤波器s3,w6=freqz(r3,q3);figure;%求数字带通滤波器的频率响应plot(w6/pi, 20*log10(abs(s3);grid;%绘制数字带通滤波器的频率响

35、应n=0:180;dt=1/Fs;t=n*dt;f1=100;f2=550;x=sin（2*pi*f1*t）+cos（2*pi*f2*t）+sin（2*pi*1*t）; figure;plot（t,x）;grid;title（'输入信号'）;xlabel（'时间 t'）;ylabel（'幅度'）;y=filter（r1,q1,x）;figure;plot（t,y）;grid;title（'低通滤波输出信号'）;xlabel（'时间 t'）;ylabel（'幅度'）y=filter（r2,q2,x）;figure;plot（t,y）;grid;title（'高通通滤波输出信号'）;xlabel（'时间 t'）;yla