引用

的

滤波是信号处理的一种基本而重要的技术，利用滤波可从复杂的信号中提取所需要的信号，抑制不需要的部分。第3章中利用FFT进行简单滤波的尝试已初步说明了这个问题。所谓滤波器是指具有一定传输特性的信号处理装置。本章首先介绍模拟滤波器设计原理，然后介绍几种常见原型滤波器的设计，尔后介绍将模拟原型滤波器转换为其他类型滤波器的频率转换方法，{zh1}介绍MATLAB中的模拟滤波器并通过应用实例进行说明。T>

5．1．1滤波原理T>

滤波器，顾名思义，其作用是对输入信号起到滤波的作用。对图5-1所示的线性时不变系统，其时域输入和输出关系为：T>

y(n)=x(n)?h(n)                 （5-1）             T>

 T>

               x(n)         y(n)T>

T>

 T>

                     图5-1 线性时不变系统示意图T>

若x(n),y(n)的Fourier变换存在，为和，则输入和输出的关系为：T>

                    （5-2）T>

式中，为系统的频率响应。再假定x(n)的Fourier变换的振幅可用图5-2的上图表示，系统幅频响应可用图5-2的中图表示，则滤波器的输出y(n)的频谱可表示为图5-2的下图。T>

T>

    图5-2   滤波器滤波示意图T>

这样，x(n)通过系统h(n)的结果是使输出y(n)中不再含有的频率成分，而使的成分“不失真”地给以通过。因此设计出不同形状的可以得到不同的滤波结果。T>

若滤波器的输入、输出都是离散时间信号，那么该滤波器的脉冲响应h(n)也必然是离散的。我们称这样的滤波器为数字滤波器（digital filter）。当用硬件实现一个数字滤波器时，所需的元件是延迟器、乘法器和加法器。当在计算机上用软件实现时，它就是一段线性卷积的程序。我们知道模拟滤波器只能用硬件来实现，其元件是电阻、电容、电感及运算放大器等。T>

5．1．2滤波器的分类T>

滤波器的种类很多，分类方法也不同，如可以从功能上分，也可以从实现方法上分，或从设计方法上来分等。但总的来说，滤波器可分为两大类，即经典滤波器和现代滤波器。经典滤波器是假定输入信号x(n)中的有用成分和希望去除的成分各自占有不同的频带，如图5-2所示，当x(n)通过一个线性系统（即滤波器）后可将欲去除的成分有效去除。如果信号和噪声的频谱相互重叠，那么经典滤波器将无能为力。现代滤波器理论研究的主要内容是从含有噪声的数据记录（又称时间序列）中估计出信号的某些特征或信号本身。一旦信号被估计出，那么估计出的信号将比原信号有更高的信噪比。本课程侧重于经典滤波器。T>

经典滤波器从功能上总的可分为四种，即低通（te="mailto:wan">lte="mailto:wan">Low pass）、高通（te="mailto:wan">hte="mailto:wan">High pass）、带通（te="mailto:wan">bte="mailto:wan">Band pass）、带阻（te="mailto:wan">bte="mailto:wan">Band stop）滤波器，当然每一种又有模拟滤波器和数字滤波器两种形式。图5-3和图5-4分别给出了模拟滤波器和数字滤波器的四种形式。图中所给的滤波器的幅频响应都是理想情况，在实际上是不可能实现的。在实际工作中，我们设计的滤波器都是在某些准则下对理想滤波器的近似，这保证了滤波器在物理上是可实现的，且是稳定的。T>

 T>

T>

图 5-3  理想模拟滤波器幅频响应T>

T>

图5-4   理想数字滤波器的幅频响应T>

图5-3和图5-4 的理想滤波器在物理上是不可实现的。不可实现的根本原因是从一个频率带到另一个频率带之间存在突变。为了物理上可实现，我们从一个带到另一个带之间设置一个过渡带，且在通带和阻带内也不应该严格为1或零，应给以较小的容限。在5.2节中，我们要介绍这种容限。T>

 T>

5.2.1  信号无失真传输的条件T>

所谓信号无失真传输是指输入信号通过系统后，输出信号的幅值和输入信号的幅值成正比。允许有一定的延时，但没有波形上的畸变。因此，系统的频率响应H(j)满足下面的特性：T>

                                          (5-3)T>

式中，k,t_d均为常数。T>

即信号无失真传输的条件是：系统的幅频响应||应为常数，相频响应应与频率成比例。或者说，滤波器应具有无限宽的定值幅频与线性相频。通常定义群延迟(Group delay)为信号系统的延迟时间，即：T>

                                                     (5-4)T>

群延迟为相频特性曲线的斜率。对于信号无失真传输，t_d为常数，即群延迟为常数；否则它是频率的非线性函数。T>

5.2.2  理想滤波器的特性T>

滤波器是一个选频装置。理想滤波器应能无失真地传输有用信号，而又能xx抑制无用信号。有用信号和无用信号往往占有不同的频带。信号能通过滤波器的频带称为通带(Passband)。信号被抑制的频带称为阻带（Stopband）。理想滤波器频率特性可写为：T>

                                        (5-5)T>

但通过后面的讲解我们会知道理想滤波器是物理不可实现系统。实际滤波器的频率特性只能“逼近”理想滤波器。图 5-5为低通滤波器的幅频特性示意图。T>

T>

图5-5低通滤波器的幅频特性示意图T>

 T>

可见，滤波器的幅频响应在通带内不是xx平直的，而是呈波纹变化；在阻带内，幅频特性也不为零，而是衰减至某个值；在通带和阻带之间存在一个过渡带，而不是突然下降。通常，实际设计要求滤波器的技术指标包括通带波纹Rp(Passband ripple)(dB)、阻带衰减Rs(Stopband attenuation)(dB)、通带边界频率、阻带边界频率、过渡带宽。T>

滤波器的通带波纹Rp为相对于频率响应{zd0}点（一般为1）的下降，因此下降越少说明通带越平直，滤波器的滤波效果越好（通常为1~5dB）。滤波器的阻带衰减Rs也是相对于频率响应{zd0}点（一般为1）的下降，因此下降越多说明信号在阻带内越不容易通过，因此滤波效果越好(通常要大于15dB)。过渡带宽越窄，滤波器的频率特性越接近于直角矩形特性，滤波效果也越好。应该指出的是，图5-5给出的是低通滤波器的情况，对于高通、带通、带阻具有相同的参数。T>

5.2.3  模拟滤波器传递函数设计原理T>

模拟滤波器的设计理论通常在Laplace域内进行讨论，本节只讨论其应用。模拟滤波器的技术指标可由平方幅值响应函数的形式给出，而和传递函数H(s)存在下面关系：T>

                                    (5-6)T>

即T>

                                             (5-7)T>

当给定模拟滤波器的技术指标后，由求出A(-s²)，再适当地分配零极点可求出H(s)。为了使滤波器稳定，H(s)的极点必须落在s平面左半平面，这是因为对于一个滤波器的极点p，te="mailto:wan">就可以写成的形式，其逆Laplace变换（对应于时间域）为，若p>0，则随着时间增大至无穷，该滤波器的输出将出现不稳定。滤波器的零点选择可任取A(-s²)的一半零点，这是因为滤波器对Laplace域表示的传递函数并无特殊要求，但如果要求H(s)具有最小相位，零点也必须选择在s左半平面。T>

本节介绍常用的模拟原型滤波器的主要特点及其MATLAB实现，包括Butterworth、Chebyshev I, ChebyshevII, Elliptical、Bessel原型低通滤波器的设计。模拟原型滤波器指的是截止频率为1的滤波器。后面所讲的各类模拟滤波器和数字滤波器可通过这些低通原型滤波器变换得到。T>

5.3.1  Butterworth滤波器T>

Butterworth模拟低通滤波器的平方幅频响应函数为：T>

                                           (5-8)T>

式中，为低通滤波器的截止频率（Cutoff frequency），N为滤波器的阶数。T>

Butterworth滤波器的特点：通带内具有{zd0}平坦的频率特性，且随着频率增大平滑单调下降；阶数愈高，特性愈接近矩形，过渡带愈窄，传递函数无零点。T>

这里的特性接近矩形，是指通带频率响应段与过渡带频率响应段的夹角接近直角。通常该角为钝角，如果该角为直角，则为理想滤波器。T>

所谓滤波器的零点就是将该点的值代入传递函数后，传递函数的值为零。所谓函数的极点就是将该点的值代入传递函数后，传递函数的值为无穷大。滤波器的增益是指传递函数表达式前的常数。若系统的传递函数表示为T>

                        (5-9)T>

则滤波器零点为[z(1),z(2),…,z(nz)]，极点为[p(1),p(2),…,p(np)]，滤波器的增益为K。T>

这里所说的零点和极点分布在一个圆上为Laplace域中的形式，感兴趣的te="mailto:wan">读者te="mailto:wan">同学可查看数学中的Laplace变换。T>

