第5章
CHAPTER 5


数字滤波器的结构












5.1〓数字滤波器的基本知识
微课视频


5.1概述
1. 数字滤波器

数字滤波器是指完成信号滤波处理功能的、用有限精度算法实现的离散时间线性非时变(DTLTI)系统,它的作用是利用DTLTI系统的特性对输入序列的波形或频谱进行加工处理,将输入数字序列通过一定的运算后转变为输出数字序列,从而达到改变信号频谱的目的。第1~4章介绍了数字信号和数字系统的主要特性,即信号频谱、系统函数、系统频率响应、稳定性和因果性,并阐述了一些基本分析方法,如DTFT、ZT、DFT(FFT)。本章将用前面所学到的基本分析方法研究数字滤波器,分析它的结构、特点及实现方法。

对于一种离散时间线性非时变系统,数字滤波器可以用时域输入输出差分方程、单位采样脉冲响应以及系统函数等表示,在研究它的实现方法,即运算结构时,一般采用方框图或信号流图表示系统。数字滤波器的运算结构对于滤波器的设计及其性能指标的实现是非常重要的,本章研究在已知系统输入输出(I/O)差分方程(或系统函数)的情况下,如何设计不同的运算结构来实现它。这里将根据系统的单位采样脉冲响应的时间特性把数字滤波器分为两类: 无限长脉冲响应(infinite impulse response,IIR)滤波器和有限长脉冲响应(finite impulse response,FIR)滤波器。

2. 表示方法

数字滤波器的数学等效描述主要有以下几种,每一种描述从不同的角度指出了数字滤波器的特性。

(1) 时域I/O差分方程:
y(n)=∑Mj=0bjx(n-j)-∑Ni=1aiy(n-i)(5.1.1)
其中,x(n)、y(n)分别为数字系统的输入输出序列; {ai,i=1,2,…,N; bj,j=0,1,2,…,M}为常数,也称为数字系统的结构参数。

(2) 卷积方程:
y(n)=∑∞i=0h(i)x(n-i)(5.1.2)

其中,x(n)、y(n)分别为数字系统的输入输出序列; h(n)是数字系统的单位采样脉冲响应。

(3) 单位采样脉冲响应: 
h(n)={h(0),h(1),…,h(M),…}(5.1.3)
(4) 系统函数: 
H(z)=∑Mj=0bjz-j1+∑Ni=1aiz-i(5.1.4)
(5) 频率响应: 
H(ω)=H(ejω)=∑∞n=-∞h(n)e-jωn(5.1.5)
(6) 零点、极点分布: 
{zj,j=1,2,…,M;pi,i=1,2,…,N}(5.1.6)
若已知数字系统的系统函数式(5.1.4),则零点zj 、极点pi 可通过求系统函数分子多项式、分母多项式的零点得到。

(7) 网络结构图和相应的(状态变量分析)样值处理算法: 这种描述表示DTLTI系统的算法结构,特别是系统内部网络结构。

3. 分类

数字滤波器按其单位采样脉冲响应h(n)的时间特性可分为无限长脉冲响应(IIR)滤波器和有限长脉冲响应(FIR)滤波器。IIR滤波器的单位采样脉冲响应h(n)包含无限个非零值,即持续时间为无限长; FIR滤波器的单位采样脉冲响应h(n)只包含有限个非零值,即持续有限长时间。

IIR滤波器的时域I/O差分方程(或系统函数的分母多项式)中至少有一个系数ai不为零,一般可表示为式(5.1.1)或式(5.1.4)的形式。因此,这类滤波器存在输出到输入的反馈,需要使用递归计算方法实现。

如果IIR滤波器的I/O差分方程(或系统函数分母多项式)中系数ai=0(i=1,2,…,N),则IIR滤波器演化为FIR滤波器,这类滤波器实现结构中没有输出到输入的反馈,主要为非递归型结构,其数学表示简化为

y(n)=∑Mi=0bix(n-i)(5.1.7)
对应的系统函数为

H(z)=∑Mi=0biz-i(5.1.8)

在实际应用中,先由给定所需的滤波器频率响应H(ω)得出滤波器设计的目标: 对IIR滤波器常求出其系统函数H(z),FIR滤波器常求出其单位采样脉冲响应h(n)。然后,由H(z)或h(n)可得到合适的网络结构图,即滤波器的算法结构实现。

5.2信号流图表示网络结构

由式(5.1.1)可以看出,实现一个数字滤波器需要几种基本的运算单元: 加法器、单位延迟器和常数乘法器。这些基本的单元有两种图形表示法: 方框图法和信号流图法,如图5.2.1所示。因而一个数字滤波器的运算结构也有这两种图形表示法,用方框图表示较明显直观,用信号流图表示则更加简单方便。



图5.2.1基本运算的方框图表示及信号流图表示


下面以信号流图表示为例,介绍并解释有关系统网络结构图描述的术语和运算规则。

(1) 节点: 在图5.2.2的信号流图中,圆点称为节点,每个节点都对应一个变量,或者说代表一个信号。

(2) 源节点: 对于一个节点,流入该节点的信号叫输入,流出该节点的信号叫输出。若一个节点只有输出支路与之相连接,则称为源节点或输入节点。如图5.2.2中的x(n)即为源节点。

(3) 吸收节点: 若一个节点只有输入支路与之相连接,则称为吸收节点或输出节点。如图5.2.2中的y(n)即为吸收节点。

(4) 支路: 连接两个节点的有向线段称为支路。两个节点间的支路相当于一个乘法器(延迟器在“z域”运算中可等效为“乘法器”),其加权系数或者说乘法器系数叫做支路传输增益。和每个节点连接的可以有输入支路和输出支路,节点变量值等于所有输入支路的输出值之和。在图5.2.2(a)中,节点变量值之间的运算关系如下:

w1(n)=w2(n-1)
w2(n)=w′2(n-1)
w′2(n)=x(n)-a1w2(n)-a2w1(n)
y(n)=b2w1(n)+b1w2(n)+b0w′2(n)(5.2.1)

显然,不同的网络信号流图代表不同的运算方法,而对于同一个系统函数可以有很多种信号流图与之相对应。从基本运算考虑,满足以下条件的信号流图,称为基本信号流图。

(1) 信号流图中所有支路都是基本的,即支路增益是常数或者是单位延迟(信号流图中的z-1物理意义上表示单位延迟,但在z域运算中,可作为符号参量进行代数运算)。

(2) 信号流图环路中必须存在延迟支路。

(3) 信号流图中节点和支路的数目是有限的。

图5.2.2(a)是基本信号流图,图中有两个环路,环路增益分别为-a1z-1和-a2z-2,且环路中都有延时支路; 而图5.2.2(b)不是基本信号流图,它不能决定一种具体算法的运算关系,不满足基本信号流图的条件。



