滤波器的研究与实现

1、1 引言随着信息技术的迅猛发展，数字信号处理已成为一个极其重要的学科和技术领域。在通信、语音、图像、自动控制和家用电器等众多领域得到了广泛的应用。数字滤波是数字信号处理的重要环节，它在数字信号处理中占有着重要的地位,它具有可靠性好、精度高、灵活性大、体积小、重量轻等优点。随着数字技术的发展，数字滤波器越来越受到人们的重视，广泛地应用于各个领域。数字滤波器的输入输出信号都是数字信号，它是通过一定的运算过程改变输入信号所含频率成分的相对比例或者滤除某些频率成分来实现滤波的，这种运算过程是由乘法器、加法器和单位延迟器组成的。数字滤波器是数字信号处理技术的重要内容，其对数字信号进行的最常见处理是保留数

2、字信号中的有用频率成分和去除信号中的无用频率成分。按照时间域的特性，数字滤波器可以分为无限冲激脉冲响应数字滤波器（IIR滤波器）和有限冲激脉冲响应数字滤波器（FIR滤波器）。从性能上来说，IIR滤波器传输函数的极点可位于单位圆内的任何地方，因此可用较低的阶数获得高的选择性，所用的存贮单元少，所以经济而效率高。但是这个高效率是以相位的非线性为代价的。选择性越好，则相位非线性越严重。相反，FIR滤波器却可以得到严格的线性相位，然而由于FIR滤波器传输函数的极点固定在原点，所以只能用较高的阶数达到高的选择性。对于同样的滤波器设计指标，FIR滤波器所要求的阶数可以比IIR滤波器高510倍。结果，成本较

3、高，信号延时也较大；如果按相同的选择性和相同的线性要求来说，则IIR滤波器就必须加全通网络进行相位较正，同样要增加滤波器的节数和复杂性。从结构上来说，IIR滤波器必须采用递归结构，极点位置必须在单位圆内，否则系统将不稳定。另外，在这种结构中，由于运算过程中对序列的舍入处理，这种有限字长效应有时会引入寄生振荡。相反，FIR滤波器主要采用非递归结构，不论在理论上还是在实际的有限精度运算中都不存在稳定性问题，运算误差也较小。此外，FIR滤波器可以采用快速傅里叶变换算法，在相同阶数的条件下，运算速度可以快得多。从设计工具上来说，IIR滤波器可以借助模拟滤波器的成果，因此一般都有有效的封闭形式的设计公式

4、可供参考，计算工作量比较小，而且对计算工具的要求不高；FIR滤波器一般没有封闭形式的设计工具。窗函数法设计FIR滤波器也仅给出了窗函数的计算公式，但是在计算通带阻带衰减时无明显表达式。一般FIR滤波器的设计只有计算程序可循，因此它对计算工具要求较高。另外，相比IIR，FIR系统只有零点，因此系统总是稳定的，所以能够获得严格的线性相位。IIR滤波器虽然设计简单，但主要用于设计具有片段常数特性的滤波器，如低通、高通、带通及带阻等，往往脱离不了模拟滤波器的格局。而FIR滤波器则要灵活得多，易于实现某些特殊功能的应用。综上所述，IIR和FIR滤波器各有所长，所以在实际应用中，选择滤波器型号时，应该从多

5、方面加以考虑。FIR数字滤波器的实现,大体可以分为软件和硬件实现两种。软件方法速度慢，难以对信号进行实时处理，虽然可以用快速傅立叶变换算法来加快计算速度 ,但很难达到实时处理的要求，因而多用于教学与科研。硬件实现主要采用以下两种方法:（1）采用DSP（Digital Signal Processing）处理器来实现。（2）采用固定功能的专用信号处理器。但这两种方法实现起来比较困难。FIR相对于IIR滤波器有许多独特的优越性，在保证满足滤波器幅频响应要求的同时还可获得严格的线性相位特性，从而保持稳定。对非线性相位FIR滤波器一般可以用 IIR滤波器来代替。由于在数据通信、语音信号处理、图像处理以