MATLAB信号处理工具箱提供Butterworth模拟低通滤波器原型设计函数buttap,函数调用形式为：T>

[z,p,k]=buttap(n)                       T>

式中，n为butterworth滤波器阶数；z,p,k分别为滤波器的零点、极点和增益。T>

Butterworth滤波器传递函数具有下面的形式：T>

                              (5-10)T>

滤波器没有零点，极点为[p(1),p(2),…,p(n)]，滤波器的增益为K。T>

在有关模拟滤波器设计的MATLAB程序中，经常遇到一些特定函数。[b,a]=zp2tf(z,p,k)为将模拟原型滤波器函数（如buttap）设计出的零点z,极点p和增益k形式转换为传递函数（Transfer function）形式；其中，b为滤波器传递函数分子多项式系数，a为滤波器传递函数分母多项式系数,见（4-4）式。[H,w]=freqs(b,a,w)求出传递函数形式（分子和分母多项式的系数为b,a）表示的滤波器的对应于频率点w的复数频率响应H（包括实部和虚部），这里w为一个矢量，表示对应的角频率。若该函数不写输出变量，则执行后绘出该滤波器的幅频响应和相频响应图。T>

【例5-1】绘制Butterworth低通模拟原型滤波器的幅频平方响应曲线，阶数分别为2，5，10，20。T>

%Samp5_1T>

n=0:0.01:2;    %频率点T>

for ii=1:4     %取4种滤波器T>

   switch iiT>

   case 1,N=2;T>

   case 2,N=5;T>

   case 3,N=10;T>

   case 4,N=20;T>

   endT>

   [z,p,k]=buttap(N);    %设计Butterworth滤波器T>

   [b,a]=zp2tf(z,p,k);   %将零点极点增益形式转换为传递函数形式T>

   [H,w]=freqs(b,a,n);   %按n指定的频率点给出频率响应T>

   magH2=(abs(H)).^2;    %给出传递函数幅度平方T>

   hold on;T>

   plot(w,magH2);     %绘制传递函数幅度平方T>

endT>

xlabel('w/wc');T>

ylabel('|H(jw)|^2');T>

title('Butterworth模拟原型滤波器');T>

text(1.5,0.18,'n=2')   %作必要的标记T>

text(1.3,0.08,'n=5')T>

text(1.16,0.08,'n=10')T>

text(0.93,0.98,'n=20')T>

grid on;T>

 T>

T>

图5-6 Butterworth滤波器原型平方幅频图T>

 T>

程序的运行结果见图5-6。可以看到，滤波器的幅频平方特性随着频率单调下降。随着滤波器阶数的增大，其幅频特性逐渐接近矩形,与前面介绍的Butterworth滤波器的特性一致。T>

5.3.2  Chebyshev I型T>

Chebvshev I型模拟低通滤波器的平方幅值响应函数为T>

                                   (5-11)T>

式中，为小于1的正数，表示通带内的幅值波纹情况；为截止频率，N为Chebyshev多项式阶数，为Chebyshev多项式，定义为：T>

                                    (5-12)T>

Chebyshev I型滤波器特点是：通带内具有等波纹起伏特性，而在阻带内则单调下降，且具有更大衰减特性；阶数愈高，特性愈接近矩形，传递函数没有零点。T>

MATLAB信号处理工具箱函数cheb1ap设计N阶Chebyshev I型模拟低通滤波器原型。因为Chebyshev I型模拟原型滤波器通带内具有等波纹特点，所以在调用时需增加通带波纹所能下降的{zd1}限度Rp（dB）的参数。由于Rp为相对于增益1的下降分贝数，所以该值越小滤波器在通带内越接近1，滤波器的性能越好，通常取值范围为1~5（dB）。此函数的调用格式为：T>

[z，p,k]=cheb1ap(N,Rp)T>

式中，N为滤波器的阶数，Rp为通带波纹，单位为dB。z,p,k分别为滤波器的零点、极点和增益。T>

滤波器传递函数具有与（5-9）式相同的形式。T>

【例5-2】绘制Chebyshev I型模拟低通滤波器的平方幅频响应曲线，阶数为2，4，6，8。T>

%Samp5_2T>

n=0:0.01:2;    %设置频率点T>

for ii=1:4T>

   switch iiT>

   case 1,N=2;T>

   case 2,N=4;T>

   case 3,N=6;T>

   case 4,N=8;T>

   endT>

   Rp=1;      %设置通带波纹为1dBT>

   [z,p,k]=cheb1ap(N,Rp);   %设计Chebyshev I型滤波器T>

   [b,a]=zp2tf(z,p,k);    %转换为传递函数形式T>

   [H,w]=freqs(b,a,n);    %求得传递函数的频率特性T>

   magH2=(abs(H)).^2;     %求得传递函数的幅频响应T>

posplot=['2,2,' num2str(ii)];%将数字ii转换为字符串，与’2,2,’合并并赋给posplot。T>

   subplot(posplot);T>

   plot(w,magH2);T>

   title(['N=' num2str(N)]);%将数字N转换为字符串与’N=’合并作为标题。T>

   xlabel('w/wc');T>

   ylabel('Chebyshev I |H(jw)|^2');T>

grid onT>

endT>

T>

图5-7 Chebyshev I型模拟原型滤波器平方幅频图T>

 T>

程序的运行结果见图5-7。与Butterworth滤波器相比（图5-6），可以看到：在相同的阶数下，Chebyshev I型模拟原型滤波器具有更窄（更陡）的过渡带。但这种特性是以牺牲了通带的单调平滑特性（而成为波纹状）为代价的。如果我们不想牺牲通带内的单调平滑特性有无办法设计阶数较小而过渡带较窄的滤波器呢？这就涉及Chebyshev II型滤波器。T>

5.3.3  Chebyshev II型滤波器T>

Chebyshev II型低通模拟滤波器的平方幅值响应函数为：T>

                               (5-13)T>

式中各项参数的意义同上。T>

Chebyshev II型模拟滤波器的特点是：阻带内具有等波纹的起伏特性，而在通带内是单调、平滑的，阶数愈高，频率特性曲线愈接近矩形，传递函数既有极点又有零点。T>

Chebyshev II型模拟滤波器传递函数具有（5-9）式的形式。滤波器零点为[z(1),z(2),…,z(nz)]，极点为[p(1),p(2),…,p(np)]，滤波器的增益为K。T>

MATLAB信号处理工具箱提供函数cheb2ap设计N阶Chebyshev II型模拟滤波器的原型。由于Chebyshev II型滤波器阻带内有波纹，通带内单调平滑，则必须给定阻带衰减Rs。前面已经讲过，通带波纹和阻带衰减都是相对于增益1的下降，因此，Rp和Rs越大则与通带增益1的差距越大。因此，阻带衰减Rs越大对应的滤波器滤波效果越好，通常Rs的取值要大于16dB。该函数通常调用格式为：T>

[z,p,k]=cheb2ap(N,Rs)T>

式中，N为滤波器的阶数；Rs为阻带波纹，单位dB; z,p,k为滤波器的零极点和增益。T>

【例5-3】 绘制Chebyshev II型原型模拟滤波器的平方幅频响应曲线，阶数分别为2,4,6,8。T>

%Samp5_3T>

n=0:0.01:2;  %设置频率点T>

for ii=1:4T>

   switch iiT>

   case 1,N=2;T>

   case 2,N=4;T>

   case 3,N=6;T>

   case 4,N=8;T>

   endT>

   Rs=16;     %阻带衰减设置为16dBT>

   [z,p,k]=cheb2ap(N,Rs);   %设计Chebyshev II型模拟原型滤波器T>

   [b,a]=zp2tf(z,p,k);    %转化为传递函数T>

   [H,w]=freqs(b,a,n);    %求出滤波器的频率响应T>

   magH2=(abs(H)).^2;     %频率响应的幅度平方T>

   posplot=['2,2,' num2str(ii)];T>

   subplot(posplot);T>

   plot(w,magH2);    %绘出幅频平方响应函数T>

   title(['N=' num2str(N)]);T>

   xlabel('w/wc');T>

   ylabel('Chebyshev II |H(jw)|^2');T>

   grid on;T>

endT>

 T>

T>

图5-8 Chebyshev II型模拟原型滤波器平方幅频图T>

 T>

程序运行结果见图5-8。可见Chebyshev II型滤波器在通带内是单调平滑的，而阻带内却出现了波纹。随着滤波器阶数的增高，其幅频特性越接近矩形。T>

由前面的讲解可知，Butterworth滤波器具有单调平滑的频率特性，但在过渡带宽相同的条件下，所需的滤波器阶数较高（即设计的滤波器滤波速度较慢，成本较高）。而Chebyshev I和II型滤波器在通带或阻带内有波纹出现，但在过渡带宽相同的条件下，所需的滤波器阶数较低（即设计的滤波器滤波速度较快，成本较低）。如果我们对通带和阻带内的平滑特性均不作限制，能否在相同的过渡带宽的情况下设计出阶数更低的滤波器呢？这就是下面要讲的椭圆滤波器。T>

