有限脉冲响应数字滤波器的设计

有限脉冲响应数字滤波器的设计第一页，共四十五页，编辑于2023年，星期三设FIR滤波器单位脉冲响应h(n)长度为N，其系统函数H(z)为：

H(z)是的N-1次多项式，它在z平面上有N-1个零点，原点z=0是N-1阶重极点。因此H(z)永远稳定。

稳定和线性相位是FIR滤波器突出的优点。第二页，共四十五页，编辑于2023年，星期三6.1线性相位FIR数字滤波器的条件和特点

本节主要介绍FIR滤波器具有线性相位的条件及幅度特性。

1.线性相位条件对于长度为N的h(n)，传输函数定义为(6.1.1)

(6.1.2)第三页，共四十五页，编辑于2023年，星期三式中，Hg(ω)称为幅度特性，θ(ω)称为相位特性。注意，这里Hg(ω)不同于|H(ejω)|,Hg(ω)为ω的实函数，可能取负值，而|H(ejω)|总是正值。H(ejω)线性相位是指θ(ω)是ω的线性函数，即

θ(ω)=-τω,τ为常数(6.1.3)

如果θ(ω)满足下式：

θ(ω)=θ0-τω,θ0是起始相位(6.1.4)

严格地说，此时θ(ω)不具有线性相位，但以上两种情况都满足群时延是一个常数，即第四页，共四十五页，编辑于2023年，星期三一般称满足(6.1.3)式是第一类线性相位；满足(6.1.4)式为第二类线性相位。下面推导与证明满足第一类线性相位的条件是：h(n)是实序列且对(N-1)/2偶对称，即

h(n)=h(N-n-1)(6.1.5)

满足第二类线性相位的条件是：h(n)是实序列且对(N-1)/2奇对称，即

h(n)=-h(N-n-1)(6.1.6)第五页，共四十五页，编辑于2023年，星期三(1)第一类线性相位条件证明：将(6.1.5)式代入上式得令m=N-n-1,则有(6.1.7)第六页，共四十五页，编辑于2023年，星期三按照上式可以将H(z)表示为

将z=ejω代入上式，得到：

按照(6.1.2)式，幅度函数Hg(ω)和相位函数分别为(6.1.8)

(6.1.9)第七页，共四十五页，编辑于2023年，星期三(2)第二类线性相位条件证明：(6.1.10)

令m=N-n-1,则有

同样可以表示为第八页，共四十五页，编辑于2023年，星期三因此，幅度函数和相位函数分别为(6.1.11)(6.1.12)该类线性相位FIR滤波器适合于在微分器及90度移相器中应用，故又称这类滤波器为正交变换网络。第九页，共四十五页，编辑于2023年，星期三6.2利用窗函数法设计FIR滤波器

设希望设计的滤波器传输函数为Hd(ejω),hd(n)是与其对应的单位脉冲响应，因此第十页，共四十五页，编辑于2023年，星期三如果已知则可求出hd(n)，经过z变换可得到滤波器的系统函数。但是一般逐段恒定，在边界频率处有不连续点，故致使hd(n)是无限时宽的，且是非因果序列。例如：理想的低通滤波器的传输函数(6.2.1)

相应的单位取样响应hd(n)为第十一页，共四十五页，编辑于2023年，星期三(6.2.2)

为了构造一个长度为N的线性相位滤波器，只有将hd(n)截取一段，并保证截取的一段对(N-1)/2对称。设截取的一段用h(n)表示，即

h(n)=hd(n)RN(n)(6.2.3)第十二页，共四十五页，编辑于2023年，星期三图6.2.1理想低通的单位脉冲响应及矩形窗第十三页，共四十五页，编辑于2023年，星期三我们知道Hd(ejω)是一个以2π为周期的函数，可以展开为傅氏级数，即对(6.2.3)式进行傅里叶变换，根据复卷积定理，得到：(6.2.4)式中，Hd(ejω)和RN(ejω)分别是hd(n)和RN(n)的傅里叶变换，即(6.2.5)