图5.2.2信号流图


根据信号流图可以求出网络(数字系统)的系统函数,方法是列出各个节点变量方程,求解该方程组,推导出输出与输入之间的关系。

例5.2.1求图5.2.2(a)信号流图决定的系统函数H(z)。

解: 对式(5.2.1)进行z变换,得到:

W1(z)=W2(z)z-1
W2(z)=W′2(z)z-1
W′2(z)=X(z)-a1W2(z)-a2W1(z)
Y(z)=b2W1(z)+b1W2(z)+b0W′2(z)

经过联立求解,消去网络内部信号参量w′2、w2λ、w1的ZT函数,得到:

H(z)=Y(z)X(z)=b0+b1z-1+b2z-21+a1z-1+a2z-2

5.3无限长脉冲响应滤波器的基本网络结构

通过前两节分析可知,IIR滤波器的系统单位采样脉冲响应h(n)为无限长序列,系统函数H(z)在有限z平面上存在极点。其网络结构的特点是含有反馈环路,即运算结构上是递归型的,但具体实现起来,其结构并不唯一。同一个系统I/O差分方程(或系统函数),可以有各种不同的结构形式,主要有5种,即直接型Ⅰ、直接型Ⅱ(典型型)、级联型、并联型和转置型。






5.3.1〓IIR直接Ⅰ型滤波器结构

微课视频


5.3.1直接型Ⅰ
1. 二阶数字滤波器



给定IIR滤波器的传输函数H(z)=N(z)D(z)=b0+b1z-1+b2z-21+a1z-1+a2z-2,相应的时域I/O差分方程为

 y(n)=-a1y(n-1)-a2y(n-2)+b0x(n)+b1x(n-1)+b2x(n-2)

直接型Ⅰ可实现上述差分方程表示的数字滤波器,如图5.3.1所示。



图5.3.1IIR滤波器直接型Ⅰ结构


由图5.3.1可以总结IIR滤波器直接型Ⅰ的结构特点如下: 

(1) 一个相加器把系统时域I/O差分方程右边的所有项都加起来了。加法器的输出是y(n)。

(2) 前向分支加权乘因子为bj,j=0,1,2,相应于H(z)的分子多项式系数,在时域I/O差分方程中是与输入x(n)有关的非递归项。

(3) 后向分支加权乘因子为-ai,i=1,2,相应于H(z)的分母多项式非零次幂项系数,在时域I/O差分方程中是与输出y(n)有关的反馈(递归)项。

注意: 网络结构图中反馈加权乘因子是传输函数分母多项式非零次幂项系数的负数。

(4) 已知直接型Ⅰ滤波器,给定输入序列x(n),求系统的响应y(n),可用样值迭代处理方法: 数字信号处理器(DSP)按设定的算法,一次处理一个采样值x(n),同时送出一个相应的输出样值y(n)。此方法常用于实时处理,如长序列的实时滤波、数字音响、数字控制系统、自适应信号处理等。与样值迭代处理算法相对应的是块处理方法: 信号处理是一段段地进行的,如有限长信号的卷积、长序列的快速卷积、语音的分析与合成、图像处理等。

图5.3.2给出了图5.3.1所示滤波器的内部状态矢量v、w及其与输入输出的关系。各内部状态值表示时间n处,图5.3.1结构图中输入样值、输出样值、各延迟寄存器的内容。当时间单元加1,各单元内容刷新。


v0(n)=x(n)w0(n)=y(n)
v1(n)=x(n-1)=v0(n-1)w1(n)=y(n-1)=w0(n-1)
v2(n)=x(n-2)=v1(n-1)w2(n)=y(n-2)=w1(n-1)
v1(n+1)=v0(n)w1(n+1)=w0(n)
v2(n+1)=v1(n)w2(n+1)=w1(n)


图5.3.2滤波器的内部状态及输入输出关系



根据图5.3.1和系统I/O差分方程,输出序列的计算可用下面的样值迭代算法: 对每个输入样值x,令

v0=x(读入当前输入样值)
w0=-a1w1-a2w2+b0v0+b1v1+b2v2(w1、w2、v0、v1、v2可按定义初始化)
y=w0(读出当前输出样值)

网络内部状态刷新要从图5.3.1所示结构的下面向上进行,只需刷新延迟寄存器的内容:

v2=v1,w2=w1
v1=v0,w1=w0

2. 任意阶数字滤波器



给定传输函数H(z)=N(z)D(z)=b0+b1z-1+b2z-2+…+bMz-M1+a1z-1+a2z-2+…+aNz-N,其阶数为max{M,N},相应的时域I/O差分方程为

y(n)=-a1y(n-1)-a2y(n-2)-…-aNy(n-N)+b0x(n)
+b1x(n-1)+…+bMx(n-M)

当M=N时,就是N阶IIR滤波器。在图5.3.1所示的二阶IIR滤波器下面加更多的延迟环节和相应乘法器即可得到其直接型Ⅰ网络结构。


例5.3.1画出下列滤波器的直接型Ⅰ实现,并写出其样值迭代算法。

H(z)=2-3z-1+4z-31+0.1z-1-0.3z-2+0.6z-4

解: 所给系统的时域I/O差分方程为

y(n)=-0.1y(n-1)+0.3y(n-2)-0.6y(n-4)
+2x(n)-3x(n-1)+4x(n-3)

按照该时域I/O差分方程画出如图5.3.3所示的IIR滤波器的直接型Ⅰ网络结构。



图5.3.3例5.3.1图


样值迭代算法: 

对每个输入样值x,令

v0=x
w0=-0.1w1+0.3w2-0.6w4+2v0-3v1+4v3
y=w0

网络内部状态刷新要从图5.3.3所示结构的下面向上进行:

w4=w3
w3=w2,v3=v2
w2=w1,v2=v1
w1=w0,v1=v0

该IIR滤波器的抽头系数:

a=[a0,a1,a2,a3,a4]=[1,0.1,-0.3,0,0.6]
b=[b0,b1,b2,b3]=[2,-3,0,4]

内部状态矢量:

w=[w0,w1,w2,w3,w4]
�瘙經=[v0,v1,v2,v3]






5.3.2〓IIR直接Ⅱ型滤波器结构

微课视频


5.3.2直接型Ⅱ(典型型)
1. 二阶数字滤波器



给定IIR滤波器的传输函数H(z)=N(z)D(z)=b0+b1z-1+b2z-21+a1z-1+a2z-2,相应的时域I/O差分方程为

y(n)={b0x(n)+b1x(n-1)+b2x(n-2)}+{-a1y(n-1)-a2y(n-2)}

(1) 把差分方程分成两部分: 递归部分和非递归部分。这种分组相当于把直接型Ⅰ的一个大相加器分裂成两个,如图5.3.4(a)所示。

