数字FIR滤波器

1、-FIR滤波器设计主要内容n线性移不变系统（线性移不变系统（LSI）；）；n数字滤波器的基本概念；数字滤波器的基本概念；nFIR数字滤波器；数字滤波器；nFIR数字滤波器的设计。数字滤波器的设计。 数字滤波器实际上就是一个离散时间系统数字滤波器实际上就是一个离散时间系统，设计，设计并实现一个数字滤波器也就是实现对于输入的离散时并实现一个数字滤波器也就是实现对于输入的离散时间信号的处理，滤波器的输出满足特定的目的。间信号的处理，滤波器的输出满足特定的目的。1. 线性移不变系统（线性移不变系统（LSI） 相对于连续时间信号，系统的分析对象是连续时间LTI系统，即线性时不变系统。在数字信号处理中，讨

2、论的是离散时间信号和系统，对于系统分析的对象是离散时间LSI系统，即线性移不变系统。LTI系统和LSI系统是等价的。 LSI系统在时域同样也是由其冲激响应决定的 y nh nx nnLSI系统的输入，冲激响应和输出之间的关系 讨论离散时间系统，手段主要有两个：n序列的卷积时域方法；nZ变换分析法频域方法。 二种方法统一于时频对应。系统函数与传输函数 n一般借助于差分方程来描述一个离散时间系统：n对该方程的时域解导致了系统分析的时域方法。相应的对差分方程进行双边Z变换（或单边Z变换，不考虑初始条件），其系统函数为： 00NMkkkka y nkb x nkNkkkMkkkzazbzH00)( 进

3、一步 单位脉冲响应为系统函数的反Z变换： 这就是分析LSI系统的频域方法。 ( ) nnH zh n z1 ( )h nZH z注注：（1）由差分方程求）由差分方程求 ，已假定系统为，已假定系统为LSI； （2）因果性和稳定性由）因果性和稳定性由 的收敛域决定的收敛域决定 （即，差（即，差 分方程不能唯一确定分方程不能唯一确定 ） ( )H z( )H z h n离散时间系统的因果性和稳定性 n对于数字滤波器（等效为LSI系统）而言，必须考虑其因果性和稳定性。LSI系统是因果系统意味着系统 时刻的输出 ，只取决于 的输入 ，LSI系统是因果系统的充要条件为：n频域： 的收敛域一定包含点。 0n

4、n 0y n0nn )(nx 0h n 0n( )H z 稳定系统意味着系统对每个有界输入，产生一个有界输出的系统。LSI系统是稳定系统的充要条件为： 系统稳定要求 的收敛域包含单位圆 | |nh n ( )H z( ) jj nz enH zh n e j nnnh n eh n 对于线性时不变、因果、稳定系统对于线性时不变、因果、稳定系统，收敛域为，收敛域为 |, 01rzr 离散时间系统的频率响应 n一个离散时间系统（滤波器）最为重要的一个考察指标就是其频率响应。滤波器的频率响应包括幅频响应和相频响应两部分，幅频响应反映信号通过系统后各频率成分衰减的情况，而相频响应反映的是信号中各频率成

5、分通过系统后在时间上发生的位移情况。相应的信号失真分别有幅度失真和相位失真。一个理想的离散时间系统，除了具有所希望的幅频特性外，还希望具有线性相位的特性，这对于实现无失真的传输具有非常重要的意义。纯粹的线性相位往往难以实现，对于窄带信号，往往采用群迟延来衡量信号的相位失真。对于滤波器的设计，往往从幅频特性和相频特性两个方面综合考虑。 2. 数字滤波器的基本概念 n数字滤波器可以分为两种：有限长冲激响应FIR滤波器和无限长冲激响应IIR滤波器。FIR系统与IIR系统是离散时间LTI和LSI系统中两类很重要的系统，它们的特性、结构以及设计方法都存在很大的差异。 滤波器有模拟滤波器和数字滤波器之分，

6、经典的滤波器在功能上分为低通、高通、带通和带阻滤波器。其中低通滤波器是最基本的滤波器，其他的滤波器都可以由低通滤波器通过适当的转换得出。我们知道一个滤波器在频域有幅频响应和相频响应之分，对于幅频响应关注其通带、阻带以及过渡带特性，对应着滤波器主要的技术指标。滤波器的设计就是依据这些指标，设计相应的系统函数，然滤波器的设计就是依据这些指标，设计相应的系统函数，然后实现对应的系统后实现对应的系统。 数字滤波器的特点 n实现的方式灵活：硬件实现（专用的数字硬件）、软件实现（编程）及软件、硬件结合的实现方式；n精度高、稳定性好（受环境条件，如温度、干扰等的影响小）；n灵活性好（改变系统参数非常容易）易

7、于实现各种功能的综合。数字滤波器的实现形式 n数字滤波器实现的具体形式有：直接型，并联型，级联型和格型（Lattice）；包含着三种基本运算：乘系数、相加、移位（延迟）。n直接型：直接型：直接直接I型和直接型和直接II型型 0101 MNkkkky nb x nka y nka0 Mkkw nb x nkn直接直接I型型 n直接直接II型型 n从系统函数的观点，将系统函数从系统函数的观点，将系统函数 分解为分解为两个部分（两个部分（ ）：）： 01a ( )H z101111( )( )1( )MMNMbb zb zH zN za za zD zn级联型和并联型的基本思想是对系统函数进行因式分

8、解（有理分解），分解为若干个子系统的并联和级联。 dx nn级联型：级联型：12( )( )( )( )KH zH z HzHzn并联型：并联型：12( )( )( )( )KH zHzHzHz值得注意的问题：值得注意的问题：n有限字长效应的存在；n级联结构具有累积误差大的缺点，但是便于准确实现滤波器的极零点，便于调整滤波器的频率响应的性能，所需要的存储器少；n并联结构的每一个子系统都是独立的，不受其他子系统的量化误差和舍入误差的影响，因此对于误差最不敏感，但是对于极零点的调整不易。滤波器的Z域极零点分析 n将系统函数n因式分解NrrMrrzdzcAzH1111)1 ()1 ()(00/Aba

9、NkkkMkkkzazbzH00)(NrrMrrMNdzczAzzH11)()()()(n设系统稳定，将 代入，得传输函数 jzeNrrjMrrjMNjjdeceAeeH11)()()()(n在Z平面上向量：n分别称为零点矢量和极点矢量，用极坐标表示为 jrrjrrc Becd Bed rrjrrjrrc Bc ed Bd e ()11()()NrjjjrNrrc BH eAH eed B n对应的幅度谱：对应的幅度谱：n相位谱为：相位谱为： 当频率 从0变化到 时，这些向量的终点B沿单位圆逆时针旋转一周。 n如果滤波器的N个零点等间隔分布在单位圆上，设N=8，极零点分布如下图。当从0变化到时

10、 ，每遇到一个零点，幅度为零。具有如下图所示幅度特性的滤波器称为梳状滤波器。梳状滤波器。 11()NrjrNrrc BH eAd B NrrNrra11)(22n梳状滤波器的极零点分布及幅度特性梳状滤波器的极零点分布及幅度特性 3. FIR数字滤波器 nFIR滤波器是全零点滤波器，即 n其冲激响应： 0 Mkky nb x nk0 Mkkh nbnknFIR数字滤波器的主要优点有：数字滤波器的主要优点有： 1、易于实现线性相位； 2、由于不存在极点，系统是稳定的。nFIR数字滤波器的缺点：数字滤波器的缺点： 需要较长的滤波器长度以满足短过渡带的要求，造成所需的计算量较大。线性相位 nFIR滤波

11、器由于是全零点滤波器，其冲激响应为有限长度，容易满足形如 的对称的形式，因此容易获得线性相位的特性n FIR滤波器的系统函数 令 ，则 1h nh Nn 1100( ) 1NNnnnnH zh n zh Nn z 1mNn (1)1( )()NH zzH z n因此， 的零点也是 的零点，反之也然。这说明线性相位的FIR滤波器其零点相对于单位圆对称分布，如图为线性相位FIR滤波器的零点分布 ( )H z1()H z4. FIR数字滤波器的设计 nFIR数字滤波器设计思想主要是建立在近似理想滤波器频率特性基础上，近似的方法有：窗函数法；频率采样法；最佳一致逼近法等。我们主要介绍窗函数法，窗函数法

12、的主要思想是对理想低通滤波器的冲激响应序列采用截断并移位的方法，使之称为一个因果的有限长的序列，将对应的理想系统转化为物理可实现的。采用不同的窗函数导致了不同的设计方法，所设计的滤波器也具有不同的特性。 理想滤波器理想的低通滤波器n低通滤波器是构建其它滤波器的基础。然而理想的低通滤波器物理上却不可实现！？窗函数法设计FIR滤波器1. 由理想的频率响应 得到理想的 ；( )dh n2. 由 得到因果、 有限长的单位抽样响应 ;( )dh n( )h n3. 对 加窗得到较好的频率响应。( )h n理想频率响应()jdHe22cc1思路与方法：思路与方法：设理想低通滤波器的幅频为1，相频为零：1(

13、 )()2jj nddh nHeed12ccjnedsin()cnn则：特点： 无限长 非因果 偶对称解决方法： 截短， 移位 保留( )(2)dh nh nM,0,1,.,nM即：隐含着使用了窗函数n时域的乘积频域的卷积0( )( )MnnH zh n z于是：注意： 是因果的，且是线性相位的，即( )2M 即事先给一线性相位2()jjMdHee为了省去每次的移位，可以令： 在通带内2()0j MjdeHe0cc这样：( )H z21sin(2)( )2(2)ccjMj ncdnMh neednM( )( )dh nh n0,1,.nM于是：使用了矩形窗上式的的表达式及设计 的思路可推广到高

14、通、带阻及带通滤波器，也可推广到其它特殊类型的滤波器。实际上，给定一个 ，只要能积分得到 ，即可由截短、移位的方法得到因果的、且具有线性相位的FIR滤波器 。( )dh n( )H z()jdHe( )H z如何设计FIR滤波器？n上述分析表明，关键在于窗函数的设计n如何由技术指标转化为窗函数的设计参数呢？n频带指标n依采样频率的归一化n频带起伏指标n简化为对于较小的passn指标转化为n窗函数和这些技术指标之间的关系如何呢？也就是说具体的窗函数如何同这些指标相对应。n进一步的工程简化 窗函数的使用在数字信号处理中是不可避免的。数据、频谱、自相关函数等都需要截短。对窗函数提出那几方面的要求？n

15、3 dB 带宽（主瓣宽度）越小越好带宽（主瓣宽度）越小越好n旁瓣起伏越小越好旁瓣起伏越小越好n旁瓣衰减速度越快越好（过渡带宽度越小越好）旁瓣衰减速度越快越好（过渡带宽度越小越好）nKaiser窗nKaiser窗的设计参数为,Nn同指标参数如何对应呢？010203000.51 boxcar-0.500.5-60-40-200010203000.51 triang-0.500.5-100-500010203000.51 hamming-0.500.5-100-500常见的几种窗函数及其对应滤波器的频率响应常见的几种窗函数及其对应滤波器的频率响应设计步骤n不同的窗函数的设计在于其参数同指标要求之间的

16、对应关系不同；n无外乎窗口的长度如何确定，窗口的形状参数如何同指标相对应；n工程设计中通常采用简化的经验设计公式，具体的内容可以查工程手册得到频率采样法（时域的设计方法）n对于所要求的一个任意的频率响应 可以采用频域采样求和的方式来近似()H j1 ()iMjnddiiMhnHeN2iiN() iMjndidnMHh n en然后对于所得到的 应用前面介绍的窗函数得到最后的冲激响应 h n dh nw n h nMclear allclffc=15; %kHz%fs=40; %kHz%f = -fs/2,-fc,-fc, 0, fc,fc,fs/2;Df=0, 0, 1, 0.5,1, 0,

17、0;axis(-fs/2,fs/2,-0.2,1.2)hold onplot(f,Df)xlabel(analog frequency f in kHz)title(desired AF frequency response D(f)gridhold offpauseclfw=2*pi*f/fs;Dw=Df;axis(-pi,pi,-0.2,1.2)hold onplot(w,Dw)xlabel(digital frequency w in rad)title(desired DF frequency response D(w)gridhold offpauseclfN=21; M=(N-1)

18、/2;delta_f=fs/N;fk=-fs/2*(1-1/N): delta_f : fs/2*(1-1/N);Dfk=(1+abs(fk)/fc)/2).*(sign(fk+fc)-sign(fk-fc)/2;axis(-fs/2,fs/2,-0.2,1.2)hold onplot(fk,Dfk,*)xlabel(analog frequency f in kHz)title(samples of D(f)gridf = -fs/2,-fc,-fc, 0, fc,fc,fs/2;plot(f,Df,r)hold offpauseclfk=-M:M;Dk=Dfk;axis(-M,M,-0.2

19、,1.2)hold onbar(k,Dk,0) plot(k,Dk,*)plot(-M,M,0,0)xlabel(k)title(D(k)gridhold offpauseclfj=sqrt(-1);wi=2*pi*(-M:M)/N;for k=-M:M d(k+M+1)=Dk*exp(-j*wi*k)/N;end dk=real(d);k=-M:M;axis(-M,M,-0.3,0.8)hold onbar(k,dk,0) plot(k,dk,*)plot(-M,M,0,0)xlabel(k)title(d(k)gridhold offpauseclfk=0:N-1;hn=dk;axis(0

20、,N-1,-0.3,0.8)hold onbar(k,hn,0) plot(k,hn,*)plot(0,N-1,0,0)xlabel(n)title(rectangular windowed h(n)gridhold offpauseclfi=1;Hf=0;for f=-fs/2:(fs/200):fs/2 Hf(i)=hn*exp(-j*2*pi*f/fs*(0:N-1); i=i+1;endHf_abs=abs(Hf);f=-fs/2:(fs/200):fs/2;axis(-fs/2,fs/2,-0.3,1.2)hold onplot(f,Hf_abs)plot(-fs/2,fs/2,0,

21、0,k)xlabel(f in kHz)title(H(f),rectangular windowed)f = -fs/2,-fc,-fc, 0, fc,fc,fs/2;Df=0, 0, 1, 0.5,1, 0, 0;plot(f,Df,k-)plot(fk,Dfk,k*)gridhold offpauseclfwhm=0.54-0.46*cos(2*pi*(0:N-1)/(N-1);axis(0,N-1,-0.3,1.2)hold onbar(0:N-1,whm,0) plot(0:N-1,whm,*)plot(0,N-1,0,0)xlabel(n)title(Hamming window)

22、gridhold offpauseclfhn_hm=hn.*whm;axis(0,N-1,-0.3,1.2)hold onbar(0:N-1,hn_hm,0) plot(0:N-1,hn_hm,*)plot(0,N-1,0,0)xlabel(n)title(Hamming windowed H(n)gridhold offpauseclfi=1;Hf=0;for f=-fs/2:(fs/200):fs/2 Hf(i)=hn_hm*exp(-j*2*pi*f/fs*(0:N-1); i=i+1;endHf_abs=abs(Hf);f=-fs/2:(fs/200):fs/2;axis(-fs/2,

23、fs/2,-0.3,1.2)hold onplot(f,Hf_abs)plot(-fs/2,fs/2,0,0,k)xlabel(f in kHz)title(H(f),Hamming windowed)f = -fs/2,-fc,-fc, 0, fc,fc,fs/2;Df=0, 0, 1, 0.5,1, 0, 0;plot(f,Df,k-)plot(fk,Dfk,k*)gridhold offpauseclfL=201;P=(L-1)/2;delta_f=fs/L;fk=-fs/2*(1-1/L): delta_f : fs/2*(1-1/L);Dfk=(1+abs(fk)/fc)/2).*(sign(fk+fc)-sign(fk-fc)/2;wi=2*pi*(-P:P)/L;for k=-P:P d(k+P+1)=Dfk*exp(-j*wi*k)/L;end dk=real(d);k=-P:P;axis(-P,P,-0.3,0.8)hold onbar(k,dk,0) plot(-P,P,0,0)xlabel(k)title(d(k),LN)gridhold offpauseclfhn=dk(P-M+1:P+M+1);k=0:N-1;axis(0,