第二章 维纳滤波和卡尔曼滤波第二章 维纳滤波和卡尔曼滤波
2.1 引言
2.2 维纳滤波器的离散形式 —— 时域解
2.3 离散维纳滤波器的 z域解
2.4 维纳预测
2.5 卡尔曼 (Kalman)滤波第二章 维纳滤波和卡尔曼滤波
2.1 引 言在生产实践中,我们所观测到的信号都是受到噪声干扰的 。 如何最大限度地抑制噪声,并将有用信号分离出来,是信号处理中经常遇到的问题 。 换句话说,信号处理的目的就是要得到不受干扰影响的真正信号 。 相应的处理系统称为滤波器 。 这里,
我们只考虑加性噪声的影响,即观测数据 x(n)是信号 s(n)与噪声
v(n)之和 ( 如图 2.1.1所示 ),即
x(n)=s(n)+v(n) ( 2.1.1)
第二章 维纳滤波和卡尔曼滤波我们的目的是为了得到不含噪声的信号 s(n),也称为期望信号,若滤波系统的单位脉冲响应为 h(n)( 如图 2.1.2所示 ),
系统的期望输出用 yd(n)表示,yd(n)应等于信号的真值 s(n);系统的实际输出用 y(n)表示,y(n)是 s(n)的逼近或估计,用公式表示为 yd(n)=s(n),y(n) = 。 因此对信号 x(n)进行处理,可以看成是对期望信号的估计,这样可以将 h(n)看作是一个估计器,
也就是说,信号处理的目的是要得到信号的一个最佳估计 。 那么,
采用不同的最佳准则,估计得到的结果可能不同 。 所得到的估计,在通信中称为波形估计 ; 在自动控制中,称为动态估计 。
)(?ns
第二章 维纳滤波和卡尔曼滤波图 2.1.1 观测信号的组成
x ( n )s ( n )
v ( n )
第二章 维纳滤波和卡尔曼滤波图 2.1.2 信号处理的一般模型
h ( n )
x ( n )
s ( n ) + v ( n )
y ( n )
第二章 维纳滤波和卡尔曼滤波假若已知 x(n-1),x(n-2),…,x(n-m),要估计当前及以后时刻的信号值 s(n+N),N≥0,这样的估计问题称为预测问题;若已知
x(n-1),x(n-2),…,x(n-m),要估计当前的信号值 s(n),称为过滤或滤波; 根据过去的观测值 x(n-1),x(n-2),…,x(n-m),估计过去的信号值 s(n-N),N≥1,称为平滑或内插 。 维纳 (Wiener)滤波与卡尔曼 (Kalman)滤波就是用来解决这样一类从噪声中提取信号的过滤或预测问题,并以估计的结果与信号真值之间的误差的均方值最小作为最佳准则 。
^
^
第二章 维纳滤波和卡尔曼滤波维纳滤波是在第二次世界大战期间,由于军事的需要由维纳提出的 。 1950年,伯特和香农给出了当信号的功率谱为有理谱时,由功率谱直接求取维纳滤波器传输函数的设计方法 。 维纳滤波器的求解,要求知道随机信号的统计分布规律 ( 自相关函数或功率谱密度 ),得到的结果是封闭公式 。 采用谱分解的方法求解,简单易行,具有一定的工程实用价值,并且物理概念清楚,但不能实时处理;维纳滤波的最大缺点是仅适用于一维平稳随机信号 。 这是由于采用频域设计法所造成的,因此人们逐渐转向在时域内直接设计最佳滤波器的方法 。
第二章 维纳滤波和卡尔曼滤波
2.2 维纳滤波器的离散形式 —— 时域解
2.2.1 维纳滤波器时域求解的方法根据线性系统的基本理论,并考虑到系统的因果性,可以得到滤波器的输出 y(n),


0
)()()()()(
m
mnxmhnhnxny
n=0,1,2,… ( 2.2.2)
设期望信号为 d(n),误差信号 e(n)及其均方值 E[ |e(n)|2]分别为
e(n)=d(n)-y(n)=s(n)-y(n) ( 2.2.3)


2
0
22 )()()(]|)()([|]|)([|
m
mnxmhndEnyndEneE
( 2.2.4)
第二章 维纳滤波和卡尔曼滤波要使均方误差为最小,须满足
0]|)([|
2

jh
neE ( 2.2.5)
这里,hj表示 h(j); 同理,可以用 aj,bj分别表示 a(j),b(j)。 由于误差的均方值是一标量,因此 ( 2.2.5) 式是一个标量对复函数的求导问题,它等价于
0]|)([|]|)([|
22

jj b
neEj
a
neE j=0,1,2,…
( 2.2.6)

jj
j bja?

j=0,1,2,…
( 2.2.7)
第二章 维纳滤波和卡尔曼滤波则( 2.2.6)式可以写为
0]|)([| 2 neEj ( 2.2.8)
将( 2.2.8)式展开





)()()()()()()()(]|)([| ****2 nje
b
nenje
b
nene
a
nene
a
neEneE
jjjj
j
( 2.2.9)
又根据( 2.2.1) ~(2.2.3)式第二章 维纳滤波和卡尔曼滤波
)(
)(
)(
)(
)(
)(
)(
)(
*
*
*
*
jnjx
b
ne
jnx
a
ne
jnjx
b
ne
jnx
a
ne
j
j
j
j




将( 2.2.10)~( 2.2.13)式代入( 2.2.9)式,得
)]()([2]|)([| *2 nejnxEneEj
( 2.2.14)
第二章 维纳滤波和卡尔曼滤波因此
E[ x*(n-j)e(n)] =0 j=0,1,2,… ( 2.2.15)
上式说明,均方误差达到最小值的充要条件是误差信号与任一进入估计的输入信号正交,这就是通常所说的正交性原理 。 它的重要意义在于提供了一个数学方法,用以判断线性滤波系统是否工作于最佳状态 。
第二章 维纳滤波和卡尔曼滤波下面计算输出信号与误差信号的互相关函数




0
*
0
**
)]()([)(
])()()([)]()([
j
j
nejnxEjh
nejnxjhEnenyE
(2.2.16)
假定滤波器工作于最佳状态,滤波器的输出 yopt(n)与期望信号 d(n)
的误差为 eopt(n),把 ( 2.2.15) 式代入上式,得到
0)]()([ *opt?nenyE opt
(2.2.17)
第二章 维纳滤波和卡尔曼滤波图 2.2.1 期望信号,估计值与误差信号的几何关系
e
o p t
( n )
d ( n )
y
o p t
( n )
第二章 维纳滤波和卡尔曼滤波图 2.2.1表明在滤波器处于最佳工作状态时,估计值加上估计偏差等于期望信号,即
)(e)()( o p to p t nnynd
注意我们所研究的是随机信号,图 2.2.1中各矢量的几何表示应理解为相应量的统计平均或者是数学期望 。 再从能量的角度来看,假定输入信号和期望信号都是零均值,应用正交性原理,则,因此在滤波器处于最佳状态时,
估计值的能量总是小于等于期望信号的能量 。
]|[| 2o p t22d o p t eEy
第二章 维纳滤波和卡尔曼滤波
2.2.2 维纳 —霍夫方程将( 2.2.15)式展开,可以得到
0)()()()(
0
***



m
mnxmhndknxE
将输入信号分配进去,得到
)()()(
0
* kmrmhkr
m
xxdx

k=0,1,2,…
对上式两边取共轭,利用相关函数的性质,ryx(-k)=r*xy(k),得到
)()()()()(
0
krkhmkrmhkr xx
m
xxxd

k=0,1,2,…
(2.2.20)
第二章 维纳滤波和卡尔曼滤波
(2.2.20)式称为维纳 -霍夫 ( Wiener- Hopf) 方程 。 当 h(n)是一个长度为 M的因果序列 (即 h(n)是一个长度为 M的 FIR滤波器 )
时,维纳 -霍夫方程表述为
)()()()()(
1
0
krkhmkrmhkr xx
M
m
xxxd
k=0,1,2,…
(2.2.21) 把 k的取值代入 (2.2.21)式,得到当 k=0时,h1rxx(0)+h2rxx(1)+…+hMrxx(M-1)=rxd(0)
当 k=1时,h1rxx(1)+ h2rxx(0)+…+ hMrxx(M-2)= rxd(+1)
当 k=M-1时,h1rxx(M-1)+ h2rxx (M-2)+…+hMrxx(0)= rxd(M-1)

(2.2.22)
第二章 维纳滤波和卡尔曼滤波定义

)0()2()1(
)2()0()0(
)1()1()0(
)1(
)1(
)0(
2
1
xxxxxx
xxxxxx
xxxxxx
xx
xd
xd
xd
xd
M
rMrMr
Mrrr
Mrrr
R
Mr
r
r
R
h
h
h
h


第二章 维纳滤波和卡尔曼滤波
( 2.2.22)式可以写成矩阵的形式,即
hRR xxxd? (2.2.23)
对上式求逆,得到
xdxx RRh 1
(2.2.24)
第二章 维纳滤波和卡尔曼滤波上式表明已知期望信号与观测数据的互相关函数及观测数据的自相关函数时,可以通过矩阵求逆运算,得到维纳滤波器的最佳解 。 同时可以看到,直接从时域求解因果的维纳滤波器,当选择的滤波器的长度 M较大时,计算工作量很大,
并且需要计算 Rxx的逆矩阵,从而要求的存贮量也很大 。 此外,
在具体实现时,滤波器的长度是由实验来确定的,如果想通过增加长度提高逼近的精度,就需要在新 M基础上重新进行计算 。 因此,从时域求解维纳滤波器,并不是一个有效的方法 。
第二章 维纳滤波和卡尔曼滤波
2.2.3
假定所研究的信号都是零均值的,滤波器为 FIR型,长度等于 M,将 (2.2.2)式和 (2.2.3)式代入 (2.2.4)式,可以得到
)]()([)1()()]()([)(
)]()([)(]|)([|]|)([|
*
1
0
1
0
**
1
0
*
1
0
*22
inxknxEhkhndknxEkh
ndknxEkhndEneE
M
k
M
i
M
k
M
k



(2.2.25)
第二章 维纳滤波和卡尔曼滤波上式可以进一步化简得到
)(])[()(
*)()()(
)()()()()()()(]|)([|
1T*11T*2
TT*T*2
1
0
1
0
*
1
0
*
1
0
*22
xdxxxxxdxxxdxxxdd
xxxdxdd
M
k
xx
M
i
M
k
xdxd
M
k
d
RRhRRRhRRR
hRhhRRh
kirihkhkrkhkrkhneE




可以看出,均方误差与滤波器的单位脉冲响应是一个二次函数关系 。 由于单位脉冲响应 h (n)为 M维向量,因此均方误差是一个超椭圆抛物形曲面,该曲面有极小点存在 。 当滤波器工作于最佳状态时,均方误差取得最小值 。
第二章 维纳滤波和卡尔曼滤波将( 2.2.24)式代入 (2.2.26)式,得到最小均方误差
optxddxdxxxdd hRRRRneE T*21T*2m i n2 )()(]|)([|
(2.2.27)
例 2.2.1 设 y(n)=x(n)+v2(n),v2(n)是一白噪声,方差 σ22=0.1。
期望信号 x1(n)的信号模型如图 2.2.2( a) 所示,其中白噪声 v1(n)
的方差 σ21=0.27,且 b0=0.8458。 x(n)的信号模型如图 2.2.2( b) 所示,b1=0.9458。 假定 v1(n)与 v2(n),x1(n)与 y(n)不相关,并都是实信号 。 设计一个维纳滤波器,得到该信号的最佳估计,要求滤波器是一长度为 2的 FIR滤波器 。
第二章 维纳滤波和卡尔曼滤波图 2.2.2 输入信号与观测数据的模型
z - 1
x
1
( n )v
1
( n )
b
0
z - 1
x ( n )x
1
( n )
b
1
y ( n )
v
2
( n )
( a ) ( b )
第二章 维纳滤波和卡尔曼滤波解 这个问题属于直接应用维纳 -霍夫方程的典型问题,其关键在于求出观测信号的自相关函数和观测信号与期望信号的互相关函数 。
图 2.2.3 维纳滤波器的框图
H
1
( z ) H
2
( z ) +
v
1
( n ) x
1
( n ) x ( n )
y ( n )
v
2
( n )


第二章 维纳滤波和卡尔曼滤波根据题意,画出这个维纳滤波器的框图,如图 2.2.3所示 。
用 H1(z)和 H2(z)分别表示 x1(n)和 x(n)的信号模型,那么滤波器的输入信号 x(n)可以看作是 v1(n)通过 H1(z)和 H2 (z)级联后的输出,
H1(z)和 H2 (z)级联后的等效系统用 H(z)表示,输出信号 y(n)就等于 x(n)和 v2(n)之和 。 因此求出输出信号的自相关函数矩阵 Ryy和输出信号与期望信号的互相关矩阵 Ryd是解决问题的关键 。 相关函数矩阵由相关函数值组成,已知 x(n)与 v2(n)不相关,那么
)()()( 22 mrmrmr vvxxyy
第二章 维纳滤波和卡尔曼滤波
(1) 求出期望信号的方差 。 根据图 2.2.2(a),期望信号的时间序列模型所对应的差分方程为
x1(n)=v1(n)-b0x1(n-1)
这里,b0=0.8458,由于 x1(n)的均值为零,其方差与自相关函数在零点的值相等 。
22
0
2
1
2
1
2
1110
2
1
2
1
2
1
11
)]1()1()(2)([)]([)0(
x
xx
b
nxbnxnvbnvEnxER



9 4 8 6.0
)8 4 5 8.0(1
27.0
1 220
2
122
1

bxd

第二章 维纳滤波和卡尔曼滤波
(2) 计算输入信号和输出信号的自相关函数矩阵 。 根据自相关函数,功率谱密度和时间序列信号模型的等价关系,已知时间序列信号模型,就可以求出自相关函数 。 这里,信号的模型 H(z)可以通过计算得到 。
)9 4 5 8.01)(8 4 5 8.01(
1)()()(
1121 zzzHzHzH
这是一个二阶系统,所对应的差分方程为
x(n)+a1x(n-1)+a2x(n-2)=v1(n)
式中,a1=-0.1,a2=-0.8。 由于 v1(n),v2(n)的均值为零,因此,
x(n)的均值为 0。 给方程两边同乘以 x*(n-m),并取数学期望,得到 rxx(m)+a1rxx(m-1)+a2rxx (m-2)=0 m> 0 (1)
rxx(0)+ a1rxx(1)+a2rxx(2)=σ21 m=0 (2)
第二章 维纳滤波和卡尔曼滤波对方程( 1)取 m=1,2,得到
rxx(1)+a1rxx(0)+a2rxx(1)=0 (3)
rxx(2)+a1rxx(1)+a2rxx(0)=0 (4)
方程( 2)、( 3)、( 4)联立求解,得
5.0
8.01
1.0
1
)1(
1
])1.0()8.01[(
27.0
8.01
8.01
])1[(1
1
)0(
2
1
222
1
2
2
2
1
2
22




a
a
r
aaa
a
r
xx
xxx
至此,输入信号的自相关矩阵 Rxx可以写出:
15.0
5.01
)0()1(
)1()0(
xxxx
xxxx
xx rr
rr
R
第二章 维纳滤波和卡尔曼滤波
v2(n)是一个零均值的白噪声,它的自相关函数矩阵呈对角形,且,2
2)0(22vvr
1.00
01.0
)0()1(
)1()0(
2222
2222
22
vvvv
vvvv
vv rr
rr
R
因此,输出信号的自相关 Ryy为


1.15.0
5.01.1
)0()0()1(
)1()1()0(
)0()1(
)1()0(
2222
2222
vvxxvvxx
vvxxvvxx
yyyy
yyyy
yy rrrr
rrrr
rr
rr
R
第二章 维纳滤波和卡尔曼滤波
(3) 计算输出信号与期望信号的互相关函数矩阵。 由于两个信号都是实信号,故
ryd(m)=E[ y(n)d(n-m)] =E[ y(n)x1(n-m)
=E[ (x(n)+v2(n))x1(n-m)] =E[ x(n)x1(n-m)] m=0,1
根据图 2.2.2系统 H2(z)的输入与输出的关系,有
x1(n)-b1x(n-1)=x(n)
推出 x
1(n)=x(n)+b1x(n-1)
这样
ryd(m) =E[ x(n)x1(n-m)] =E[ x(n)(x(n-m)+b1x(n-1-m))]
=rxx(m)+b1rxx(m-1)
第二章 维纳滤波和卡尔曼滤波将 m=0,m=1代入上式,得
ryd(0)=rxx(0)+b1rxx(-1)=1-0.9458× 0.5=0.5272
ryd(1)=rxx(1)+ b1rxx(0)=0.5-0.9458× 1=-0.4458
因此,输出信号与期望信号的互相关 Ryd为




4458.0
5272.0
)1(
)0(
yd
yd
yd r
r
R
求出输出信号自相关的逆矩阵,并乘以 Ryd,就可以得到维纳滤波器的最佳解 Wopt:

1 4 5 6.15 2 0 8.0
5 2 0 8.01 4 5 6.1
)0()1(
)1()0(
)1()0(
1
)0()1(
)1()0(
22
1
1
rr
rr
rrrr
rr
R
yy
第二章 维纳滤波和卡尔曼滤波






7 8 5 3.0
8 3 6 0.0
4 5 8.0
5 2 7 2.0
1 4 5 6.15 4 2 8.0
5 2 0 8.01 4 5 6.11
ydyyopt RRW
把 Wopt代入 ( 2.2.27) 式,可以计算出该维纳滤波达到最佳状态时均方误差,即取得了最小值 E[ |e(n)|2] min,
1 5 7 9.0
7 5 8 3.0
8 3 6 0.0
]4 4 5 8.05 2 7 2.0[9 4 8 6.0)(]|)([| o p tT*2m i n2
WRneE ydd?
第二章 维纳滤波和卡尔曼滤波
2.3 离散维纳滤波器的 z域解若不考虑滤波器的因果性,( 2.2.20)式可以写为
)()()()()( krkhmkrmhkr xx
m
xxxd


设定 d(n)=s(n),对上式两边做 Z变换,得到
Sxs(z)=Hopt(z)Sxx(z)
)(
)()(
zS
zSzH
xx
xs
o p t?
第二章 维纳滤波和卡尔曼滤波假设信号和噪声不相关,即 rsv(m)=0,则
Sxs(z)=Sss(z)
Sxx(z)=Sss(z)+Svv(z)
( 2.3.2)式可以写成
)()(
)(
)(
)()(
zSzS
zS
zS
zSzH
vvss
ss
xx
xs
o p t
( 2.3.5) 式表示,当噪声为 0时,信号全部通过;当信号为 0时,
噪声全部被抑制掉,因此维纳滤波确有滤除噪声的能力 。 把信号的频谱用 Pss(ejω)表示,噪声的频谱用 Pvv(ejω)表示,那么非因果的维纳滤波器的传输函数 Hopt(ejω)的幅频特性如图 2.3.1所示 。
第二章 维纳滤波和卡尔曼滤波

0
1
1
)e( j ωo ptH
Pss(ejω)≠0,Pvv(ejω)=0
Pss(ejω)≠0,Pvv(ejω) ≠ 0
Pss(ejω)=0,Pvv(ejω) ≠ 0
然而实际的系统都是因果的 。 对于一个因果系统,不能直接转入频域求解的原因是由于输入信号与期望信号的互相关序列是一个因果序列,如果能够把因果维纳滤波器的求解问题转化为非因果问题,求解方法将大大简化 。 那么怎样把一个因果序列转化为一个非因果序列呢?
第二章 维纳滤波和卡尔曼滤波图 2.3.1 非因果维纳滤波器的传输函数的幅频特性
H
o p t
(e
j?
)
P
S S
(e
j?
)
P
v v
(e
j?
)
0?
第二章 维纳滤波和卡尔曼滤波回顾前面讲到的时间序列信号模型,假设 x(n)的信号模型
B(z)已知 ( 如图 2.3.2(a)所示 ),求出信号模型的逆系统 B-1(z),
并将 x(n)作为输入,那么逆系统 B-1(z)的输出 ω(n)为白噪声 。 一般把信号转化为白噪声的过程称为白化,对应的滤波器称为白化滤波器 ( 如图 2.3.2(b)所示 ) 。
图 2.3.2 x(n)的时间序列信号模型及其白化滤波器
B ( z )? ( n ) x ( n ) B - 1 ( z )? ( n )x ( n )
第二章 维纳滤波和卡尔曼滤波具体思路如图 2.3.3所示 。 用白噪声作为待求的维纳滤波器的输入,设定 1/B(z)为信号 x(n)的白化滤波器的传输函数,那么维纳滤波器的传输函数 G(z)的关系为
)(
)()(
zB
zGzH? (2.3.7)
因此,维纳滤波器的传输函数 H(z)的求解转化为 G(z)的求解。
图 2.3.3 维纳滤波解题思路
( n )
x ( n ) G ( z ) y ( n ) = s ( n ))(1 zB ^
第二章 维纳滤波和卡尔曼滤波
2.3.1
假设待求维纳滤波器的单位脉冲响应为 ω(n),期望信号
d(n)=s(n),系统的输出信号 y(n)=s(n),g(n)是 G(z)的逆 Z变换,
如图 2.3.3所示 。



k
knkgngnnsny )()()()()(?)(
第二章 维纳滤波和卡尔曼滤波





























k
s
k
s
k k
s
k
s
k k
kk
k
rkr
kgr
krkgkrkgkgr
nsknkgnsknkgE
rnknrgkgEnsE
knkgnsEneE
2
2
2
ss
**22
ss
***
**2
2
2
||)(
)()0(
)()()()(|)(|)0(
)()()()()()(
)]()()()([]|)([|
)()()(]|)([|




(2.3.9)
第二章 维纳滤波和卡尔曼滤波可以看出,均方误差的第一项和第三项都是非负数,要使均方误差为最小,当且仅当
0)()(

krkg s -∞< k< ∞
(2.3.10)
因此 g(n)的最佳值为
2
)()(
krkg s
o p t?
-∞< k< ∞
(2.3.11)
对上式两边同时做 Z变换,得到
2
)()(
zSzG s
opt?
(2.3.12)
第二章 维纳滤波和卡尔曼滤波这样,非因果维纳滤波器的最佳解为
)(
)(1
)(
)()(
2
o p t
o p t zB
zS
zB
zGzH s?


(2.3.13)
因为 s(n)=s(n)*δ(n),且 x(n)=ω(n)*b(n),根据相关卷积定理
(1.4.15)式,得到
rxs(m)=rωs(m)*b(-m) (2.3.14)
对上式两边做 Z变换,得到
Sxs (z)=Sωs(z)B(z-1)
因此
)(
)()(
1 zB
zSzS xs
s?
(2.3.15)
第二章 维纳滤波和卡尔曼滤波将上式代入 (2.3.13)式,并根据 x(n)的信号模型,得到非因果的维纳滤波器的复频域最佳解的一般表达式
)(
)(
)(
)(
)(
11)(
12opt zS
zS
zB
zS
zBzH xx
xsxs

(2.3.16)
假定信号与噪声不相关,即当 E[ s(n)v(n)] =0时,有
rxs(m)=E[ (s(n)+v(n))*s(n+m)] =rss(m)
rxx(m)=E[ (s(n)+v(n))*(s(n+m)+v(n+m))] =rss(m)+rvv(m)
对上边两式做 Z变换,得到
Sxs(z)=Sss(z)
Sxx(z)=Sss(z)+Svv(z)
(2.3.17)
(2.3.18)
第二章 维纳滤波和卡尔曼滤波把 (2.3.17)式代入 (2.3.15)式,得到
)(
)()(
1 zB
zSzS ss
s?
( 2.3.19)
将 ( 2.3.18) 式和 ( 2.3.19) 式代入 ( 2.3.16) 式,得到信号和噪声不相关时,非因果维纳滤波器的复频域最佳解和频率响应分别为
)()(
)(
)(
)()(
o p t zSzS
zS
zS
zSzH
vvss
ss
xx
xs

)()(
)(
)e()e(
)e()e(
j ωj ω
j ω
j ω
o p t
vvss
ss
vvss
ss
PP
P
SS
SH

( 2.3.20)
( 2.3.21)
第二章 维纳滤波和卡尔曼滤波下面我们推出该滤波器的最小均方误差 E[ |e(n)|2] min的计算,重新写出 (2.3.9)式的最佳解



k
s
ss
krrneE
2
2
mi n
2 |)(|)0(]|)([|
根据围线积分法求逆 Z变换的公式,rss(m)用下式表示,
C mssss zzzSmr d)(πj2 1)( 1
(2.3.22)
得出
C ssss zzzSr d)(πj2 1)0( 1
(2.3.23)
第二章 维纳滤波和卡尔曼滤波由复卷积定理
z
z
zYzXnynx Cn
d1)(
πj2
1)()(
*
**



(2.3.24)
取 y(n)=x(n),有
zzzXzXnx
Cn
d)(
πj2
1|)(| 12

(2.3.25)
因此
z
zzSzSkr
C ssn s
d)()(
πj2
1|)(| 12


(2.3.26)
把 (2.3.23)式和 (2.3.26)式代入 (2.3.9)式,得到
z
zzSzSzSneE
sssC
d)()(1)(
πj2
1]|)([| 1
2m i n
2?



(2.3.27)
第二章 维纳滤波和卡尔曼滤波将( 2.3.19)式代入上式,得到

z
z
zSzHzS
z
z
zB
zS
zB
zS
zSneE
sso p tss
C
sss
ss
C
d
)()()(
πj2
1
d
)(
)(
)(
)(1
)(
πj2
1
]|)([|
1
1
12m i n
2


(2.3.28)
因为实信号的自相关函数是偶函数,即 rss(m)=rss(-m),因此
Sss(z)=Sss(z-1)
假定信号与噪声不相关,E[ s(n)v(n)] =0,则

C
xx
vvss
C
ss
xx
ss
ss
z
z
zS
zSzS
z
z
zS
zS
zS
zSneE
d
)(
)()(
πj2
1
d
)(
)(
)(
)(
πj2
1
]|)([|
1
m i n
2
(2.3.30)
第二章 维纳滤波和卡尔曼滤波
2.3.2
若维纳滤波器是一个因果滤波器,要求
g(n)=0 n< 0
则滤波器的输出信号


0
)()()()()(?)(
k
knkgngnnsny
估计误差的均方值
E[ |e(n)|2] =E[ |s(n)-y(n)|2]
(2.3.32)
(2.3.31)
类似于( 2.3.9)式的推导,得到




0
2
2
2
0
2 |)(|1)()()0(]|)([|
k
s
k
s
ss kr
krkgrneE


(2.3.33)
第二章 维纳滤波和卡尔曼滤波要使均方误差取得最小值,当且仅当
)(
)(
00
0
)(
)(
2
2
opt
nu
nr
n
n
nr
ng
s
s

(2.3.34)






)]([
1
)]([ZT)(
])([])()([)]([
2opt
0
zSngzG
znrznunrzS
sopt
n
n
s
n
n
ss

(2.3.35)
(2.3.36)
第二章 维纳滤波和卡尔曼滤波又由( 2.3.15)式得到



)(
)(1)(
12o p t zB
zSzG xs

(2.3.37)
所以因果维纳滤波器的复频域最佳解为



)(
)(
)(
11
)(
)()(
12
o p t
zB
zS
zBzB
zGzH xs
o p t

(2.3.38)
第二章 维纳滤波和卡尔曼滤波维纳滤波的最小均方误差为

z
z
zSzHzS
z
z
zB
zS
zB
zS
zS
z
z
zSzSr
krkukrr
kr
rneE
xsss
C
xsxs
ss
C
s
C
sss
k
ssss
k
s
ss
d
)()()(
πj2
1
d
)(
)(
)(
)(1
)(
πj2
1
d
)()]([
1
πj2
1
)0(
)()]()([
1
)0(
|)(|
)0(]|)([|
1
opt
1
12
1
2
*
2
0
2
2
m i n
2










(2.3.39)
第二章 维纳滤波和卡尔曼滤波比较 ( 2.3.28) 式和 ( 2.3.39) 式,可以看出因果维纳滤波器的最小均方误差与非因果维纳滤波器的最小均方误差的形式相同,但公式中的 Hopt(z)的表达式不同,(2.3.16)
式和 (2.3.38)式 。 前面已经导出,对于非因果情况,



k
s
ss
krrneE
2
2
mi n
2 |)(|)0(]|)([|
对于因果情况,


0
2
2
mi n
2 |)(|)0(]|)([|
k
s
ss
krrneE
比较两式,它们的第二项求和域不同,因为因果情况下,k=0~+∞,
因此可以说明非因果情况的 E[ |e(n)|2] min一定小于等于因果情况 E[ |e(n)|2] min。 在具体计算时,可以选择单位圆作为积分曲线,
应用留数定理,计算积分函数在单位圆内的极点的留数来得到 。
第二章 维纳滤波和卡尔曼滤波通过前面的分析,因果维纳滤波器设计的一般方法可以按下面的步骤进行:
(1) 根据观测信号 x(n)的功率谱求出它所对应的信号模型的传输函数,即采用谱分解的方法得到 B(z)。 具体方法为
Sxx(z)=σ 2ωB(z)B(z-1),把单位圆内的零极点分配给 B(z),单位圆外的零极点分配给 B(z-1),系数分配给 σ 2ω。
(2)求 的 Z反变换,取其因果部分再做 Z变换,即舍掉单位圆外的极点,得
(3) 积分曲线取单位圆,应用 ( 2.3.38) 式和 ( 2.3.39)
式,计算 Hopt(z),E[ |e(n)|2] min。


)(
)(
1zB
zS xs


)(
)(
1zB
zS xs 。
第二章 维纳滤波和卡尔曼滤波例 2.3.1 已知
)8.01)(8.01(
36.0)(
1 zzzS ss
信号和噪声不相关,即 rsv(m)=0,噪声 v(n)是零均值,单位功率的白噪声 ( σ 2v=1,mv=0),求 Hopt(z)和 E[ |e(n)|2] min。
解 根据白噪声的特点得出 Svv(z)=1,由噪声和信号不相关,
得到 rxx(m)=rss(m)+rvv(m)。
对上式两边做 Z变换,并代入已知条件,对 x(n)进行功率谱分解,
)()(
)8.01)(8.01(
)5.01)(5.01(6.1
1
)8.01)(8.01(
36.0
)()()(
12
1
1
1




zBzB
zz
zz
zz
zSzSzS
vvssxx
第二章 维纳滤波和卡尔曼滤波考虑到 B(z)必须是因果稳定的系统,得到
6.1,8.01 5.01)( 21
1

z
zzB
(1) 首先分析物理可实现情况,应用公式( 2.3.38),








)5.01)(8.01(
36.0
)5.01(6.1
8.01
)(
)(
)(
11)(
11
1
12opt zzz
z
zB
zS
zBzH
xs





)5.01)(8.01(
36.0)(,
)5.01)(8.01(
36.0)(
11
zz
zF
zz
zF
F(z)的极点为 0.8和 2,考虑到因果性,稳定性,仅取单位圆内的极点 zi=0.8,f(n)为 F(z)的 Z反变换 。 用 Res表示留数,应用留数定理,有第二章 维纳滤波和卡尔曼滤波
n
z
n
n
z
zz
z
z
zz
snf
8.06.0)8.0(
)5.01)(8.01(
36.0
8.0,
)5.01)(8.01(
36.0
Re)(
8.01
1
1
1




取因果部分,
f+(n)=0.6× 0.8n× u(n)
111
1
o p t
1
5.01
1
8
3
8.01
6.0
)5.01(6.1
8.01
)(
8.01
6.0
)](8.06.0[)]([)(





zzz
z
zH
z
nuZTnfZTzF
n
第二章 维纳滤波和卡尔曼滤波
z
z
zz
z
z
zzzzz
z
z
zSzHzS
neE
C
C
C
xsss
d
)8.01)(8.01(
36.0
πj2
1
d
)8.01)(8.01(
36.0
5.01
8
3
)8.01)(8.01(
36.0
πj2
1
d
))()()((
πj2
1
]|)([|
1
111
1
opt
m i n
2






zzz
zE 1
)5.01)(8.01(
8
536.0
)( 1?


第二章 维纳滤波和卡尔曼滤波单位圆内只有极点 Zi=0.5,
8
3
)5.0(
1
)5.01)(8.01(
36.0
]5.0),([R e s]|)([|
5.01m i n
2



z
z
zzz
zEneE
未经滤波器的均方误差
1]|)([||)()([|]|)([| 2222 vnvEnsnxEneE?
第二章 维纳滤波和卡尔曼滤波
(2) 对于非物理可实现情况,应用公式 (2.3.20)和 (2.3.28),有
)5.01)(5.01(
2 2 5.0
)()(
)(
)(
)()(
1opt zzzSzS
zS
zS
zSzH
vvss
ss
xx
xs

z
z
zzzz
zz
z
z
zzzzzz
z
z
zSzHzSneE
C
C
xsss
d
)5.01)(5.01)(8.01)(8.01(
)5.05.0025.1(36.0
πj2
1
d
)8.01)(8.01(
36.0
)5.01)(5.01(
225.0
)8.01)(8.01(
36.0
πj2
1
d
))()()((
πj2
1
]|)([|
11
1
111
1
optm i n
2








第二章 维纳滤波和卡尔曼滤波令
)5.01)(5.01)(8.01)(8.01(
)5.05.00 2 5.1(36.0)(
11
1
zzzz
zzzF



单位圆内有两个极点 0.8和 0.5,应用留数定理,有
10
3]5.0),([R e s]8.0),([R e s)]([
m i n
2 zFzFneE
比较两种情况下的最小均方误差,可以看出非物理可实现情况的最小均方误差小于物理可实现情况的均方误差 。
第二章 维纳滤波和卡尔曼滤波
2.4 维 纳 预 测
2.4.1 维纳预测的计算在维纳滤波中,期望的输出信号 yd(n)=s(n),实际的输出为
y(n)=s(n)。 在维纳预测中,期望的输出信号 yd(n)=s(n+N),
实际的输出 y(n)=s(n+N)。 前面已经推导得到维纳滤波的最佳解为
^
^
)(
)(
)(
)()(
o p t zS
zS
zS
zSzH
xx
xy
xx
xs d
( 2.4.1)
其中,Sxx(z)是观测数据的功率谱 ;Sxyd(z)是观测数据与期望信号的互功率谱,即互相关函数 rxyd(k)的傅里叶变换
)]()([ *)( knynxEr dkxy d ( 2.4.2)
第二章 维纳滤波和卡尔曼滤波对应于维纳预测器,其输出信号 y(n)和预测误差信号
e(n+N)分别为
)(?)()(
)()()(?)(
0
NnsNnsNne
mNnxmhNnsny
m



( 2.4.3)
( 2.4.4)
同理,要使预测误差的均方值为最小,须满足
0]|)([|
2

kh
NneE ( 2.4.5)
其中,hk表示 h(k)。
第二章 维纳滤波和卡尔曼滤波观测数据与期望的输出的互相关函数 rxyd(k)和互谱密度
Sxyd(z)分别为
N
xsxy
xsdxy
zzSzS
kNrkNnsnxEknynxEkr
d
d
)()(
)()]()([)]()([)( **
( 2.4.6)
( 2.4.7)
这样,非因果维纳预测器的最佳解为
)(
)(
)(
)()(
opt zS
zSz
zS
zSzH
xx
xs
N
xx
xy d ( 2.4.8)
因果维纳预测器的最佳解为






)(
)(
)(
11
)(
)(
)(
11)(
1212opt zB
zSz
zBzB
zS
zBzH
xs
N
xy d

( 2.4.9)
第二章 维纳滤波和卡尔曼滤波维纳预测的最小均方误差为


C
xsss
C
xyss
z
dz
zSzHzS
z
dz
zSzHzSNneE
d
)]()()([
πj2
1
)]()()([
πj2
1
]|)([|
1
opt
1
optm i n
2
从上面分析可以看出,维纳预测的求解和维纳滤波器的求解方法是一致的。
第二章 维纳滤波和卡尔曼滤波
2.4.2 纯预测假设 x(n)=s(n)+v(n),式中 v(n)是噪声,且 v(n)=0,期望信号为 s(n+N),N> 0,此种情况称为纯预测 。
假定维纳预测器是因果的,仍设 s(n)与 v(n)不相关,纯预测情况下的输入信号的功率谱及维纳预测器的最佳解分别为


)]([
)(
1
)(
)(
)(
11
)(
)()()()()(
12opt
12
zBz
zBzB
zSz
zB
zH
zBzBzSzSzS
Nxs
N
ssxsxx
( 2.4.11)
( 2.4.12)
第二章 维纳滤波和卡尔曼滤波纯预测器的最小均方误差为







C
NN
C
N
N
C
N
xsss
z
dz
zBzzzBzBzB
z
dz
zzBzB
zB
zBz
zBzB
z
dz
zzSzHzS
NneE
)]([)()()(
πj2
)()(
)(
)]([
)()(
πj2
1
))()()((
πj2
1
]|)([|
11
2
1212
1
opt
m i n
2


( 2.4.13)
应用复卷积定理
z
z
zYzXnynx Cn
d1)(
πj2
1)()(
*
**


( 2.4.14)
第二章 维纳滤波和卡尔曼滤波取 y(n)=x(n)
z
zzXzXnx
Cn
d)()(
πj2
1)( 12

( 2.4.15)
将上式代入 (2.4.13)式,并考虑到 b(n)是因果系统,得到
)(
})()({
)()]()([)(]|)([|
1
0
22
0 0
222
22
m i n
2
nb
Nnbnb
NnbnuNnbnbNneE
N
n
n n
n n










可以看到,随着 N增加,E[ |e(n+N)|2] min也增加 。 这一点也容易理解,当预测的距离越远,预测的效果越差,偏差越大,因而 E
[ |e(n+N)|2] min越大 。
第二章 维纳滤波和卡尔曼滤波例 2.4.1 已知
)1)(1(
1)(),()(
1
2
azaz
azSnsnx
xx

其中 -1< a< 1,求,
(1) 最小均方误差下的 s(n+N);
(2) E[ |e(n+N)|2] min。
^
解 首先对 Sxx(z)进行功率谱分解。因为
1
22
1
2
1
1
)(,1
)1)(1(
1
)(


az
zBa
azaz
a
zS xx
所以第二章 维纳滤波和卡尔曼滤波其次,求出 B(z)的 Z反变换
)(1 1 T1 1 nuaaz nZ
然后,应用 Z变换的性质,得到
1
ZTT
1 1)()(1
1


az
anuaNnua
az
z NNnNnZN 取因果部分
N
N
N a
az
aazzBz
zBzH
1
1
o p t 1)1()]([)(
1)(
( 2.4.17)
第二章 维纳滤波和卡尔曼滤波图 2.4.1 纯预测维纳滤波器
a Nx ( n ) y ( n )
第二章 维纳滤波和卡尔曼滤波由 Hopt(z)=aN,此时可以把纯预测的维纳滤波器看作是一个线性比例放大器 ( 如图 2.4.1所示 ) 。 根据 x(n)的信号模型
11
1)(
azzB
可以写出 x(n)的时间序列模型所对应的输入输出方程
x(n)=ω (n)+ax(n-1)
将信号 x(n)通过纯预测维纳滤波器,随着时间的递增,可以得到当 N=1时,x(n+1)=ax(n)=as(n)
当 N=2时,x(n+2)=ax(n+1)=a2s(n)
当 N=N时,x(n+N)=ax(N+n-1)=aNs(n)
… ( 2.4.19)
第二章 维纳滤波和卡尔曼滤波以上推导结果相当于在 n+N时刻,ω(n+N)=0,即去掉噪声时的结果 。 设 N> 0时,ω (n+N)=0,则
x(n+N)=ax(n+N-1)
此时,从统计意义上讲,当 N> 0时,白噪声信号 ω(n+N)对 x(n)
无影响 。
这一结论还可以推广,对于任何均值为零的 x(n),要估计 s
(n+N)时,只需要考虑 B(z)的惯性,即可认为 ω(n+N)=0,N> 0,
这样估计出来的结果将有最小均方误差 。
终值定理
^
2
1 )(l im)()1(l im xxxmxxz mmrzSz
第二章 维纳滤波和卡尔曼滤波表明一个信号的功率谱在单位圆上没有极点与信号均值等于 0等价,因此对于功率谱在单位圆上没有极点的信号,要估计 s(n+N)时,可认为 ω(n+N)=0,N> 0,即仅需要考虑 B(z)的惯性,这样估计出来的结果将有最小均方误差 。
^
第二章 维纳滤波和卡尔曼滤波
2.4.3 一步线性预测的时域解已知 x(n-1),x(n-2),…,x(n-p),预测 x(n),假设噪声
v(n)=0,这样的预测称为一步线性预测 。 设定系统的单位脉冲响应为 h(n),根据线性系统的基本理论,输出信号

p
k
knxkhnxny
1
)()()(?)(
令 apk=-h(k),则

p
k
pk knxanx
1
)()(?
( 2.4.21)
预测误差


p
k
pk
p
k
pk knxaknxanxnxnxne
01
)()()()(?)()(
( 2.4.22)
第二章 维纳滤波和卡尔曼滤波其中,ap0=1,

2
1
2 )()(]|)([|
p
k
pk knxanxEneE
( 2.4.23)
要使均方误差为最小值,要求
pla neE
pl
,,2,10]|)([|
2

同维纳滤波的推导过程一样,可以得到
E[ e*(n)x(n-l)] =0 l=1,2,…,p ( 2.4.24)
把 (2.4.22)式代入 (2.4.24)式,得到

p
k
xxpkxx lkralr
1
0)()(
l=1,2,…,p
( 2.4.25)
第二章 维纳滤波和卡尔曼滤波由于预测器的输出 是输入信号的线性组合,参见
( 2.4.21)式,得到
)(? nx
0)](?)([ *?nxneE ( 2.4.26)
(2.4.24)式说明误差信号与输入信号满足正交性原理,
(2.4.26)式说明预测误差与预测的信号值同样满足正交性原理 。
预测误差的最小均方值







p
k
xxpkxx
p
k
pk krarnxknxanxE
nxneEnxnxneEneE
11
**
**
m i n
2
)()0()()()(
)]()([) ) ](?)()(([]|)([|
(2.4.27)
第二章 维纳滤波和卡尔曼滤波将 (2.4.25)式和 (2.4.27)式联立,得到下面的方程组:


p
k
xxpkxx
p
k
xxpkxx
pllkralr
neEkrar
1
1
m i n
2
,,2,10)()(
]|)([|)()0(
( 2.4.28)
将方程组写成矩阵形式
0
0
]|)([|1
)0()1()(
)1()0()1(
)()1()0(
m i n
2
1


neE
a
a
rprpr
prrr
prrr
pp
p
xxxxxx
xxxxxx
xxxxxx
( 2.4.29)
第二章 维纳滤波和卡尔曼滤波这就是有名的 Yule-Walker方程,Yule-Walker
方程具有以下特点:
(1) 除了第一个方程外,其余都是齐次方程 ;
(2) 与维纳 -霍夫方程相比,不需要知道观测数据 x(n)与期望信号 s(n)的互相关函数 。
该方程组有 p+1个方程,对应地,可以确定 apk,k=1,2,…,
p和 E[ e2(n)] min,共计 p+1个未知数,因此可用来求解 AR模型参数 。 这就是后面要介绍的 AR模型法进行功率谱估计的原理,
它再一次揭示了时间序列信号模型,功率谱和自相关函数描述一个随机信号的等价性 。
第二章 维纳滤波和卡尔曼滤波例 2.4.2 已知
)8.01)(8.01(
36.0)(
1 zzzS xx
x(n)为 AR模型,求 AR模型参数 。
解 求解 AR模 型参 数包括 确定 AR模型 的阶 数 p及 系数
ap1,ap2,…,app。
首先对 Sxx(z)做傅里叶反变换,得到 x(n)的自相关函数 rxx(m),
rxx(m)=0.8|m|
采用试验的方法确定模型阶数 p。 首先取 p=2,各相关函数值由上式计算,并代入 (2.4.29)
0
0
1
18.064.0
8.018.0
64.08.01 2
2
1

a
a
第二章 维纳滤波和卡尔曼滤波计算得到
a1=-0.8,a2=0,σ 2ω =0.36
如果取 p=3,可计算出 a1=-0.8,a2=a3=0,σ 2ω=0.36,说明 AR
模型的阶数只能是一阶的 。 采用谱分解的方法,即对 Sxx(z)进行谱分解,得到的模型也是一阶的,其时间序列模型和差分方程为
)1(8.0)()(
8.01
1
)(
1

nxnnx
z
zB
第二章 维纳滤波和卡尔曼滤波
2.5 卡尔曼 (Kalman)滤波卡尔曼滤波是用状态空间法描述系统的,由状态方程和量测方程所组成 。 卡尔曼滤波用前一个状态的估计值和最近一个观测数据来估计状态变量的当前值,并以状态变量的估计值的形式给出 。 卡尔曼滤波具有以下的特点:
(1) 算法是递推的,且状态空间法采用在时域内设计滤波器的方法,因而适用于多维随机过程的估计;离散型卡尔曼算法适用于计算机处理 。
(2) 用递推法计算,不需要知道全部过去的值,用状态方程描述状态变量的动态变化规律,因此信号可以是平稳的,也可以是非平稳的,即卡尔曼滤波适用于非平稳过程 。
(3) 卡尔曼滤波采取的误差准则仍为估计误差的均方值最小 。
第二章 维纳滤波和卡尔曼滤波
2.5.1 卡尔曼滤波的状态方程和量测方程假设某系统 k时刻的状态变量为 xk,状态方程和量测方程
( 也称为输出方程 ) 表示为
kkkk xAx 1
(2.5.1a)
kkkk vxCy
(2.5.1b)
其中,k表示时间,这里指第 k步迭代时,相应信号的取值;输入信号 ωk是一白噪声,输出信号的观测噪声 vk也是一个白噪声,
输入信号到状态变量的支路增益等于 1,即 B=1; A表示状态变量之间的增益矩阵,可以随时间发生变化,用 Ak表示第 k步迭代时,
增益矩阵 A的取值;
第二章 维纳滤波和卡尔曼滤波
C表示状态变量与输出信号之间的增益矩阵,可以随时间变化,
第 k步迭代时,取值用 Ck表示,其信号模型如图 2.5.1所示 。
将状态方程中时间变量 k用 k-1代替,得到的状态方程和量测方程如下所示:
xk=Ak-1xk-1+ω k-1
yk=Ckxk+vk
其中,xk是状态变量 ;ω k-1表示输入信号是白噪声 ; vk是观测噪声 ; yk是观测数据 。
第二章 维纳滤波和卡尔曼滤波图 2.5.1 卡尔曼滤波器的信号模型
z
- 1
A
k - 1
C
k
k - 1
x
k - 1
x
k
v
k
y
k
第二章 维纳滤波和卡尔曼滤波为了后面的推导简单起见,假设状态变量的增益矩阵 A不随时间发生变化,ω k,vk都是均值为零的正态白噪声,方差分别是 Qk和 Rk,并且初始状态与 ω k,vk都不相关,γ 表示相关系数 。

kjkvvkvkk
kjkkkk
RRvEv
QQE
jk
jk




,
2
,
2
,,0][:
,,0][:
其中

jk
jk
kj
0
1
第二章 维纳滤波和卡尔曼滤波
2.5.2 卡尔曼滤波的递推算法卡尔曼滤波是采用递推的算法实现的,其基本思想是先不考虑输入信号 ω k和观测噪声 vk的影响,得到状态变量和输出信号 ( 即观测数据 ) 的估计值,再用输出信号的估计误差加权后校正状态变量的估计值,使状态变量估计误差的均方值最小 。
因此,卡尔曼滤波的关键是计算出加权矩阵的最佳值 。
,状态方程和量测方程为
1
'
1
'
'?

kkkkkk
kkk
xACxCy
xAx
(2.5.4)
(2.5.5)
第二章 维纳滤波和卡尔曼滤波显然,由于不考虑观测噪声的影响,输出信号的估计值与实际值是有误差的,用 表示
kky~
'?~ kkk yyy (2.5.6)
为了提高状态估计的质量,用输出信号的估计误差 来校正状态变量
ky~
)?(?
)?(
11
'
1



kkkkkkk
kkkkkk
xACyHxA
yyHxAx (2.5.7)
其中,Hk为增益矩阵,实质是一加权矩阵 。 经过校正后的状态变量的估计误差及其均方值分别用 和 Pk表示,把未经校正的状态变量的估计误差的均方值用 表示k
x~
'kP
第二章 维纳滤波和卡尔曼滤波
])?)(?[(
])?)(?[(
~~
T'''
T
kkkkk
kkkkk
kkk
xxxxEP
xxxxEP
xxx



(2.5.8)
(2.5.9)
(2.5.10)
卡尔曼滤波要求状态变量的估计误差的均方值 Pk为最小,
因此卡尔曼滤波的关键就是要得到 Pk与 Hk的关系式,即通过选择合适的 Hk,使 Pk取得最小值 。
首先推导状态变量的估计值 和状态变量的估计误差,
然后计算 的均方值 Pk,并通过化简 Pk,得到一组卡尔曼滤波的递推公式 。
kx? kx?
kx?
第二章 维纳滤波和卡尔曼滤波将( 2.5.3),(2.5.5)式代入( 2.5.7)式
kkkkkkkkkkk
kkkkkkkkk
kkkkkkkkk
kkkkkk
vHxACHxACHI
vHxCHxACHI
xACvxCHxA
yyHxAx






)(?)(
)(
)?(
~
)?(
~~
111
1
11
'
1
(2.5.11)
同理,状态变量的估计误差 为x~

kkkkkkkk
kkkkkkkkkk
kkkkkkkkkkkkkkk
kkkkkkkkkkkkkk
kkk
vHxxACHI
vHCHIxxACHI
vHACHxxACHxxA
vHxACHxACHIxA
xxx









111
111
111111
11111
)?()(
)()?()(
)?(?(
})(?){(
~~


(2.5.12)
第二章 维纳滤波和卡尔曼滤波由上式可以看出,状态变量的估计误差 由三部分组成,可记为
x~
cbax~
其中
kk
kkk
kkkkk
vHc
CHIb
xxACHIa



1
11
)(
)?()(
(2.5.13b)
(2.5.13c)
(2.5.13d)
那么,状态变量的估计误差的均方值 Pk就由 9项组成,
][
]))([(]~,~[
TTTTTTTTT
TT
cbcabcbaacabccbbaaE
cbacbaExxEP kkk


(2.5.14a)
第二章 维纳滤波和卡尔曼滤波其中
T
k
T
k
kkk
kkkkk
vHc
CHIb
CHIAxxa



T
TT
1
T
TTT
11
T
)(
)()?(
(2.5.14b)
(2.5.14d)
(2.5.14c)
下面化简 Pk的表达式,根据假设的条件,状态变量的增益矩阵 A不随时间发生变化,起始时刻为 k0,则 ( 2.5.2) 式经过迭代,得到

0
0
0
1
1 )(
kk
l
l
k
kk
k lkAxAx?
令 l=k-k0-j,得到


1
0
0
1
0
0
0
0 )(
kk
j
jkk
k
kk
k jkAxAx?
第二章 维纳滤波和卡尔曼滤波取 k0=0,k=k-1,得到
j
k
j
jkk
k AxAx


2
0
2
0
1
1
(2.5.15)
所以 xk-1仅依赖于 x0,ω 0,ω 1,…,ω k-2,与 ω k-1不相关,即
0][][ T 11T 11 kkkk xexE (2.5.16)
又据 (2.5.7)式和 (2.5.3)式,得
)?( 2111111211 kkkkkkkkkk xACvxCHxAx
(2.5.17)
第二章 维纳滤波和卡尔曼滤波所以 仅依赖于 xk-1,vk-1,而与 vk不相关,即
1kx
0])?([])?[(
0])?([])?[(
T
111
T
111
T
11
T
11




kkkkkk
kkkkkk
xxExxE
xxvEvxxE

(2.5.18)
(2.5.19)
把( 2.5.15)~( 2.5.19)式代入 (2.5.14)式,Pk中的 9项可以分别化简为
TT
1
TTT
1111
T
)()(
)()?)(?()[(][
kkkkkkk
kkkkkkkkkk
CHIAPACHI
CHIAxxxxACHIEaaE



T
1
TT
11
T
)()(
])()[(][
kkkkk
kkkkkk
CHIQCHI
CHICHIEbbE



(2.5.20a)
(2.5.20b)
第二章 维纳滤波和卡尔曼滤波
0])((][
0])()?((][
0])(][
0])()?()(][
0])?()(][
0])?()?()[(][
][][
TT
1
TTT
TTT
11
T
TT
1
T
TTT
111
T
TT
11
T
T
11
T
111
T
TT











kkkkk
kkkkkkk
kkkkk
kkkkkkkk
kkkkkkk
kkkkkkkk
kkkk
T
kkk
CHIvHEcbE
CHIAxxvHEcaE
HvCHIEbcE
CHIAxxCHIEbaE
HvxxACHIEacE
xxxxACHIEabE
HRHHvvHEccE
(2.5.20c)
(2.5.20d)
(2.5.20e)
(2.5.20f)
(2.5.20g)
(2.5.20h)
(2.5.20j)
第二章 维纳滤波和卡尔曼滤波也就是说,Pk仅有其中的三项不为零,化简成
T
kkkkkk
T
kkkkk
T
kkkkkkkk
kk
T
kkkkk
k
HRHCHIQAPACHI
HRHCHIQCHI
CHIAPACHI
ccEbbEaaEP





T
11
T
1
T
1
TTT
)()(
)()(
)()(
][][][
(2.5.21)
第二章 维纳滤波和卡尔曼滤波为了进一步化简 Pk,推导未经误差校正的状态估计误差的均方值 Pk′,由下面推导结果可以看出,Pk′ 是一对称矩阵,满足
Pk′ =(Pk′ )T。

1
T
1
T
11
TT
1111
T
111111
T
111111
T'''
][])?)(?[(
])?(] [ ()?([(
])?)(?[(
])?)(?[(d e f









kkkk
kkkkkkkk
kkkkkkkk
kkkkkkkkkk
kkkkk
QAPA
EAxxxxEA
xxAxxAE
xAxAxAxAE
xxxxEP



(2.5.22)
第二章 维纳滤波和卡尔曼滤波将 (2.5.22)式代入 (2.5.21)式,即把 Pk′ 代入 Pk,
TT'T'''
TT'T'''
TT'
)(
)()(
kkkkkkk
T
kkkkkk
kkkkkkkk
T
kkkkkk
kkkkkkkkk
HRHPCHHCPPCHP
HRHHPCHHCPPCHP
HRHCHIPCHIP



(2.5.23) 其中,是正定阵,记
kTkkk RCPC?'
T' SSRCPC kTkkk
(2.5.24)

'T'TT'T )( kkkkkk PCPCCPU
(2.5.25)
第二章 维纳滤波和卡尔曼滤波
TTTT' )( kkTkkkkk HSSHUHUUUHPP
将上式代入( 2.5.23)式,得
(2.5.26)
将 (2.5.26)式后三项配对
'1T'T'
'T1T1T
T1T'T1T1T
)(
])(][)([
)(])(][)([
kkkkkkkk
kkk
kkkk
PCPCPCCP
PSUSHSUSH
USSUPSUSHSUSHP





(2.5.27)
第二项和第三项均与 Hk无关,第一项为一半正定阵,因此使 Pk最小的 Hk应满足
1T''
1T11T
opt
T
)(
)()(
01)(




kkkk
T
kk
k
RCPCCP
SSUSSUH
SUSH
(2.5.28)
(2.5.29)
第二章 维纳滤波和卡尔曼滤波将 Hopt代入 Pk,得到最小均方误差阵
'''
'1T'T''
)(
)(
kkkkkkk
kkkkkkkkkk
PCHIPCHP
PCRCPCCPPP


将 (2.5.7),(2.5.22),(2.5.29)式和 (2.5.30)式联立,得到一组卡尔曼递推公式






'
1
T
1
'
1T'T'
11
)(
)(
)?(
kkkk
kkkkk
kkkkkkk
kkkkkkkk
PCHIP
QAPAP
RCPCCPH
xACyHxAx
(2.5.30)
(2.5.31a)
(2.5.31b)
(2.5.31c)
(2.5.31d)
第二章 维纳滤波和卡尔曼滤波假设初始条件 Ak,Ck,Qk,Rk,yk,xk-1,Pk-1已知,其中 x0=E[ x0],
P0=var[ x0],那么,递推流程见图 2.5.2。
^ ^
图 2.5.2 卡尔曼滤波递推流程
( 2,5,3 1 c ) 式 ( 2,5,3 1 b ) 式
( 2,5,3 1 a ) 式
( 2,5,3 1 d ) 式P k H k
x
k
P
k
x k - 1,P k - 1
^
^ ′

第二章 维纳滤波和卡尔曼滤波例 2.5.1 已知
0,1)(,)8.01)(8.01( 36.0)( 1 vvvxx mzSzzzS
信号与噪声不相关,yk=xk+vk,求卡尔曼信号模型中的 Ak和 Ck。
解 由 yk=xk+vk知道,Ck=1。
对 Sxx(z)进行谱分解,确定 x(n)的信号模型 B(z),从而确定 Ak。
根据 Sxx(z) =σ 2ωB(z)B(z-1),得出
)()1(8.0)(
8.01
1
)(,36.0 12
nnxnx
z
zB



第二章 维纳滤波和卡尔曼滤波上式与卡尔曼状态方程相比,不同之处在于输入信号 ω(n)
的时间不同,因此将 Sxx(z)改写为
z
z
z
zzS
xx 8.018.0136.0)( 1
1

再对 Sxx(z)进行谱分解,得到
kkk
xx
nnxnx
z
z
zB


8.0
)1()1(8.0)(
8.01
)(
1
1
1
(解毕 )
第二章 维纳滤波和卡尔曼滤波卡尔曼滤波和维纳滤波都是采用均方误差最小的准则来实现信号滤波的,但维纳滤波是在信号进入了稳态后的分析,
卡尔曼滤波是从初始状态采用递推的方法进行滤波 。 对于平稳随机信号,当过渡过程结束以后,卡尔曼滤波与维纳滤波的结果间存在什么关系呢? 下面举一例说明 。
第二章 维纳滤波和卡尔曼滤波例 2.5.2 已知
1]v a r [,0?
,0)(,1)(,
)8.01)(8.01(
36.0
)(
001
1



xPx
zSzS
zz
zS vxvvxx
在 k=0时开始观察 yk,yk=xk+vk,用卡尔曼过滤的计算公式求 xk,
并与维纳过滤的方法进行比较 。
解 (1) 由 x(n)功率谱及量测方程,确定卡尔曼递推算法 。
首先对 Sxx(z)进行功率谱分解,由例 2.5.1的结果,得到卡尔曼滤波的状态方程为
xk=0.8xk-1+ω k-1,确定 Ak=0.8
第二章 维纳滤波和卡尔曼滤波由量测方程 yk=xk+vk,确定 Ck=1,
1)0(2 vvvvk SP?
将参数矩阵 Ak,Ck,Rk代入卡尔曼递推公式( 2.5.30),得到





'
1
'
1''
11
)1(
36.064.0
)1(
)?8.0(?8.0?
kkk
kk
kkk
kkkkk
PHP
PP
PPH
xyHxx (2.5.32a)
(2.5.32b)
(2.5.32c)
(2.5.32d)
第二章 维纳滤波和卡尔曼滤波
(2) 求出卡尔曼滤波的输出 。 由卡尔曼递推公式,以及
,P0=var[ x0] =1,可得到 Pk′,Hk,Pk及 xk( k表示迭代次数 ),迭代流程为:
由具体迭代结果可以看出,原先的增益矩阵 Ak,由于只选择了一个状态变量,变成了加权系数 。 见表 2.5.1。
01x
kkkkk PdxaHbPcP 32.5.2?32.5.232.5.232.5.2 '1?
第二章 维纳滤波和卡尔曼滤波表 2.5.1 Kalman滤波迭代结果第二章 维纳滤波和卡尔曼滤波
(3) 求出卡尔曼滤波的稳态解 。 将 (2.5.32b)式代入方程
(2.5.32d),得到第 5个方程
kkk
k
kkk
kkkkkk
HPP
P
PPP
PPPPHP



)1(
1
)()1(
])1(1[)1(
''
'
2'''
'1'''
(2.5.32e)
将方程 (2.5.32c),(2.5.32e) 代入方程 (2.5.32d),消去
Pk′,可以得到 Pk的递推关系,
Pk=(1-Pk)[ 0.64 Pk-1+0.36] =0.64 Pk-1-0.64 Pk -1 Pk +0.36-0.36 Pk
第二章 维纳滤波和卡尔曼滤波化简上式,得到
1.36Pk+0.64Pk-1Pk=0.64Pk-1+0.36
36.164.0
36.064.0
1
1

k
k
k P
PP
要求的是稳态解,因此将 Pk,Pk-1都用 P∞ 代替,得到
036.072.064.0
036.064.036.164.0
2
2




PP
PPP
8
3
64.02
36.064.0472.072.0 2?

P
第二章 维纳滤波和卡尔曼滤波根据 P∞,可以确定达到稳态后的卡尔曼滤波的状态方程,
kkkkkk
k
k
k
kk
yxxyxx
P
P
P
H
PP
8
3
5.0)?8.0(
8
3
8.0?
8
3
6.1
6.0
1
6.036.064.0
111
'
'
1
'




(2.5.33)
第二章 维纳滤波和卡尔曼滤波
(4) 用维纳滤波的方法分析 。 采用功率谱分解的方法,得到 x(n)的时间序列信号模型的传输函数 H(z):
15.01
8
3
)(
)()(
zzY
zXzH
第二章 维纳滤波和卡尔曼滤波上式说明 x是一阶 AR模型,对 H(z)做 Z反变换得到
)
2
1
(
8
3
8
75.0
8
5.1
8
3
)2()1()0()2(?
2
1
8
3
8
5.1
8
3
)1()0()1(?
8
3
)0()0()0(?
)()()(?)(?
)()5.0(
8
3
)(
12
012012
010101
0
*
yy
yyyyhyhyhx
yyyyyhyhx
yyhx
mnymhnxns
nunh
m
n






第二章 维纳滤波和卡尔曼滤波




n
m
nn
yymnymhnx
yy
yyyy
yhyhyhyhx
0
1
23
0123
0123
2
1
8
3
)()()(?
2
1
8
3
8
375.0
8
75.0
8
5.1
8
3
)3()2()1()0()3(?
(2.5.34)
比较 ( 2.5.33) 式和 ( 2.5.34) 式,可以看出卡尔曼滤波的稳态解与维纳解是相等的 。
(解毕 )
第二章 维纳滤波和卡尔曼滤波通过上面的例题,可以看出维纳滤波是已知前 p个观测数据及信号与噪声的相关函数,通过建立模型的方法分析的 。 卡尔曼滤波要求已知前一个时刻的状态估计值 x(k-1)和当前的观测值 yk,由状态方程和量测方程递推得到结果 。 维纳滤波的解以 H(z)的形式给出,卡尔曼滤波是以状态变量的估计值给出解的形式 。 它们都采用均方误差最小的准则,但卡尔曼滤波有一个过渡过程,其结果与维纳滤波不完全相同,但到达稳态后,
结果相同 。
^
第二章 维纳滤波和卡尔曼滤波
2.5.3 应用举例下面举一个雷达跟踪目标物的例子说明卡尔曼滤波的应用 。
雷达跟踪目标的基本原理是通过发射脉冲,根据接收到的脉冲与发射脉冲的时间间隔,来确定目标物的距离和速度 。
由于干扰的影响,接收到的脉冲波形变化很大,那么一次的测量结果可能存在很大的误差 。 为了减小误差,往往采取发射一串脉冲的方法进行测量 。
第二章 维纳滤波和卡尔曼滤波我们知道,空间中的一点需要由径向距离和方位角来确定 。
假设雷达跟踪的目标为飞行器,发射的脉冲时间间隔为 T,在时间 k,径向距离为 R+ρ (k),在时间 k+1,距离为 R+ρ (k+1),
两者之间有秒的延时,因此 T表示空间一次扫描的时间间隔 。
平均距离用 R表示,ρ (k)和 ρ (k+1)表示对平均值的偏差 。 假定偏差是统计随机的,均值为零,那么可以写出距离方程
)()()1( kTkk
式中,表示速度。用 u表示加速度,则可以得到加速度方程)(k
)()1()1( kkkTu
第二章 维纳滤波和卡尔曼滤波假定加速度 u(k)是零均值的平稳白噪声,即满足
22 )]([
0)]()1([
ukuE
kukuE


为了后面叙述方便,定义 x(k)表示第 k个雷达回波脉冲获得的目标距离,z(k)表示在第 k个雷达脉冲进行数据处理之后的目标距离估计,z(k)表示在第 k个雷达脉冲进行数据处理之后的目标速度估计 。 设定状态变量为 x(k),选择的状态变量有 4个,
分别表示径向距离,径向速度,方位角和方位角速度,即
.
T4321 )](),(),(),([)](),(),(),([)( kkkkkxkxkxkxkX
( 2.5.35)
第二章 维纳滤波和卡尔曼滤波根据状态变量的物理含义,得到以下方程:
)()()1(
)()()1(
)()()1(
)()()1(
244
433
122
211
kukxkx
kTxkxkx
kukxkx
kTxkxkx




式中,u1(k)和 u2(k)表示在区间 T径向速度和方向的变化。 将上式写成矩阵形式
)(
0
)(
0
)(
)(
)(
)(
1000
100
0010
001
)1(
)1(
)1(
)1(
2
1
4
3
2
1
4
3
2
1
ku
ku
kx
kx
kx
kx
T
T
kx
kx
kx
kx
第二章 维纳滤波和卡尔曼滤波由此得到卡尔曼滤波信号模型的状态方程
)()(),1()1(? kkXkkAkX ( 2.5.36)
再来看量测方程,距离和方向的估计值为
)()()(
)()()(
232
111
kvkxkz
kvkxkz


这里 v1(k),v2(k)为观测偏差,将上两式分别写成向量形式和矩阵形式,Z(k)=CX(k)+V(k)

)(
)(
)(
)(
)(
)(
0100
0001
)(
)(
2
1
4
3
2
1
2
1
kv
kv
kx
kx
kx
kx
kz
kz
( 2.5.37a)
( 2.5.37b)
第二章 维纳滤波和卡尔曼滤波观测噪声 V(k)假定为高斯噪声,均值为 0,方差为 σ ρ 2和 σ θ 2 。
状态方程激励信号的协方差阵为

2
2
2
1T
T
000
0000
000
0000
)]()([)(
)()]()([


kkEkQ
kQjkE
kj
(2.5.38)
其中,,为径向加速度在 T时刻的方差;,
为角加速度在 T时刻的方差 。
][ 2121 uE ][ 2222 uE
第二章 维纳滤波和卡尔曼滤波量测方程的噪声协方差阵

)(0
0)(
)]([)]()([
)]()([)]([
)]()([)(
)()]()([
2
2
22
22
2
212
21
2
1T
T
2,21,2
2,11,1
k
k
kvEkvkvE
kvkvEkvE
kvkvEkR
kRjvkvE
vv
vv
kj


(2.5.39)
第二章 维纳滤波和卡尔曼滤波为了简单起见,假定在各个方向,加速度服从均匀分布,
其概率密度函数,并将 u的值限制在 ± M之间,那么,加速度的方差
Muf 21)(?
22
3
1 M
u
根据误差信号协方差阵 P(k)的定义
P(k)=E[ e(k)eT(k)]
可以计算出,单个信号的均方误差和两个信号的协方差矩阵分别为

)()(
)()(
)]([)]()([
)]()([)]([
)(
)]([)()(
2,21,2
2,11,1
2
212
21
2
1
2
111
kpkp
kpkp
keEkekeE
kekeEkeE
EkP
keEkpkP
第二章 维纳滤波和卡尔曼滤波根据接收到的相邻两个回波脉冲,可以测量出飞行器的距离
z1(1)和 z2(2),方位角 z2(1)和 z2(2)。 根据这四个数据,用
(2.5.35)式计算状态变量
)]2()2([
1
)2()(?
)2()2(?)(?
)]2()2([
1
)2()(?
)2()2(?)(?
124
23
111
11
zz
T
kx
zkx
zz
T
kx
zkx




T
22
2
11
1
)1()2(),2(,)1()2(),2()(?



T
zzz
T
zzzkX
因此,k时刻的状态向量 可以写为
)(? kX
第二章 维纳滤波和卡尔曼滤波取 k=2,将上式中的观测信号 z1(1),z1(2),z2(1),z2(2)用 (2.5.37)
式代入,得到状态变量


T
vv
x
vx
T
vv
x
vx
T
vxvx
vx
T
vxvx
vx
X
)1()2(
)1(
)2()2(
)1()2(
)1(
)2()2(
))1()1(())2()2((
)2()2(
))1()1(())2()2((
)2()2(
)2(?
33
4
33
11
2
11
3333
33
1111
11
( 2.5.40)
根据 (2.5.35)式,k=2时的状态向量为
)1()1(
)2(
)1()1(
)2(
)2(?
24
3
12
1
ux
x
ux
x
X ( 2.5.41)
第二章 维纳滤波和卡尔曼滤波计算 k=2时刻的协方差阵
})]2(?)2()][2(?)2({[)2( TXXXXEP
式中

T
vv
u
v
T
vv
u
v
XX
)1()2(
)1(
)2(0
)1()2(
)1(
)2(0
)2(?)2(
11
2
2
11
1
1
第二章 维纳滤波和卡尔曼滤波因此,误差协方差阵是 4× 4阶矩阵,假设噪声源 u和 v是独立的,
则协方差阵为
4443
3433
2221
1211
00
00
00
00
)2(
pp
pp
pp
pp
P
其中
2
22
2
44
2
4334
2
33
2
12
2
22
2
2112
2
11
2
,,
2
,,





T
p
T
ppp
T
p
T
ppp
pp
p
第二章 维纳滤波和卡尔曼滤波
2.5.4
从理论上讲,卡尔曼滤波的递推算法可以无限地继续下去 。
然而在实际问题中的某些条件下,可能产生发散问题 。 也就是说,实际应用中发现估计误差大大地超过了理论误差的预测值,
而且误差不但不减小,反而越来越大,即不收敛 。
导致发散的一个原因是舍入误差的影响以及递推算法使得舍入误差积累的影响 。 计算机存贮单元的长度有限,使得舍入误差不可避免地存在,它相当于在状态方程和量测方程中又加入了噪声,带来的后果是有可能改变某些矩阵的性质,引起误差矩阵失去正定性和对称性 。 如果均方误差阵受到扰动而离开稳定解,只要它没有失去正定性,那么仍可能返回稳定解 。
第二章 维纳滤波和卡尔曼滤波舍入误差引起的发散现象可以采用双精度运算得以改善,
但运算量要增加许多 。 目前多采用平方根法,即把递推公式中的均方误差阵 P改用其平方根 P1/2实现 。 具体的分析参见文献
[ 16] 。
另一种类型的发散问题是由于待估计过程模型的不精确引起的 。 人们在设计卡尔曼滤波时,认为分析过程是按某一规律发展的,但实际上是按另一规律演变的 。 如假定待分析过程的模型是一随机数,而实际过程是一个随机斜面,这样滤波器将连续地试着用错误曲线去拟合观测数据,结果导致发散 。
第二章 维纳滤波和卡尔曼滤波当选择系统模型不准确时,由于新观测值对估计值的修正作用下降,陈旧观测值的修正作用相对上升,是引发滤波发散的一个重要因素 。 因此逐渐减小陈旧观测值的权重,相应地增大新观测值的权重,是抑制这类发散的一个可行途径 。 常用的方法有衰减记忆法,限定记忆法,限定下界法等 。 另外,通过人为地增加8模型输入噪声方差,用扩大了的系统噪声来补偿模型误差,抑制模型不准确所造成的发散现象,也是一种常见的策略,常用的方法有伪随机噪声法等 。
还存在第三种发散问题,它是由于系统不可观察引起的 。
所谓不可观察,是指系统有一个或几个状态变量是隐含的,现有的观测数据不能提供足够的信息来估计所有的状态变量 。 这种发散问题表现为估计值误差不稳定或者均方误差阵的主对角线上有一项或几项无限增长 。