(2) 从传输函数的等效性,整个系统可看成是两个传输函数分别为N(z)、1D(z)的滤波器级联: 

H(z)=N(z)1D(z)



图5.3.4IIR滤波器直接型Ⅱ结构



(3) 理论上,级联的次序对整个系统是没关系的,可前后交换,H(z)=1D(z)N(z),如图5.3.4(b)所示。

(4) 设图5.3.4(b)中第一个滤波器1D(z)的输出为w(n),可见前后两个滤波器的两套延时器的内容完全一样。所以两套延时器可合并成一套,这就构成了直接型Ⅱ(canonical form,典型型),方框图表示法见图5.3.4(c),信号流图表示见图5.3.4(d)。由于直接型Ⅱ比直接型Ⅰ的延迟单元数量少,有些文献中把直接型Ⅱ就叫直接型。

(5) I/O差分方程描述的时域运算可通过下列迭代算法实现:

w(n)=x(n)-a1w(n-1)-a2w(n-2)
y(n)=b0w(n)+b1w(n-1)+b2w(n-2)

(6) 对每一时刻n,w(n-1)、w(n-2)是两共享延时寄存器w1(n)和w2(n)的内容。它们也是滤波器的内部状态。

(7) 样值迭代算法: 定义图5.3.4(c)所示IIR滤波器的内部状态矢量为w=[w0,w1,w2]。对每一输入样值x,令

w0=x-a1w1-a2w2
y=+b0w0+b1w1+b2w2

内部状态刷新要从图5.3.4(c)所示结构的下面向上进行:

w2=w1
w1=w0

图5.3.4(c)所示IIR滤波器的内部状态变量间的时序关系如下:

w0(n)=w(n)
w1(n)=w(n-1)=w0(n-1)w1(n+1)=w0(n)
w2(n)=w(n-2)=w1(n-1)w2(n+1)=w1(n)

2. 任意阶数字滤波器



给定传输函数H(z)=N(z)D(z)=b0+b1z-1+b2z-2+…+bMz-M1+a1z-1+a2z-2+…+aNz-N,其阶数=max{M,N},相应的时域I/O差分方程为

y(n)=-a1y(n-1)-a2y(n-2)-…-aNy(n-N)
+b0x(n)+b1x(n-1)+…+bMx(n-M)

当M=N时,就是N阶IIR滤波器。其典型型网络结构是在图5.3.4(c)所示的二阶滤波器典型型下面加更多的延迟环节和相应乘法器。

当M≠N时,公共延时器的数目K=max{M,N}。不失一般地,设N>M,此时: 

(1) 任意阶IIR滤波器时域I/O差分方程可等效为

w(n)=x(n)-a1w(n-1)-…-aNw(n-N)
y(n)=b0w(n)-b1w(n-1)-…-bMw(n-M)

(2) 样值迭代算法: 定义IIR滤波器典型型网络的内部状态矢量w=[w0,w1,w2,…,wN]。

对每一输入样值x,令

w0=x-a1w1-a2w2-…-aNwN
y=b0w0+b1w1+…+bMwM

内部状态变量刷新要从典型型网络结构的下面向上进行: 

wi=wi-1,i=K,K-1,…,1





5.3.3〓级联型、并连型网络结构

微课视频


5.3.3级联型

直接型Ⅰ和直接型Ⅱ(典型型)结构简单,但IIR滤波器各抽头系数对其性能控制不直接,且相互影响,为此可采用级联型结构。用实抽头系数的一阶、二阶滤波器的级联可实现任意阶的数字滤波器。

1. 二阶节



实际物理可容易实现的滤波器系统函数H(z)是实系数的多项式的分式。相应的系统函数的零点和极点是实数,或是成对的共轭复数。共轭成对的零点、极点对应实系数二次多项式。级联型的滤波器单级基本结构是“二阶滤波器”,称二阶节(second order section,SOS)。二阶节可用直接型Ⅰ或典型型实现。式(5.3.1)给出了L级SOS级联型滤波器总系统函数H(z)的一般表示式:

H(z)=A∏Li=1Hi(z)(5.3.1)

其中,Hi(z)=N(z)D(z)=bi0+bi1z-1+bi2z-21+ai1z-1+ai2z-2,i=1,2,…,L是第i级二阶节的系统函数; A是归一化常数。

因为可把实际一阶系统看作二阶节的特例,所以总级联型实现的系统函数H(z)分子、分母多项式的实际阶数≤2L。级联型实现的数字滤波器方框图如图5.3.5所示。



图5.3.5采用级联型的H(z)方框图


2. 级联型滤波器实现特点

硬件实现时,可用一个二阶节“时分复用”单元,以节省存储单元。

因为每一个二阶节都对应一对零点和一对极点,所以可单独调整滤波器的零点、极点,从而方便地调整滤波器的频率响应。

级联次序不同,运算误差不同,且极间电平移动不同(极间输出电平太大,则滤波器运算可能溢出; 极间输出电平太小,则输出信号的信噪比会很小)。

例5.3.2画出给定滤波器的级联型和典型型实现。并写出相应的差分方程和样值迭代算法。

H(z)=3-4z-1+2z-21-0.6z-1+0.5z-23+4z-1+2z-21+0.6z-1+0.5z-2=H0(z)H1(z)
=9-4z-2+4z-41+0.64z-2+0.25z-4
 
解: 

(1) 图5.3.6(a)是滤波器的级联型实现; 图5.3.6(b)是滤波器的典型型实现。可见,级联型的抽头系数(5×2=10)比典型型的(5)多,硬件实现上需要的乘法器数量就比典型型的多。实时要求高的场合,常用典型型。级联型滤波器的模块化结构规范,便于根据需要实时配置抽头系数,实现SOS的分时复用,应用更加灵活。



图5.3.6例5.3.2图


(2) 级联型。

时域I/O差分方程为

w0(n)=x(n)+0.6w0(n-1)-0.5w0(n-2)
x1(n)=3w0(n)-4w0(n-1)+2w0(n-2)
w1(n)=x1(n)-0.6w1(n-1)-0.5w1(n-2)
y(n)=3w1(n)+4w1(n-1)+2w1(n-2) 

样值迭代算法: 对每一级定义内部状态矢量w0=[w00,w01,w02],w1=[w10,w11,w12]。

对每一输入样值x,做下列循环:

w00=x+0.6w01-0.5w02
x1=3w00-4w01+2w02 
w02=w01w01=w00(SOS内部状态刷新)
w10=x1-0.6w11-0.5w12
y(n)=3w10+4w11+2w12
w12=w11w11=w10(SOS内部状态刷新)

(3) 典型型。

时域I/O差分方程为