5.3.4  椭圆滤波器T>

椭圆(Elliptic)模拟低通原型滤波器的平方幅值响应函数为T>

                                  (5-14)T>

式中，为小于1的正数，表示波纹情况；为低通滤波器的截止频率（Cutoff frequency），N为滤波器的阶数，为椭圆函数，其定义已超出本课程的范围，我们直接利用。T>

Elliptic滤波器的特点：在通带和阻带内均具有等波纹起伏特性，与以上滤波器原型相比，相同的性能指标所需的阶数最小，但相频响应具有明显的非线性。T>

MATLAB信号处理工具箱提供Elliptic模拟低通滤波器原型设计函数ellipap。由于椭圆滤波器通带和阻带均有波纹，因此函数输入参数中通带波纹和阻带衰减均需给出。该函数调用形式为：T>

[z,p,k]=ellipap(n,Rp,Rs)T>

式中，n为椭圆滤波器阶数；Rp,Rs分别为通带波纹和阻带衰减，单位dB,通常滤波器的通带波纹的范围为1~5dB，阻带衰减的范围大于15dB。z,p,k分别为滤波器的零点、极点和增益。T>

滤波器传递函数具有与（5-9）式相同的形式。    T>

【例5-4】 绘制Elliptic低通模拟原型滤波器的幅频平方响应曲线，阶数分别为2，te="mailto:wan">3te="mailto:wan">4，te="mailto:wan">6te="mailto:wan">4，te="mailto:wan">8te="mailto:wan">5。T>

%Samp5_4T>

n=0:0.01:2;   %设置频率点T>

for ii=1:4T>

   switch iiT>

   case 1,N=2;T>

   case 2,N=3;T>

   case 3,N=4;T>

   case 4,N=5;T>

   endT>

   Rp=1; Rs=15;   %设置通带波纹为1dB，阻带衰减为15dBT>

   [z,p,k]=ellipap(N,Rp,Rs);   %设计椭圆滤波器T>

   [b,a]=zp2tf(z,p,k);    %将零点极点增益形式转化为传递函数形式T>

   [H,w]=freqs(b,a,n);     %求得传递函数的复数频率响应T>

   magH2=(abs(H)).^2;     %给出幅度平方函数T>

   posplot=['2,2,' num2str(ii)];  %给出绘图区域T>

   subplot(posplot);T>

   plot(w,magH2); ylim([0 1])       %绘出幅度平方函数T>

   title(['N=' num2str(N)]);T>

   xlabel('w/wc');T>

   ylabel('椭圆 |H(jw)|^2');T>

   grid on;T>

endT>

 T>

T>

图5-9 椭圆模拟原型滤波器平方幅频图T>

 T>

程序运行结果见图5-9。可见阶数为4的椭圆滤波器的过渡带已相当窄（陡），但这种特性的获得是以牺牲通带和阻带的单调平滑特性为代价的。可以看到滤波器的阶数越高平方幅频响应越接近于矩形。T>

5.3.5  Bessel滤波器T>

前面讲过的各类原型滤波器均没有绘出其相位随频率的变化特性（相频特性）。在后面的数字信号处理学习中将会看到它们的相位特性是非线性的。本节所介绍的Bessel滤波器就能{zd0}限度地减少相频特性的非线性，使得通带内通过的信号形状不变(拷贝不走样)。T>

Bessel模拟低通滤波器的特点是在零频时具有最平坦的群延迟，并在整个通带内群延迟几乎不变。在零频时的群延迟为。由于这一特点，Bessel模拟滤波器通带内保持信号形状不变。但数字Bessel滤波器没有平坦特性，因此MATLAB信号处理工具箱只有模拟Bessel滤波器设计函数。T>

函数besselap用于设计Bessel模拟低通滤波器原型，调用格式为：T>

[z,p,k]=besselap(N)T>

式中，N为滤波器的阶数，应小于25。z,p,k为滤波器的零点、极点和增益。T>

滤波器的传递函数无零点，具有与（5-10）式相同的形式。下面用实例观看Bessel滤波器的幅频和相频特性。T>

【例5-5】绘制5阶和10阶Bessel低通滤波器原型的平方幅频和相频图。T>

%Samp5_5T>

clfT>

n=0:0.01:2;   %设置频率点T>

for ii=1:2T>

   switch iiT>

   case 1,pos=1;N=5;T>

   case 2,pos=3;N=10;T>

   endT>

   [z,p,k]=besselap(N);   %设计Bessel模拟滤波器T>

   [b,a]=zp2tf(z,p,k);   %将零点极点增益形式转化为传递函数形式T>

   [H,w]=freqs(b,a,n);   %求得传递函数的复数频率响应T>

   magH2=(abs(H)).^2;T>

   phaH=unwrap(angle(H));    %求得函数的相位角并进行解缠绕运算T>

   phaH=phaH*180/pi;       %将相位角由弧度转化为度T>

   posplot=['2,2,' num2str(pos)];   %设置绘图位置字符串T>

   subplot(posplot);T>

   plot(w,magH2); grid on       %绘出平方幅频响应T>

   title(['N=' num2str(N)]);T>

   xlabel('w/wc');   ylabel('Bessel |H(jw)|^2');T>

   grid onT>

   subplot(['2,2,' num2str(pos+1)]);T>

   plot(w,phaH);T>

   xlabel('w/wc');ylabel('Bessel 相位/^o'); T>

 title(['N=' num2str(N)]);T>

   grid onT>

endT>

图5-10 Bessel模拟原型滤波器相频图T>

 T>

可见，Bessel滤波器具有{zy}线性相频的特点，但这个特点的获得是以牺牲窄过渡带为代价的,即滤波器的幅频平方特性与矩形特性相差甚远。T>

对所有的模拟原型滤波器做一总结可知：Butterworth滤波器在通带和阻带内均具有平滑单调的特点，但在相同过渡带宽的条件下，该滤波器所需的阶数最多。Chebyshev I和II型滤波器在通带或阻带内具有波纹，但在相同过渡带宽的条件下，该滤波器所需的阶数比Butterworth滤波器要少。椭圆滤波器在通带和阻带内均有波纹出现，但在相同过渡带宽的条件下，该滤波器所需的阶数最少。Bessel滤波器具有最宽的过渡带，但具有{zy}的线性相频特性。因此没有{jd1}“好”的滤波器，要根据解决问题的不同选择不同的滤波器，因此，每一种滤波器的设计方法我们都要熟练掌握。T>

前面所讲的模拟原型滤波器均是截止频率为1的滤波器，在实际设计中是很难遇到的，然而它是设计其他各类滤波器的基础。我们通常遇到的是截止频率任意的低通滤波器、高通滤波器、带通滤波器和带阻滤波器。如何te="mailto:wan">以te="mailto:wan">由低通原型模拟滤波器te="mailto:wan">为基础设计这些滤波器呢？这就要用到我们今天要讲的频率变换。所谓频率变换是指各类滤波器（低通、高通、带通、带阻）和低通滤波器原型的传递函数中频率自变量之间的变换关系。通过频率变换，我们可以从模拟低通滤波器原型获得模拟的低通滤波器、高通滤波器、带通滤波器和带阻滤波器，再借助于s域至z域的变换关系又可以设计各类后面所讲的无限冲激响应数字滤波器，这是滤波器设计的重要方法之一。T>

MATLAB信号处理工具箱有lp2lp,lp2hp,lp2bp,lp2bs四个频率变换函数。下面分别叙述其使用方法及各参量的意义。T>

（1）       函数lp2lp用于实现由低通模拟原型滤波器至低通滤波器的频率变换，调用格式为T>

[bt,at]=lp2lp(b,a,)T>

其中，a,b为模拟原型滤波器的分母和分子多项式的系数，为低通滤波器所期望的截止频率（rad/s），若给定的单位为Hz，应乘以2。bt,at为返回的低通滤波器的分母和分子多项式的系数。该函数将模拟原型滤波器传递函数执行下面变换：T>

                                                    (5-15)T>

式中，H(p)为低通原型滤波器传递函数，H(s)为低通滤波器传递函数。该项操作可以执行模拟原型滤波器由截止频率为1到指定截止频率的变换，其原理讨论已超出本课程的范围，可参看其他信号处理参考书。下面的例子说明如何进行模拟原型低通滤波器变换为截止频率不为1的模拟低通滤波器。T>

【例5-6】将4阶椭圆模拟原型滤波器变换为截止频率为0.5的椭圆模拟低通滤波器，其中通带波纹Rp=2dB,阻带衰减Rs=30dB。T>

%Samp5_6T>

Rp=2;Rs=30;   %模拟原型滤波器的通带波纹为2dB,阻带衰减为30dB。T>