6、及自适应处理等领域往往要求信号在传输过程中不能有明显的相位失真 ,而IIR存在频率色散的问题，所以 FIR滤波器获得了更广泛的应用。随着通信与信息技术的发展，数字信号处理在该领域显得越来越重要。同时数字信号处理在语音、自动控制、航空航天和家用电器等领域也得到了广泛应用，它已成为当今一门极其重要的学科和技术。在数字信号处理中起重要作用并获得广泛应用的是数字滤波器，数字滤波器是数字信号处理的基础。Matlab Matrix laborator是美国MathWorks公司推出的具有强大数值分析、矩阵运算、图形绘制和数据处理等功能的软件，现已广泛应用到教学、科研、工程设计等领域。随着 Matlab软件

7、信号处理工具箱的推出，Matlab已成为信息处理，特别是数字信号处理DSP应用中分析和设计的主要工具。就Matlab信号处理中的滤波器设计而言，在很大程度上能快速有效地实现滤波器的分析、设计及仿真，大大节约了设计时间，相对传统设计而言，简化了滤波器设计难度。2 FIR数字滤波器的特性FIR数字滤波器的冲激响应的Z变换为： （2-1）其中滤波器最重要的两个特性为线性相位特性和幅度特性。2.1FIR滤波器的线性相位特性由式（2-1）可以看出，H（z）是的N-1次多项式，它在z平面内有N-1个零点，同时在原点有N-1个重极点。因为FIR数字滤波器的单位冲激响应是有限长的，所以它永远都是稳定的。在数字

8、信号处理的许多领域中，常常要求滤波器具有线性相位，达到这一要求，仅需要对FIR数字滤波器的冲激响应施加一定的约束。令，就可由得到FIR数字滤波器的频率响应：（2-2）式中是的幅频特性，是的相频特性（2-3）当要求滤波器具有严格的线性相位，或者说具有相位不失真时，应有：（2-4）即（2-5）也即（2-6）将式（2-6）展开可得：（2-7）式（2-7）的左边具有傅里叶级数的形式。由于FIR数字滤波器的冲激响应为实数，若呈偶对称，即（2-8a）其对称中心在（2-9）此时有（2-10）因此：（2-11）根据式（2-2），有(2-12a)（2-12b）若呈奇对称，（2-8b）其对称中心也在处。此时有：（