w(n)=x(n)-0.64w(n-2)-0.25w(n-4)
y(n)=9w(n)-4w(n-2)+4w(n-4)

样值迭代算法: 定义内部状态矢量w=[w0,w1,w2,w3,w4]。

对每一输入样值x,做下列循环:

w0=x-0.64w2-0.25w4
y(n)=9w0-4w2+4w4
w4=w3w3=w2w2=w1w1=w0(内部状态刷新)

讨论: 两种滤波器实现形式,若仅考虑延迟寄存器的状态为独立内部状态,则总的内部状态变量数一样,均为4个。

5.3.4并联型

若把数字滤波器系统函数表示为实系数二阶多项式部分因式之和的形式,则总系统函数是各二阶节表示的子系统的并联:


H(z)=H1(z)+H2(z)+…+HL(z)(5.3.2)

并联型滤波器总输出y(n)的z变换Y(z)可以写为

Y(z)=H1(z)X(z)+H2(z)X(z)+…+HL(z)X(z)(5.3.3)

这意味着输入序列x(n)通过并联的L个子滤波器后,在输出端把它们累加起来就可得到总输出y(n)。每个子滤波器通常选用式(5.3.4)所示的SOS:


图5.3.7采用并联型H(z)

方框图


Hi(z)=bi0+bi1z-1+bi2z-21+ai1z-1+ai2z-2,i=1,2,…,L(5.3.4)



每个Hi(z)由一对共轭复极点产生,并且可用直接型Ⅰ或典型型的形式实现。并联型实现的数字滤波器方框图如图5.3.7所示。

并联型滤波器的实现特点如下:  

(1) 运算速度快,各子网络对输入同时计算。它较之直接型、级联型结构,运算速度最快。

(2) 可单独调整滤波器极点。

(3) 各基本节的运算误差互不影响。它较之直接型、级联型结构,累积运算误差最小。

5.3.5转置型


根据线性信号流图中的转置理论: 把基本信号流图网络中所有支路的方向颠倒,且将输入x(n)和输出y(n)位置互换,则转置后的网络结构的输入输出关系和原来网络的一样(系统传输函数保持不变)。

所有上述各种滤波器的实现都可得到其转置型。两者传输函数一样,但有限字长的影响不同。例如,对图5.3.4(d)的直接型Ⅱ结构,转置后的网络如图5.3.8所示,画成输入在左方,输出在右方的习惯形式,则如图5.3.9所示。




图5.3.8直接型Ⅱ结构的转置




图5.3.9输入在左,输出在右








5.4有限长脉冲响应滤波器的基本结构

FIR滤波器的系统单位采样脉冲响应h(n)为有限长序列,系统函数H(z)的极点在z=0处,其网络结构的特点是没有反馈环路,即运算结构是非递归型的。






5.4 1. FIR滤波器的基本结构——直接型、级联型
微课视频


1. 直接型FIR滤波器

从卷积形式的滤波器描述知道,IIR滤波器的单位采样脉冲响应h(n)无限长; 而FIR滤波器的h(n)是有限长,故FIR滤波器的直接型实现可直接用h(n)作为滤波器的系数(延时抽头的权值),即横向滤波器(直接卷积型结构)。其时域I/O差分方程为

y(n)=∑N-1i=0h(i)x(n-i)(5.4.1)

按照式(5.4.1)的时域I/O差分方程直接画出FIR滤波器直接型结构图如图5.4.1所示。



图5.4.1FIR滤波器直接型结构图


2. 级联型FIR滤波器

将滤波器系统函数H(z)分解成实系数二阶因子的乘积形式:

H(z)=∑N-1n=0h(n)z-n=∏N2k=1(β0k+β1kz-1+β2kz-2)(5.4.2)

式中,N2表示取N2的整数部分。若N为偶数,则N-1为奇数,故系数β2kk=1,2,…,N2中至少有一个为零。这是因为,这时H(z)在z域上有奇数个根,其中复数根成共轭对,为偶数个,综合起来,H(z)应有奇数个非零实根,而H(z)的实数根对应一阶因子。图5.4.2画出了N为奇数时,FIR滤波器的级联结构,其中每个二阶因子用图5.4.1所示的FIR滤波器直接型结构。



图5.4.2FIR滤波器的级联型结构(N为奇数)







5.4 3. FIR滤波器的基本结构——线性相位型
微课视频


3. 线性相位型FIR滤波器

FIR滤波器系统的最主要特性之一就是它可以构成具有线性相位特性的滤波器。

(1) 线性相位特性: 滤波器对不同频率ωi的输入正弦波信号,输出响应所产生的相移φi和该输入正弦波的频率呈线性关系,即

φi=kωi

其中,k是滤波器线性相频特性的斜率。

(2) 信号无频率失真: 对线性相位的滤波器,当滤波器通带内通过多个频率信号时,输出各频率分量的瞬时相位均满足

φT(t)=ωit+kωi=ωi(t+k)

其中,k可认为是所有频率分量的延迟。因为滤波器线性相频特性的斜率决定了输入信号的所有频率分量均产生延迟T=k(k<0表示右移——延迟,k>0表示左移——超前),所以滤波器通带内的信号输入无频率失真。

(3) 线性相位型FIR滤波器结构特点: 可证明单位采样脉冲响应长度为N的线性相位的因果FIR滤波器,其单位采样脉冲响应h(n)是实的,且相对N-12一定是奇(或偶)对称的(N=奇数时具有实际的对称中心点,N=偶数时是虚的对称点)。线性相位型FIR滤波器h(n)的对称性使其网络结构只需N2(N为偶数)或N-12(N为奇数)个乘法器进行相应的乘法运算。线性相位型FIR滤波器的结构如图5.4.3所示。


例5.4.1FIR滤波器传递函数为H(z)=1+16+116z-4+z-8。

(1) 确定并画出滤波器的直接型、线性相位型和级联型网络结构。

(2) 画出用这3种形式表示时,滤波器的单位采样脉冲响应波形。

解: 例5.4.1涉及的MATLAB程序参考实现见5.6节。

(1) 直接型时域I/O差分方程为

y(n)=x(n)+16.0625x(n-4)+x(n-8)

线性相位型时域I/O差分方程为

y(n)=[x(n)+x(n-8)]+16.0625x(n-4)

级联型的FIR滤波器系统函数的系数可编写以下MATLAB程序计算得到。

b0 =1 ; (FIR滤波器级联型系统函数归一化系数)

B =1.00002.82844.0000 ; (输出数组B的每行对应式(5.4.2)中一级,k=1)

1.00000.70710.2500 ; (顺序地对应式(5.4.2)中的系数β0k、β1k、β2k,k=2)

1.0000-0.70710.2500 ; (k=3)

1.0000-2.82844.0000 ; (k=4)

A =100

100

100  

100


