第5章 CHAPTER 5 数字滤波器基本结构 及状态变量分析法 数字信号处理的目的之一,是设计某种设备或建立某种算法分析处理序列,使序列具有某些确定的性质,这种设备或是算法结构就是数字滤波器。与FFT一样,网络结构是数字滤波器设计中的一个非常重要的内容,也是数字信号处理的重要内容,因为数字滤波器的稳定性、运算速度以及系统的成本和体积等许多重要性能都取决于其网络结构。本章的主要内容就是理解数字滤波器结构的表示方法; 无限长脉冲响应(IIR)基本网络结构,掌握IIR滤波器的直接型、级联型和并联型结构; 有限长脉冲响应(FIR)基本网络结构,掌握FIR滤波器的直接型、级联型、频率采样型结构; 滤波器网络结构的状态变量分析方法。 视频讲解 5.1引言 1. 描述数字滤波器的方法 一般,数字滤波器可以采用下面4种方法描述。 (1) 系统单位脉冲响应h(n)(系统的时域特性)。 (2) 系统频率响应(变换域特性): H(ejω)=∑∞n=-∞h(n)·e-jωn,Y(ejω)=X(ejω)·H(ejω)(51) (3) 系统函数H(z)(变换域特性): H(z)=∑∞n=-∞h(n)·z-n(52) (4) 差分方程(输入输出序列间的关系): y(n)=∑Mi=0bix(n-i)+∑Ni=1aiy(n-i)(53) 2. 实现方法 硬件实现: 根据描述数字滤波器的数学模型或信号流图,用数字硬件设计成一台专门的设备,构成专用的信号处理机。 软件实现: 直接利用通用计算机,将所需要的运算编成程序让计算机执行。 为了用计算机或专用硬件完成对输入信号的处理(运算),必须把式(52)或者式(53)变换成一种算法,按照这种算法对输入信号进行运算。其实,式(53)就是对输入信号的一种直接算法,如果已知输入信号x(n)以及ai、bi和n时刻以前的y(n-i),则可以递推出y(n)值。但给定一个差分方程,不同的算法有多种,例如: H(z)=11-3z-1+2z-2=21-2z-1-11-z-1=11-2z-1·11-z-1(54) 不同的算法将直接影响系统运算误差、运算速度以及系统的复杂程度和成本等,因此研究实现信号处理的算法是一个很重要的问题。我们用网络结构表示具体的算法,因此网络结构实际表示的是一种运算结构。本章是第6、7章数字滤波器设计的必要基础。在介绍数字系统的基本网络结构之前,先介绍网络结构的表示方法。 5.2用信号流图表示网络结构 观察式(53)可知,数字信号处理中有三种基本算法,即乘法、加法和单位延迟。三种基本运算框图及其流图如图51所示。 图51三种基本运算的流图表示 例如,二阶数字滤波器为 y(n)=a1y(n-1)+a2y(n-2)+b0x(n) 其方框图及信号流图结构如图52所示。 图52二阶网络方框图及信号流图 图52(b)中,1、2、3、4、5为网络节点。x(n): 输入节点或源节点(没有输入支路),y(n): 输出节点或阱节点(没有输出支路)。节点之间用有向支路连接,每个节点可以有几条输入支路和几条输出支路,任意节点的节点值等于它所有输入支路的信号和。而输入支路的信号值等于这一支路起点处节点信号值乘以支路上的传输系数。如果支路上不标传输系数值,则认为其传输系数为1。节点2处的节点值可以用w2(n)表示,其他节点类似,则有如下方程组: w2(n)=y(n)=w1(n) w3(n)=w2(n-1)=y(n-1) w4(n)=w3(n-1)=y(n-2) w5(n)=a1w3+a2w4=a1y(n-1)+a2y(n-2) w1(n)=b0x(n)+w5(n)=b0x(n)+a1y(n-1)+a2y(n-2) 不同的信号流图代表不同的运算方法,而对于同一个系统函数,可以有很多种信号流图与其对应。从基本运算考虑,满足以下三个条件,称为基本信号流图(Primitive Signal Flow Graghs)。 (1) 信号流图中所有支路都是基本的,即支路增益是常数或者是z-1。 (2) 流图环路中必须存在延时支路。 (3) 节点和支路的数目是有限的。 从该例中可以看出,用信号流图表示系统的运算情况(网络结构)是比较简明的。以下我们均用信号流图表示网络结构。一般将网络结构分成两类,一类称为有限长单位脉冲响应(Finite Impulse Response,FIR)网络,另一类称为无限长单位脉冲响应(Infinite Impulse Response,IIR)网络。FIR网络中一般不存在输出对输入的反馈支路,因此差分方程为 y(n)=∑Mi=0bix(n-i)(55) 其单位脉冲响应h(n)是有限长的,按照式(55),则 h(n)=bn,0≤n≤M 0,其他 系统函数为 H(z)=∑Mn=0bnz-n IIR网络结构存在输出对输入的反馈支路,也就是说,信号流图中存在反馈环路。这类网络的单位脉冲响应是无限长的: y(n)=∑Mi=0bix(n-i)+∑Ni=1aiy(n-i) 系统函数为 H(z)=∑Mi=0biz-i1-∑Ni=1aiz-i(56) 例如,一个简单的一阶IIR网络的差分方程为 y(n)=ay(n-1)+x(n) 其单位脉冲响应h(n)=anu(n)。 综上所述,这两类不同的网络结构各有不同的特点,下面分类叙述其网络结构。 5.3无限长单位冲激响应滤波器的基本结构 视频讲解 无限长单位冲激响应(IIR)数字滤波器的结构特点: ①系统的单位冲激响应h(n)是无限长的; ②系统函数H(z)在有限Z平面上有极点存在; ③结构上存在输出到输入的反馈,即结构是递归的。其基本结构有直接Ⅰ型、直接Ⅱ型、级联型和并联型。 5.3.1直接型 1. 直接Ⅰ型 系统输入输出关系的N阶差分方程为 y(n)=∑Mi=0bix(n-i)+∑Ni=1aiy(n-i) 对应的系统函数为 H(z)=N(z)·1D(z)=∑Mi=0biz-i·11-∑Ni=1aiz-i=H1(z)·H2(z)(57) 直接Ⅰ型IIR滤波器结构如图53所示。 图53直接Ⅰ型IIR滤波器结构 直接Ⅰ型网络结构的特点: (1) ∑Mi=0bix(n-i)表示将输入及延时后的输入组成M阶的延时网络,即横向延时网络,实现零点。 (2) ∑Ni=1aiy(n-i)表示输出及其延时组成N阶延时网络,实现极点。 (3) 总的网络由上面两个网络级联而成。 (4) 直接Ⅰ型需要N+M级延时单元。 2. 直接Ⅱ型 按照差分方程可以直接画出网络结构如图53所示。图中第一部分系统函数用H1(z)表示,第二部分用H2(z)表示,那么H(z)=H1(z)·Hz(z),当然也可以写成H(z)=H2(z)·H1(z),相当于将图53中两部分流图交换位置,如图54(a)所示。该图中节点变量w1=w2,因此前后两部分的延时支路可以合并,形成如图54(b)所示的网络结构流图,图54(b)所示的这类流图称为IIR直接Ⅱ型网络结构。 图54直接Ⅰ型向直接Ⅱ型IIR滤波器结构的转变 直接Ⅱ型网络结构的优点: 结构简单、清晰,延时支路比直接Ⅰ型减少1/2。 直接Ⅱ型网络结构的缺点: (1) ai、bi常数对滤波器的性能控制作用不明显; (2) 零、极点关系不明显,调整困难; (3) 系数量化效应敏感度高。 【例51】IIR数字滤波器的系统函数H(z)为 H(z)=8-4z-1+11z-2-2z-31-54z-1+34z-2-18z-3 画出该滤波器的直接Ⅱ型结构。 【解】由H(z)写出差分方程如下: y(n)=54y(n-1)-34y(n-2)+18y(n-3)+8x(n) -4x(n-1)+11x(n-2)-2x(n-3) 该滤波器的直接Ⅱ型结构如图55所示。 图55例51直接Ⅱ型结构 【例52】一个LTI系统单位脉冲响应为h(n)=δ(n)+212nu(n)-14nu(n),求解如下问题: (1) 确定系统函数(包括收敛域)、系统频率响应。 (2) 给出系统零极点并确定系统的稳定性。 (3) 求解系统的差分方程。 (4) 绘制出在使用最少存储单元及单位延时数情况下的系统结构图。 【解】(1)根据Z(anu(n))=11-az-1可以直接给出系统函数 H(z)=1+21-12z-1-11-14z-1=2-34z-1+18z-21-34z-1+18z-2,|z|>12 频率响应为 H(jω)=H(z)|z=ejω=2-34e-jω+18e-j2ω1-34e-jω+18e-j2ω (2) 系统函数可以写成 H(z)=2-34z-1+18z-21-34z-1+18z-2=2z2-34z+18z2-34z+18 系统零点为 z1,2=3±7j16 系统极点为 p1=12,p2=14 由于两个极点都在单位圆中,所以系统是稳定的。 (3) 由于 H(z)=2-34z-1+18z-21-34z-1+18z-2 因此 Y(z)1-34z-1+18z-2=X(z)2-34z-1+18z-2 这里,X(z)和Y(z)分别是输入输出信号的Z变换。 图56例52直接Ⅱ型结构 因此,差分方程可表示成 y(n)-34y(n-1)+18y(n-2) =2x(n)-34x(n-1)+18x(n-2) (4) 题目要求最少存储单元及单位延时数,系统结构应该是直接Ⅱ型结构(图56)。 5.3.2级联型 对于系统函数,有 H(z)=Y(z)X(z)=∑Mi=0biz-i1-∑Ni=1aiz-i 分子分母均为多项式,且多项式的系数一般为实数,现将分子分母多项式分别进行因式分解,得到 H(z)=A∏Mr=1(1-crz-1)∏Nr=1(1-drz-1)(58) 式中,A为常数; cr,dr分别表示系统的零点、极点,为实数或共轭成对的复数。将共轭成对的零点(极点)放在一起,形成一个二阶多项式,系数仍为实数,将分子、分母均为实数的二阶多项式放在一起,形成一个二阶网络。二阶网络系统函数为 Hi(z)=β0i+β1iz-1+β2iz-21-α1iz-1-α2iz-2(59) 式中,β0i、β1i、β2i、α1i和α2i均为实数。这样,H(z)就分解成一阶或二阶数字网络的级联形式,如 H(z)=H1(z)H2(z)…Hk(z) 级联型示意图见图57。 图57级联型IIR滤波器结构示意图 图58式(59)对应的直接Ⅱ型结构 级联型结构不是唯一的,式中Hi(z)表示一个一阶或二阶的数字网络的系统函数,每个Hi(z)的网络结构均采用前面介绍的直接型网络结构表示,式(59)对应的直接Ⅱ型结构见图58。 级联型IIR滤波器结构的优点: (1) 所需存储器最少,系统结构组成灵活。 (2) 每个一阶网络决定一个零点、一个极点,每个二阶网络决定一对零点、一对极点。调整一阶网络和二阶网络系数可以改变零极点位置,所以零极点调整方便,便于调整频响。 级联型IIR滤波器结构的缺点: (1) 存在误差积累,级联结构中前面网络误差会输出到后面的网络。 (2) 零、极点配合关系着网络最优化的问题,而最佳配合关系不易确定。 【例53】已知IIR数字滤波器的系统函数,画出该滤波器的级联型结构。 H(z)=8-4z-1+11z-2-2z-31-1.25z-1+0.75z-2-0.125z-3 【解】为了减少单位延迟的数量,将一阶的分子、分母多项式组成一个一阶网络,二阶的分子、分母多项式组成一个二阶网络。将H(z)的分子、分母进行因式分解,得 H(z)=(2-0.379z-1)(4-1.24z-1+5.264z-2)(1-0.25z-1)(1-z-1+0.5z-2) 则H(z)的级联型结构见图59。 图59例53级联型结构图 5.3.3并联型 图510并联型结构图 将H(z)展成部分分式形式得到IIR并联型结构,如图510所示,即 H(z)=H1(z)+H2(z)+…+Hk(z) 式中,Hi(z)通常为一阶网络或二阶网络。二阶网络的系统函数一般为 Hi(z)=β0i+β1iz-11-α1iz-1-α2iz-2(510) 式中,β0i,β1i,α1i,α2i都是实数。如果α2i=0,则构成一阶网络。 其输出Y(z)表示为 Y(z)=H1(z)X(z)+H2(z)X(z)+…+Hk(z)X(z) 表明: 将x(n)送入每个二阶(或一阶)网络后,将所有输出相加得到输出y(n)。 并联型IIR滤波器结构的优点: (1) 所需存储器最小。 (2) 无误差积累,各级误差互不影响,仅极点调整方便。所以,在要求准确传输极点的场合,宜采用这种结构。 并联型IIR滤波器结构的缺点: 零点调整不方便,当H(z)有多阶极点时,部分分式不易展开。 【例54】若系统函数 H(z)=8-4z-1+11z-2-2z-31-1.25z-1+0.75z-2-0.125z-3 求H(z)的并联型结构。 【解】对H(z)展开成部分分式 H(z)=16+81-0.25z-1+-16+20z-11-z-1+0.5z-2 将每一部分用直接型结构实现,其并联型网络结构如图511所示。 图511例54并联型结构图 视频讲解 5.4有限长单位冲激响应滤波器的基本结构 有限长单位冲激响应(FIR)数字滤波器的结构特点: ①系统的单位冲激响应h(n)是有限长序列; ②系统函数|H(z)|在|z|>0处收敛,有限Z平面上只有零点,极点全部在z=0处(即FIR一定为稳定系统); ③结构上主要是非递归结构,没有输出到输入反馈。但有些结构中(如频率抽样结构)也包含有反馈的递归部分。 5.4.1直接型 设单位脉冲响应h(n)长度为N,其系统函数H(z)为 H(z)=∑N-1n=0h(n)z-n 差分方程表示为 y(n)=∑N-1m=0h(m)x(n-m)=h(0)x(n)+h(1)x(n-1)+…(511) 式(511)是线性移不变系统的卷积和公式,是x(n)延时链的横向结构,称为横截型结构或卷积型结构,也可称为直接型结构。直接按H(z)或者差分方程画出没有反馈支路的结构图(图512)。 图512FIR直接型网络结构 5.4.2级联型 将H(z)进行因式分解,并将共轭成对的零点放在一起,形成一个系数为实数的二阶网络,形式如下: H(z)=∑N-1n=0h(n)z-n=∏Mi=1(a0i+a1iz-1+a2iz-2)(512) 这样级联型网络结构就是由一阶或二阶因子构成的级联结构,如图513所示,其中每一个因式都用直接型实现。 Hz=H1zH2z…HMz 图513级联型FIR滤波器结构 级联型FIR滤波器结构的特点: (1) 由于这种结构所需的系数比直接型多,所需乘法运算也比直接型多,很少用。 (2) 由于这种结构的每一节控制一对零点,因而只能在需要控制传输零点时用。 【例55】设FIR网络系统函数H(z)如下: H(z)=0.96+2.0z-1+2.8z-2+1.5z-3 画出H(z)的直接型结构和级联型结构。 【解】将H(z)进行因式分解,得到 H(z)=(0.6+0.5z-1)(1.6+2z-1+3z-2) 其级联型结构和直接型结构如图514所示。 图514例55级联型和直接型结构 5.4.3频率采样型 由频域采样定理可知,对有限长序列h(n)的Z变换H(z)在单位圆上做N点的等间隔采样,N个频率采样值的离散傅里叶反变换所对应的时域信号hN(n)是原序列h(n)以采样点数N为周期进行周期延拓的结果,当N大于或等于原序列h(n)长度M时,hN(n)=h(n),不会发生信号失真,此时H(z)可以用频域采样序列H(k)内插得到,各个过程间关系如图515所示,H(k)、H(z)关系用内插公式(513)表示。 图515频率采样与系统函数、离散傅里叶变换间关系 H(z)=(1-z-N)1N∑N-1k=0H(k)1-W-kNz-1(513) 式中 H(k)=H(z)z=ej2πNk,k=0, 1, 2, …,N-1 1. 频率采样型滤波器结构 由式(513)得到FIR滤波器的另一种结构——频率采样型结构,它由两部分级联而成。H(z)也可以重写为 H(z)=1NHc(z)∑N-1k=0Hk(z)(514) 式中 Hc(z)=1-z-N(515) Hk(z)=H(k)1-W-kNz-1(516) 2. 梳状滤波器 Hc(z)是一个由N阶延时单元组成的梳状滤波器,信号流图及幅频特性如图516所示。它在单位圆上有N个等间隔的零点: zi=ej2πNi=W-iN,i=0, 1,2, …,N-1 将z=ejω代入式(515)中,可得到幅频特性表达式为 Hc(ejω)=1-e-jNω=1-cosNω+jsinNω Hc(ejω)=(1-cosNω)2+sin2Nω =2(1-cosNω)=2sinNω2 图516梳状滤波器信号流图及幅频特性 3. 谐振器 图517谐振器信号流图 式(516)表示的网络称为谐振器,谐振器的极点为 z=W-kN=ej2πNk=ejω 谐振器信号流图(IIR)如图517所示。 此一阶网络在频率ω=2πNk处响应为无穷大,即Hk(z)→∞。故等效为谐振频率为2πNk的无损耗谐振器。 4. 谐振矩 ∑N-1k=0Hk(z)=∑N-1k=0H(k)1-W-kNz-1(517) 它是由N个谐振器并联而成的。这个谐振矩的极点正好与梳状滤波器的一个零点(i=k)相抵消,从而使这个频率(ω=2πk/N)上的频率响应等于H(k) Hc(z)·Hk(z)=zk-ej2πkNH(k)(zk-ej2πkN)=H(k)(518) 把这两部分级联起来就可以构成FIR滤波器的频率采样型结构,如图518所示。 图518频率采样型FIR滤波器结构 频率采样型FIR滤波器结构的优点: (1) 频响特性调整方便,在频率采样点ωk,H(ejωk)=H(k),只要调整H(k)(即一阶网络Hk(z)中乘法器的系数H(k)),可有效地调整频响特性。 (2) 易于标准化、模块化: 只要h(n)长度N相同,对于任何频响形状,其梳状滤波器的部分完全相同,N个一阶网络部分的结构也完全相同,只是各支路增益H(k)不同。这样,相同部分便于标准化、模块化。 频率采样型FIR滤波器结构的缺点: (1) 系统稳定是靠位于单位圆上的N个零极点对消来保证的,由于寄存器的长度有限,有限字长效应可能使零极点不能完全抵消,影响系统的稳定性。 (2) 由于H(k)和W-kN一般为复数,要求乘法器完成复数乘法运算,这对硬件实现是不方便的。 为了克服以上缺点,采取下面修正措施。 将单位圆上的零极点向单位圆内收缩一点,收缩到半径r<1 且r≈1,这样,以z/r代替式(513)H(z)表示式中z,得 H(z)=(1-rNz-N)1N∑N-1k=0Hr(k)1-rW-kNz-1(519) 此时,Hrk≈Hk,零极点均为rej2πNk,如果由于量化原因零、极点不能相互抵消时,极点也都在单位圆内,系统仍然是稳定的。若h(n)是实序列,根据其DFT变换对称性,H(k)=H*(N-k),旋转因子W-kN=W-N-kN,将Hk(z)和HN-k(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-rWkNz-1 =a0k+a1kz-11-2rcos2πNkz-1+r2z-2(520) 其中 a0k=2Re[H(k)] a1k=-2Re[rH(k)WkN],k=1,2,…,N2-1(521) 当N为偶数时,H(0)和H(N/2)为实数,H(z)可用式(522)表示,此时,修正的Z平面采样点如图519所示,FIR频率采样修正型子网络结构如图520所示,修正的频率采样结构如图521所示。 H(z)=(1-rNz-N)1NH(0)1-rz-1+HN21+rz-1+∑N2-1k=1a0k+a1kz-11-2cos2πNkz-1+r2z-2(522) 图519修正的Z平面采样点图 图520FIR频率采样修正型子网络结构 图521频率采样修正结构(N为偶数) 当N为奇数时,只有H(0)为实数,H(z)可用式(523)表示,此时,修正的频率采样结构如图522所示。 H(z)=(1-rNz-N)NH(0)1-rz-1+∑N-12k=1a0k+a1kz-11-2rcos2πkNz-1+r2z-2(523) 图522频率采样修正结构(N为奇数) 可见,当采样点数N很大时,网络结构很复杂,需要的乘法器和延时单元很多,对于窄带滤波器,大部分采样值为零,使二级网络个数大大减少,所以频率采样结构适合窄带滤波器的设计。 5.5状态变量分析法 本书其他章节都把离散时间系统看作一个具有输入—输出两个端口的黑盒,而不关心系统内部状态。这种两端口系统的时域特性用冲激响应或差分方程来表征,而频域特性则用频率响应或系统函数来表征。当然,这些特性存在着一定的数学关系。随着系统理论的发展,人们不再满足只对系统输出变化的了解,对系统内部某些变量的变化规律也感兴趣。此外,实际应用中还会遇到多输入多输出的情况,而且将系统看成是时不变系统也往往不符合客观实际。这样,仅仅用两端口分析法就不够了,所以要引入另一种系统分析方法,即状态变量分析法。利用矩阵线性代数等数学工具对系统进行分析研究。 “状态”指系统内一组变量,它包含了系统全部过去的信息,由这一组变量和现在与将来的输入,可求出系统现在和将来的全部输出; 一般来说,描述一个系统有输入输出法(差分方程、系统函数、单位脉冲响应)和状态变量法。状态变量分析法具有如下优越性: (1) 更好地了解系统的内部结构。 (2) 既能处理线性系统,又能处理非线性系统。 (3) 既能处理单输入输出变量的系统,又能处理多输入输出的系统。 5.5.1由信号流图建立状态方程 状态变量分析法有两个基本方程: 状态方程和输出方程。状态方程把系统内部一些称为状态变量的节点变量和输入联系起来; 输出方程则把输出信号和那些状态变量联系起来。具体规则为 (1) 状态变量选在基本信号流图中单位延迟时支路输出节点处。 (2) 围绕加法器写出状态方程和输出方程。 图523是二阶网络基本信号流图,信号流图有两个延时支路,因此,建立两个状态变量w1n和w2n。 图523二阶网络基本信号流图 图523流图中输入x(n)、输出y(n)与状态变量之间的关系为 w′2=w2(n+1) w2(n+1)=-a2w1(n)-a1w2(n)+x(n) w1(n+1)=w2(n) y(n)=b2w1(n)+b1w2(n)+b0w′2=(b2-a2b0)w1(n)+(b1-a1b0)w2(n)+b0x(n) 将以上w1(n+1)、w2(n+1)和y(n)写成矩阵形式: w1(n+1) w2(n+1)=0,1 -a2,-a1w1(n) w2(n)+0 1x(n) y(n)=[b2-a2b0b1-a1b0][w1(n)w2(n)]T+b0x(n) 再分析一个信号流图,信号流图中两个延时支路输出节点定为状态变量w1n和w2n,如图524所示。 图524一般二阶网络基本信号流图 这个一般的二阶网络基本信号流图,按照信号流图写出以下方程组,即 w′1(n)=w1(n+1) w1(n+1)=a11w1(n)+a12w2(n)+b1x(n) w′2(n)=w2(n+1) w2(n+1)=a21w1(n)+a22w2(n)+b2x(n) y(n)=c1w1(n)+c2w2(n)+dx(n) 将以上写成矩阵形式,即 w1(n+1) w2(n+1)=a11,a12 a21,a22w1(n) w2(n)+b1 b2x(n) y(n)=[c1 c2][ω1(n) ω2(n)]T+dx(n) 用矩阵符号表示为 状态方程:W(n+1)=AW(n)+Bx(n) 输出方程:Y(n)=CW(n)+Dx(n) 其中 A=a11a12 a21a22,B=b1b2T C=c1c2,D=d 若系统中有N个单位延时支路,M个输入信号: x1(n),x2(n),…,xM(n),L个输出信号y1(n),y2(n),…,yL(n),状态变量分析法信号流图如图525所示,用矩阵符号表示为 状态方程: W(n+1)=AW(n)+Bx(n)(524) 输出方程: Y(n)=CW(n)+Dx(n)(525) 图525状态变量分析法 其中 A=a11a12 a21a22,B=[b1b2], C=[c1c2],D=[d](526) 【例56】建立图526的状态方程和输出方程。 图526例56网络信号流图 【解】(1) 按顺序在z-1支路输出端建立状态变量wi(n),z-1支路的输入端为wi(n+1)。 (2) 列出所有节点变量方程,找出状态变量wi(n+1)与wi(n)和输入x(n)之间的关系,用矩阵方程表示。 (3) 找出输出信号和状态变量wi(n)以及输入信号关系用矩阵方程表示。 w1(n+1)=a1w1(n)+a2w2(n)+x(n) w2(n+1)=w1(n) 将上式写成矩阵方程,即 w1(n+1) w2(n+1)=a1a2 10w1(n) w2(n)+1 0x(n) 输出信号y(n)的方程推导如下: y(n)=b0w1(n+1)+b1w1(n)+b2w2(n) 将w1(n+1)代入上式,则 y(n)=a1b0w1(n)+b0a2w2(n)+b0x(n)+b1w1(n)+b2w2(n) =(a1b0+b1)w1(n)+(b0a2+b2)w2(n)+b0x(n) y(n)=[a1b0+b1,a2b0+b2]w1(n) w2(n)+b0x(n) 【例57】IIR直接型系统流图及其参数如图527所示。给出系统函数,并写出状态方程和输出方程。 图527例57网络信号流图 【解】根据图527得系统函数为 H(z)=8-4z-1+11z-2-2z-31-1.25z-1+0.75z-2-0.125z-3 本题与例56相似,只是阶数不同,参考例56可用直接写出状态方程和输出方程 w1(n+1) w2(n+1) w3(n+1)=-a1-a2-a3 100 010w1(n) w2(n) w3(n)+1 0 0x(n) =1.25-0.750.125 100 010w1(n) w2(n) w3(n)+1 0 0x(n) y(n)=[b1+a1b0b2+a2b0b3+a3b0][w1(n)w2(n)w3(n)]T+b0x(n) =[65-1][w1(n)w2(n)w3(n)]T+8x(n) 5.5.2由系统函数建立状态方程 由系统函数建立状态方程的步骤为 (1) 画出系统网络结构。 (2) 根据信号流图写出其状态方程和输出方程。 【例58】系统函数 H(z)=30z2-10z+90z3-6z2+11z-6 写出其状态方程和输出方程。 【解】根据部分分式展开法,得到 H(z)=5z-1+10z-2+15z-3 Y(z)=5z-1X(z)+10z-2X(z)+15z-3X(z) =5z-11-z-1X(z)+10z-11-2z-1X(z)+15z-11-3z-1X(z) 并联型结构图如图528所示。 图528例58网络信号流图 状态方程和输出方程为 w1(n+1)=w1(n)+5x(n) w2(n+1)=2w2(n)+10x(n) w3(n+1)=3w3(n)+15x(n) y(n)=w1(n)+w2(n)+w3(n) 写成矩阵形式为 w1(n+1) w2(n+1) w3(n+1)=100 020 003w1(n) w2(n) w3(n)+5 10 15x(n) y(n)=111w1(n) w2(n) w3(n) 【例59】已知FIR滤波网络系统函数H(z)为 H(z)=∑3i=0aiz-i 写出其状态方程和输出方程。 【解】画出直接型结构如图529所示。 图529例59网络信号流图 则状态方程为 w1(n+1) w2(n+1) w3(n+1)=000 100 010w1(n) w2(n) w3(n)+1 1 0x(n) 输出方程为 y(n)=[a1a2a3][w1(n)w2(n)w3(n)]T+a0x(n) 5.5.3由状态变量分析法转换到输入输出分析法 设系统的状态方程和输出方程为 W(n+1)=AW(n)+Bx(n)(527) y(n)=CW(n)+dx(n)(528) 将上面两式进行Z变换: zW(z)=AW(z)+BX(z) (529) Y(z)=CW(z)+dX(z)(530) 式中 W(z)=[W1(z),W2(z),…,Wm(z)]T Wi(z)=ZT[wi(n)] X(z)=ZT[x(n)] Y(z)=ZT[y(n)] 由式(529)得 [zI-A]W(z)=BX(z) W(z)=[zI-A]-1BX(z)(531) 由式(530)和式(531)得 Y(z)=C[zI-A]-1BX(z)+dX(z) H(z)=Y(z)X(z)=C[zI-A]-1B+d(532) 【例510】已知二阶网络的4个参数矩阵如下: A=01 -a2-a1,B=0 1 C=[b2-a2b0b1-a1b0],d=b0 求该网络的系统函数。 【解】 [zI-A]=z-1 a2z+a1 [zI-A]-1B=1z(z+a1)+a2z+a11 -a2z0 1 =1z2+a1z+a21 z 由式(533)可得 H(z)=C[zI-A]-1B+d =b0z2+b1z+b2z2+a1z+a2=b0+b1z-1+b2z-21+a1z-1+a2z-2 5.6MATLAB应用实例 【例511】求下列直接型系统函数的零、极点,并将它转换成二阶节形式。 H(z)=1-0.1z-1-0.3z-2-0.3z-3-0.2z-41+0.1z-1+0.2z-2+0.2z-3+0.5z-4 MATLAB代码如下: num=[1 -0.1 -0.3 -0.3 -0.2]; den=[1 0.1 0.2 0.2 0.5]; [z,p,k]=tf2zp(num,den); m=abs(p); disp('零点');disp(z); disp('极点');disp(p); disp('增益系数');disp(k); sos=zp2sos(z,p,k); disp('二阶节');disp(real(sos)); zplane(num,den) 输入到“num”和“den”的分别为分子和分母多项式的系数。计算求得零、极点增益系数和二阶节的系数。 零点: 0.9615-0.5730-0.1443+0.5850i-0.1443-0.5850i 极点: 0.5276+0.6997i0.5276-0.6997i-0.5776+0.5635i-0.5776-0.5635i 增益系数: 1 二阶节: 1.0000-0.3885-0.55091.00001.15520.6511 1.00000.28850.36301.0000-1.05520.7679 系统函数的二阶节形式为 Hz=1-0.3885z-1-0.5509z-21+0.2885z-1+0.3630z-2·1+1.1552z-1+0.6511z-21+1.0552z-1+0.7679z-2 系统零极点如图530所示。 图530系统函数的零极点图 【例512】周期序列发生器。 在用状态变量法分析数字系统时,N阶系统被视为具有N个抽头的延迟线结构。这种结构既可以用来构建梳状滤波器,也可以用作周期序列发生器。在作为周期序列发生器时,应周期地刷新延迟线各抽头状态,并且必须有反馈支路才能实现自激振荡。例如,若要产生图531所示的6点周期序列,则需要具备图531所示的延迟线。图中延迟线的6个状态变量,即w1(n),w2(n),…,w6(n)。由图531知x(i)=w(1),y(i)=x(i)。x(i)是系统内部的输入点,而系统外部并无信号加到这一点。状态变量w(1)~w(6)被初始化如图531所示。每次得到一个输出点y(i)后,都要将延迟线按图中所示方向刷新,即先使临时变量temp=w(1),然后依次使w(1)=w(2) ,w(2)=w(3),…,w(5)=w(6) ,w(6)=temp。运行下面代码,得到如图532所示的序列发生器输出波形。 图5316点周期序列 图532序列发生器输出波形 MATLAB代码如下: set(gcf, 'color', 'w') %设定图面背景色是白色 N=34;%序列总点数 w=[1,2,3,zeros(1,3)]; %初始化延迟线的状态 [1,2,3,0,0,0,0] x=[zeros(1,N)]; %初始化输入序列 x(1)=w(1); %初始化第一个输入样点 for i=1 : N y(i)=x(i); %得到第 i 个输出样点 temp=w(1); %暂存 for j=1 : 5 w(j)=w(j + 1); %延迟线状态刷新 end w(6)=temp %反馈 x(i+1)=w(1); %下一个输入样点 end H=stem(0:N-1,y(1:N));set(H, 'markersize', 3); grid; axis([0,N-1,0,4]); 【例513】从例题57中可以看出,对下列系统函数用状态分析法求出了的状态方程和输出方程,从而得到参数矩阵。这个过程是烦琐的。 H(z)=8-4z-1+11z-2-2z-31-1.25z-1+0.75z-2-0.125z-3 现在,根据系统函数可以利用转换函数tf2ss从系统函数的分子分母多项式系数组容易求出 系统的状态方程和输出方程的参数矩阵。 MATLAB代码如下: b = [8 -4 11 -2]; %系统函数的分子多项式系数组 a = [1 -1.25 0.75 -0.125]; %系统函数的分母多项式系数组 [A B C D] = tf2ss(b, a) %求出状态方程和输出方程的参数矩阵 程序运行结果为: A =1.2500-0.75001.2500 1.000000 01.00000 B =100 C =65-1 D =8 答案与例57结果是一致的。 【本章习题】 5.1填空题 (1) 无限长单位冲激响应滤波器的基本结构有直接Ⅰ型、直接Ⅱ型、和4种。 (2) 数字信号处理中有3种基本算法,即乘法、加法和。 (3) 状态变量分析法有2个基本方程: 和输出方程。 5.2选择题 (1) 下面信号流图(题5.2(1)图)表示的系统函数为()。 A. H(z)=1+12z-11-23z-1+14z-2B. H(z)=1-23z-1+14z-21+12z-1 C. H(z)=1+12z-11+23z-1-14z-2D. H(z)=1+23z-1-14z-21+12z-1 题5.2(1)图 (2) 下列各种滤波器的结构中()不是IIR滤波器的基本结构。 A. 直接型B. 级联型 C. 并联型D. 频率抽样型 (3) 下列结构中不属于FIR滤波器基本结构的是()。 A.横截型B. 级联型 C. 并联型D. 频率抽样型 5.3画出H(z)=(2-0.379z-1)(4-1.24z-1+5.264z-2)(1-0.25z-1)(1-z-1+0.5z-2)级联型网络结构。 5.4已知某三阶数字滤波器的系统函数为H(z)=3+53z-1+23z-21-13z-11+12z-1+12z-2,试画出其并联型网络结构。 5.5已知系统函数H(z)为 H(z)=2(1-z-1)(1-1.414z-1+0.7z-2)(1+0.5z-1)(1-0.9z-1+0.81z-2) (1) 画出H(z)的级联型网络结构; (2) 根据已画出的流图写出其状态方程和输出方程。 5.6已知某数字系统的系统函数为 H(z)=z3(z-0.4)(z2-0.6z+0.25) 试分别画出直接型、级联型、并联型结构框图。 5.7一线性时不变系统用题5.7图的流图实现。 题5.7图 (1) 写出该系统的差分方程和系统函数; (2) 计算每个输出样本需要多少次实数乘法和实数加法。 5.8求题5.8图各系统的单位脉冲响应。 题5.8图 5.9一线性时不变系统的单位脉冲响应 h(k)=ak,0≤k≤7 0,其他 (1) 画出该系统的直接型FIR结构流图。 (2) 证明该系统的系统函数为 H(z)=1-a8z-81-az-1 并由该系统函数画出由FIR系统和IIR系统级联而成的结构流图。