[z,p,k]=ellipap(4,Rp,Rs);   %设计椭圆滤波器T>

[b,a]=zp2tf(z,p,k);      %由零点极点增益形式转换为传递函数形式T>

[H,w]=freqs(b,a,0:0.01:2);   %给出复数频率响应T>

subplot(2,1,1),plot(w,abs(H).^2);   %绘出平方幅频函数T>

xlabel('w/wc');ylabel('椭圆 |H(jw)|^2');T>

title('原型低通椭圆滤波器(wc=1)')T>

[bt,at]=lp2lp(b,a,0.5);  %将模拟原型低通滤波器的截止频率变换为0.5T>

[Ht,wt]=freqs(bt,at,0:0.01:2);  %给出复数频率响应T>

subplot(2,1,2),plot(wt,abs(Ht).^2); %绘出平方幅频函数T>

xlabel('w/wc');ylabel('椭圆 |H(jw)|^2');T>

title('低通椭圆滤波器(wc=0.5)')T>

T>

         图5-11 将4阶椭圆模拟原型滤波器变换为截止频率为0.5的椭圆模拟低通滤波器的结果T>

 T>

程序的运行结果见图5-11。清楚地表明将截止频率由1变换到0.5。T>

（2）       函数lp2hp用于实现由低通模拟滤波器至高通滤波器的频率变换。调用格式：T>

[bt,at]=lp2hp(b,a,)T>

式中，为高通模拟滤波器所期望的截止频率（rad/s），若给定的频率单位为Hz，应乘以2。该函数将模拟原型滤波器传递函数执行下面变换： T>

                                                    (5-16)T>

【例5-7】将6阶Chebyshev I型原型滤波器变换为截止频率为0.8的模拟高通滤波器，其中通带波纹Rp=0.5dB。T>

%Samp5_7T>

Rp=0.5;   %设置滤波器的通带波纹为0.5dBT>

[z,p,k]=cheb1ap(6,0.5);   %设计Chebyshev I型模拟原型滤波器T>

[b,a]=zp2tf(z,p,k);    %由零点极点增益形式转化为传递函数形式T>

[H,w]=freqs(b,a,0:0.01:2);   %计算传递函数的复数频率响应T>

subplot(2,1,1),plot(w,abs(H).^2); %绘制传递函数的平方幅频响应T>

xlabel('w/wc');ylabel('Chebyshev I |H(jw)|^2');T>

title('Chebyshev I 低通原型滤波器(wc=1)')T>

[bt,at]=lp2hp(b,a,0.8);   %由低通原型滤波器转换为截止频率为0.8的高通滤波器T>

[Ht,wt]=freqs(bt,at,0:0.01:2);   %计算滤波器的复数频率响应T>

subplot(2,1,2),plot(wt,abs(Ht).^2);  %绘制传递函数的平方幅频响应T>

xlabel('w/wc');ylabel('Chebyshev I |H(jw)|^2');T>

title('Chebyshev I 高通滤波器(wc=0.8)')T>

T>

图5-12 将6阶Chebyshev I型原型滤波器变换为截止频率为0.8的模拟高通滤波器的结果图T>

 T>

程序运行结果见图5-12。可见实现了从低通滤波器到高通滤波器的转换，并且截止频率也符合例题要求。T>

（3）       函数lp2bp用于实现由低通模拟原型滤波器至带通滤波器的频率变换。调用格式：T>

[bt,at]=lp2bp(b,a,,Bw)T>

式中，为带通滤波器的中心频率(rad/s)，Bw为带通滤波器带宽(rad/s)。而T>

                                          (5-17)T>

式中, 为带通滤波器的下边界频率，为带通滤波器上边界频率。若给定的边界频率为Hz需乘以2。T>

该函数将模拟原型滤波器传递函数执行下面变换运算：T>

                                              (5-18)T>

下面用例题说明截止频率和滤波器类型的转换。这里要注意，输出的带通滤波器阶数为模拟原型滤波器阶数的2倍。T>

【例5-8】将6阶Chebyshev II型原型滤波器变换为模拟带通滤波器，其中上边界截止频率为0.8 rad/s,下边界截止频率为1.4 rad/s,阻带衰减Rs=20dB。T>

%Samp5_8T>

Rs=20;    %滤波器的阻带衰减为20dBT>

[z,p,k]=cheb2ap(6,Rs);   %设计6阶阻带衰减为20dB的Chebyshev II型原型滤波器T>

[b,a]=zp2tf(z,p,k);    %将零点极点增益形式转换为传递函数形式T>

[H,w]=freqs(b,a,0:0.01:2);   %计算模拟原型滤波器的复数频率响应T>

subplot(2,1,1),plot(w,abs(H).^2);  %绘制平方幅频响应T>

xlabel('w/wc');ylabel('Chebyshev II |H(jw)|^2');T>

title('Chebyshev II 型低通原型滤波器(wc=1)')T>

w1=0.8;w2=1.4;   %给定欲设计滤波器通带的下限和上限频率T>

w0=sqrt(w1*w2);bw=w2-w1; %计算中心点频率和频带宽度T>

[bt,at]=lp2bp(b,a,w0,bw);  %频率转换T>

[Ht,wt]=freqs(bt,at,0:0.01:2);  %计算滤波器的复数频率响应T>

subplot(2,1,2),plot(wt,abs(Ht).^2); %绘制平方幅频响应T>

xlabel('w/wc');ylabel('Chebyshev II |H(jw)|^2');T>

title('Chebyshev II型带通滤波器(wc=0.8~1.4)')T>

T>

图5-13 将6阶Chebyshev II型原型滤波器变换为模拟带通滤波器的结果图T>

 T>

程序运行结果为图5-13。可见该程序实现了由低通原型滤波器到带通滤波器的转换。T>

（4）       函数lp2bs用于实现由低通模拟原型滤波器至带阻滤波器的频率变换。调用格式T>

[bt,at]=lp2bs(b,a,,Bw)T>

式中，为带阻滤波器的中心频率(rad/s)，Bw为带阻滤波器带宽(rad/s)。而T>

                                          (5-19)T>

式中, 为带阻滤波器的下边界频率，为带阻滤波器上边界频率。若给定的边界频率为Hz需乘以2。T>

该函数将模拟原型滤波器传递函数执行下面变换运算：T>

                     (5-20)T>

 T>

注意：输出的带阻滤波器和带通滤波器是滤波器原型阶数的2倍。T>

【例5-9】 将6阶Butterworth原型滤波器变换为模拟带阻滤波器，其中上边界频率为0.7 rad/s,下边界频率为1.5 rad/s。T>

%Samp5_9T>

[z,p,k]=buttap(6);   %设计阶数为6的Butterworth模拟原型滤波器T>

[b,a]=zp2tf(z,p,k);  %将零点极点增益形式转换为传递函数形式T>

[H,w]=freqs(b,a,0:0.01:2); %计算指定频率点的复数频率响应T>

subplot(2,1,1),plot(w,abs(H).^2); %绘制平方幅频响应T>

xlabel('w/wc');ylabel('Butterworth |H(jw)|^2');T>

title('Btterworth低通原型滤波器(wc=1)')T>

w1=0.7;w2=1.5;   %带阻滤波器的下限及上限频率T>

w0=sqrt(w1*w2);bw=w2-w1;   %计算中心点频率和阻带宽度T>

[bt,at]=lp2bs(b,a,w0,bw);   %频率转换T>

[Ht,wt]=freqs(bt,at,0:0.01:2);  %计算带阻滤波器指定点的复数频率响应T>

subplot(2,1,2),plot(wt,abs(Ht).^2);%绘制平方幅频响应T>

xlabel('w/wc');ylabel('Butterworth |H(jw)|^2');T>

title('Butterworth带阻滤波器(wc=0.7~1.5) ')T>

T>

图5-14 将6阶Butterworth原型滤波器变换为模拟带阻滤波器的结果图T>

 T>

程序运行结果为图5-10。可见，该程序实现了由低通原型滤波器到带阻滤波器之间的转换，并且频带参数与所要求的xx一致。T>

前面所述的模拟滤波器设计中，滤波器阶数是我们在编程时任意指定的。其实它是决定滤波器品质的主要参数之一。通常在满足性能指标的前提下，阶数应该尽可能小，以满足易于实现、提高运算速度的要求。而在滤波器阶数和滤波器性能之间存在一定的函数关系，我们通过这一函数关系可以求出满足滤波性能指标的{zd1}阶数。下面介绍常用低通滤波器最小阶数的确定原理及MATLAB实现；接着介绍MATLAB信号处理工具箱中用来计算最小阶数和截止频率的工具函数。T>

5.5.1  滤波器最小阶数选择原理：以Butterworth低通模拟滤波器为例T>

模拟低通滤波器的设计指标为：通带边界频率，阻带边界频率，通带波纹Rp(dB)、阻带衰减Rs(dB)。T>