此FIR滤波器的直接型、线性相位型和级联型结构分别如图5.4.4(a)、图5.4.4(b)和图5.4.4(c)所示。

(2) 滤波器单位采样脉冲响应h(n)波形可通过运行本书编写的MATLAB程序实现,如图5.4.5所示。



图5.4.3线性相位型FIR滤波器的结构




图5.4.4例5.4.1中FIR滤波器结构图





图5.4.5例5.4.1中FIR滤波器的单位采样脉冲响应图






5.4 4. FIR滤波器的基本结构——频率采样型
微课视频


4. 频率采样型FIR滤波器

设函数H(k)是FIR滤波器单位采样脉冲响应h(n)的N点DFT。根据DFT与ZT的函数关系,FIR滤波器的系统函数可写为

H(z)=1N(1-z-N)∑N-1k=0H(k)1-W-kNz-1(5.4.3)



其中,H(k)=H(z)z=ej2πNk,k=0,1,2,…,N-1。

根据式(5.4.3)构成的FIR滤波器网络结构,就是其频率采样型结构,即用系统函数为Hc(z)=(1-z-N)的FIR(梳状)滤波器与IIR滤波器级联而成。这里的IIR滤波器是N个一阶网络Hk(z)的并联: 

∑N-1k=0Hk(z)=∑N-1k=0H(k)1-W-kNz-1

每个系统函数为Hk(z)的一阶子网络,均是一个无损耗的IIR谐振器; 谐振频率分别是数字域频率2πNk。

频率采样型FIR的滤波器结构如图5.4.6所示。



图5.4.6频率采样型FIR的滤波器结构


1) 梳状陷波滤波器(notch filter)

梳状滤波器的系统函数一般表达式为

H(z)=1-rNz-N1-aNz-N(5.4.4)

其中,当a=0,r=1时,梳状滤波器称为“FIR(梳状)陷波滤波器”,系统函数退化为

H(z)=1-z-N(5.4.5)

式(5.4.5)所表示的FIR梳状滤波器网络结构如图5.4.7(a)所示。

FIR梳状滤波器系统函数的零点是令1-z-N=0,在z域求得

zk=ej2πNk,k=0,1,2,…,N-1

即z平面单位圆上的N个等分点就是FIR梳状滤波器系统函数的N个零点。

FIR梳状滤波器的极点是z平面上z=0处的N阶极点。综合起来,FIR梳状滤波器零点、极点分布如图5.4.7(b)所示。

令z=ejω代入式(5.4.5),可得其频率响应:

H(ejω)=1-e -jωN(5.4.6)

由式(5.4.6)可得FIR梳状滤波器幅频响应:

H(ω)=|H(ejω)|=2sinN2ω(5.4.7)

式(5.4.7)所对应的幅频响应曲线如图5.4.7(c)所示,是状如梳子的齿。由图5.4.7(c)可知,在频率ωk=2πNk(k=0,1,…,N-1)处的FIR梳状滤波器幅频响应为零,即陷波。



图5.4.7FIR梳状滤波器的结构,零点、极点分布和幅频响应


式(5.4.4)中,当a<1,且接近于1,r=1时,梳状滤波器为“IIR(梳状)陷波滤波器”,其系统函数具体为:

H(z)= 1-z-N1-aNz-N(5.4.8)


式(5.4.8)的零点、极点分布如图5.4.8(a)所示。由图5.4.8(a)可知,IIR梳状滤波器零点是z平面单位圆上的N个等分点:

zk=ej2πNk,k=0,1,2,…,N-1

N个极点在z平面单位圆内,接近单位圆:

azk=aej2πNk,k=0,1,2,…,N-1


z平面单位圆上旋转点z=ejω对零点、极点对(zi, azi)的距离近似相等, 其中i=0,1,…,N-1,但当单位圆上的点z旋转到频率ωi=i2πN弧度处,与该处零点zi的距离等于零,故IIR梳状滤波器幅频响应|H(ω)|近似是平坦的,只在频率点i2πN 上有极小值,其幅频响应如图5.4.8(b)所示,形状如梳子的齿。



图5.4.8IIR梳状滤波器零点、极点分布和幅频响应



例5.4.2梳状滤波器有两种形式: 

FIR梳状滤波器H1(z)=1-z-N

IIR梳状滤波器H2(z)=1-z-N1-aNz-N

分别对N=8,a=0.8和0.98计算并图示两种形式的梳状滤波器零点、极点分布及幅频特性曲线。最后讨论FIR、IIR梳状滤波器的特点。

解: 例5.4.2的MATLAB程序参考实现见5.6节。

MATLAB程序运行得零点、极点分布及幅频特性曲线如图5.4.9所示。


讨论: 

(1) 梳状滤波器可滤除输入信号中ω=k2πN(k=0,1,2,…,N-1)的频率分量,如用于消除电网工频谐波干扰,在彩色电视接收机中进行亮度、色度信号分离、三色信号分离等。

(2)  IIR梳状滤波器具有更平坦的通带特性和更窄的过渡带。在滤波器极点接近z平面上的单位圆时幅频响应通带的平坦特性更加明显。

2) 无损耗IIR并联谐振器

(1) 无损耗IIR谐振器。

前面学习频率采样型FIR滤波器网络结构时,已经介绍了无损耗IIR谐振器的系统函数为Hk(z)=H(k)1-W-kNz-1,是一阶数字系统,H(k)为式(5.4.3)所示FIR滤波器频率响应函

数的频域采样值,是常数。Hk(z)在z域单位圆上zk=ej2πNk处有一极点,物理上可解释为: 系统函数为Hk(z)的IIR滤波器对频率为ωk=2πNk的输入序列的响应是无穷大(∞),说明IIR滤波器Hk(z)对频率为ωk的输入信号不稳定,工程实践上表现为系统对该频率的信号“谐振”,且无损耗。故此一阶系统是谐振频率为ωk=2πNk的无损耗IIR谐振器。如果谐振器的频率响应不是无穷大,则是有损耗的。



图5.4.9例5.4.2梳状滤波器零点、极点分布和幅频特性曲线


(2) 并联IIR谐振器。

由多个谐振频率不同的一阶无损耗IIR谐振器并联,就构成了并联型IIR滤波器(并联IIR谐振器),其系统函数为

H(z)=∑N-1k=0H(k)1-W-kNz-1(5.4.9)

式(5.4.9)所示的系统函数在z域单位圆上共有N个极点zk=ej2πNk(k=0,1,2,…,N-1),说明分别对N个频率分量ωk=2πNk(k=0,1,2,…,N-1)是谐振的。

3) 修正的频率采样型FIR滤波器

(1) 频率采样型FIR滤波器的特点。