第十四页，共四十五页，编辑于2023年，星期三RN(ω)称为矩形窗的幅度函数；将Hd(ejω)写成下式：按照(6.2.1)式，理想低通滤波器的幅度特性Hd(ω)为将Hd(ejω)和RN(ejω)代入(6.2.4)式，得到：第十五页，共四十五页，编辑于2023年，星期三将H(ejω)写成下式：(6.2.6)该式表明了滤波器的幅度特性等于理想低通滤波器的幅度特性与矩形窗幅度特性的卷积。第十六页，共四十五页，编辑于2023年，星期三

图6.2.2矩形窗对理想低通幅度特性的影响第十七页，共四十五页，编辑于2023年，星期三通过以上分析可知，对hd(n)加矩形窗处理后，H(ω)和原理想低通Hd(ω)差别有以下两点：

(1)在理想特性不连续点ω=ωc附近形成过渡带。过渡带的宽度，近似等于RN(ω)主瓣宽度，即4π/N。

(2)通带内增加了波动，最大的峰值在ωc-2π/N处。阻带内产生了余振，最大的负峰在ωc+2π/N处。通带与阻带中的波动与窗函数的幅度谱有关系。RN(ω)波动越快(N加大时)，通带、阻带波动越快，RN(ω)旁瓣的大小直接影响H(ω)波动的大小。

以上两点就是对hd(n)用矩形窗截断后，在频域的反映，称为吉布斯效应。

第十八页，共四十五页，编辑于2023年，星期三通过实验分析：调整窗口的长度N可以有效的控制过渡带的宽度。减少带内波动以及加大阻带的衰减只能从窗函数的形状上找解决办法。思路：能否找到窗函数的形状，使其谱函数的主瓣包含更多的能量，相应旁瓣幅度就减小了，旁瓣的减小可使通带、阻带波动减小，从而加大阻带衰减。第十九页，共四十五页，编辑于2023年，星期三下面介绍几种常用的窗函数。设

h(n)=hd(n)w(n)

式中w(n)表示窗函数。

1.矩形窗(RectangleWindow)

wR(n)=RN(n)

前面已分析过，按照(6.2.5)式，其频率响应为第二十页，共四十五页，编辑于2023年，星期三2.三角形窗(BartlettWindow)(6.2.8)

其频率响应为

(6.2.9)

第二十一页，共四十五页，编辑于2023年，星期三3.汉宁(Hanning)窗——升余弦窗当N1时，N-1≈N第二十二页，共四十五页，编辑于2023年，星期三图6.2.3汉宁窗的幅度特性第二十三页，共四十五页，编辑于2023年，星期三4.哈明(Hamming)窗——改进的升余弦窗(6.2.11)其频域函数WHm(ejω)为其幅度函数WHm(ω)为当N>>1时，可近似表示为第二十四页，共四十五页，编辑于2023年，星期三5.布莱克曼(Blackman)窗(6.2.13)

其频域函数为其幅度函数为(6.2.14)

第二十五页，共四十五页，编辑于2023年，星期三6.凯塞－贝塞尔窗(Kaiser-BaselWindow)

式中

I0(x)是零阶第一类修正贝塞尔函数，可用下面级数计算：第二十六页，共四十五页，编辑于2023年，星期三图6.2.4常用的窗函数第二十七页，共四十五页，编辑于2023年，星期三

图6.2.5常用窗函数的幅度特性(a)矩形窗；(b)巴特利特窗(三角形窗)；(c)汉宁窗；(d)哈明窗；(e)布莱克曼窗第二十八页，共四十五页，编辑于2023年，星期三

图6.2.6理想低通加窗后的幅度特性(N=51,ωc=0.5π)(a)矩形窗；(b)巴特利特窗(三角形窗)；(c)汉宁窗；

(d)哈明窗；(e)布莱克曼窗第二十九页，共四十五页，编辑于2023年，星期三表6.2.2六种窗函数的基本参数第三十页，共四十五页，编辑于2023年，星期三MATLAB提供了窗函数的子程序：1)w=boxcar(M)在数组w中产生M点的矩形窗函数；

2)w=triang(M)在数组w中产生M点的三角形窗函数；

3)w=hanning(M)在数组w中产生M点的hanning窗函数；

4)w=hamming(M)在数组w中产生M点的hamming窗函数；

5)w=blackman(M)在数组w中产生M点的blackman窗函数；