当时T>

                                 (5-21)T>

以截止频率（幅值下降3dB）为1，化为相对的相对频率，则上式可写成：T>

                       (5-22)T>

同理，当时，T>

                                             (5-23)T>

综合上面两式可得：T>

                     (5-24)T>

式中，N应向上取整，则N为该滤波器的最小阶数。表示为：T>

                         (5-25)T>

或T>

                                                   (5-26)T>

下面用MATLAB编程可以计算滤波器的最小阶数N和截止频率。T>

【例5-10】设计一个模拟Butterworth滤波器。设计指标为通带边界频率，阻带边界频率，通带波纹1dB，阻带衰减16dB,试确定最小阶数N和截止频率。T>

%Samp5_10T>

wp=200*pi;ws=300*pi;Rp=1;Rs=16;%滤波器的通带阻带边界频率、通带波纹和阻带衰减T>

N=ceil(log10((10^(Rp/10)-1)/(10^(Rs/10)-1))/(2*log10(wp/ws))) T>

% （5-22）式，这里ceil函数表示向上取整T>

Wcp=wp/((10^(Rp/10)-1)^(1/(2*N)))   %（5-23）式T>

Wcs=ws/((10^(Rs/10)-1)^(1/(2*N)))   %（5-24）式T>

给程序的输出为：T>

N =   7T>

Wcp =  691.9839T>

Wcs =  725.7292T>

程序输出的两个截止频率分别为运用通带波纹和阻带衰减得到的。为使信号保留较多的有用信息，对于低通滤波器可选用较大的截止频率，即725.7292。T>

其他几种滤波器的最小阶数选择思路与此相同，由于篇幅所限，我们略去。T>

5.5.2  滤波器最小阶数选择函数T>

上面给出了Butterworth滤波器的最小阶数和截止频率的选择公式及程序。其实MATLAB工具箱中运用滤波器的最小阶数选择公式给出了滤波器最小阶数选择函数。几种滤波器最小阶数的选择函数如下：T>

[n,wc]=buttord(wp,ws,Rp,Rs, 's');     Butterworth 滤波器T>

[n,wc]=cheb1ord(wp,ws,Rp,Rs, 's');     Chebyshev I 滤波器T>

[n,wc]=cheb2ord(wp,ws,Rp,Rs, 's');     Chebyshev II 滤波器T>

[n,wc]=ellipord(wp,ws,Rp,Rs, 's');     Elliptic 滤波器T>

式中，wp为通带边界频率，ws为阻带边界频率，单位为rad/s。Rp,Rs分别为通带波纹和阻带衰减，单位为dB。二者分别表示通带内的{zd0}允许幅值损失和阻带下降的分贝数。's'表示模拟滤波器（缺省时，该函数适用于数字滤波器）；函数返回值n为模拟滤波器的最小阶数；wc为模拟滤波器的截止频率，单位为rad/s。这四个函数适用于低通、高通、带通、带阻滤波器。T>

若wp<ws，对应于低通模拟滤波器，当wp>ws时对应于高通模拟滤波器，对于带通和带阻滤波器存在两个过渡带，wp和ws均应该为含有两个元素的向量，分别表示两个过渡带的边界频率。即：wp=[通带下界频率, 通带上界频率]，ws=[阻带下界频率, 阻带上界频率]。对于带通滤波器，这四个频带界线的大小排列为:阻带下界频率<通带下界频率<通带上界频率<阻带上界频率；对于带阻滤波器，这四个频带界线的大小排列为:通带下界频率<阻带下界频率<阻带上界频率<通带上界频率。这时返回值wc包括两个元素（{dy}个元素小于第二个元素），分别为通带和阻带之间的界线频率。函数自动判断是带通还是带阻滤波器T>

【例5-11】设计一个模拟低通Butterworth滤波器和Chebyshev II型滤波器，满足通带边界频率，阻带边界频率，通带波纹Rp=1dB，阻带衰减Rs=16dB,试确定最小阶数N和通带下降3dB的截止频率。T>

%Samp5_11T>

Wp=200*pi;Ws=300*pi;Rp=1;Rs=16;%滤波器的通带阻带截止频率、通带波纹和阻带衰减T>

disp('Butterworth滤波器:'); %在命令窗口显示括号内的文本T>

[N,Wc]=buttord(Wp,Ws,Rp,Rs,'s') %计算Butterworth滤波器的最小阶数和截止频率T>

disp('Chebyshev II型滤波器:');T>

[N,Wc]=cheb2ord(Wp,Ws,Rp,Rs,'s') %计算ChebyshevII滤波器的最小阶数和截止频率T>

程序的输出为：T>

Butterworth滤波器:T>

N =7T>

Wc =725.7292T>

Chebyshev II型滤波器:T>

N =4T>

Wc =  839.8376T>

程序输出的Butterworth滤波器的最小阶数和截止频率与例5-10的程序输出相同。证明调用求取最小阶数和截止频率的程序就是采用了（5-24）~（5-26）的公式。运用Chebyshev II型滤波器所需的滤波器的最小阶数远小于Butterworth滤波器，说明运用Chebyshev II型滤波器比Butterworth滤波器更容易实现、计算速度更快。大家可以自己设计一些例子运用本节给出的求取滤波器最小阶数的函数求解滤波器最小阶数和截止频率。T>

滤波器设计好之后，一般要进行各方面的测试。在正式设计滤波器之前，我们先介绍如何测试滤波器的性能。对于模拟滤波器，在本章第二节我们已采用函数freqs来求模拟滤波器的频率响应，这里我们详细介绍该函数，若其调用格式为：T>

[h[,w]]=freqs(b,a[,n])T>

式中，b,a分别为模拟滤波器传递函数分子和分母多项式系数； h为对应频率点的传递函数值。上面的表示中，[,w]和[,n]表示可有可无的参数，本书中的该类表示均为可有可无的参数, 由用户根据需要确定。w为频率点的值,n为频率点数。若n=128，则用128个频率点来给出此模拟滤波器的频率特性（给定频率点的传递函数值），缺省时为200。若该函数不写输出变量，则执行后绘出该滤波器的幅频响应和相位响应图。此函数模拟滤波器的传递函数形式为：T>

                       (5-27)T>

MATLAB工具箱还提供了两个函数abs和angle，由频率响应求幅频响应和相频响应。其中angle的输出单位为rad。可采用rad2deg函数转化为度。另外注意函数的幅频响应经常用分贝(dB)来表示。te="mailto:wan">T>

求出的幅频响应te="mailto:wan">后，通过下式转换为分贝:(dB)。T>

【例5-12】已知模拟滤波器的传递函数为：，绘制滤波器的幅频响应和相频响应。T>

%Samp5_12T>

b=[0.2 0.3 1];   %滤波器传递函数分子多项式系数T>

a=[1 0.4 1]; %滤波器传递函数分母多项式系数T>

figure(1),freqs(b,a);  %{dy}种输出方法T>

[h,w]=freqs(b,a); %计算滤波器的复数频率响应T>

mag=abs(h);pha=angle(h);  %得到滤波器的幅频和相频响应T>

figure(2),subplot(2,1,1),loglog(w,mag);  %运用双对数坐标绘制幅频响应T>

grid on;xlabel('角频率/rad\cdots^-^1');ylabel('振幅');T>

subplot(2,1,2),semilogx(w,pha*180/pi) %运用半对数坐标绘相频响应T>

grid on;xlabel('角频率/rad\cdots^-^1');T>

ylabel('相位/^o');T>

 T>

T>

图5-15 例5-12中滤波器的幅频响应和相频响应 T>

 T>

程序运行结果的两幅图xx一致（图5-15），说明freqs的绘图程序与本例绘图程序xx一致。T>

我们知道，除了用传递函数描述滤波器特性外，还可用脉冲(冲激)响应来描述滤波器，因为在模拟滤波器中，脉冲响应与传递函数是Laplace变换对。此外还可以用阶跃响应（输入一个阶跃时系统的输出）来描述滤波器特性。下面我们介绍在MATLAB中如何得到模拟滤波器的脉冲响应和阶跃响应。T>

将滤波器的传递函数表示成分子和分母多项式系数的形式，如分子和分母多项式的系数为[b,a]，则在MATLAB中用H=[tf(b,a)]来表示此模拟滤波器，采用[[y,t]=]impulse(H)给出该系统的模拟脉冲响应。采用[[y,t]=]step(H)来得到该系统的阶跃响应。这两个函数与freqs一样，若没有输出则程序自动绘图模拟该滤波器的脉冲响应或阶跃响应。输出值y,t分别为该滤波器的脉冲响应或阶跃响应及其对应的时间序列。T>

另外，我们还可以运用一个输入信号模拟该滤波器的输出。若给定滤波器的输入值序列和对应的时间序列为u,t，则系统的输出可用y=lsim(H,u,t)来模拟, y为对应t的输出。若该函数没有输出变量则程序自动绘图显示。T>