由上面对梳状滤波器和并联IIR谐振器的零点、极点分析可知,FIR滤波器频率采样型网络中梳状滤波器在z平面单位圆上的N个零点恰好可以和并联IIR谐振器在z平面单位圆上的N个极点一一对应相消。总体上,FIR滤波器对频率分量ωk=2πNk的响应正好是该


滤波器频率响应的频域采样值H(k),频率采样型网络结构控制FIR滤波器的响应很直接。由于上述零点、极点对消机理,FIR滤波器频响函数只在频率ωk=2πNk上取H(k),在其他频率点上,系统频率响应近似为常数,恰好保证FIR滤波器是稳定系统。

(2) 频率采样型FIR滤波器的最大问题。

上面分析得到频率采样型FIR滤波器中并联IIR谐振器的N个极点在z平面单位圆上,极点位置是由FIR滤波器理论上的系数决定的。在实际数字实现中,由于计算机或数字器件的有限字长效应等,滤波器的系数要量化,量化中一定有误差存在,故实现的并联IIR滤波器谐振点会移动(并联IIR谐振器的极点偏离理论位置),这可能使梳状滤波器的零点和并联IIR谐振器的极点不能完全对消,导致FIR滤波器系统不稳定。

(3) 修正频率采样型FIR滤波器。

为了改善因梳状滤波器零点和并联IIR谐振器极点可能不完全对消,导致FIR滤波器可能不稳定的问题,实际应用中可以如下构造“修正频率采样型FIR滤波器”,即把并联IIR谐振器的极点设置在z平面上半径r小于1但近似于1的圆周上; 相应地把梳状滤波器的零点也移到同一个半径小于1的圆周上,实现与并联IIR谐振器的极点一一对消。这时可实现的“修正”频率采样型FIR滤波器传输函数:

H(z)≈1N (1-rNz-N) ∑N-1k=0H(k)1-rW-kNz-1,r≤1(5.4.10)

将式(5.4.10)加以简化,求得并联IIR谐振器的各个极点(即H(z)的极点):

zk=rej2πNk,k=0,1,…,N-1(5.4.11)

(4) 根据各个极点的情况讨论“修正”频率采样型FIR滤波器网络结构。

式(5.4.11)中的共轭复极点在半径为r的圆周上,且以实轴为轴成对称分布,满足

zN-k=z*k(5.4.12)

即

rW-(N-k)N=rej2πN(N-k)=rej2πNk=(rW-kN)(5.4.13)

式中,WN=e-j2πN。

由于FIR滤波器单位采样脉冲响应h(n)是实数,故其N点DFT,H(k)=DFT[h(n)]也是共轭对称的。因此,可以将式(5.4.10)中的第k个与第(N-k)个一阶IIR谐振器合并为一个实系数的二阶网络,以Hk(z)表示:

Hk(z)=H(k)1-rW-kNz-1+H(N-k)1-rW-(N-k)Nz-1
=H(k)1-rW-kNz-1+H(k)1-r(W-kN)z-1
=a0k+a1kz-11-z-12rcos2πNk+r2z-2,
k=1,2,…,N-12,N为奇数
k=1,2,…,N2-1,N为偶数
(5.4.14)


其中,

a0k=2Re[H(k)]
a1k=-2rRe[H(k)WkN](5.4.15)

由于式(5.4.14)所示的二阶网络的极点在z平面单位圆内,而不在单位圆上,因而从滤波器频率响应的几何解释可知,它相当于一个有限Q值的IIR谐振器,谐振频率为ωk=2πNk,其网络结构如图5.4.10所示。



图5.4.10式(5.4.14)所示的二阶IIR谐振器Hk(z)


式(5.4.11)中的实极点对应的是一阶IIR网络,其系统函数为

H0(z)=H(0)1-rz-1(5.4.16)
HN2(z)=HN21+rz-1,当N为偶数时(5.4.17)

当然,由式(5.4.11)可知,当N为偶数时,有一对实极点(即k=0及k=N/2时的极点),对应有式(5.4.16)和式(5.4.17)所示两个一阶网络; 当N为奇数时,只有一个实极点(即k=0的极点),对应只有式(5.4.16)所示一个一阶网络H0(z)。

最后,将并联IIR谐振器的实极点、复极点以及梳状滤波器合起来,得到修正频率采样型FIR滤波器的实系数系统函数H(z)及其对应的总网络结构。

当N为偶数时,修正频率采样型FIR滤波器的系统函数H(z):

H(z)=(1-rNz-N)1NH0(z)+HN2(z)+∑N/2-1k=1Hk(z)(5.4.18)

对应的网络结构图如图5.4.11所示,图中由一阶网络H0(z),HN/2(z)和各二阶网络Hk(z)构成,k=1,2,…,N2-1。

当N为奇数时,修正频率采样型FIR滤波器的系统函数H(z):

H(z)=(1-rNz-N)1NH0(z)+∑(N-1)/2k=1Hk(z)(5.4.19)

对应的网络结构图,只需将图5.4.11所示的网络结构中去掉一阶网络HN/2(z),只由一阶网络H0(z)和各二阶网络Hk(z)构成,k=1,2,…,N-12。

图5.4.11中已经画出了一阶网络H0(z)和HN/2(z)的信号流图,而其中二阶网络Hk(z)的结构如图5.4.10所示。



图5.4.11修正频率采样型FIR滤波器结构



5.5格型结构

在数字滤波器网络结构中,格型网络是一类十分重要并且得到广泛应用的网络结构。本节分别讨论全零点(FIR)格型滤波器和全极点(IIR)格型滤波器。

1. 全零点(FIR)格型滤波器

一个N阶的FIR滤波器的系统函数H(z),不失一般地可写成式(5.5.1)的形式:

H(z)=B(z)=∑Ni=0biz-i=1+∑Ni=1b(i)Nz-i(5.5.1)

其中,b(i)N表示N阶FIR滤波器的第i个系数,并假设首项(零次幂)系数b0=1。H(z)对应的格型结构如图5.5.1所示。



图5.5.1全零点(FIR)格型滤波器网络结构




图5.5.2全零点(FIR)格型

结构基本单元


图5.5.1所示的网络结构可看成由N个如图5.5.2所示的格型网络单元级联而成,其中kn(n=1,2,…,N)又称为反射系数,每个格型网络单元有两个输入端和两个输出端。由图5.5.1可见,


输入信号x(n)同时送到第一级格型网络单元的两个输入端,而在输出端仅取最后一级格型网络单元上面的一个输出端作为整个FIR格型滤波器的输出信号y(n)。

下面以二阶FIR滤波器为例,掌握如何从已知的滤波器系统函数求出其格型网络结构系数,并画出其格型网络结构图。

例5.5.1设二阶FIR滤波器系统函数为

H(z)=1+14z-1-12z-2

求其格型网络结构系数,并画出格型网络结构图。