6)w=kaiser(M,beta)在数组w中产生beta值M点的矩形窗函数；第三十一页，共四十五页，编辑于2023年，星期三下面介绍用窗函数设计FIR滤波器的步骤。

(1)根据技术要求确定待求滤波器的单位取样响应hd(n)。如果给出待求滤波器的频响为Hd(ejω)，那么单位取样响应用下式求出：如果Hd(ejω)较复杂，或者不能用封闭的公式表示时，则不能用上式求出hd(n)。可以对Hd(ejω)从0至采样M点，并用代替(6.2.17)中的，(6.2.17)近似写为：(6.2.17)

(6.2.18)

第三十二页，共四十五页，编辑于2023年，星期三根据频率采样定理，hM(n)与hd(n)应满足如下关系：因此，如果M取得比较大，可以保证在窗口内hM(n)有效的逼近hd(n)。实际计算时，可以用Hd(ejω)的M点采样值，进行M点IDFT得到。

如果给出通带阻带衰减和边界频率的要求，可选用理想滤波器作为逼近函数，从而用理想滤波器的特性作傅立叶逆变换，求出hd(n)。第三十三页，共四十五页，编辑于2023年，星期三

(2)根据对过渡带及阻带衰减的要求，选择窗函数的形式，并估计窗口长度N。设待求滤波器的过渡带用Δω表示，它近似等于窗函数主瓣宽度。因过渡带Δω近似与窗口长度N成反比，，A取决于窗口形式。

(3)计算滤波器的单位取样响应h(n)

h(n)=hd(n)w(n)第三十四页，共四十五页，编辑于2023年，星期三其中是已选好的窗函数。如果要求线性相位，则要求hd(n)和均对(N-1)/2对称。

(4)验算技术指标是否满足要求。设计出的滤波器频率响应用下式计算：计算时可用FFT算法。如果H(ejω)不满足要求，根据具体情况重复(2),(3),(4)步，直到满足要求。第三十五页，共四十五页，编辑于2023年，星期三例6.2.1设计线性相位FIR低通滤波器，给定采样频率，通带截止频率，阻带截止频率，阻带衰减不少于-50dB，幅度特性如图所示。解：1求各对应的数字频率。第三十六页，共四十五页，编辑于2023年，星期三2理想的低通线性相位滤波器求故取

满足线性相位。第三十七页，共四十五页，编辑于2023年，星期三3.由确定窗函数形状，确定N。，查表可选用哈明窗。故4.由确定h(n).窗函数过渡带宽阻带最小衰减/dB矩形窗4π/N－21三角形窗8π/N－25汉宁窗8π/N－44哈明窗8π/N－53布来克曼窗12π/N－74第三十八页，共四十五页，编辑于2023年，星期三

故

5.检验各项指标是否满足要求，不满足则需改变N或窗形状。第三十九页，共四十五页，编辑于2023年，星期三例设计数字FIR低通滤波器，技术指标如下：

解：使用hamming和blackman窗函数都能提供大于50dB阻带衰减。这里选用hamming窗函数，它给出比较小的过渡带，因此有较低的阶。wp=0.2*pi;ws=0.3*pi;tr_width=ws-wpM=ceil(8*pi/tr_width)n=[0:1:M-1];wc=(ws+wp)/2hd=ideal_lp(wc,M);w_ham=(hamming(M))';第四十页，共四十五页，编辑于2023年，星期三

h=hd.*w_ham;[db,mag,pha,grd,w]=freqz(h,[1]);%plotssubplot(1,1,1)subplot(2,2,1);stem(n,hd);title('IdealImpulseResponse')axis([0M-1-0.10.3]);xlabel('n');ylabel('hd(n)')subplot(2,2,2);stem(n,w_ham);title('HammingWindow')axis([0M-101.1]);xlabel('n');ylabel('w(n)')subplot(2,2,3);stem(n,h);title('ActualImpulseResponse')axis([0M-1-0.10.3]);xlabel('n');ylabel('h(n)')subplot(2,2,4);plot(w/pi,db);title('MagnitudeResponseindB');gridaxis([01-10010]);xlabel('frequencyinpiunits');ylabel('Decibels')set(gca,'XTickMode','manual','XTick',[0,0.2,0.3,1])set(gca,'YTickMode','manual','YTick',[-50,0])set(gca