DSP研究性学习报告数字滤波器设计.doc

数字信号处理课程研究性学习报告数字滤波器设计专题研讨【目的】(1) 掌握IIR和FIR数字滤波器的设计方法及各自的特点。(2) 掌握各种窗函数的时频特性及对滤波器设计的影响。(3) 培养学生自主学习能力，以及发现问题、分析问题和解决问题的能力。【研讨题目】 基本题 1IIR 数字滤波器设计设计一个IIR数字低通滤波器，其能取代下列指标的模拟低通滤波器(系统的抽样频率为44.1kHz) fp=2kHz， fs=10kHz , Ap=0.5dB, As=50dB(1) 分别用双线性变换和冲激响应不变法设计一个BW型数字低通滤波器，并进行比较。(2) 用双线性变换分别设计Chebyshev I型Chebyshev I I型和椭圆型数字低通滤波器，并进行比较。【温磬提示】在数字滤波器的设计中，不管是用双线性变换法还是冲激响应不变法，其中的参数T的取值对设计结果没有影响。但若所设计的数字滤波器要取代指定的模拟滤波器时，则抽样频率（或抽样间隔T）将对设计结果有影响。【仿真】（1）脉冲响应不变法clear;Fs=44100;fs=10000;fp=2000;ws=fs*2*pi;wp=fp*2*pi;%确定模拟滤波器指标Ws=ws/Fs;Wp=wp/Fs;Ap=0.5;As=50;N=buttord(wp,ws,Ap,As,s);%确定滤波器阶数fprintf(N=%.0fn,N);wc=wp/10(0.1*Ap-1)(1/N/2);%计算3DB截频numa,dena=butter(N,wc,s);%确定BW AF，求出滤波器系统函数numd,dend=impinvar(numa,dena,Fs); w=linspace(0,pi,512);h=freqz(numd,dend,w); %模拟滤波器到数字滤波器的转换norm=max(abs(h);numd=numd/norm;%计算频谱plot(w/pi,20*log10(abs(h)/norm);xlabel();ylabel(幅值);w=Wp Ws;h=freqz(numd,dend,w);fprintf(Ap=%.4fn,-20*log10(abs(h(1);fprintf(As=%.4fn,-20*log10(abs(h(2); %计算通带衰减、阻带衰减grid on;Ap=0.4618As=60.4050N=5原本As应该为50.000，但是由于有混叠，所以As大于50.000双线性变换法clear;FS=44100;fp=2000;fs=10000;Ap=0.5;As=50;wp=fp*2*pi; ws=fs*2*pi;wp1=wp/FS;ws1=ws/FS;OmegaP=2*FS*tan(wp1/2);OmegaS=2*FS*tan(ws1/2); N,wc=buttord(OmegaP,OmegaS,Ap,As,s); bt,at=butter(N,wc,s); bz,az=bilinear(bt,at,FS);w=linspace(0,pi,512);h=freqz(bz,az,w); norm=max(abs(h);bz=bz/norm;plot(w/pi),20*log10(abs(h)/norm);xlabel();ylabel(幅值);w=wp1 ws1;h=freqz(bz,az,w);fprintf(Ap= %.4fn,-20*log10(abs(h(1);fprintf(As= %.4fn,-20*log10(abs(h(2); grid on;Ap= 0.2445As= 50.0000利用双线性变换法，使得混叠现象被克服（2）Chebyshev I型clear;FS=44100;fp=2000;fs=10000;Ap=0.5;As=50;wp=fp*2*pi; ws=fs*2*pi;wp1=wp/FS;ws1=ws/FS;OmegaP=2*FS*tan(wp1/2);OmegaS=2*FS*tan(ws1/2); N,wc=cheb1ord(OmegaP,OmegaS,Ap,As,s); bt,at=cheby1(N,Ap,wc,s); bz,az=bilinear(bt,at,FS);w=linspace(0,pi,512);h=freqz(bz,az,w); norm=max(abs(h);bz=bz/norm;subplot(2,1,1);plot(w/pi),20*log10(abs(h)/norm);xlabel();ylabel(幅值);r,p,k=residuez(bz,az);subplot(2,1,2);zplane(bz,az);w=wp1 ws1;h=freqz(bz,az,w);fprintf(Ap= %.4fn,-20*log10(abs(h(1);fprintf(As= %.4fn,-20*log10(abs(h(2);Ap= 0.4999As= 71.0563Chebyshev I I型clear;FS=44100;fp=2000;fs=10000;Ap=0.5;As=50;wp=fp*2*pi; ws=fs*2*pi;wp1=wp/FS;ws1=ws/FS;OmegaP=2*FS*tan(wp1/2);OmegaS=2*FS*tan(ws1/2); N,wc=cheb2ord(OmegaP,OmegaS,Ap,As,s); bt,at=cheby2(N,As,wc,s); bz,az=bilinear(bt,at,FS);w=linspace(0,pi,512);h=freqz(bz,az,w); norm=max(abs(h);bz=bz/norm;subplot(2,1,1);plot(w/pi),20*log10(abs(h)/norm);xlabel();ylabel(幅值);r,p,k=residuez(bz,az);subplot(2,1,2);zplane(bz,az);w=wp1 ws1;h=freqz(bz,az,w);fprintf(Ap= %.4fn,-20*log10(abs(h(1);fprintf(As= %.4fn,-20*log10(abs(h(2);Ap= 0.5000As= 53.0577椭圆型滤波器clear;FS=44100;fp=2000;fs=10000;Ap=0.5;As=50;wp=fp*2*pi; ws=fs*2*pi;wp1=wp/FS;ws1=ws/FS;OmegaP=2*FS*tan(wp1/2);OmegaS=2*FS*tan(ws1/2); N,wc=ellipord(OmegaP,OmegaS,Ap,As,s); bt,at=ellip(N,Ap,As,wc,s); bz,az=bilinear(bt,at,FS);w=linspace(0,pi,512);h=freqz(bz,az,w); norm=max(abs(h);bz=bz/norm;subplot(2,1,1);plot(w/pi),20*log10(abs(h)/norm);xlabel();ylabel(幅值);r,p,k=residuez(bz,az);subplot(2,1,2);zplane(bz,az);w=wp1 ws1;h=freqz(bz,az,w);fprintf(Ap= %.4fn,-20*log10(abs(h(1);fprintf(As= %.4fn,-20*log10(abs(h(2);Ap= 0.5000As= 51.2355Chebyshev I型 Chebyshev I I型 椭圆型滤波器【结果分析】双线性变换和冲激响应不变法所设计的滤波器的性能有什么不同。脉冲响应不变法的主要优点是模拟频率与数字频率之间的关系是线性的，其主要缺点是存在频谱混叠，使得阻带衰减不满足条件。双线性变换法主要优点是避免了频谱混叠，依靠频率的非线性关系得到s平面与z平面的单值对应关系，整个 轴单值对应于单位圆一周，消除了脉冲响应不变法固有的频谱混叠效应，其缺点是模拟频率与数字频率之间的关系是非线性的。Chebyshev I型、Chebyshev I I型和椭圆型滤波器的零极点分布各有什么特点。在相同的设计指标下，BW型滤波器的阶数最高，椭圆滤波器的阶数最低。即使阶数相等，它们的裕量也不同。零极点分布都是关于虚轴对称且都在单位园内。在滤波器的实现过程中， BW型滤波器最容易实现，而椭圆滤波器不易实现。【自主学习内容】利用脉冲响应不变法和双线性变换法设计滤波器【阅读文献】数字信号处理【问题探究】脉冲响应不变法的优缺点 优点：数字滤波器和模拟滤波器的频率关系为线性缺点：存在频谱混叠，故不能用脉冲响应不变法设计高通、带阻等滤波器。 双线性变换法的优缺点优点：当模拟滤波器的幅度响应为分段常数时，双线性变换后的数字滤波器能保持模拟滤波器的幅度响应，但分段边界点将产生畸变，这种畸变可以利用式 在数字滤波器指标转换成相应的模拟滤波器指标时进行预畸变校正。双线性变换法一般适合于设计幅度相应为分段常数的数字滤波器，缺点：幅度响应不是常数时会产生幅度失真，不适合设计像数字微分器等幅度相应为非常数的数字滤波器。BW型、Chebyshev I型、Chebyshev I I型和椭圆型滤波器的特点;1. 在相同的设计指标下，一般来说， BW型滤波器的阶数最高，椭圆滤波器的阶数最低。即使阶数相等，它们的裕量也不同。2在滤波器的实现过程中， BW型滤波器最容易实现，而椭圆滤波器不易实现（因为它的系统函数H(s)的极点离jw轴最近）。3 .由图可以看出，BW型设计的的滤波器没有震荡，而CB I型在通带会有波动，CB I I型在阻带有较大波动，椭圆型在通阻带都有波动。且由图可以看出不同滤波器的过渡带的宽度也有差距。【研讨题目】 基本题 2窗函数研究分析矩形窗、汉纳窗、哈明窗、布莱克曼窗、凯泽窗的频域特性，并进行比较。【仿真结果】设计一满足下列指标的线性相位 FIR高通滤波器：Wp=0.67,Ws=0.53,Ap=0.3dB,As=50dB。采用Hamming窗加窗截断，设计过程如下：Wp=0.67*pi;Ws=0.53*pi;Ap=0.3;As=50;N=ceil(7*pi/(Wp-Ws); %确定滤波器阶数，并使滤波器为I型 N=mod(N+1,2)+N; % 进行取余运算M=N-1;fprintf(N=%.0fn,N); %输出Nw=hamming(N);Wc=(Wp+Ws)/2; %理想低通截频k=0:M;hd=-(Wc/pi)*sinc(Wc*(k-0.5*M)/pi); hd(0.5*M+1)=hd(0.5*M+1)+1;h=hd.*w;omega=linspace(0,pi,512); %0-取512个点mag=freqz(h,1,omega);magdb=20*log10(abs(mag);plot(omega/pi,magdb);N=51哈明窗：通带衰减Ap0dB 阻带衰减As54 dB设计一满足下列指标的线性相位 FIR高通滤波器：Wp=0.67,Ws=0.53,Ap=0.3dB,As=50dB。采用blackman窗加窗截断，设计过程如下：Wp=0.67*pi;Ws=0.53*pi;Ap=0.3;As=50;N=ceil(7*pi/(Wp-Ws);N=mod(N+1,2)+N;M=N-1;fprintf(N=%.0fn,N);w=blackman(N);Wc=(Wp+Ws)/2;k=0:M;hd=-(Wc/pi)*sinc(Wc*(k-0.5*M)/pi);hd(0.5*M+1)=hd(0.5*M+1)+1;h=hd.*w;omega=linspace(0,pi,512);mag=freqz(h,1,omega);magdb=20*log10(abs(mag);plot(omega/pi,magdb);布莱克曼窗通带衰Ap0dB阻带衰减 As75dB 设计一满足下列指标的线性相位 FIR高通滤波器：Wp=0.67,Ws=0.53,Ap=0.3dB,As=50dB。采用hanning窗加窗截断，设计过程如下：Wp=0.67*pi;Ws=0.53*pi;Ap=0.3;As=50;N=ceil(7*pi/(Wp-Ws);N=mod(N+1,2)+N;M=N-1;fprintf(N=%.0fn,N);w=hanning(N);Wc=(Wp+Ws)/2;k=0:M;hd=-(Wc/pi)*sinc(Wc*(k-0.5*M)/pi);hd(0.5*M+1)=hd(0.5*M+1)+1;h=hd.*w;omega=linspace(0,pi,512);mag=freqz(h,1,omega);magdb=20*log10(abs(mag);plot(omega/pi,magdb);汉纳窗 通带衰减 Ap0dB阻带衰减 As44dB设计一满足下列指标的线性相位 FIR高通滤波器：Wp=0.67,Ws=0.53,Ap=0.3dB,As=50dB。采用矩形窗加窗截断，设计过程如下：Wp=0.67*pi;Ws=0.53*pi;Ap=0.3;As=50;N=ceil(7*pi/(Wp-Ws);N=mod(N+1,2)+N;M=N-1;fprintf(N=%.0fn,N);w=boxcar(N);Wc=(Wp+Ws)/2;k=0:M;hd=-(Wc/pi)*sinc(Wc*(k-0.5*M)/pi);hd(0.5*M+1)=hd(0.5*M+1)+1;h=hd.*w;omega=linspace(0,pi,512);mag=freqz(h,1,omega);magdb=20*log10(abs(mag);plot(omega/pi,magdb);矩形窗通带衰减 Ap0dB阻带衰减 As21dB设计一满足下列指标的线性相位 FIR高通滤波器：Wp=0.67,Ws=0.53,Ap=0.3dB,As=50dB。采用凯泽窗加窗截断，设计过程如下：Ap=0.3;As=50;Rp=1-10.(-0.05*Ap);Rs=10.(-0.05*As);f=0.53,0.67;a=0,1;dev=Rp,Rs;M,Wc,beta,ftype = kaiserord(f,a,dev);M=mod(M,2)+M;h = fir1(M,Wc,ftype,kaiser(M+1,beta)omega=linspace(0,pi,512);mag=freqz(h,1,omega);plot(omega/pi,20*log10(abs(mag);凯泽窗通带衰减Ap0dB阻带衰减As52dB【结果分析】1.哈明窗 通带衰减 Ap0dB 阻带衰减 As54 dB2.汉纳窗 通带衰减 Ap0dB 阻带衰减 As44dB3.矩形窗 通带衰减 Ap0dB 阻带衰减 As21dB4.凯泽窗 通带衰减 Ap0dB 阻带衰减 As52dB5.布莱克曼窗 通带衰减 Ap0dB 阻带衰减 As75dB【结果分析】各种窗有何特点?矩形窗设计出的FIR滤波器阻带衰减最小,布莱克曼窗阻带衰减最大。1.Hamming(哈明)窗(w=hamming(M+1)哈明窗设计的FIR滤波器在不连续点附近的幅度函数Ap 0.019dB, As 53dB2.Blackman窗(w=blackman(M+1)布莱克曼窗设计的FIR滤波器在不连续点附近的幅度函数Ap 0.0017dB，As 74dB3.矩形窗矩形窗设计的FIR滤波器在不连续点附近的幅度函数A（）Ap =-20lg(1-dp) 0.82dB, As = -20lg(ds ) 21dB4.汉纳窗图为汉纳窗设计的FIR滤波器在不连续点附近的幅度函数A（）Ap 0.056dB, As 44dB5.凯泽窗是一可调参数，通过改变的值可以调节窗函数的形状。式中A= -20lg (minp，s )滤波器阶数M或凯泽窗的长度则可由下式估计矩形窗、汉纳窗、哈明窗、布莱克曼窗比较【自主学习内容】各种窗函数的设计方法【阅读文献】数字信号处理【发现问题】 (专题研讨或相关知识点学习中发现的问题)：在谱分析中如何选择窗函数，在滤波器设计中如何选择窗函数？【问题探究】在谱分析中如何选择窗函数，在滤波器设计中如何选择窗函数？由结果分析可知，在矩形窗、汉纳窗、哈明窗、布莱克曼窗中，矩形窗的过渡带最窄，但利用它设计出的FIR滤波器的阻带衰减最小。利用布莱克曼窗设计出的FIR滤波器阻带衰减最大，但其过渡带也最宽。显然，减小了窗函数旁瓣的相对幅度却增加了其主瓣的宽度，即提高FIR滤波器阻带衰减是以增加过渡带宽度为代价的。在工程应用中，在满足阻带衰减的前提下，尽可能地选择主瓣宽度较小的窗函数。而在实际设计中，可由待设计的FIR数字滤波器阻带衰减或通带波动来确定窗函数的类型，有过渡带宽度估计窗函数的长度N(N=M+1)。 而凯泽窗是一种应用广泛的可调窗，它可以通过改变窗函数的形状来控制窗函数旁瓣的大小，而在设计中可根据滤波器的衰减指标来确定窗函数的形状。【研讨题目】 基本题 3 窗函数法设计FIR 数字滤波器(1)分别用Blackman窗和Kaiser窗法设计一个满足下列指标的线性相位的FIR低通滤波器Wp=0.4p rad, Ap=0.5 dB, Ws=0.6p rad, As=55dB (2)（M5-5）在用窗口法设计FIR滤波器时，由于理想滤波器的幅度响应在截频处发生突变，使得设计出的滤波器的幅度响应发生振荡，这个现象被称为Gibbs现象。解决这个问题的一个方案是本书中介绍的用逐步衰减的窗函数。另一个方案是使理想滤波器过渡带为渐变的，如下图所示具有线性过渡带的理想低通滤波器的频率响应，试用窗口法设计逼近该频率响应的FIR滤波器。 题3图【(2)单位脉冲响应证明】试证该滤波器的单位脉冲响应为其中：，【设计步骤】（1）确定线性相位滤波器的类型（I型）（2）确定理想滤波器的幅度函数Ad()（3）确定滤波器相位（4）计算hdk（5）利用窗函数截断hdk【仿真结果】所设计滤波器的幅度响应和相位响应（1）wp=0.4*pi;ws=0.6*pi;Ap=0.5;As=55;N=ceil(11.4*pi/(ws-wp);M=N-1;w=blackman(N);wc=(wp+ws)/2;k=0:M;hd=(wc/pi)*sinc(wc*(k-0.5*M)/pi);h=hd.*w;omega=linspace(0,pi,512);mag=freqz(h,1,omega);plot(omega/pi,20*log10(abs(mag)Rs=0.01;f=0.4 0.6 ;a=0.975,0;dev=Rs*ones(1,length(a);M Wc beta ftype=kaiserord(f,a,dev);h=fir1(M,Wc,ftype,kaiser(M+1,beta);omega=linspace(0,pi,512);mag=freqz(h,1,omega);plot(omega/pi,20*log10(abs(mag)（2）渐变的窗函数选择hamming窗。为了简便研究过程，设Wp=0.55pi、Ws=0.45pi、As=25dB、Ap=1dB。设hamming窗阶数为M，矩形窗的长度为M1，图中蓝线为第一种方案涉及到滤波器，红线为第二种方案设计的滤波器。M1=5时 M1=8时M1=15时 M1=30M1=80Wp=0.55*pi; Ws=0.45*pi; Ap=1; As=25;N=ceil(7*pi/(Wp-Ws);N=mod(N+1,2)+N;M=N-1;w=hamming(N);Wc=(Wp+Ws)/2;k=0:M;hd=-(Wc/pi)*sinc(Wc*(k-0.5*M)/pi); h=hd.*w;omega=linspace(0,pi,512);mag=freqz(h,1,omega);magdb=abs(mag);plot(omega/pi,magdb);grid;W=Ws-Wp;M1=15;k2=-M1:M1;Wc=(Wp+Ws)/2;hd=sinc(W*k2/2).*(sin(Wc*k2)./(k2.*pi);hd(M1+1)=Wc/pi;omega2=linspace(0,pi,512);mag2=freqz(hd,1,omega2);magdb2=abs(mag2);hold on;plot(omega2/pi,magdb2,r);【结果分析】由仿真结果可知，第一种方案的过渡带明显短于第二种方案；当M=7时，第一种方案几乎可以完全消除Gibbs现象，但M1=8时，第二种方案仍然可以看到明显的通阻带波动，故消除Gibbs现象方案二需要的阶数更高。【研讨题目】 中等题 4频率取样法FIR 数字滤波器（1）（M5-6）利用频率取样法设计某I型线性相位带通FIR滤波器，其通带截频分别为Wp1=0.3p rad, Wp2=0.5p rad（2）（M5-7）在通带和阻带间增加1个过渡点，重新设计该滤波器。过渡点的最佳幅度由实验确定。【设计步骤】由题意我们确定设计滤波器的方法：频率取样法。我们设计滤波器的类型为一型滤波器。采用以下步骤：确定理想滤波器的幅度函数Ad( )确定理想滤波器的相位 ，对于类型一 确定 在M+1个取样点上的值Hdm对Hdm做M+1点IDFT，即可得到hk【仿真结果】(1)滤波器阶数M与阻带衰减的关系（2）确定过渡点的最佳幅度（以M=44为例）【结果分析】滤波器的阻带衰减，滤波器的阶数与设计结果的关系。增加滤波器的阶数M对阻带的波动几乎没有改善。出现这种现象的原因是从通带到阻带所给定的样本点发生了从0到1的跳变。改善的方法是增加非0的过渡点具体方法是：在通、阻带交界处人为地安排一到几个过渡点，其值介于零和1之间，这样可减小样点间的落差，使过渡平缓，阻带最小衰耗增大。本题中经试验表明，两个过渡点分别取值0.32， 0.38。获得最佳阻带衰减。结论：在过渡带插入过渡点后，改善了滤波器的幅度特性，在滤波器的阶数不变的情况下，得到较好的滤波效果。在进行滤波器设计时候，是一种优化的好方法。【仿真程序】（1）不加过渡点M=input(Input a number to M=);m=0:M/2;Wp1=0.3*pi;Wp2=0.5*pi;Wm=2*pi*m./(M+1);Ad=double(Wm=wp1)&(Wm=Wp1)&(Wm=Wp);Ad(mtr)=0.28;Hd=Ad.*exp(-j*0.5*M*Wm);Hd=Hd conj(fliplr(Hd(2:M/2+1);h=real(ifft(Hd);w=linspace(0,pi,1000);H=freqz(h,1,w);plot(w/pi,20*log10(abs(H)Parks-McClellan算法Fp=0.5;Fs=0.4;ds=0.0017;dp=0.0017;f=Fs Fp;a=0 1;dev=ds dp;M,fo,ao,w=remezord(f,a,dev);h=remez(M,fo,ao,w);w=linspace(0,pi,1000);mag=freqz(h,1,w);plot(w/pi,20*log10(abs(mag)【研讨题目】 中等题 5设计幅度响应逼近下图所示的数字高通滤波器，其中Wc=0.5p。要求：(1) 将该数字高通滤波器设计成IIR数字滤波器，具体要求：分别设计成BW型、CBI型、CBII型和椭圆型滤波器，并比较设计结果。(2) 将该数字高通滤波器设计成FIR数字滤波器，具体要求：a. 采用窗函数法，分别利用矩形窗、汉纳窗、哈明窗、布莱克曼窗、凯泽窗截断，并将设计结果进行分析比较。b. 采用频率取样法，讨论过渡点对滤波器阻带衰减的影响。c. 采用Parks-McClellan算法。d. 比较窗函数法、频率取样法和Parks-McClellan算法所设计的FIR滤波器。(3) 所设计的IIR和FIR数字滤波器应具有基本相同的幅度响应。根据设计结果，从幅度响应、相位响应、滤波器阶数等对两类滤波器进行分析比较。【数字高通滤波器设计指标】因为Wc=0.5p，且该滤波器为高通滤波器，为了方便研究过程，设Wp=0.55p,Ws=0.45p,As=30dB, AP=1dB.【设计步骤】（1） IIR数字滤波器设计步骤：双线性法数字滤波器H(Z)模拟滤波器H(S)模拟滤波器指标频率转换待设计的数字滤波器指标模拟设计转换高通滤波器HH(S)低通滤波器HL(S)WP ，Ws ，Ap, AsWP ，Ws ，Ap, As频率转换确定过渡带和窗函数以及函数类型(2) a.FIR数字滤波器设计步骤窗函数法加窗单位脉冲响应hdk滤波器频率响应Wp，Ws，Ap, Ashk= hdk*wkIDTFT确定过渡带和窗函数以及函数类型bFIR数字滤波器设计步骤频率抽样法hk抽样后的频率HmDTFT抽样滤波器频率响应Wp，Ws，Ap, As【仿真结果】(1) Wp=0.55p,Ws=0.45p,As=30dB, AP=1dB.转化为模拟滤波器指标，双线性法取T=2Wp=0.55*pi; Ws=0.45*pi; Ap=1; As=30;T=2;Fs=1/T; wp=2*tan(Wp/2)/T;ws=2*tan(Ws/2)/T;wp1=1/wp;ws1=1/ws;N,wc=buttord(wp1,ws1,Ap,As,s);num,den=butter(N,wc,s);numa,dena=lp2hp(num,den,1);numd,dend=bilinear(numa,dena,Fs);w=linspace(0,pi,1024);h=freqz(numd,dend,w);plot(w/pi,20*log10(abs(h);axis(0 1 -50 0);grid;xlabel(Normalized frequency);ylabel(Gain,dB);BW型:以下的CBI型、CBII型和椭圆型滤波器型只需该上面相应的BW型函数既可，在此就不赘述程序。CBI型：Wp=0.55*pi; Ws=0.45*pi; Ap=1; As=30;T=2;Fs=1/T; wp=2*tan(Wp/2)/T;ws=2*tan(Ws/2)/T;wp1=1/wp;ws1=1/ws;N,wc=cheb1ord(wp1,ws1,Ap,As,s);num,den=cheby1(N,Ap,wc,s);numa,dena=lp2hp(num,den,1);numd,dend=bilinear(numa,dena,Fs);w=linspace(0,pi,1024);h=freqz(numd,dend,w);plot(w/pi,20*log10(abs(h);axis(0 1 -50 0);grid;xlabel(Normalized frequency);ylabel(Gain,dB);明显CBI型的幅度响应在通带有等波纹特性。CBII型：Wp=0.55*pi; Ws=0.45*pi; Ap=1; As=30;T=2;Fs=1/T; wp=2*tan(Wp/2)/T;ws=2*tan(Ws/2)/T;wp1=1/wp;ws1=1/ws;N,wc=cheb2ord(wp1,ws1,Ap,As,s);num,den=cheby2(N,As,wc,s);numa,dena=lp2hp(num,den,1);numd,dend=bilinear(numa,dena,Fs);w=linspace(0,pi,1024);h=freqz(numd,dend,w);plot(w/pi,20*log10(abs(h);axis(0 1 -50 0);grid;xlabel(Normalized frequency);ylabel(Gain,dB);明显CBII型的幅度响应阻带有等波纹特性。椭圆型滤波器型：Wp=0.55*pi; Ws=0.45*pi; Ap=1; As=30;T=2;Fs=1/T; wp=2*tan(Wp/2)/T;ws=2*tan(Ws/2)/T;wp1=1/wp;ws1=1/ws;N,wc=ellipord(wp1,ws1,Ap,As,s);num,den=ellip(N,Ap,As,wc,s);numa,dena=lp2hp(num,den,1