【例5-13】设计一个5阶的Chebyshev I型带通滤波器，通带波纹3dB，下边界频率100Hz，上边界频率500Hz，绘制幅频响应图。给出该滤波器的脉冲响应、阶跃响应。假定输入sin(2*pi*30*t)+0.5*cos(2*pi*300*t)+2*sin(2*pi*800*t)的信号，求模拟滤波器的输出并给出模拟输入信号和模拟输出信号的Fourier振幅谱。T>

%Samp5_13T>

N=5;Rp=3;f1=100;f2=500;  %滤波器阶数、边界频率（Hz）T>

w1=2*pi*f1;w2=2*pi*f2;    %边界频率（rad/s）T>

[z,p,k]=cheb1ap(N,Rp);   %设计Chebyshev I型原型低通滤波器T>

[b,a]=zp2tf(z,p,k);  %转换为传递函数形式T>

Wo=sqrt(w1*w2);    %中心频率T>

Bw=w2-w1;    %频带宽度T>

[bt,at]=lp2bp(b,a,Wo,Bw);   %频率转换T>

[h,w]=freqs(bt,at);    %计算复数频率响应T>

figure(1)T>

subplot(2,2,1),semilogy(w/2/pi,abs(h));   %绘制幅频响应T>

xlabel('频率/Hz');grid on;title('幅频图');T>

subplot(2,2,2),plot(w/2/pi,angle(h)*180/pi);%绘制相频响应T>

xlabel('频率/Hz');T>

ylabel('相位/^o');grid on;title('相频图')T>

H=[tf(bt,at)];  %在MATLAB中表示此滤波器T>

[h1,t1]=impulse(H);  %te="mailto:wan">输出te="mailto:wan">绘te="mailto:wan">出系统的脉冲响应te="mailto:wan">图T>

subplot(2,2,3),plot(t1,h1);xlabel('时间/s');title('脉冲响应')T>

[h2,t2]=step(H);     %te="mailto:wan">输te="mailto:wan">绘出系统的阶跃响应te="mailto:wan">图T>

subplot(2,2,4),plot(t2,h2);xlabel('时间/s');title('阶跃响应')T>

figure(2)T>

dt=1/2000;T>

t=0:dt:0.1;    %给出模拟滤波器输出的时间范围T>

u=sin(2*pi*30*t)+0.5*cos(2*pi*300*t)+2*sin(2*pi*800*t);  %模拟输入信号T>

subplot(2,2,1),plot(t,u)  %绘制模拟输入信号T>

xlabel('时间/s');title('输入信号')   T>

[ys,ts]=lsim(H,u,t);   %模拟系统的输入u时的输出T>

subplot(2,2,2),plot(ts,ys)  %绘制模拟输出信号T>

xlabel('时间/s'),title('输出信号');  T>

subplot(2,2,3),plot((0:length(u)-1)/(length(u)*dt),abs(fft(u))*2/ length(u));%绘输入信号振幅谱T>

title('输入信号振幅谱'),xlabel('频率/Hz')T>

Subplot(2,2,4),T>

Y=fft(ys);T>

plot((0:length(Y)-1)/(length(Y)*dt),abs(Y)*2/length(Y)); %绘制输出信号振幅谱T>

title('输出信号振幅谱')T>

xlabel('频率/Hz')T>

T>

图5-16 例5-13设计的Chebyshev I型滤波器的幅频响应、相频响应、脉冲响应和阶跃响应T>

 T>

T>

图5-17 例5-13所设计滤波器模拟输入和输出信号的时间域图形和振幅谱T>

 T>

程序输出结果见图5-16和图5-17。图5-16给出了该程序得到的幅频图、相频图、脉冲响应和阶跃响应。幅频图清楚地给出了te="mailto:wan">的通带范围和阻带范围。图5-17给出了模拟输入、输出信号的时间域波形及振幅谱。输入信号中含有30Hz,300Hz和500Hz的波，由滤波器的幅频特性（图5-16左上图）可知，30Hz和500Hz全在阻带内，只有300Hz的波可通过滤波器。从输入信号和模拟输出信号的时间域波形和振幅谱，可以看到达到了预期的效果。注意，由该滤波器的相频特性（图5-16的相频图）可知，该滤波器并不是线性相位。T>

5.7.1  模拟滤波器设计步骤T>

用户对设计的滤波器提出设计要求，我们可以针对滤波器的设计要求设计滤波器。通常用户对模拟滤波器提出的要求有：（1）滤波器的性能指标，包括截止频率(对于低通和高通)或上下边界频率、，通带波纹、阻带衰减等；（2）滤波器的类型，通常为Butterworth、Chebyshev I、 Chebyshev II、 Elliptic或Bessel滤波器。我们根据滤波器的类型通常按下列步骤设计滤波器。T>

（1）       给定模拟滤波器的性能指标，如截止频率(对于低通和高通)或上下边界频率、；通带波纹、阻带衰减以及滤波器类型等（用户给定）。T>

（2）       确定滤波器阶数（根据本章第5节介绍的方法和函数）。T>

（3）       设计模拟低通原型滤波器。MATLAB信号处理工具箱的滤波器原型设计函数有butterap、cheb1ap、cheb2ap、ellipap、besselap。（本章第5.3节介绍的方法和函数）T>

（4）       按频率变换设计模拟滤波器（低通、高通、带通、带阻）。MATLAB信号处理工具箱的频率变换函数有lp2lp、lp2hp、lp2bp、lp2bs。（本章第5.4节介绍的方法和函数）。T>

下面我们用例子说明如何根据用户的要求设计滤波器。T>

【例5-14】设计带通Chebyshev I型模拟滤波器Ws1=0.2*pi（rad/s）; Ws2=0.8*pi（rad/s）; Rs=60dB; Wp1=0.35*pi（rad/s）; Wp2=0.65*pi（rad/s）; Rp=1dB。 T>

%Samp5_14T>

wp=[0.35 0.65]*pi;   %给出滤波器的通带边界频率，用弧度表示T>

ws=[0.2 0.8]*pi;   %给出滤波器的阻带边界频率，te="mailto:wan">，用弧度表示T>

Rp=1;Rs=60;   %给出滤波器的通带波纹和阻带衰减T>

[n,wn]=cheb1ord(wp,ws,Rp,Rs,'s');  %根据滤波器的要求求滤波器的最小阶数T>

[z,p,k]=cheb1ap(n,Rp);  %用求得的滤波器的最小阶数设计Chebyshev I型滤波器T>

[b,a]=zp2tf(z,p,k); %将零点极点增益形式转换为传递函数形式T>

Wo=sqrt(wn(1)*wn(2)); %计算中心点频率T>

Bw=wn(2)-wn(1);  %计算频带宽度T>

[bt,at]=lp2bp(b,a,Wo,Bw);  %模拟原型低通滤波器转换为带通滤波器T>

[h,w]=freqs(bt,at,128);  %计算滤波器的复数频率响应T>

plot(w/pi,20*log10(abs(h)));   %绘制滤波器的振幅响应T>

xlabel('角频率/pi');T>

ylabel('|H(jw)|/dB')T>

grid on;hold on;T>

plot([0.2 0.2],ylim, 'r'); plot([0.8 0.8],ylim, 'r');%绘出阻带界限以显示达到的阻带衰减T>

%plot([0.35 0.35],ylim); plot([0.65 0.65],ylim);%可绘出通带界限以显示达到的通带衰减T>

 T>

T>

图5-18  例5-14设计的Chebyshev I型带通滤波器的幅频响应T>

 T>

程序的运行结果见图5-18。可以看到，设计的滤波器的通带范围确实为0.35*pi~0.65*pi,并且在阻带边界处下降分贝数大于60dB，满足了用户的设计要求。程序中ylim函数是提取现行坐标轴y轴的下限和上限。程序用它绘出了阻带边界。T>

5.7.2  模拟滤波器设计函数T>

上面滤波器的设计步骤比较麻烦，根据设计要求求解滤波器的最小阶数和边界频率之后需要设计模拟原型滤波器并进行频率转换。其实MATLAB将这一系列的过程组合成了更为方便的设计函数：butter,cheby1,cheby2,ellip,besself。这些函数称为模拟滤波器xx设计函数。用户在求得滤波器的最小阶数和截止频率之后只需调用一次xx设计函数就可以自动完成所有设计过程，编程十分简单。这些工具函数适用于模拟滤波器的设计，但同样也适用于数字滤波器。本节只讨论这些函数在模拟滤波器设计中的应用。但要注意，MATLAB是将上述一系列的步骤复合而已，并不是一种新的设计模拟滤波器的方法。T>

[b,a]=butter(n,wn[,'ftype'], 's')T>

[z,p,k]=butter(n,wn[,'ftype'], 's')T>

[b,a]=cheby1(n,Rp,wn[,'ftype'], 's')T>