解: 由图5.5.1所示的一般FIR滤波器格型网络结构单元的输入输出关系可知

r1(n)=k1e0(n)+r0(n-1)
e0(n)=r0(n)=x(n)
e1(n)=e0(n)+k1r0(n-1)

e2(n)=e1(n)+k2r1(n-1)
=e0(n)+k1r0(n-1)+k2k1e0(n-1)+k2r0(n-2)
=x(n)+(k1+k2k1)x(n-1)+k2x(n-2)

因此,图5.5.1所示的FIR滤波器“二级格型网络”的输入输出关系可用二阶差分方程描述:

y(n)=x(n)+(k1+k2k1)x(n-1)+k2x(n-2)


题目已知二阶FIR滤波器系统函数为H(z)=1+14z-1-12z-2,对应的时域I/O差分方程为

y(n)=x(n)+14x(n-1) -12x(n-2)

所以令

k2=-12k1+k2k1=14

解上述方程组,求出所给FIR滤波器格型网络结构系数为

k1=12,k2=-12

最后画出所给全零点(FIR)滤波器格型网络结构图,如图5.5.3所示。



图5.5.3例5.5.1中的FIR滤波器格型网络结构


2. 全极点(IIR)格型滤波器

给定全极点滤波器的系统函数,可以根据FIR滤波器格型网络结构开发IIR滤波器的格型网络结构。设一个全极点滤波器的系统函数由式(5.5.2)给定:

H(z)=11+∑Ni=1a(i)Nz-i=1A(z)(5.5.2)



其中,a(i)N为分母多项式中z-i项的系数。与式(5.5.1)比较可知,H(z)=1A(z)是FIR滤波器系统B(z)=A(z)的逆系统,所以可以按照系统求逆准则得到H(z)=1A(z)的格型网络结构如图5.5.4所示。



图5.5.4全极点(IIR)格型滤波器网络结构


比较图5.5.1和图5.5.4,可总结格型网络结构系统求逆准则如下: 

(1) 将输入至输出的无延迟通路全部反向,并把该通路的常数值支路增益改变为原来增益值的倒数,图5.5.1中上面e0-e1-e2-…-eN的支路,增益全为1; 


(2) 再把指向上面这条新通路的各节点的其他支路的增益(指向上面节点的各支路增益: 反射系数)改变为原来增益值的负值; 

(3) 将输入和输出交换位置;

(4) 按照常规,将输入画在网络左端,输出画在网络右端,最终所得的格型网络(如图5.5.4所示)的系统函数正好就是原来FIR滤波器系统函数的倒数。

5.6MATLAB实现

本节给出的MATLAB实现参考程序,均在MATLAB 7.4.0环境下调试通过。

(1) 例5.4.1参考的MATLAB实现(即FIR直接型、线性相位型和级联型的单位采样脉冲响应程序)

b=[1,0,0,0,16+1/16,0,0,0,1];

[b0,B,A]=dir2cas(b,1); %级联型结构系数

N=24;n=1:N+1;

delta=impseq(0,0,N); h1=filter(b,1,delta); h2=casfiltr(b0,B,A,delta);

fprintf('\n***由图可知这三种形式的单位冲激响应是一致的\n')

subplot(3,1,1); stem(n,h1);gtext('Response to \delta(n), Direct Form')

subplot(3,1,2); stem(n,h2);gtext(' Response to \delta(n),Linear Phase')

axis([0,25,0,20])

subplot(3,1,3); stem(n,h2);gtext('Response to \delta(n), Cascade Form')

axis([0,25,0,20])

(2) 例5.4.2参考的MATLAB实现(即画出梳状滤波器的零点、极点及幅频特性曲线)

% zero/pole plotH(z)=sum[b(i)z^(-i)]/sum[1-a(i)z^(-i)], let N=8

%clear all;close all;

N=8;      %2*pi/N=0.25*pi;

r1=0.8;

r2=0.95;

b=[1 zeros(1,N-1) -1];

a0=1;

a1=[1 zeros(1,N-1) -r1^N];

a2=[1 zeros(1,N-1) -r2^N];

subplot(2,1,1)

zplane(b,a0)

st1=sprintf('zero-pole locations for FIR comb filter,  N=%d',N);

title(st1);

 [h,w]=freqz(b,a0);

%[z,p]=tf2zp(b,a)

subplot(2,1,2);

plot(w/pi,abs(h));

axis([0,1, min(abs(h))-0.05, max(abs(h))+0.05]); 

ylabel('|H(e^j^w)|');

xlabel('frequency unit  in /pi rad');

title('Frequency response to magnitude of FIR')



figure

[h1,w]=freqz(b,a1);

subplot (2,1,1); zplane(b,a1)

st1=sprintf('zero-pole locations for IIR comb filter,  N=%d',N);

st2=sprintf(',     a=%d',r1);

st=[st1 st2];

title(st);

subplot(2,1,2);plot(w/pi,abs(h1));

axis([0,1, min(abs(h1))-0.05, max(abs(h1))+0.05]); 

ylabel('|H(e^j^w)|');

%subplot(2,1,2);plot(w/pi,angle(h));

%axis([0,1, min(angle(h))-0.05, max(angle(h))+0.05]); 

xlabel('frequency unit  in /pi rad');

%ylabel('Phase of H(jw)    in rad');

title('Frequency response to magnitude of IIR');



figure

[h2,w]=freqz(b,a2);

subplot (2,1,1); zplane(b,a2)

st1=sprintf('zero-pole locations for IIR comb filter,  N=%d',N);

st2=sprintf(',     a=%d',r2);

st=[st1 st2];

title(st);

subplot(2,1,2);plot(w/pi,abs(h2));

axis([0,1, min(abs(h2))-0.05, max(abs(h2))+0.05]); 

ylabel('|H(e^j^w)|');

%subplot(2,1,2);plot(w/pi,angle(h));

%axis([0,1, min(angle(h))-0.05, max(angle(h))+0.05]); 

xlabel('frequency unit  in /pi rad');

%ylabel('Phase of H(jw)    in rad');

title('Frequency response to magnitude of IIR');


5.7习题

5.1设数字系统用下面时域I/O差分方程描述: 

y(n)-34y(n-1)+18y(n-2)=x(n)+13x(n-1)

试分别画出系统的直接型、级联型和并联型网络结构。差分方程中x(n)和y(n)分别表示系统的输入和输出信号。

5.2设数字滤波器的时域I/O差分方程为

y(n)=(a+b)y(n-1)-aby(n-2)+x(n-2)+(a+b)x(n-1)+abx(n)

式中|a|<1,|b|<1,试画出系统的直接型、级联型结构,x(n)和y(n)分别表示系统的输入和输出信号。

5.3用直接型Ⅰ及典型型结构实现以下系统函数描述的数字系统:

