10.6 卡尔曼滤波器简介本节讨论如何从带噪声的测量数据把有用信号提取出来的问题。
通常,信号的频谱处于有限的频率范围内,而噪声的频谱则散布在很广的频率范围内。
如前所述,为了消除噪声,可以把 FIR滤波器或IIR滤波器设计成合适的频带滤波器,进行频域滤波。
但在许多应用场合,需要进行时域滤波,从带噪声的信号中提取有用信号。
虽然这样的过程其实也算是对信号的滤波,但所依据的理论,即针对随机信号的估计理论,是自成体系的。
人们对随机信号干扰下的有用信号不能“确知”,只能“估计”。
为了“估计”,要事先确定某种准则以评定估计的好坏程度。
最小均方误差是一种常用的比较简单的经典准则。
典型的线性估计器是离散时间维纳滤波器与卡尔曼滤波器。
对于平稳时间序列的最小均方误差估计的第一个明确解是维纳在1942年2月首先给出的。
当时美国的一个战争研究团体发表了一个秘密文件,其中就包括维纳关于滤波问题的研究工作。
这项研究是用于防空火力控制系统的。
维纳滤波器是基于最小均方误差准则的估计器。
为了寻求维纳滤波器的冲激响应,需要求解著名的维纳-霍夫方程。
这种滤波理论所追求的是使均方误差最小的系统最佳冲激响应的明确表达式。
这与卡尔曼滤波(Kalman filtering)是很不相同的。
卡尔曼滤波所追求的则是使均方误差最小的递推算法。
在维纳进行滤波理论研究并导出维纳-霍夫方程的十年以前,在1931年,维纳和霍夫在数学上就已经得到了这个方程的解。
对于维纳-霍夫方程的研究,20世纪五十年代涌现了大量文章,特别是将维纳滤波推广到非平稳过程的文章甚多,但实用结果却很少。
这时正处于卡尔曼滤波问世的前夜。
维纳滤波的困难问题,首先在上世纪五十年代中期确定卫星轨道的问题上遇到了。
1958年斯韦尔林(Swerling)首先提出了处理这个问题的递推算法,并且立刻被承认和应用。
1960年卡尔曼进行了比斯韦尔林更有意义的工作。
他严格地把状态变量的概念引入到最小均方误差估计中来,建立了卡尔曼滤波理论。
空间时代的到来推动了这种滤波理论的发展。
维纳滤波与卡尔曼滤波所研究的都是基于最小均方误差准则的估计问题。
维纳滤波理论的不足之处是明显的。
在运用的过程中,它必须把用到的全部数据存储起来,而且每一时刻都要通过对这些数据的运算才能得到所需要的各种量的估值。
按照这种滤波方法设置的专用计算机的存储量与计算量必然很大,很难进行实时处理。
虽经许多科技工作者的努力,在解决非平稳过程的滤波问题时,给出能用的方法为数甚少。
到五十年代中期,随着空间技术的发展,这种方法越来越不能满足实际应用的需要,面临了新的挑战。
尽管如此,维纳滤波理论在滤波理论中的开拓工作是不容置疑的,维纳在方法论上的创见,仍然影响着后人。
五十年代中期,空间技术飞速发展,要求对卫星轨道进行精确的测量。
为此,人们将滤波问题以微分方程表示,提出了一系列适应空间技术应用的精练算法。
1960年和1961年,卡尔曼(R.E.Kalman)和布西(R.S.Bucy)提出了递推滤波算法,成功地将状态变量法引入到滤波理论中来,用消息与干扰的状态空间模型代替了通常用来描述它们的协方差函数,将状态空间描述与离散时间更新连系起来,适于计算机直接进行运算,而不是去寻求滤波器冲激响应的明确公式。
这种方法得出的是表征状态估计值及其均方误差的微分方程,给出的是递推算法。
这就是著名的卡尔曼理论,或称卡尔曼-布西滤波。
从维纳-霍夫方程来看,维纳滤波算法是十分低效的。
这种算法要求设置大量存储器保存过去的测量数据,新的一个数据来到后,要进行刷新,重新计算自相关和互相关序列。
再者,求解这个方程需要耗费大量机时对高阶矩阵求逆。
因此,维纳滤波算法难以用于实时处理中,尤其是无法用于军事、航空航天等领域。
在维纳滤波中,必须首先求解维纳-霍夫方程,得到有限长度的冲激响应序列,在这个意义上,维纳滤波器属于FIR滤波器的范畴。
卡尔曼滤波不要求保存过去的测量数据,当新的数据到来时,根据新的数据和前一时刻的诸量的估值,借助于系统本身的状态转移方程,按照一套递推公式,即可算出新的诸量的估值。
这一点说明卡尔曼滤波器属于IIR滤波器范畴。
这就是说,与维纳滤波器不同,卡尔曼滤波器能够利用先前的运算结果,再从当前数据提供的最新信息,即可得到当前的估值。
卡尔曼递推算法大大减少了滤波装置的存储量和计算量,并且突破了平稳随机过程的限制,使卡尔曼滤波器适用于对时变信号的实时处理。
因此,卡尔曼滤波器在应用上有更加广泛的可能性和更加美好的前景。
请读者回顾第2章的两个例子,即例2.5.1 和例2.6.5。
这两个例子有助于理解维纳滤波和卡尔曼滤波的主要区别。
例2.5.1用FIR滤波器实现5点滤波。
这有点像维纳滤波。
例2.5.2则用IIR滤波器实现同样的功能,而且滤波效果还稍好。
在这个例子中,利用前一个时刻的滤波输出,再添上当前的输入(新息),即可递归地求出当前的滤波输出。
这有点像卡尔曼滤波。
为了导出维纳滤波器和卡尔曼滤波器,以下各小节先介绍两种估计器:非递归估计器和递归估计器。
然后依次介绍标量维纳滤波、向量维纳滤波、标量卡尔曼滤波和向量卡尔曼滤波。
重点是对标量卡尔曼滤波算法的推导。
整个过程不涉及繁琐的数学推导。
10.7 非递归估计器y与设x表示慢变化(相对于叠加于其上的噪声)的时变信号,信号的测量值)(ix成线性关系,并加上均值为零而方差为2σ的白色噪声成分)(i v。
这个噪声成分源v自测量或其他原因引起的随机误差。
这样就有=i y+)()(i vx(10.7.1)这里把信号看作是均值为x i y E =)]([,方差为2v σ的随机变量。
人们能够测到的信号是)(i y 。
由于噪声的存在,人们只能够估计x 值,而永远无法直接测到有用信号x 。
)(i y 的平均值)]([ˆi y E x=是一个最好的估值。
如果用如下的N 阶FIR 滤波器对信号)(i y 进行滤波,∑==Nn i y i h Nx1)()(1ˆ(10.7.2A ) 式中,)(,),2(),1(N y y y 是N 个数据。
则滤波值xˆ可以被作为x 的估计值。
由于)(i v 的均值为零,故随着平均次数N的增大,xˆ中的随机噪声成分趋于零。
若1)()2()1(===N h h h ,则x i v Nx N i y NxN i N i ≈+==∑∑-=-=])([1)(1ˆ110 (10.7.3)所以,FIR 滤波器可以作为一个估计器(estimator )。
随着平均次数N 的增大,估计值无限地接近真值。
这样的估计称为无偏估计(unbiased estimation)。
估计器的均方误差为2110}])([1{x i v x N E p N i N i e -+=∑∑-=-=N i v NE v N i /])(1[2210σ==∑-=(10.7.4)由此可见,经滤波后,估计量的方差是原信号方差除以N 。
增大N 值可以提高估计器质量。
但N 值不能过大,因为估计器本质上是一个低通滤波器,N 值增大时,有用信号的细节(高频分量)也将被滤除了。
上述估计器只是简单地进行求和-平均运算,所依据的准则是使平均误差最小,而不是均方误差最小。
10.3节已导出基于LME 准则的最佳系统方程,即W-H 方程(10.3.6)。
基于这个方程,曾导出10.5节的维纳滤波器。
以这个滤波器作为估计器,就得到均方误差最小的非递归估计器。
在目前的场合,输入信号已按式(10.7.1)定义为)()(i v x i y +=,即变量y 代表输入信号,而x 则代表期望信号。
所以,在对非递归估计器进行分析前,应将式(10.3.6)修改为下式:∑-===-10110 ),()()(N k xy yy ,N-,,l l r k l r k h(10.7.5A )相应的矩阵方程为:⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡-=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡---- )1( )1( )0( )1( )1( )0( )0( )2( )1( )2( )0( )1( )1( )1( )0( N r r r N h h h r N r N r N r r r N r r r xy xy xy yy yy yy yy yy yy yy yy yy (10.7.6)这里,假定信号是平稳过程,故自相关矩阵是对称的。
应该记住,被估计量是式(10.7.1A )中的x ,而估计量xˆ则由式(10.7.2)算出。
估计误差为∑=-=-=Ni i y i h x xx e 1)()(ˆ (10.7.7) 方差为][2e E p e =(10.7.8)在按LME 准则导出W-H 方程时,曾令0)(])()([2)(1=--=∂∂∑=j y i y i h x E j h p Ni eMj ,,2,1 = (10.7.9)由此得N j j ey E ,,2,1 0)]([ ==(10.7.10) 式中,xx eˆ-=是估计误差。
在估计理论中,把 0)]([=j ey E 称为正交原理。
它的意思是:误差xx e ˆ-=与被测采样值)(j y 的乘积在期望(即平均值)意义上等于零。
根据这个式子,由)(]})()([{][12ex E i y i h x e E e E p Ni e =-==∑=得到对应于上述最优解的最小均方误差。
因此有 ])([)()(12∑=-=Ni e i xy E i h x E p或 ∑=-=Ni xy e i r i h x E p 12)()()((10.7.11)其中,)(i r xy 是互相关函数 )]([)(i xy E i r xy =(10.7.12A)这样,用方程组(10.7.5A)、估计方程(10.7.2A) 以及由式(10.7.11)给出的最小均方误差就得到了估计问题的完整解。
这三个式子的矩阵形式为xy yy r h r = (10.7.5B) 式中,yy r 是)(M M⨯自相关矩阵,h 和xy r 是)1(⨯M 列向量。
方程(10.7.5B)的形式解为xyyy r r h 1-=(10.7.13) 估计方程可写为y h xT =ˆ(10.7.2B)这里yy r 和y 是)1(⨯M列向量,T h 是行向量。
把式(10.7.13)代入式(10.7.2B),由于矩阵yy r 是对称的,故x 的估计为y r r x 1yy T xy -=ˆ (10.7.14)用同样的方法可以得到最小均方误差为=e p )(2x E xy 1yy T xyr r r -- (10.7.15) 上面分析了基于维纳滤波器的标量估计器。
应该指出,在推导中还没有使用关系式(10.7.1),即)()(k v x k y +=。
因此所得结论更具有一般性。