[z,p,k]=cheby1(n,Rp,wn[,'ftype'], 's')T>

[b,a]=cheby2(n,Rs,wn[,'ftype'], 's')T>

[z,p,k]=cheby2(n,Rs,wn[,'ftype'], 's')T>

[b,a]=ellip(n,Rp,Rs,wn[,'ftype'], 's')T>

[z,p,k]=cheby2(n,Rp,Rs,wn[,'ftype'], 's')T>

[b,a]=besself(n,wn[,'ftype'], 's')T>

[z,p,k]=besself(n,wn[,'ftype'], 's')T>

在上面的调用方式中，n为滤波器的阶数，wn为滤波器的截止频率，单位rad/s（wn>0）；’s’为模拟滤波器，缺省时为数字滤波。T>

   ‘ftype’滤波器的类型可取为：T>

    ‘high’高通滤波器，截止频率为wn。T>

    ‘stop’带阻滤波器，截止频率为wn=[w1,w2] (wte="mailto:wan"><te="mailto:wan">1>w2)。T>

‘ftype’缺省时为低通或带通滤波器。T>

当设计带通滤波器时，截止频率也为wn=[w1,w2] (w1te="mailto:wan"><te="mailto:wan">>w2)。a,b分别为滤波器的传递函数分子和分母多项式系数向量；z,p,k分别为滤波器的零点、极点和增益。Rp,Rs分别为所设计滤波器的通带波纹和阻带衰减, 单位为dB。T>

滤波器的传递函数具有下面的形式：T>

                        (5-26)T>

若滤波器为带通或带阻型，则滤波器的阶数为2n，否则阶数为n。下面用一些例子说明滤波器的设计。T>

【例5-15】设计一个Butterworth模拟带通滤波器，设计指标为：通带频率：1000~2000Hz，两侧过渡带宽500Hz，通带波纹1dB，阻带衰减100dB。假设一个信号，其中f1=100Hz,f2=1500Hz,f3=2900Hz。信号的采样频率为10000Hz。试将原信号与通过该滤波器的模拟信号进行比较。T>

%Samp5_15T>

wp=[1000 2000]*2*pi;ws=[500 2500]*2*pi;Rp=1;Rs=100; %滤波器设计参数，对于给定Hz应乘以2T>

[N,Wn]=buttord(wp,ws,Rp,Rs,'s'); %求得滤波器的最小阶数和截止频率  T>

w=linspace(1,3000,1000)*2*pi;   %设置绘制频率响应的频率点T>

[b,a]=butter(N,Wn,'s');    %设计模拟Butterworth滤波器T>

H=freqs(b,a,w);    %计算给定频率点的复数频率响应T>

magH=abs(H);phaH=unwrap(angle(H)); %计算幅频响应和相频响应T>

plot(w/(2*pi),20*log10(magH));  %以频率为横坐标绘制幅频响应T>

xlabel('频率/Hz');ylabel('振幅/dB');T>

title('Butterworth模拟带通滤波器');T>

hold on;plot([1000 1000],ylim,'r');plot([2000 2000],ylim,'r');%绘te="mailto:wan">通te="mailto:wan">阻带边界grid onT>

figure(2)T>

dt=1/10000;  %模拟信号采样间隔T>

f1=100;f2=1500;f3=2900;%输入信号的三个频率成分T>

t=0:dt:0.04;  %给定模拟时间段T>

x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号T>

H=[tf(b,a)];     %滤波器在MATLAB系统中的表示T>

[y,t1]=lsim(H,x,t);   %模拟输出T>

subplot(2,1,1),plot(t,x),title('输入信号') %绘出输入信号T>

subplot(2,1,2),plot(t1,y)   %绘制输出信号T>

title('输出信号'),xlabel('时间/s')T>

T>

图5-19   例5-15设计的模拟滤波器幅频响应，两条竖线表示阻带边界T>

 T>

T>

图5-20  例5-15中滤波器的输入信号和输出信号T>

 T>

程序的运行结果为图5-19和图5-20。可见，滤波器的幅频响应xx符合要求。在阻带边界确实下降到100dB。大家修改程序可以观看通带边界处下降的分贝数不超过1个dB。模拟的输入信号中含有三个频率成分，根据图5-19可知，2900Hz和100Hz的频率均在阻带内，只有1000Hz的频率在通带内，通过模拟系统后，将2900Hz的高频成分和100Hz的低频成分滤除，输出只含有1500Hz频率成分的振动。T>

【例5-16】 设计一个Chebyshev I型模拟带阻滤波器，设计指标为：te="mailto:wan">阻te="mailto:wan">通带频率：1000～2000Hz，两侧过渡带宽500Hz，通带波纹1dB，阻带衰减大于50dB。假设一个信号，其中f1=100Hz,f2=1500Hz,f3=2900Hz。信号的采样频率为20000Hz。试将原信号与通过该滤波器的模拟信号进行比较。T>

%Samp5_16T>

close all;figure(1)T>

ws=[1000 2000]*2*pi;wp=[500 2500]*2*pi;Rp=1;Rs=50;  %滤波器设计参数，乘以2转换为弧度/sT>

[N,Wn]=cheb1ord(wp,ws,Rp,Rs,'s');   %求模拟滤波器的最小阶数和截止频率T>

w=linspace(1,3000,1000)*2*pi;   %设置计算复数频率响应的频率点T>

[b,a]=cheby1(N,Rp,Wn,'stop','s');  %设计带阻滤波器T>

H=freqs(b,a,w);    %计算给定频率点的复数频率响应T>

magH=abs(H);phaH=unwrap(angle(H));  %求滤波器的幅频响应和相频响应T>

plot(w/(2*pi),20*log10(magH));  %以频率为横轴绘制幅频响应T>

xlabel('频率/Hz');ylabel('振幅/dB');T>

title('Butterworth模拟带阻滤波器');T>

hold on;plot([1000 1000],ylim, 'r');plot([2000 2000],ylim, 'r');%绘阻带边界T>

grid onT>

figure(2)T>

dt=1/20000;T>

f1=100;f2=1500;f3=2900;   %输入信号含有的频率T>

t=0:dt:0.04;%时间序列T>

x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号T>

H=[tf(b,a)];   %滤波器在MATLAB系统中的表示T>

[y,t1]=lsim(H,x,t);   %模拟输出信号T>

subplot(2,1,1),plot(t,x),title('输入信号')  %绘输入信号T>

subplot(2,1,2),plot(t1,y)   %绘制输出信号T>

title('输出信号'),xlabel('时间/s')T>

T>

图5-21 例5-16设计的带阻滤波器的幅频响应T>

图中，两条竖线表示阻带边界T>

 T>

T>

图5-22 例5-16检测滤波器的输入输出信号T>

 T>

程序的输出结果见图5-21和图5-22。可见设计的带阻滤波器xx符合要求，阻带边界达到50dB的衰减。模拟的输出信号xx滤除了1500Hz的信号，模拟输出信号中只含用100Hz和2900Hz的信号。T>

【例5-17】设计一个Chebyshev II型模拟带阻滤波器，设计指标与上例相同。T>

%Samp5_17T>

ws=[1000 2000]*2*pi;wp=[500 2500]*2*pi;Rp=1;Rs=50;  %滤波器设计参数T>

[N,Wn]=cheb2ord(wp,ws,Rp,Rs,'s');   %求滤波器的最小阶数和截止频率T>

w=linspace(1,3000,1000)*2*pi;   %设置计算复数频率响应的频率点T>

[b,a]=cheby2(N,Rs,Wn,'stop','s');  %设计带阻滤波器T>

H=freqs(b,a,w);   %计算给定频率点的复数频率响应T>

magH=abs(H);phaH=unwrap(angle(H));  %计算幅频响应和相频响应T>

plot(w/(2*pi),20*log10(magH));  %以频率为横轴绘制幅频响应T>

xlabel('频率/Hz');ylabel('振幅/dB');T>

title('Chebyshev II 模拟带阻滤波器')T>

hold on; plot([1000 1000],ylim, 'r');plot([2000 2000],ylim, 'r');%阻带边界T>

grid onT>

T>

图5-23  例5-17设计滤波器的幅频响应图T>

 T>

程序运行结果为图5-23。可以看出滤波器在阻带边界处达到了下降50dB的要求，此滤波器在阻带内有振荡，这正是Chebyshev II型滤波器的特点。T>

【例5-18】设计一个高通椭圆滤波器，设计性能指标为：通带边界频率wp=1500Hz,阻带边界频率ws=1000Hz,通带波纹Rp=1dB,阻带衰减Rs=30dB。假设一个信号，其中f1=400Hz,f2=1600Hz。信号的采样频率为10000Hz。试将原信号与通过该滤波器的模拟信号进行比较。T>

%Samp5_18T>

wp=1500*2*pi;ws=1000*2*pi;Rp=1;Rs=30;   %滤波器设计参数T>