H(z)=3+4.2z-1+0.8z-22+0.6z-1-0.4z-2

5.4用级联型结构实现以下系统函数描述的数字系统,试问一共能构成几种级联型网络?

H(z)=4(z+1)(z2-1.4z+1)(z-0.5)(z2+0.9z+0.8)

5.5图5.7.1中画出了4个系统,试用各子系统的单位采样脉冲响应分别表示各总系统的单位采样脉冲响应h(n)。并求其总系统函数H(z)。

5.6写出图5.7.2中各信号流图表示的数字系统的系统函数及时域I/O差分方程。

5.7给出以下系统函数所描述数字系统的并联型实现:

H(z)=5.2+1.58z-1+1.41z-2-1.6z-3(1-0.5z-1)(1+0.9z-1+0.8z-2)

5.8已知FIR滤波器的单位采样脉冲响应为

h(n)=δ(n)+0.3δ(n-1)+0.72δ(n-2)+0.11δ(n-3)+0.12δ(n-4)

试画出其级联型结构实现。

5.9写出图5.7.3中各信号流图描述的数字系统的系统函数。

5.10已知滤波器的单位采样脉冲响应为h(n)=0.9nR5(n),求出该滤波器的系统函数,并画出其直接型网络结构。

5.11已知滤波器的单位采样脉冲响应为h(n)=δ(n)-δ(n-1)+δ(n-4),试用频率采样结构来实现该滤波器。设频域采样点数N=5,要求画出频率采样型网络结构图,写出滤波器的参数计算公式。




图5.7.1习题5.5用图




图5.7.2习题5.6用图








图5.7.3习题5.9用图




5.12设某FIR滤波器的系统函数为

H(z)=15(1+3z-1+5z-2+3z-3+z-4)

试画出此滤波器的线性相位型网络结构。

5.13设滤波器的时域I/O差分方程为

y(n)=x(n)+x(n-1)+13y(n-1)+14y(n-2)

(1) 试用直接型Ⅰ,典型型及一阶节的级联型、一阶节的并联型网络结构实现此差分方程描述的滤波器。

(2) 求系统的频率响应(幅度响应及相位响应)。

(3) 设采样频率为10kHz,输入正弦波幅度为5,频率为1kHz,试求该滤波器的稳态输出。

5.14分别给定3个数字滤波器的系统函数为

H1(z)=1-0.6z-1-1.4145z-2+0.864z-3
H2(z)=1-0.98z-1+0.9z-2-0.898z-3
H3(z)=H1(z)/H2(z)

分别画出它们的直接型网络结构。

思政小课堂


思政主题:  基于FIR和IIR的位置姿态系统的低延时降噪方法

思政方向:  课程知识在我国先进技术领域的应用、培养学生科研精神与爱国情操




通过前几章专业知识的学习,同学们已经掌握了DTFT、z变换、DFT、FFT。在学习的过程中可能会疑惑课程中的专业知识以及复杂的公式究竟是如何应用到我们的国家发展和日常生活中的。在课程的绪论中我们曾经学习过,数字信号处理广泛地应用于数字通信、雷达、遥感、图像处理、航空航天等领域,随着2022年神舟十二号载人航天飞船的成功发射,以及航天站的投入使用,航空航天也成为近期大众关注的热门话题。数字信号处理中提及的FIR和IIR滤波器的设计,在解决北斗卫星的发射和执行通信任务过程中,信号实时传输过程遇到的延时和精确度的问题发挥着重要作用。

我国的航天事业的建设起始于1956年,从1970年4月24日成功发射的第一颗人造地球卫星东方红一号,到近十年自主创新技术的显著增强,目前我国的航天成就已包括神舟、嫦娥、空间站、天问、长征系列运载火箭、北斗卫星导航系统、通信卫星、气象卫星、地球资源卫星、反卫星武器等多个领域(图E.1)。其中北斗卫星导航系统是中国自主研制的全球卫星导航系统,可在全球范围内全天候、全天时为各类用户提供高精度、高可靠定位、导航、授时服务。卫星在发射和运行的过程中,对其位置和姿态信息的实时获取,是检测卫星是否处于正常工作状态的重要手段。地面控制中心将利用位置姿态测量系统(position and orientation system,POS)对卫星进行监测。它是一种航空遥感专用的惯性/卫星组合测量系统,主要用于测量与其固连的航空遥感载荷相位中心的位置、速度和姿态信息,并对成像进行运动补偿。



图E.1我国航天卫星发展历程



受距离和带宽等因素的限制,实际中对卫星信息的获取无法实现极小间隔的连续采样。同时,获取到的信息通常都带有一定的因为机械抖动偏频而产生的噪声,这些都会影响获取到的信号的精度,为了确保信息的精度,往往需要进行精密的滤波降噪。

FIR 滤波器能够实现严格的线性相位,有利于POS 三轴惯性信息测量的同步和POS 时延补偿。但随着航空遥感载荷成像精度的不断提高,要求POS 能够测量更宽频率范围内的载荷运动,导致FIR 滤波器阶数提高,运算量增加,使POS 无法同时完成降噪解调、数据平滑和补偿等所有信号预处理功能,而且滤波器阶数的提高导致滤波时延大幅提高,限制了航空遥感载荷成像的实时性。
IIR滤波器则可以弥补FIR滤波器的这些缺陷。它具有噪声抑制效果好、阶数低和运算量小的优点,但常规的IIR 滤波器相位非线性严重,不利于POS惯性信息测量同步和POS 时延补偿,而且常规IIR滤波器时延也很大,不利于POS 实时测量。

为了解决这一问题,科研人员提出一种线性相位IIR滤波降噪解调方法(图E.2)。首先在时域内通过卡尔曼滤波参数估计使IIR 低通滤波器相位线性化,解决常规IIR滤波器相位非线性严重的问题;  然后在频域内建立IIR低通滤波器目标函数,通过拟牛顿法进行优化求解,以降低时延,获得最终在线使用的IIR 低通滤波器系数;  最后设计低阶线性相位IIR 带阻滤波器,降低线性相位IIR 低通滤器的阶数和运算量,获得最终在线使用的IIR 带阻滤波器系数。




实验结果表明,从幅频特性来看,同等阶数的IIR滤波器可以逼近更高阶的FIR 滤波器的幅频特性以获得更好的噪声抑制效果和通频带增益稳定性。22阶的IIR的幅频性能可以很好地等效为90阶的FIR滤波器。而且从相频特性来看,同等阶数的IIR逼近更低阶的线性相位FIR滤波器的相频特性,22阶的IIR的相频性能可以很好地等效为38阶的FIR滤波器。这样信息获取过程中POS的时延就可以大幅度地降低。



图E.2利用线性IIR滤波降噪解决方案