中科院研究生院—现代数字信号处理课件

1、第二章 维纳滤波和卡尔曼滤波 n 2.1 引言引言 n 2.2 离散维纳滤波器的时域解离散维纳滤波器的时域解n 2.3 离离散维纳滤波器的散维纳滤波器的z z域解域解 n 2.4 维纳预测维纳预测n 2.5 卡尔曼卡尔曼(Kalman)(Kalman)滤波滤波 2.1 引引 言言n最优滤波n维纳滤波和卡尔曼滤波简介n本章讨论的主要内容1、最优滤波n信号处理的目的是从噪声中提取信号，得到不受干扰影响的真正信号。采用的处理系统称为滤波器。n滤波器的分类：线性滤波器、非线性滤波器；FIR滤波器、IIR滤波器；时域滤波器、频域滤波器；图 2.1.1 信号处理的一般模型 x(n)=s(n)+v(n) (

2、 )( )( )( )( ) ()my ns nx nh nh m x nm( )( )( )( )( )e ns ny ns ns nn最优准则最优准则：最大输出信噪比准则匹配滤波器最小均方误差准则误差绝对值的期望值最小 误差绝对值的三次或高次幂的期望值最小2min| ( )| E e nmin| ( )|E e nmin| ( )| kE e nx(n)=s(n)+v(n) ( )( )( ) ()my ns nh m x nm( )( )( )e ns ny nWiener滤波器的一般结构滤波器的一般结构 2min| ( )| E e n2、维纳滤波和卡尔曼滤波简介n维纳(Wiener)