[N,Wn]=ellipord(wp,ws,Rp,Rs,'s');   %求滤波器的最小阶数和截止频率T>

w=linspace(1,3000,1000)*2*pi;   %给出计算复数频率响应的频率点T>

[b,a]=ellip(N,Rp,Rs,Wn,'high','s');  %设计高通椭圆滤波器T>

H=freqs(b,a,w);   %计算给定频率的复数频率响应T>

magH=abs(H);phaH=unwrap(angle(H));  %求幅频响应和相频响应T>

figure(1)T>

plot(w/(2*pi),20*log10(magH));  %以频率为横轴绘幅频响应T>

xlabel('频率/Hz');ylabel('振幅/dB');T>

title('椭圆模拟高通滤波器');T>

hold on;plot([1000 1000],ylim, 'r'); %阻带边界T>

grid onT>

figure(2)   %给出另一个图形窗口T>

dt=1/10000;   %采样间隔T>

f1=400;f2=1600;%信号中所含频率成分T>

t=0:dt:0.04;%时间序列T>

x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t);  %输入信号T>

H=[tf(b,a)];   %滤波器在MATLAB系统中的表示T>

[y,t1]=lsim(H,x,t);  %模拟输出T>

subplot(2,1,1),plot(t,x),title('输入信号') %绘输入信号T>

subplot(2,1,2),plot(t1,y)  %绘制输出信号T>

title('输出信号'),xlabel('时间/s')T>

T>

图5-24  例5-18设计滤波器的幅频响应图T>

 T>

T>

图5-25例5-18滤波器的输入和输出信号T>

 T>

程序运行结果见图5-24和图5-25。可见所设计的高通滤波器xx符合要求，在阻带边界下降的分贝数达30dB。并且通带和阻带均有波纹，这正是椭圆滤波器的特点。将该滤波器模拟输入400Hz和1600Hz的信号后，输出结果xx滤除了400Hz低频成分，只含有1600Hz的信号。T>

【例5-19】设计一个五阶Bessel低通模拟滤波器，截止频率为te="mailto:wan">200Hzte="mailto:wan">1000rad/s,绘制滤波器的频率特性图。假设一个信号，其中f1=100Hz,f2=1000Hz。信号的采样频率为10000Hz。试将原信号与通过该滤波器的模拟信号进行比较。T>

%Samp5_19T>

N=5;  %滤波器阶数T>

Wn=te="mailto:wan">200*2*pite="mailto:wan">1000;  %边界频率T>

[b,a]=besself(N,Wn); %设计Bessel滤波器T>

figure(1)T>

[H,w]=freqs(b,a,512);   %绘出滤波器的复数频率特性T>

magH=abs(H);phaH=unwrap(angle(H));  %求幅频响应和相频响应T>

subplot(2,1,1),plot(w/(2*pi),20*log10(magH));  %以频率为横轴绘幅频响应T>

grid on;xlabel('频率/Hz');ylabel('振幅/dB');T>

subplot(2,1,2),plot(w/2/pi,angle(H)*180/pi);%绘制相频响应T>

grid on;xlabel('频率/Hz');ylabel('相位/^o')T>

figure(2)    %新的绘图框T>

dt=1/10000;   %信号采样间隔T>

f1=100;f2=1000;  %输入两个频率成分T>

t=0:dt:0.1;T>

x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t); %输入信号T>

H=[tf(b,a)];   %滤波器在MATLAB系统中的表示T>

[y,t1]=lsim(H,x,t);  %模拟滤波器的输出T>

subplot(2,1,1),plot(t,x),title('输入信号')  %绘出输入信号T>

subplot(2,1,2),plot(t1,y)  %绘出输出信号T>

title('输出信号'),xlabel('时间/s') T>

te="mailto:wan">te="mailto:wan">T>

图5-26 例5-19设计Bessel滤波器的频率特性T>

 T>

T>

图5-27 例5-19中Bessel滤波器的输入和输出T>

 T>

 T>

程序运行结果为图5-26和图5-27。虽然模拟Bessel滤波器的幅频特性不如其他滤波器好，但其相频特性为线性的。因此图5-27中的输出信号与输入的低频信号形状不变，只是有延迟。T>

第五章  习题T>

1.     模拟滤波器可分为哪几种类型？设计时要求的主要性能指标有哪些？T>

2.     设计程序研究本章第二节设计Butterworth、Chebyshev I,II、椭圆滤波器例题中滤波器的非线性相位特性，并与Bessel滤波器的线性相位特性进行比较。T>

3.     设计一个Butterworth模拟低通滤波器，要求满足下面指标：截止频率wc=1000Hz,通带波纹小于3dB,阻带始点频率ws=1500Hz，阻带衰减大于35dB。（1）确定滤波器的最小阶数。（2）滤波器的传递函数；（3）绘制滤波器的平方幅频和相频特性图。（4）当输入信号为x(t)=sin(2*pi*f1*t)+cos(2*pi*f2*t),f1=500Hz,f2=2000Hz时，绘制滤波器输出信号曲线。（5）绘制滤波器的单位脉冲响应曲线。T>

4.     确定三阶Chebyshev I型低通模拟滤波器的传递函数，性能指标为：通带波纹为1dB，截止频率为wc=2rad/s。T>

5.     低通模拟滤波器满足下列条件：边界频率为500Hz和750Hz,通带波纹为2dB,阻带衰减为30dB。试求满足此要求的Butterworth,Chebyshev I, ChebyshevII和te="mailto:wan">椭圆te="mailto:wan">Elliptic滤波器的最小阶数。T>

6.     设计一个高通滤波器，边界频率为100Hz和150Hz，通带波纹不大于2dB,阻带衰减不小于30dB。试设计满足此要求的Butterworth, Chebyshev I, ChebyshevII和te="mailto:wan">Ellipticte="mailto:wan">椭圆滤波器并绘制它们的平方幅频响应图。绘制此系统的脉冲响应和阶跃响应。T>

7.     设计满足下列条件的带通Butterworth滤波器：中心频率1kHz,带宽为200Hz,通带波纹在1dB以下，阻带衰减在40dB以上。确定滤波器的阶数、系统函数并绘制幅频特性、相频特性和零极点图。T>

8.     FSS-3地震计零极点为：零点：[0.0,0.0,1.127E3, -6.25E+01+j*6.307089E+02, -6.25E+01-j*6.307089E+02, -1.827E+02]。极点：[-4.442212-j*4.443554, -4.442212+j*4.443554,-1.3645E+02-j*1.137163E+02,-1.3645E+02+j*1.137163E+02,-5.51E+01-j*4.409807E+02, -5.51E+01+j*4.409807E+02]。归一化因子：-7.342E-2。增益：4.405E+2。将增益与归一化因子相乘作为计算滤波器的增益。试绘出该地震计的频率特性。T>

9.     FSS-3型地震仪系统传递函数为：零点：[1.355083E-4+j*1.043789E-2, 1.355083E-4-j*1.043789E-2, 1.684441E2+j*7.147554E1, 1.684441E2-j*7.147554E1, 6.572567+j*2.38467E2, 6.572567-j*2.38467E2, 8.87029E1+j*1.476857E+2, 8.87029E1-j*1.476857E+2, 3.297725E1+j*1.584205E2, 3.297725E1-j*1.584205E2, 1.018952E1+j*1.653736E2, 1.018952E1-j*1.653736E2, 2.924976+j*1.530131E2, 2.924976-j*1.530131E2, 9.947939+j*1.550421E2, 9.947939-j*1.550421E2]; 极点：[-4.442212+j*4.443554, -4.442212-j*4.443554, -8.49856E3, -7.043488E1+j*3.404871E2, -7.043488E1-j*3.404871E2, -4.773638+j*1.452003E2, -4.773638-j*1.452003E2, -6.337485+j*1.39173E2, -6.337485-j*1.39173E2, -7.902612+j*1.330283E2, -7.902612-j*1.330283E2, -2.083218E1+j*1.247569E2, -2.083218E1-j*1.247569E2, -2.792995E1+j*1.203601E2, -2.792995E1-j*1.203601E2, -3.405637E01+j*9.91384E1, -3.405637E01-j*9.91384E1, -4.302115E1+j*7.238421e1, -4.302115E1-j*7.238421e1, -4.992166e1+j*3.912072e1, -4.992166e1-j*3.912072e1, -5.206292e1];归一化因子te="mailto:wan">：te="mailto:wan">为1.264053E12; 增益te="mailto:wan">为：4.405E2, 试绘出该地震系统的频率特性。T>

10.  SLJ-100型强震加速度计传递函数及其零极点：零点：[-3.601E02+j*2.532257488882203E3, -3.601E02+j*2.532257488882203E3]；极点：[-2.299999999999999E2+ j*4.516635916254485E02, -2.299999999999999E2-j*4.516635916254485E02,  -5.636E02]增益：k=21.955。试绘制其频率特性。T>