9、2-13a）（2-13b）即呈奇对称的滤波器相位与呈偶对称的相位产生了的相移。式（2-8）和式（2-9）是FIR数字滤波器具有线性相位的充要条件。即它要求FIR数字滤波器的单位冲激响应的序列必须满足式（2-8）所表示的特定的对称性，其相位延迟等于长度的一半，即个采样周期。具有线性相位的FIR数字滤波器既有恒定的延迟群，又有恒定的相延迟。2.2 FIR滤波器的幅频特性FIR滤波器的幅频特性可以分为以下4种情况：（1）为偶对称，且N为奇数(型滤波器）当为偶对称，且N为奇数时，根据式（2-12）滤波器的幅频函数可以表示为：（2-14）其中，。此时，对呈偶对称。型滤波器的幅频函数可以通过例程2-1的M

10、ATLAB语言实现，曲线图如图2.1所示。例程2-1 型滤波器幅频函数FunctionHr,w,a,L=hr_type1(h);%计算所设计的型滤波器的振幅响应%Hr=振幅响应%a=型滤波器的系数%L=Hr的阶次%h=型滤波器的单位冲激响应M=length(h);L=(M-1)/2;a=h(L+1) 2*h(L:-1:1);n=0:1:L;w=0:1:500'*2*pi/500;Hr=cos(w*n)*a't=-3.14:0.314:3.14;h=cos(t);Hr,w,a,L=hr_type1(h);plot(Hr)图2.1型滤波器幅频函数（2）为偶对称，且N为偶数（型滤波器

11、）当为偶对称，且N为偶数时，根据式（2-12）滤波器的幅频函数可以表示为：（2-15）其中：。此时对呈奇对称。但是，故高通滤波器不能用这种方法实现。型滤波器的幅频函数可以通过例程2-2的MATLAB语言实现，曲线图如图2.2所示。 例程2-2 型滤波器幅频函数functionHr,w,b,L=hr_type2(h);%计算所设计的型滤波器的振幅响应%Hr=振幅响应%a=型滤波器的系数%L=Hr的阶次%h=型滤波器的单位冲激响应M=length(h);L=M/2;b=2*h(L:-1:1);n=1:1:L;n=n-0.5;w=0:1:500*2*pi/500;Hr=cos(w*n)*b;t=-3

12、.14:0.33:3.14;h=cos(t);Hr,w,b,L=hr_type2(h);plot(Hr)图2.2 型滤波器幅频函数（3）为奇对称，且N为奇数（型滤波器）当为奇对称，且为奇数时，根据式（2-13）滤波器的幅频函数可以表示为：（2-16）其中。此时对呈奇对称。但是当时，所以低通、高通滤波器不能采用这种形式。型滤波器的幅频函数可以通过例程2-3的MATLAB语言实现，曲线图如图2.3所示。 例程2-3 型滤波器幅频函数functionHr,w,c,L=hr_type3(h);%计算所设计的型滤波器的振幅响应%Hr=振幅响应%b=型滤波器的系数%L=Hr的阶次%h=型滤波器的单位冲激响

13、应M=length(h);L= (M-1)/2;c=2*h(L+1:-1:1);n=0:1:L;w=0:1:500*2*pi/500;Hr=sin(w*n)*c;t=-3.14:0.314:3.14;h=sin(t);Hr,w,c,L=hr_type3(h);plot(Hr)图2.3 型滤波器幅频函数4、为奇对称，且N为偶数（型滤波器）当为奇对称，且N为偶数时，根据式(2-13)滤波器的幅频函数可表示为：（2-17）其中：。此时对呈奇对称，对呈偶对称。但是当时，所以低通滤波器不能采用这种形式。型滤波器的幅频函数可以通过例程2-4的MATLAB语言实现，曲线图如图2.4所示。例程2-4 型滤波器

14、的幅频函数functionHr,w,d,L=hr_type4(h);%计算所设计的型滤波器的振幅响应%Hr=振幅响应%d=型滤波器的系数%L=Hr的阶次%h=型滤波器的单位冲激响应M= length (h);L= M/2;d=2*h(L:-1:1);n=1:1:L;n=n-0.5;w=0:1:500*2*pi/500;Hr=sin(w*n)*d;t=-3.14:0.33:3.14;h=sin(t);Hr,w,d,L=hr_type4(h);plot(Hr)图2.4 型滤波器的幅频函数3基于窗函数的FIR数字滤波器设计3.1常用窗函数及MATLAB实现矩形窗（Rectangular window

15、）矩形窗函数的时域形式可以表示为: （3-1-1）它的频域特性为： （3-1-2）矩形窗函数的时域幅度与频域幅度特性曲线的MATLAB实现如例程3-1所示，曲线图如图3.1所示。 例程3-1 矩形窗函数的时域幅度与频域幅度特性曲线%draw trectangle window time and frequency amplitudeclear all;N=64;w=rectwin(N);wvtool(w);图3.1 矩形窗函数的时域幅度与频域幅度特性曲线三角窗（Triangular window）三角窗函数的时域形式可以表示为：当窗长N为奇数时：（3-1-3）当窗长N为偶数时：（3-1-4）它

16、的频域特性为：（3-1-5）三角窗函数的主瓣宽度为，比矩形窗函数的主瓣宽度增加一倍，但是它的旁瓣却小得多。三角窗函数的时域幅度与频域幅度特性曲线的MATLAB实现如例程3-2所示，曲线图如图3.2所示。 例程3-2 三角窗函数的时域幅度与频域幅度特性曲线%draw triangle window time and frequency amplitudeclear all;N=64;w=triang(N);wvtool(w);图3.2三角窗函数的时域幅度与频域幅度特性曲线汉宁窗（Hanning window）汉宁窗函数是余弦平方函数，又称之为升余弦函数，它的时域形式可以表为：（3-1-6）其中。

17、它的频域幅度特性函数为：（3-1-7）其中为矩形窗函数的幅度频率特性函数。汉宁窗函数的最大旁瓣值比主瓣值低31dB，但是主瓣宽度比矩形窗函数的主瓣宽度增加了一倍，为。汉宁窗函数的时域幅度与频域幅度特性曲线的MATLAB实现如例程3-3所示，曲线图如图3.3所示。 例程3-3 汉宁窗函数的时域幅度与频域幅度特性曲线%draw hanning window time and frequency amplitudeclear all;N=64;w=hann(N);wvtool(w); 图3.3 汉宁窗函数的时域幅度与频域幅度特性曲线海明窗（Hamming window）海明窗函数是一种改进的升余弦函

18、数，定义为：（3-1-8）其中。它的幅度频率特性为：（3-1-9）和汉宁窗（Hanning window）函数相比，海明窗函数的主瓣宽度和汉宁窗相同，但是它的旁瓣又被进一步压低，其最大旁瓣值比主瓣值低41dB。海明窗函数的时域幅度与频域幅度特性曲线的MATLAB实现如例程3-4所示，曲线图如图3.4所示。 例程3-4 海明窗函数的时域幅度与频域幅度特性曲线%draw hanning window time and frequency amplitudeclear all;N=64;w=hanning(N);wvtool(w);图3.4 海明窗函数的时域幅度与频域幅度特性曲线布拉克曼窗（Blac

19、kman window）为了进一步抑制旁瓣，对升余弦函数再加上一个二次谐波的余弦分量，便得到了布拉克曼窗函数，也称之为二阶升余弦函数。它的时域形式可以表示为：（3-1-10）其中。其幅度频率特性为：（3-1-11）其中为矩形窗函数的幅度频率特性函数。布拉克曼窗函数主瓣的宽度是矩形窗函数主瓣宽度的三倍，为。它的最大旁瓣值比主瓣值低57dB。布拉克曼窗函数的时域幅度与频域幅度特性曲线的MATLAB实现如例程3-5所示，曲线图如图3.5所示。 例程3-5 布拉克曼窗函数的时域幅度与频域幅度特性曲线%draw blackman window time and frequency amplitudecl

20、ear all;N=64;w=blackman(N);wvtool(w);图3.5 布拉克曼窗函数的时域幅度与频域幅度特性曲线切比雪夫窗（Chebyshev window）切比雪夫窗又称Dolph-Chebyshev窗，它是由一个切比雪夫多项式在单位圆上作N点等间隔抽样，然后再作DFT反变换得到的。当给定窗的长度N，最大旁瓣衰减r（切比雪夫窗函数傅里叶变换旁瓣幅度比主瓣低rdB），则其3dB带宽B的表达式为：（3-1-12）切比雪夫窗函数的时域幅度与频域幅度特性曲线的MATLAB实现如例程3-6所示，曲线图如图3.6所示。 例程3-6 切比雪夫窗函数的时域幅度与频域幅度特性曲线%draw ch

21、ebyshev window time and frequency amplitudeclear all;N=64;w=chebwin(N);wvtool(w);图3.6 切比雪夫窗函数的时域幅度与频域幅度特性曲线3.2数字低通滤波器的窗函数设计假设理想低通数字滤波器的频率响应如图3.7所示，其幅频特性为。那么该滤波器的：（3-2-1）式中表示截止频率（rad），表示采样延迟。0低通1图3.7 理想低通数字滤波器幅频特性由此可得理想数字低通滤波器的单位冲激响应为：（3-2-2）根据式（3-2-2），可知为无限非因果序列，关于对称。为了从得到一个因果线性相位的FIR滤波器，必须利用有限长度N的窗

22、函数对进行截取，截取后的冲激响应函数可表示为：（3-2-3）此时为关于偶对称的有限因果序列。当N为奇数时，所设计的FIR数字高通滤波器为型滤波器；当N为偶数时，为型滤波器。而表示的滤波器的频率特性为：（3-2-4）是否能够很好地逼近取决于窗函数的频率特性。若将理想滤波器的频率响应写成：（3-2-5）其中幅度频率特性：（3-2-6）此时，滤波器的频率特性可表述为：（3-2-7）由此可以得到所设计的滤波器的幅度频率特性为：（3-2-8）由式（3-2-8）可见，对实际FIR滤波器有影响的只是窗函数的幅度频率特性。实际中的FIR滤波器的幅度频率特性，是理想低通滤波器的幅度频率特性和窗函数的幅度频率特性

23、的复卷积。复卷积给带来过冲和波动，所以加窗函数后，对滤波器的理想特性的影响有以下几点：（1）在截止频率的间断点变成了连续的曲线，使得出现了一个过渡带，它的宽度等于窗函数的主瓣宽度。由此可知，如果窗函数的主瓣越宽，过渡带就越宽。（2）由于窗函数旁瓣的影响，使得滤波器的幅度频率特性出现了波动，波动的幅度取决于旁瓣的相对幅度。旁瓣范围的面积越大，通带波动和阻带的波动就越大，也就是说阻带的衰减越小。而波动的多少，取决于旁瓣的多少。（3）增加窗函数的长度，只能减少窗函数的幅度频率特性的主瓣宽度，而不能减少主板和旁瓣的相对值，该值取决于窗函数的形状，即增加截取函数的长度N只能相应的减少过渡带，而不能改变滤

24、波器的波动程度。为了满足工程上的需要，可以通过改变窗函数的形状来改善滤波器的幅度频率特性，而窗函数的选择原则是：（1） 具有较低的旁瓣幅度，尤其是第一旁瓣的幅度；（2） 旁瓣的幅度下降的速率要快，以利用增加阻带的衰减；（3） 主板的宽度要窄，这样可以得到比较窄的过渡带。通常上述的几点难以同时满足。当选用主板宽度较窄时，虽然能够得到比较陡峭，但是通带和阻带的波动明显增加。当选用比较小的旁瓣幅度时，虽然能够得到比较平坦和匀滑的幅度频率响应，但是过渡带将加宽，因此实际中选用的窗函数往往是它们的折中。在保证主瓣的宽度达到一定要求的条件下，适当牺牲主瓣的宽度来换取旁瓣的波动减少。以上是从幅度频率特性设计

25、方面对窗函数提出的要求，实际中设计FIR数字滤波器往往要求是线性相位的，因此要求满足线性相位的条件，即要求满足：（3-2-9）所以窗函数不仅有截短的作用，而且能够起到平滑的作用，在很多领域得到了应用。型理想低通滤波器的单位冲激响应计算的MATLAB实现如例程3-9，曲线图如图3.9所示。例程3-9 型理想低通滤波器的单位冲激响应计算function hd=ideal_lp(Wc,N)%compute the ideal lowpass fiter unit pulse respondence hd(n)%wc:cutoff frequency%N:window length%hd:unit p

26、ulse respondencealpha=(N-1)/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(Wc*m)./(pi*m);Wc=100;N=128; hd= ideal_lp(Wc,N);plot(hd);图3.9 型理想低通滤波器的单位冲激响应4、基于频率抽样法的FIR滤波器设计窗函数设计FIR数字滤波器是从时域出发的，把理想滤波器的单位取样响应用合适的窗函数截短成为有限长度的，使得逼近理想的。以实现所设计的滤波器的频率响应逼近于理想滤波器的频率响应。一个有限长的序列，如果满足频率采样定理的条件，可以通过频谱的有限个采样点的值准确地恢复。频率抽样法正是采用这种思想

27、来设计FIR数字滤波器的。设理想滤波器的频率响应是，它在连续频率的周期函数。现对其抽样，使每一个周期有N个抽样值，即：（4-1）对作IDFT,可得到N点的单位抽样序列，即：（4-2）将作为所设计滤波器的单位冲激响应，这样就可以求出该滤波器的转移函数，即：（4-3）当然也可以用来表示，将（4-2）代入式（4-3），得（4-4）经推导得：（4-5）该系统的频率响应为：（4-6）经推导有：（4-7）这样，由连续的抽样得到，由的反变换得到，由做DTFT又得到连续谱。若再对进行抽样，令抽样点数，为大于1的整数，得抽样结果，那么：（4-8）如果，那么：（4-9）这说明由式（4-5）求出的滤波器，其频率响应

28、在的抽样点上严格地等于所希望的值，而在的点上，则由内插函数的插值决定。所以这种滤波器的设计方法被称之为频率抽样法。该内插函数是：（4-10）由此可得：（4-11）由此可以看出，连续函数是由N个离散值作为权重和插值函数线性组合的结果。显然对抽样点N取的越大，对的近似程度越好。N的选取要视对在通带和阻带的技术要求而定。是移位后的sinc函数，它正是一个离散矩形窗函数的频谱。在频率抽样法中指定要比在窗函数法中指定复杂。指定的原则是：（1）在通带内可以令，阻带内，且在通带内赋给一相位函数；（2）指定的应保证由式求出的是实的；（3）由求出的应具有线性相位。 为保证所涉及的滤波器具有线性相位，必须对频率采

29、样值进行约束。假设可以表示为：（4-12）那么根据式（4-2）有以下4种情况：（1）偶对称，且N为奇数当偶对称，且N为奇数时，式（4-12）的约束条件为：（4-13）（4-14） （2）偶对称，且N为偶数当偶对称，且N为偶数时，式（4-12）的约束条件为：（4-15）（4-16） （3）奇对称，且N为奇数当奇对称，且N为奇数时，式（4-12）的约束条件为：（4-17）（4-18）（4）奇对称，且N为偶数当奇对称，且N为偶数时，式（4-12）的约束条件为：（4-19）（4-20）对于理想滤波器的频率响应，根据第一个采样点的不同，可分为两种采样方法：方法1：第一个采样点在处，根据式（4-1）、式（

30、4-2）、式（4-5）、式（4-6）、式（4-7）可以推导出：（4-21）（4-22）（4-23）（4-24）其中根据式（4-3）（4-22），可以表示为：（4-25）所以当为实数时，是周期为N的函数，且满足：（4-26）（4-27）1、N为奇数（1）若理想滤波器的频率响应不存在相移时，即设计型FIR滤波器时，则线性相位约束条件为：（4-28）（2）若理想滤波器的频率响应存在的相移时，即设计型滤波器时，则线性相位约束条件为：（4-29）2、N为偶数（1）若理想滤波器的频率响应不存在相移时，即设计型FIR滤波器时，则线性相位约束条件为：（4-30）（2）若理想滤波器的频率响应存在的相移时，即设计

31、型滤波器时，则线性相位约束条件为：（4-31）方法2：第一个采样点在处，根据式（4-1）、式（4-2）、式（4-5）、式（4-6）、式（4-7）可以推导出：（4-32）（4-33）（4-34）（4-35）其中。根据式（4-3）和式（4-33），可以表示：（4-36）所以当为实数时，是周期为N的函数，且满足：（4-37）（4-38）1、 N为奇数（1）若理想滤波器的频率响应不存在相移时，即设计型FIR滤波器时，则线性相位约束条件为：（4-39）（2）若理想滤波器的频率响应存在的相移时，即设计型FIR滤波器时，则线性相位约束条件为：（4-40）2、N为偶数（1）若理想滤波器的频率响应不存在相移时，

32、即设计型FIR滤波器时，则线性相位约束条件为：（4-41）（2）若理想滤波器的频率响应存在的相移时，即设计型FIR滤波器时，则线性相位约束条件为：（4-42）例程4-1：根据下列指标，利用频率抽样法，设计一个型FIR高通滤波器：wp=0.8pi,ws=0.7pi,Ap=1db,As=40db程序：%use frequency sampling method 1 to design highpass digital filterclear all;N=61;T1=0.1095;T2=0.598;alpha=(N-1)/2;l=0:N-1;wl=(2*pi/N)*l;Hrs=zeros(1,22)

33、,T1,T2,ones(1,14),T2,T1,zeros(1,21); %理想振幅响应采样Hdr=0,0,1,1;wdl=0,0.75,0.75,1;k1=0:floor(N-1)/2);k2=floor(N-1)/2)+1:N-1;angH=-alpha*(2*pi)/N*k1,alpha*(2*pi)/N*(N-k2); %相位约束条件Hdk=Hrs.*exp(j*angH); %构成Hd(k)h=real(ifft(Hdk,N); %实际单位冲激响应db,mag,pha,w=freqz_m2(h,l);Hr,ww,a,L=hr_type1(h); %实际振幅响应subplot(221)

34、plot(wl/pi,Hrs,'.',wdl,Hdr)title('频率样本Hd(k):N=61')axis(0 1 -0.1 1.2)subplot(222)stem(l,h)title('实际单位脉冲响应h(n)')subplot(223)plot(ww/pi,Hr,wl/pi,Hrs,'.')title('实际振幅响应H(w)')axis(0 1 -0.1 1.2)subplot(224)plot(w/pi,db)title('幅度响应(db)')axis(0 1 -80 10)图4.1利用频

35、率抽样法，设计一个型FIR高通滤波器：5、基于切比雪夫逼近法直接设计FIR数字滤波器基于切比雪夫逼近法直接设计FIR数字滤波器的MATLAB函数主要有firls、remez和remezord三种函数。函数firls是函数fir1和fir2的扩展，函数firls与remez和remezord在使理想滤波器的频率响应与实际滤波器的频率响应之间的整体误差最小时，均采用切比雪夫原理，只是误差方程求解的算法不一样，firls函数采用了最小二乘算法，remez和remezord函数采用了Parks-McClellan算法。Parks-McClellan算法利用Remez交换算法和切比雪夫逼近理论来设计滤波

36、器，使实际滤波器的幅频响应最优地拟合理想滤波器的幅频响应，使它们之间最大误差最小化，故函数remez设计的滤波器是最优的，有时称之为最优滤波器。又因为其在频域内滤波器呈现等波纹特点，因此也称为等纹滤波器。Parks-McClellan滤波器设计方法是FIR滤波设计中最流行的和应用最广的一种方法。remezord函数主要利用Parks-McClellan算法来优化估计FIR滤波器的阶次N。在MATLAB中，firls函数的调用格式为：remez函数的调用格式为：remezord函数的调用格式为：例程5-1 基于firls函数的FIR带通滤波器设计%to test firls and to des

37、ign bandpass FIR filter;clear all;F=0 0.3 0.4 0.6 0.7 0.9;A=0 1 0 0 0.5 0.5;b=firls(24,F,A,'hilbert');H,f=freqz(b,1,512,2);plot(f,abs(H),'k.')hold onfor i=1:2:6,plot(F(i) F(i+1),A(i) A(i+1),'k-')hold onendgrid onhold offlegend('实际滤波器','理想滤波器')图5.1 基于firls函数设计的

38、FIR带通滤波器与理想滤波器的幅频特性曲线例程5-2 基于remez函数的FIR带通滤波器设计%to test remez and to design bandpass FIR filter;clear all;f=0 0.3 0.4 0.6 0.7 1;a=0 0 1 1 0 0;b=remez(17,f,a);h,w=freqz(b,1,512);plot(f,a,'k-',w/pi,abs(h),'k.')legend('理想滤波器','实际滤波器')图5.2基于remez函数设计的FIR带通滤波器与理想滤波器的幅频特性曲线例程5-3 基于remezord函数的FIR低通滤波器设计%to test remezord and to design lowpass FIR filter;clear all;rp=3; %Passband ripplers=40; %