3、滤波与卡尔曼(Kalman)滤波以估计的结果与信号真值之间的误差的均方值最小作为最优准则。)( )()(nsnsne22( )E e nEss假设信号的真值与估计值间的误差为： 均方误差最小即误差的平方的统计平均值最小： 最小3、本章讨论的主要内容n主要内容：维纳滤波器（FIR维纳滤波器和IIR维纳滤波器）、维纳预测器、卡尔曼滤波。n分析思路：在均方误差最小的前提下，求得系统的单位脉冲响应h(n)或传递函数H(z)，进而计算滤波器的最小均方误差2min| ( )| E e n22min( )( )( )optmin E e nhnE e n2.2 离散维纳滤波器的时域解离散维纳滤波器的时域解n

4、本节要解决的主要问题及方法本节要解决的主要问题及方法n正交性原理正交性原理n维纳维纳霍夫方程霍夫方程nFIR维纳滤波器的时域解维纳滤波器的时域解1、本节要解决的主要问题及方法本节要解决的主要问题及方法n要解决的问题：要解决的问题：寻求在均方误差最小情况下的单位脉冲响应h(n)或传递函数H(z)的表达式，这一过程称为设计维纳滤波器的过程。n解决方法：解决方法：实质是求解维纳霍夫(Wiener-Hopf)方程，即2( )( )optmin E e nhn 本节讨论维纳滤波器的时域求解方法，即在时域本节讨论维纳滤波器的时域求解方法，即在时域求最小均方误差下的求最小均方误差下的 。( )opthn2、

5、 维纳滤波器时域求解的方法维纳滤波器时域求解的方法0( )( )( )( )()kky ns nx nh nh x nkn 因果维纳滤波器的输出y(n) ：n=0,1, 2, 设期望信号为d(n)，误差信号e(n)及其均方值E|e(n)|2分别为 e(n)=d(n)-y(n)=s(n)-y(n) 2*( )| ( )| ( )( )J nE e nE e n e n代价函数为( ),0,1,2,.kkh kajb k22| ( )| | ( )| 0kkE e nE e njabk=0, 1, 2, 记梯度算子为 kkkjab k=0, 1, 2, ( )khmin J nn 要使均方误差为最

6、小，须满足 0kJ nJ nh *kkkE e n enE e n enJ njab上式展开为 *2*( )( )( )( )| ( )| ( )( )( )( )kkkkke ne ne ne nE e nEe ne nje nje naabb又00( )( )( )( )()( )( )( )()kkke ns ny ns nh x nks na kjb kx nk将上述4式代入得 2*| ( )| 2 () ( )kkJ nE e nE x nk e n *( )()( )()( )()( )()kkkke nx nkae njx nkbe nx nkae njx nkb 分析：分析：上

7、式说明，若使滤波器的均方误差达到，则误差信号与输入信号正交，这就是通常所说的正交性原理。 *0()( )0,0,1,2,.koptJ nE x nk enkn正交性原理：正交性原理：n正交性原理的重要意义：提供了一个数学方法，用以判正交性原理的重要意义：提供了一个数学方法，用以判断线性滤波系统是否工作于最佳状态。断线性滤波系统是否工作于最佳状态。0 s x1 sse s w2x2 x2 w1x1 *0*0 ( )( )()( ) ()( )kkkkE y n e nEh x nk e nh E x nk e n 正交性原理的引理：最佳状态时，由滤波器输出定义的期望正交性原理的引理：最佳状态时，

8、由滤波器输出定义的期望响应的估计响应的估计yopt(n)与估计误差与估计误差eopt(n)正交：正交：0)()(*optnenyEopt3、 维纳维纳霍夫方程霍夫方程*,0()( )0,0,1,2,.()( )()0optopt iiE x nk enkE x nkd nhx ni将输入信号分配进去， 得到 *,0()()dxopt ixxirkhrikk=0, 1, 2, 维纳维纳-霍夫（霍夫（WienerHopf）方程：）方程：k=0, 1, 2, ,0( )()xdopt i xxirkhrki4、FIR维纳滤波器的时域解n FIR维纳滤波器的维纳维纳滤波器的维纳-霍夫方程霍夫方程 当h

9、(n)是一个长度为M的因果序列时，FIR维纳滤波器的维纳-霍夫方程表述为 10( )( )()Mxdxxirkh i rkik=0, 1, 2, ,M-1 (2.2.21)把k的取值代入(2.2.21)式， 得到 当k=0时，h0rxx(0)+h1rxx(-1)+hM-1rxx(-M+1)=rxd (0)当k=1时，h0rxx(1)+ h1rxx(0)+ hM-1rxx(-M+2)= rxd (1) 当k=M-1时，h0rxx(M-1)+ h1rxx (M-2)+hM-1rxx(0)= rxd (M-1) (2.2.22) 定义 011(0)(1)(1)(0)(1)(1)(1)(0)(2)(1

10、)(2)(0)xdxdxdxdMxxxxxxxxxxxxxxxxxxxxhrhrhRrMhrrrMrrrMRrMrMr（2.2.22）式可以写成矩阵形式矩阵形式， 即 xdxxRR h 对上式求逆，得到 1optxxxdhhR R 这里涉及到计算相关矩阵和逆矩阵，计算量可能较大。n FIR维纳滤波器的估计误差的均方值维纳滤波器的估计误差的均方值 假定所研究的信号都是零均值的，滤波器为FIR型，长度等于M， 12*01*012*min012*20| ( )| ( ) ( )( ) ()( )( )( )()| ( )| ( )( ) ( )( ) ()( )( )( )MmMmMmMdxddmE

11、 e nE e n d nh m x nmE e n dnhm E e n x nmE e nE e n dnEd nh m x nm dnh m rmHxdR h 结论：在所有在所有N阶阶FIR滤波器中，最优滤波器的均方误差值是滤波器中，最优滤波器的均方误差值是最小的，从这个意义上说，它是最优的。最小的，从这个意义上说，它是最优的。其阶数越高，采用的已知信息就越多，它的最小均方误差就越小，但相应的计算量也越大。例例1 假设有一实的广义平稳随机信号s(n)的自相关函数（序列）为 ,且伴随有实的噪声w(n)，方差为 ，与s(n)无关。试设计一个M=2的FIR维纳滤波器来估计s(n) ，并计算最小

12、均方误差。 0.6mssrm 21w1optxxxshhR R22Hmin| ( )| ()sxsoptE e nRhn解：已知 0.6 ,msswwrmrmm 0.6mxxxsxwsswwrmrmrmrmrmm 0120.6100.62xxxxxxxxxxrrRrr 01010.610.60.6xsxsxsrRr由此，M=2最佳FIR维纳滤波器如下：10.4510.165optxxxshR Rn或者，利用下式求解10( )( )()Mxdxxirkh i rkik=0, 1当k=0时，2h0+0.6h11当k=1时，0.6h0+2 h10.6010.451;0.165hh估计该滤波器的输出误

13、差的最小均方值： 122*Tmin0| ( )| ( )( )0()0.451110.60.1650.45MdxdssxsoptmE e nh m rmrRh 222( )( )( )( )(0)1wwE e nEx ns nE w nr经过此滤波器以前的均方误差为 2.3 离散维纳滤波器的离散维纳滤波器的z域解域解n本节要解决的主要问题及方法本节要解决的主要问题及方法n白化滤波器白化滤波器n非因果非因果IIR维纳滤波器的维纳滤波器的Z域解域解n因果因果IIR维纳滤波器的维纳滤波器的Z域解域解1、本节要解决的主要问题及方法本节要解决的主要问题及方法n待解决的问题：待解决的问题：当h(n)是物理

14、可实现的因果序列时，所得到的Wiener-Hopf方程 将存在k0的约束，不能直接转到Z域求解。这使得在要求满足物理可实现条件下，求解维纳-霍夫方程成为一个十分困难的问题。 n解决方法：解决方法：采用将观测序列x(n)白化的方法，求解Wiener-Hopf方程的Z域解。n 若不考虑滤波器的因果性，维纳霍夫方程可以改写为 ( )( )()( )*xdxxxxmrkh m rkmh krk 设定d(n)=s(n)，对上式两边做Z变换，得到 Sxs (z)=Hopt(z)Sxx(z) ( )( )( )xsoptxxSzHzSzx(n)=s(n)+v(n) 假设信号和噪声不相关，即rsv(m)=0，

15、则 Sxs (z)=Sss(z)Sxx(z)=Sss(z)+Svv(z) ( )( )( )( )( )( )xsssoptxxssvvSzSzHzSzSzSzn 对于因果IIR维纳滤波器，其维纳霍夫方程为 0( )( )()( )( )xdxxxxmrkh m rkmh krkk=0, 1, 2, 因为存在k0的约束，使得上式不能直接转到Z域求解。如有可能将其转化为非因果问题，则求解会大大简化。 如果滤波器的输入是白噪声，即x(n)= w(n)，w(n)是方差为2w的白噪声，由于 2( )( )xxwwwrkrkk 则因果IIR维纳滤波器的维纳霍夫方程变为：220( )( )( )( )xd

16、wdwwmrkrkh mkmh k k=0, 1, 2, 2( )wdwrkh kk=0, 1, 2, 由此可见，只要将输入信号转化为白噪声，就可以解得因果IIR维纳滤波器的单位脉冲响应。2、白化滤波器、白化滤波器n任何具有有理功率谱密度的随机信号都可以看成是由一白色噪声w(n)激励某个物理网络所形成。x(n)的时间序列信号模型 21( )( ) ()xxwPzB z B z)()()(zWzBzX一般把信号转化为白噪声的过程称为白化白化，对应的滤波器称为白化滤波器白化滤波器。 x(n)的白化滤波器 如果B(z)是一个最小相移网络函数，那么1/B(z)显然也是一个物理可实现的最小相移网络，因此

17、可以利用上式白化x(n)。)()(1)(zXzBzW利用白化x(n)的方法求解维纳-霍夫方程 n 利用白化x(n)的方法求解维纳-霍夫方程:)()()(zBzGzH于是，在最小均方误差准则下，求最佳于是，在最小均方误差准则下，求最佳Hopt(z)的问题就归结的问题就归结为求最佳为求最佳G(z)的问题了。的问题了。G(z)当然也需分因果性或非因果性的当然也需分因果性或非因果性的约束情况加以讨论。约束情况加以讨论。21( )( ) ()xxPzB z B z 如果已知信号的Pxx(z)，即可由下式求得B(z) 。 n 计算计算Hopt (z)： ( )( )( )( )( ) ()ky ns nw

18、 ng ng k w nk3、 非因果非因果IIR维纳滤波器的求解维纳滤波器的求解22*22*2ss( )() ( )( ) ()(| ( )| ( )( )( ) ()(0)|( )|( ) ( )() ()| kkwkrkkEg kE e nEs ng kEg k g r wE s nwnk w nrnw nkrk s ng k w nk s ng kg k r 第一项第二项*22ss2( )( )( )( )|(0)( )wswskkwswswkkwwkg k rkrkrrg k 第三项(2.3.9) 求满足最小均方误差条件下的求满足最小均方误差条件下的g(k)：为求得相对于g(k)的最

19、小均方误差值，令( )( )0wswwrkg k -k 2| ( )| 0E e ng k2( )( )wsoptwrkgk -k 2( )( )wsoptwSzGzZ变换后 非因果非因果IIR维纳滤波器的最佳解：维纳滤波器的最佳解： optopt2( )( )1( )( )( )wswGzSzHzB zB zs(n)=s(n)*(n)，x(n)=w(n)*b(n)rxs(m)=rws(m)*b(-m) Sxs (z)=Sws(z)B(z-1) 1( )( )()xswsSzSzB z 非因果IIR维纳滤波器的复频域最佳解的一般表达式 opt21( )( )11( )( )()( )xsxsw

20、xxSzSzHzB z B zSz假定信号与噪声不相关，即当Es(n)v(n)=0时可以得到： Sxs(z)=Sss(z) Sxx(z)=Sss(z)+Svv(z) 1( )( )()sswsSzSzB z 信号和噪声不相关时，非因果IIR维纳滤波器的复频域最佳解和频率响应分别为 )()()()()()(optzSzSzSzSzSzHvvssssxxxs)()()()e ()e ()e ()e (jjjjoptvvssssvvssssPPPSSSHn 由上式可知： 当噪声为0时，Hopt=1，信号全部通过； 当信号为0时， Hopt=0，噪声全部被抑制掉； 当即有信号又有噪声时， Hopt1，

21、大小随Pvv的增加而减小，从而达到降低噪声影响的目的。011)e (joptHPss(ej)0, Pvv(ej)=0 Pss(ej)0, Pvv(ej) 0 Pss(ej)=0, Pvv(ej) 0 图 2.3.6 非因果维纳滤波器的传输函数的幅频特性 n 计算最小均方误差计算最小均方误差E|e(n)|2min： 22min2|( )| ( )| (0)wssskwrkE e nr 第一项根据围线积分法求逆Z变换的公式, rss(m)用下式表示: CmsssszzzSmrd)(j21)(1得出 CsssszzzSrd)(j21)0(1 第二项由帕塞伐尔定理：zzzYzXnynxCnd1)(j2

22、1)()(*取y(n)=x(n), 有 zzzXzXnxCnd)(j21| )(|12因此 211d|( )|( )()2jwswswsCnzrkSz Szz得到 21min211d| ( )| ( )( )()2jsswswsCwzE e nSzSz Szz112min121()( )( )()( )( )()11d| ( )| ( )2j1d( )()2jwsxsopxsssCsstxsCwsszszSzzE e nSzzzSSzB zB zHzzSzz 1( )( )()xswsSzSzB z21min211d| ( )| ( )( )()2jsswswsCwzE e nSzSz Szz

23、 假定信号与噪声不相关，Es(n)v(n)=0，又因为实信号的自相关函数是偶函数，即rss(m)=rss(-m)，则2m1(1)in( )1d| ( )| ( )()2j( )( )( )( )()1d2j( )( )( )( )1d2j( )( )( )1d2j( )optHssssssCxxssxxssssCxxssxxssCxxssvvCxzxSzzE e nSzSzSzzSz SzSz SzzSzzSzSzSzzSzzSz SzzSzz Sxs(z)=Sss(z) ，Sxx(z)=Sss(z)+Svv(z)； Sss(z)=Sss(z-1) 4 4、 因果因果IIRIIR维纳滤波器的求

24、解维纳滤波器的求解n 若维纳滤波器是一个因果滤波器， 要求 g(n)=0 n0 则滤波器的输出信号 0( )( )( )( )( ) ()ky ns nw ng ng k w nk估计误差的均方值 E|e(n)|2=E|s(n)-y(n)|2 类似于（2.3.9）式的推导，得到 222200( )1| ( )| (0)( )|( )|wssswwskkwwrkE e nrg krk要使均方误差取得最小值， 当且仅当 2opt2( )0( )00( )( )wswwswrnngnnrnu n令 0opt221( )( ) ( )( )( )11( )ZT( )( )()nnwswswsnnxso

25、ptwswwSzrn u n zrn zSzGzgnSzB z因果维纳滤波器的复频域最佳解为 )()()(11)()()(12optzBzSzBzBzGzHxsopt维纳滤波的最小均方误差为 22min20*212121|( )| ( )| (0)1(0)( ) ( )( )11d(0)( )()2j( )()11d( )2j()( )1( )2jwssskwsswswsksswswsCxsxsssCwssCrkE e nrrrk u k rkzrSzSzzSzSzzSzB zB zzSz 1optd( )()xszHz Szz 非因果情况时，滤波器的最小均方误差为22min2|( )| (

26、)| (0)wssskwrkE e nr 对于因果情况， 22min20|( )| ( )| (0)wssskwrkE e nr 比较两式，可以看出非因果情况的E|e(n)|2min一定小于等于因果情况E|e(n)|2min。 因果维纳滤波器设计的一般方法： (1) 根据观测信号x(n)的功率谱求出它所对应信号模型的传输函数，即采用谱分解的方法得到B(z)，Sxx(z)=2wB(z)B(z-1)。 (2) 求的Z反变换，取其因果部分再做Z变换，即舍掉单位圆外的极点，得 (3) 积分曲线取单位圆，计算Hopt(z), E|e(n)|2min。 )()(1zBzSxs)()(1zBzSxs例例 2

27、.3.1 已知 )8 . 01)(8 . 01 (36. 0)(1zzzSss信号和噪声不相关，即rsv(m)=0，噪声v(n)是零均值、单位功率的白噪声（2v=1,mv=0)，求Hopt(z)和E|e(n)|2min。 解解 根据白噪声的特点得出Svv(z)=1, 由噪声和信号不相关， 得到rxx(m)=rss(m)+rvv(m)。 11211( )( )( )0.361(1 0.8)(1 0.8 )1.6 (1 0.5)(1 0.5 )( ) ()(1 0.8)(1 0.8 )xxssvvwSzSzSzzzzzB z B zzz考虑到B(z)必须是因果稳定的系统，得到 1211 0.5(

28、),1.61 0.8wzB zz (1)、 首先分析物理可实现情况:1opt2111( )111 0.80.36( )( )()1.6 (1 0.5)(1 0.8)(1 0.5 )xswSzzHzB zB zzzz因为 111110.810.360.36Re,0.8(1 0.8)(1 0.5 )(1 0.8)(1 0.5 )0.360.8(1 0.8)(1 0.5 )3(0.8)5nnznZszzzzzzzzz1opt1111 0.80.631( )1.6 (1 0.5)1 0.881 0.5zHzzzz取其因果部分 110.3633/5(0.8)(1 0.8)(1 0.5 )51 0.8nZ

29、Tu nzzz 110.363(0.8)(1 0.8)(1 0.5 )5nZu nzz2min1opt111| ( )| 1d( )( )()2j310.360.36d82j(1 0.8)(1 0.8 )1 0.5(1 0.8)(1 0.8 )50.45(0.5)182(0.8)(1/0.8)(0.5)ssxsCCCE e nzSzHz Szzzzzzzzzzdzjzzz取单位圆为积分围线，上式等于单位圆内的极点 )5 . 08 . 0(zz及的留数之和，即 未经滤波器的均方误差 1| )(| )()(| )(|2222vnvEnsnxEneE20.80.5min550.450.50.450.

30、588( )(1/0.8)(0.5)(0.8)(1/0.8)3/8zzzzE e nzzzz 所以通过因果维纳滤波器后均方误差下降8/3(2.7)倍。 (2)、 对于非物理可实现情况有 opt111( )( )( )( )( )( )0.36(1 0.8)(1 0.8 )0.361(1 0.8)(1 0.8 )0.225(1 0.5)(1 0.5 )xsssxxssvvSzSzHzSzSzSzzzzzzz21minopt111111d| ( )| ( )( )()2j10.360.2250.36d2j (1 0.8)(1 0.8 )(1 0.5)(1 0.5 )(1 0.8)(1 0.8 )1

31、0.360.22512(1 0.8)(1 0.8 )(1 0.5)(1 0.5 )ssxsCCzE e nSzHz Szzzzzzzzzzdzjzzzzz11110.36 (1.0250.50.5 )d2j(1 0.8)(1 0.8 )(1 0.5)(1 0.5 )Czzzzzzzz令 1110.36 (1.0250.50.5 )( )(1 0.8)(1 0.8 )(1 0.5)(1 0.5 )zzF zzzzz z单位圆内有两个极点0.8和0.5, 应用留数定理，有 1035 . 0),(Res8 . 0),(Res)(min2zFzFneE结论：比较两种情况下的最小均方误差,可以看出非物理

32、可实现情况的最小均方误差小于物理可实现情况的均方误差。 维纳滤波部分的总结：n主要内容：主要内容：FIR维纳滤波求解、非因果IIR维纳滤波求解、因果IIR维纳滤波求解；n知识点：知识点：最小均方误差准则、正交性原理、维纳霍夫方程、白化滤波器；n结论：结论：在所有N阶FIR滤波器中，最优滤波器的均方误差值是最小的，从这个意义上说它是最优的；与非最优滤波相比，最优滤波的优势在于能对滤波的质量（逼近的好坏）做出评价；E|e(n)|2min与Sss(z)和Svv(z)重叠部分大小有关；最小均方误差比较：非因果IIR因果IIRFIR维纳滤波的最小均方误差2.4 维维 纳纳 预预 测测n本节讨论的主要问题

33、及方法本节讨论的主要问题及方法n预测的可能性预测的可能性n维纳预测的计算维纳预测的计算n纯预测纯预测n一步线性预测的时域解一步线性预测的时域解1、本节讨论的主要问题及方法、本节讨论的主要问题及方法n讨论的主要问题：本节将讨论维纳预测器，以观测到的全部过去数据来估计当前或将来的值n解决方法：以均方误差最小为估计原则 0),( NNnsmin)( )(2NnsNnsE H(z) x(n)=s(n)+(n) )( )(nsny( H(z) x(n)=s(n)+(n) ( )()y ns nN 图2.4.1(b) 维纳预测器图2.4.1(a) 维纳滤波器2、预测的可能性、预测的可能性n信号可以预测是由

34、于信号内部存在关联性。数据间关联越密切，预测越准确；完全不关联，则无法预测。 2( )wwwrmm 0( )0wwmrm时， 21( ),0 xxwxxSzB z B zrm输入：输出： x(n)在各不同时间点上的值的相关性是起因于B(z)的惯性。由观测到的x(n)的所有过去值和当前值来估计将来值 时：如果 ，则 仅由B(z)的惯性决定，如果 ，则 将由B(z)的惯性和 共同决定；()0 x nNN()00w nNN()0 x nNN()00w nNN()0 x nNN()0w nNN()x nNx nNn随机信号预测的特点：随机信号预测的特点：以信号的统计特性作为预测的主要依据；不可能作预测

35、误差为零的绝对精确的预测；实际信号通常带有噪声干扰，使得预测和滤波联系在一起，成为带滤波的预测。3 3、 维纳预测的计算维纳预测的计算0( )()( ) ()()()()my ns nNh m x nNme nNs nNs nN( )()dyns nN )()(Nnsnyd )()()(nnsnx )( )(Nnsny H（z） 同理，要使预测误差的均方值为最小，须满足 2| ()| 0kE e nNh其中，hk表示h(k)。 02jiiixxhNnsE0( )( )(), 0dxyoptxxmrkhm rkmk即 2()jhmin E e nNNxsxyxsdxyzzSzSkNrkNnsnx

36、EknynxEkrdd)()()()()()()()(*n 非因果维纳预测器的最佳解为 )()()()()(optzSzSzzSzSzHxxxsNxxxydn 因果维纳预测器的最佳解为 )()()(11)()()(11)(1212optzBzSzzBzBzSzBzHxsNxyd维纳预测的最小均方误差为 21minopt1opt1| ()| ( )( )()2j1( )( )()2jdssxyCNssxsCdzE e nNSzHz SzzdzSzHz Szzz 维纳预测的求解和维纳滤波器的求解方法是一致的。 4 4、 纯预测纯预测n 假设假设x(n)=s(n)+v(n)，纯预测问题是在，纯预测问

37、题是在v(n)=0情况下对情况下对s(n+N), N0的预测，此时的预测，此时x(n)=s(n)。 因果情况下，假设s(n)与v(n)不相关，纯预测情况下的输入信号的功率谱及维纳预测器的最佳解分别为 )()(1)()()(11)()()()()()(12opt12zBzzBzBzSzzBzHzBzBzSzSzSNxsNssxsxx纯预测器的最小均方误差为 CNNCNNCNxssszdzzBzzzBzBzBzdzzzBzBzBzBzzBzBzdzzzSzHzSNneE)()()()(j2)()()()()()(j21)()()(j21| )(|11212121optmin2应用复卷积定理 zzz

38、YzXnynxCnd1)(j21)()(*取y(n)=x(n) zzzXzXnxCnd)()(j21)(12得到 )( )()()()()()(| )(|10220022222min2nbNnbnbNnbnuNnbnbNneENnnnnn 可以看到，随着N增加，E|e(n+N)|2min也增加。这一点也容易理解，当预测的距离越远，预测的效果越差，偏差越大，当预测的距离越远，预测的效果越差，偏差越大，因而因而E|e(n+N)|2min越大。越大。 例例2.4.1 已知 )1)(1 (1)(),()(12azazazSnsnxxx其中-1a1， 求: (1) 最小均方误差下的s(n+N)； (2)

39、 E|e(n+N)|2min。 解解 （1）、求B(z)：首先对Sxx(z)进行功率谱分解。因为 1221211)(,1)1)(1 (1)(azzBaazazazSxx所以 （2）、求H(z)： 求出B(z)的Z反变换 )(11T11nuaaznZ 对于因果维纳预测器有： 1ZTT11)()(11 azanuaNnuaazzNNnNnZN取因果部分NNNaazaazzBzzBzH11opt1)1 ()()(1)(11111)()(1)(azzazzBzzBzHNNopt图 2.4.1 纯预测维纳滤波器 由Hopt(z)=aN，此时可以把纯预测的维纳滤波器看作是一个线性比例放大器（如图2.4.1

40、所示）。)()( nsaNnsN x(n)=s(n) )()( )(nsaNnsnyN a N u 为什么预测值只与当前值s(n)有关，而和s(n+1)、 s(n+2)无关呢？ 将x(n)看成由白噪声(n)通过B(z)产生的，根据x(n)的信号模型 111)(azzB可以写出x(n)的时间序列模型所对应的输入输出方程 x(n)=(n)+ax(n-1) 以上推导结果从统计意义上讲，当N0时，白噪声信号(n+N)对x(n)无影响。 23(1)(1)( )( )( ) (1) 0(2)(2)(1)(1)( ) (2) 0(3)(3)(2)(2)( ) x nnax nax nas nnx nnax

41、nax na s nnx nnax nax nas n 当当(3) 0 ()()(1)(1)( ) () 0Nnx n Nn Nax n Nax n Na s nn N 当当将信号x(n)通过纯预测维纳滤波器，随着时间的递增，基于 和 ) 1()()(naxnnx)()(nsnx x的均值等于零，正说明 的影响就统计平均来讲等于零。当当mx=0，估计，估计 时，只需要考虑系统时，只需要考虑系统B(z)的惯性而的惯性而可认为可认为 ，这样估计出来的结果将有最小均方误，这样估计出来的结果将有最小均方误差。差。 ( )( )nx n对()s nN()0(0)nNN 211( )0(1)(1)xxxa

42、SzmE x nazaz实际上，由可知，21lim(1)( )lim( )0 xxxxxzmzSzrmm0 xm 这是因为，如果，由终值定理可知那么Sxx(z)在z=1处有一个极点，而现在在z1处无极点，故有 0 xm （3）、 求E|e(n+N)|2min1222min0122012202()( )( )(1)1 11NnNnNnnNE e nNbna u naaaa 结论：由上式可知，N越大，误差越大，如果N=0则没有误差。 ()( )( )( )0Ns nNa s ns ns ne n5、 一步线性预测的时域解一步线性预测的时域解n 一步线性预测一步线性预测：采用p个最近的采样值来预测时

43、间序列下一时刻的值，包括前向预测和后向预测两种。 前向预测前向预测：在噪声v(n)=0的情况下，已知x(n-1), x(n-2)，,x(n-p), 预测当前时刻x(n)； 后向预测后向预测：在噪声v(n)=0的情况下，已知x(n)，x(n-1)，x(n-p+1)基础上，估计x(n-p)。 图 2.4.2 前后向预测数据之间的关系 x(n p) , x(n p1) , , x(n2) , x(n1) , x(n)后向预测前向预测（1）、前向预测）、前向预测n设定系统的单位脉冲响应为h(n)，其输出信号为1( )( )( )( ) ()pky ns nx nh k x nk令apk=-h(k)，则

44、 pkpkknxanx1)()( n 前向预测误差为 pkpkpkpkknxaknxanxnxnxne01)()()()( )()(其中, ap0=1, 一步前向预测器结构图 n前向预测误差的均方值为： 212)()(| )(|pkpkknxanxEneEplaneEpl, 2 , 10| )(|2或Ee (n)x* (n-l)=0 l=1, 2, , p 即*12( )()()01,2,PpkkEx na x nkx nllppkxxpkxxlkralr10)()( 由于预测器的输出 是输入信号的线性组合，故预测误差与预测的信号值同样满足正交性原理：)( nx0)( )(*nxneEn 前向

45、预测误差的最小均方值为： 2*min*11| ( )| ( )( ( )( )( ) ( )( )()( )(0)( )ppkkpxxpk xxkE e nE e n x nx nE e n x nEx na x nkx nra rkpkxxpkxxpkxxpkxxpllkralrneEkrar11min2, 2 , 10)()(| )(|)()0(将方程组写成矩阵形式 （Yule-Walker方程）方程）00| )(|1)0() 1()() 1()0() 1 ()() 1 ()0(min21neEaarprprprrrprrrpppxxxxxxxxxxxxxxxxxx1optxxxdhhR

46、R12(0)(1)(1)(0)(1)(0)(2)(1)(1)(1)(2)(0)xxxxxxxdxxxxxxxdxdMxxxxxxrrrMhrhrrrMrrMhrMrMr 维纳霍夫方程00| )(|1)0() 1()() 1()0() 1 ()() 1 ()0(min21neEaarprprprrrprrrpppxxxxxxxxxxxxxxxxxx Yule-Walker方程 前向预测误差为 1( )( )( )( )()ppkke nx nx nx na x nk AR信号模型为 1( )()pkkx na x nkw n 222min,pkkwaa e nw nE enE wn对比两式可知，

47、（2）、后向预测）、后向预测n假设前、后向预测器具有相同的系数，即 1( )()()()ppkky ns npx npa x npk n 后向预测误差为 pkpkkpnxapnxpnxpnxnb1)()()( )()(n后向预测误差的均方值为： 2122)()()( )()(PkpkkpnxapnxEpnxpnxEnbE2| ( )| 01,2,plE b nlpa或Eb (n)x* (n-p+l)=0 l=1, 2, , p 即*12()()()01,2,ppkkEx npa x npkx npllppkxxpkxxlkralr10)()( 由于预测器的输出 是输入信号的线性组合，故预测误差

48、与预测的信号值同样满足正交性原理：()x np*( ) ()0E b n x npn 后向预测误差的最小均方值为： 2*min*1( )( )()()( ) ()( ) ()( ) ()()()()ppkkE b nE b nx npx npE b n x npb n x npE b n x npEx npa x npkx np同理，可以得到下面方程组： 2min11(0)( )| ( )| ( )()01,2,pxxpk xxkpxxpk xxkra rkE b nrla rkllp将方程组写成Yule-Walker方程方程形式2min1(0)(1)( )1| ( )| (1)(0)(1)0

49、0( )(1)(0)xxxxxxpxxxxxxppxxxxxxrrrpE b narrrparprpr Yule-Walker方程具有以下特点： (1) 除了第一个方程外，其余都是齐次方程; (2) 与维纳-霍夫方程相比，不需要知道rxs(m)。 (3) 由方程组的p+1个方程，可以确定apk,k=1, 2, , p和Ee2(n)min，共计p+1个未知数。n Levinson-Durbin算法算法 Levinson-Durbin算法首先由一阶AR模型开始，一阶AR模型(p=1)的Yule-Walker为 01)0() 1 ()2()0(211 , 1arrrrxxxxxxxx由该方程解出：

50、)0()1 ()0() 1 (21 , 1211 , 1xxxxxxrarra然后增加一阶，即令p=2，得到 001)0() 1 ()2() 1 ()0() 1 ()2() 1 ()0(222, 21 , 2aarrrrrrrrrxxxxxxxxxxxxxxxxxx 由上面方程解出： 2122, 2221 , 12, 21 , 1221 , 2211 , 12222, 2)1 ()1 ()0(/)2() 1 () 1 ()0(/)1 ()2()1 ()0(/)1 ()2()0(aaaarrrrrrararrrrrraxxxxxxxxxxxxxxxxxxxxxxxxxx 然后令p=3, 4, ,

51、 以此类推， 可以得到Levinson-Durbin的一般递推公式如下： )()0()1 (1, 3 , 2 , 1 )()(2202122, 1, 1,21111nxErkpkakaaakkpraprkxxpppkpppkpkppppppkxxpxxp例例2.4.2 已知 )8 . 01)(8 . 01 (36. 0)(1zzzSxxx(n)为AR模型，求AR模型参数（包括模型阶数和系数）。rxx(m)=0.8|m| 解解 首先对Sxx(z)做傅里叶反变换，得到x(n)的自相关函数rxx(m), （1）、采用试验的方法确定模型阶数p。首先取p=2，各相关函数值由上式计算 00118 . 06

52、4. 08 . 018 . 064. 08 . 01221aa计算得到 a1=-0.8， a2=0, 2=0.36 （2）、如果取p=3，可计算出a1=-0.8, a2=a3=0, 2=0.36，说明AR模型的阶数只能是一阶的。（3）、采用谱分解的方法，即对Sxx(z)进行谱分解，得到的模型也是一阶的，其时间序列模型和差分方程为 ) 1(8 . 0)()(8 . 011)(1nxnnxzzB2.5 卡尔曼卡尔曼(Kalman)滤波滤波 n本节讨论的主要问题及方法本节讨论的主要问题及方法n卡尔曼滤波的状态方程和量测方程卡尔曼滤波的状态方程和量测方程n卡尔曼滤波的递推算法卡尔曼滤波的递推算法n发散

53、问题及其抑制发散问题及其抑制1、本节讨论的主要问题及方法、本节讨论的主要问题及方法n讨论的主要问题：本节将主要讨论卡尔曼滤波的状态方程和量测方程，及其递推算法。n解决方法：利用状态方程和递推方法寻找最小均方误差下状态变量 的估计值 ，即kxkxkkkxxxminTkkkE x xx 2、卡尔曼滤波的状态方程和量测方程卡尔曼滤波的状态方程和量测方程 假设某系统k时刻的状态变量为xk，状态方程和量测方程（也称为输出方程）表示为 11kkkkxA xwkkkkvxCy Ak为状态转移矩阵，描述系统状态由时间k-1的状态到时间k的状态之间的转移； Ck为量测矩阵，描述状态经其作用，变成可量测或可观测的

54、； xk为状态向量，是不可观测的；yk为观测向量； wk为过程噪声；vk为量测噪声。图 2.5.1 卡尔曼滤波器的信号模型 kskx 假设状态变量的增益矩阵A不随时间发生变化，wk,vk都是零均值白噪声，方差分别是Qk和Rk，并且初始状态x0与wk,vk都不相关，且噪声向量wk,vk也互不相关，即2200:0,:0,0;,00TkkwkkjkkjTkkvkkjkkjkkTkkwE wQE w wQvE vRE v vRCov x wCov x vE w v其中 jkjkkj013、 卡尔曼滤波的递推算法卡尔曼滤波的递推算法 n 基本思想： 先不考虑输入信号k和观测噪声vk的影响，得到状态变量和

55、输出信号（即观测数据）的估计值 和 再用输出信号的估计误差 加权后校正状态变量的估计值 ，使状态变量估计误差 的均方值最小。 kxkyky kxkxminkkkkkkkkkTkkkxyyyyxxxxE x xx 量测方程校正 当不考虑观测噪声和输入信号时，状态方程和量测方程为：11 kkkkkkkkkxACxCyxAx 输出信号的估计误差（新息）为：kkkyyy 为了提高状态估计的质量，用输出信号的估计误差 来校正状态变量 ky)()(111kkkkkkkkkkkkkxACyHxAyyHxAx其中，Hk为增益矩阵，实质是一加权矩阵。 校正后状态变量的估计误差及其均方值分别为：kkkxxxT()

56、() TkkkkkkkPE x xE xxxx 未经校正的状态变量估计误差的均方值为：T()() kkkkkPE xxxx 卡尔曼滤波要求状态变量的估计误差的均方值Pk为最小， 因此卡尔曼滤波的关键就是要得到卡尔曼滤波的关键就是要得到Pk与与Hk的关系式，即通过选择的关系式，即通过选择合适的合适的Hk,使使Pk取得最小值。取得最小值。 minkkkPHx1111111()()()()()kkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkxA xHyyA xHC xvC A xIH CA xH C xH vIH CA xH CA xH vn 递推步骤：1. 计算状态变量的估计

57、值 ： kx2. 计算状态变量的估计误差 ： x11111111111111111()()()()()()()()()kkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkxxxA xIH CA xH CA xH vA xxH C A xxH C AH vIH CA xxIH CH vIH CA xxH v11kkkkxA x 由上式可以看出，状态变量的估计误差 由三部分组成， 可记为 xcbax其中 kkkkkkkkkkvHcCHIbxxACHIa111)()()(TTTT11TTT1T()()()kkkkkkkkTTkkaxxAIH CbIH

58、Ccv H3. 计算状态变量估计误差的均方值Pk ： TTTTTTTTTTT()() kkkPE x xE abc abcE aabbccabacbabccacb 其中，T1111TTTTT1()()()()()kkkkkkkkkkkkkkkkkE aaE IH CAAIH CIH CAAIHxxxPCxTTTT111()() ()()kkkkkkkkkkkE bbE IH CIH CQIH CIH CTTkkkTkkkkv vE ccEHH RHHT1TTTTTTTTT1TT111TTTT11T111()0() 0()()()()() )()(kkkkkkkkkkkkkkkkkkkkkkkk

59、kkkkkkkkkkvvxxIH CE bcE IH CHE cbE HIH CE abE IH CAE acE IH CAHE baE IH CAIH CE caEvHxxxxvT11TT) )(kkkkkkAIH Cxx 根据假设的条件，状态变量的增益矩阵A不随时间发生变化，起始时刻为0，则10022110012121001100 kkkjkjjkkkjkjjxAxxAxA xAxAxAxA xA xk-1仅依赖于x0,0, 1,k-2,与k-1不相关，即 TT11110kkkkE xEx 又112111121111112()()kkkkkkkkkkkkkkkxAxHyyAxHCxvCAx

60、 仅依赖于xk-1,vk-1，而与vk不相关，即 1kx0)()(0)()(T111T111T11T11kkkkkkkkkkkkxxExxExxvEvxxE由此可知，TTTT0;0;0;0E abE acE baE ca 也就是说，Pk仅有其中的三项不为零， 化简成 TTTT1T11T1()()()()()()kTkkkkkkkTkkkkkkkkTkkkkkkTkkkkkPE aaE bbE ccIH CA PAIH CIH CQIH CH R HA PAIH CICH R HQH4. 计算未经误差校正的状态变量估计误差的均方值Pk：TT111111T111111TTT111T11111def