CN108681621A - 基于Chebyshev正交多项式扩展RTS Kalman平滑方法 - Google Patents

基于Chebyshev正交多项式扩展RTS Kalman平滑方法 Download PDF

Info

Publication number
CN108681621A
CN108681621A CN201810309568.XA CN201810309568A CN108681621A CN 108681621 A CN108681621 A CN 108681621A CN 201810309568 A CN201810309568 A CN 201810309568A CN 108681621 A CN108681621 A CN 108681621A
Authority
CN
China
Prior art keywords
chebyshev
matrix
rts
polynomials
smoothing
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201810309568.XA
Other languages
English (en)
Other versions
CN108681621B (zh
Inventor
丁国强
张焕龙
娄泰山
杨存祥
张铎
王晓雷
方洁
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhengzhou University of Light Industry
Original Assignee
Zhengzhou University of Light Industry
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhengzhou University of Light Industry filed Critical Zhengzhou University of Light Industry
Priority to CN201810309568.XA priority Critical patent/CN108681621B/zh
Publication of CN108681621A publication Critical patent/CN108681621A/zh
Application granted granted Critical
Publication of CN108681621B publication Critical patent/CN108681621B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Feedback Control In General (AREA)

Abstract

本发明提出了一种基于Chebyshev正交多项式扩展RTS Kalman平滑方法,用以解决传统的平滑算法无法对非线性系统状态变量开展滤波后平滑操作的问题。本发明建立SLAM系统的非线性的状态模型;对非线性SLAM系统的状态变量参数Kalman滤波;基于Chebyshev多项式拟合逼近SLAM系统实施Chebyshev多项式逼近计算操作,计算平滑算法的预测均值、预测方差矩阵和协方差矩阵;获取非线性系统方程的Chebyshev多项式拟合逼近计算的平滑均值及平滑方差矩阵;根据估计数据开展Chebyshev多项式RTS平滑计算。本发明利用Chebyshev多项式拟合SLAM系统的模型方程,实现状态向量的滤波平滑计算,具有较好的计算优势和计算效能。

Description

基于Chebyshev正交多项式扩展RTS Kalman平滑方法
技术领域
本发明涉及导航制导与控制中的航空航天系统信息处理的技术领域,尤其涉及一种基于Chebyshev正交多项式扩展RTS(Rauch-Tung-Striebel)Kalman平滑方法,对运动载体的即时定位与地图构建(Simultaneous Localization And Mapping,SLAM)SLAM问题进行平滑处理。
背景技术
根据估计时刻用到的观测量信息情况,最优估计理论分为预测、滤波和平滑三种类型,其中滤波计算是利用当前时刻以及以前时刻的所有观测信息对当前系统状态变量进行估计计算,而平滑方法除了利用滤波计算用到的观测信息,还要用到当前时刻以后的部分或者所有观测信息。因此,理论上平滑方法是一种离线处理方法,在滤波计算基础上对系统状态变量进一步做出改善,从而获得更加精确的计算数据结果。
平滑算法和滤波算法一样,其理论基础也是基于Bayesians最优滤波理论,假设系统状态变量满足高斯分布。平滑估计算法一般可分为固定点平滑、固定滞后平滑和固定区间平滑等,其中固定区间平滑是利用某一时间区间内的所有观测信息对所有状态变量进行估计计算的一种方法,其应用范围最为广泛,而Rauch-Tung-Striebel平滑算法就是一种固定区间平滑计算方法,但是传统的平滑算法都是针对线性系统状态变量开展滤波后的平滑操作,但是近年来也有很多种非线性平滑算法,如其中应用最广泛的基于Taylor级数扩展表达式的一阶或者二阶RTS非线性平滑算法,被称为扩展RTS平滑算法,以及基于Sigma点非线性RTS平滑算法、中心差分RTS平滑算法、Gauss-Hermite数值积分逼近的RTS平滑算法和容积RTS平滑算法等。
发明内容
针对传统的平滑算法都是针对线性系统状态变量开展滤波后平滑操作的技术问题,本发明提出一种基于Chebyshev正交多项式扩展RTS Kalman平滑方法,将Chebyshev正交多项式的优异性质应用于SLAM问题系统设计中,实现SLAM系统非线性模型状态参数最优滤波平滑计算,具有较好的计算优势和计算效能。
为了达到上述目的,本发明的技术方案是这样实现的:一种基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其步骤如下:
步骤一:建立SLAM系统的非线性的状态模型,包括状态方程和非线性观测方程;
步骤二:对非线性SLAM系统状态模型的状态变量参数的Kalman滤波计算,经由k=T步迭代计算获得第T步的最优滤波结果,并存储各个时刻的估计数据;
步骤三:从k=T开始,对步骤二获取的各个时刻滤波的估计数据开展逆向平滑操作,基于Chebyshev多项式拟合逼近SLAM系统的非线性状态方程,实施Chebyshev多项式逼近计算操作,根据第T步滤波数据,计算平滑算法的预测均值、预测方差矩阵和协方差矩阵;
步骤四:计算Chebyshev多项式拟合逼近计算的平滑增益,获取非线性系统方程的Chebyshev多项式拟合逼近计算的平滑均值及平滑方差矩阵;
步骤五:令k=T-1,根据T-1步的估计数据再次开展Chebyshev多项式RTS平滑计算,获得T-1步的Chebyshev多项式拟合逼近的平滑均值和方差;进而计算直至k=0时的Chebyshev多项式拟合逼近平滑数据,从而完成SLAM问题系统状态参数变量的平滑计算任务。
所述步骤一中SLAM系统的非线性的状态模型为
其中,xk∈Rn是第k步的系统状态变量,yk∈Rm是第k步的系统观测变量,Rn和Rm分别表示n和m维的实数空间,qk-1~N(0,Qk-1)和rk~N(0,Rk)分别表示高斯过程噪声和观测噪声,Qk-1表示系统状态变量第k-1步的过程噪声方差,Rk表示系统观测变量第k步的过程噪声方差,f(·)和h(·)分别表示系统模型的动态模型函数和观测模型函数。
所述步骤三中从k=T-1开始,利用第k步的估计数据mk|k和Pk|k,利用Chebyshev多项式拟合逼近计算第k+1步组合导航系统状态变量的预测均值为:
mk+1|k=E[f(xk)|mk|k,Pk|k]
利用有限N项的Chebyshev多项式逼近整理可得到:
进一步整理可以获得预测均值表达式为,
其中,表示第k步已知均值mk|k表达的矩阵矩阵表示为:
利用Chebyshev多项式拟合逼近组合导航系统状态变量的平滑预测方差矩阵,由第一算法的卷积计算获得系数矩阵V0:nN,系统状态变量的预测方差矩阵表达为:
其中,Qk表示系统状态变量的过程噪声方差;V0:2N=Conv(b0:N,b0:N),b0:N为Chebyshev多项式第k步的系数:b0:N=[a0-mk+1|k,a1,a2,…,aN];
协方差矩阵Ck+1|k计算为:
所述第一算法为:
首先利用Chebyshev多项式的系数赋值V0:nN←b0:N
从i=1到n-1,计算卷积V0:nN=Conv(V0:nN,b0:N);其中,Conv(·)表示卷积计算;
经过n-1步迭代计算获得系数向量V0:nN
所述步骤三中利用Chebyshev多项式拟合逼近SLAM系统非线性状态方程开展多项式拟合逼近计算操作:利用有限N项Chebyshev多项式开展SLAM系统状态方程拟合逼近计算为:
其中,y=f(x)表示SLAM系统状态方程的非线性项,g(x)表示有限N项Chebyshev多项式函数;且在自变量x的取值区间[-1,1]上,N+1个Chebyshev多项式系数cj为:
那么,
所述有限N项Chebyshev多项式整理获得:
其中,c0:N=[c0,c1,…,cN]是Chebyshev级数多项式的系数项,AN矩阵由Chebyshev系数组成,表示矩阵AN的第n列向量:
其中,列向量是由N阶Chebyshev多项式中的包括第N个单项式的所有系数组成。
所述步骤三中计算SLAM系统状态方程拟合逼近的预测均值:非线性函数y的一阶矩期望均值可写为
其中,mi表示随机变量Δx的第i阶矩,上式可整理为:
系数 是一个下三角矩阵,将其整理为
其中,
其中,⊙表示Hadamard乘积符号,是一个常值矩阵,通过构建第N+1行向量
将RN行向量截断来获取矩阵的其他行向量。
所述步骤四中平滑增益Gk为:
计算第k+1步的平滑均值为:
mk+1|k+1=mk|k+Gk(xk+1-mk+1|k);
相应的平滑方差矩阵为:
最终经由T步迭代平滑计算获得系统状态变量的平滑均值和方差矩阵分别为:
mk|T=mk|k+Gk(xk+1|T-mk+1|k)
本发明的有益效果:利用Chebyshev多项式扩展逼近SLAM问题系统的非线性多项式模型方程,利用Chebyshev多项式拟合SLAM问题系统的模型方程,实现SLAM问题系统状态向量的滤波平滑计算,从而达到进一步改善SLAM问题系统状态变量参数最优估计的精度要求。经由机器人SLAM系统仿真实验,并与传统的扩展Kalman平滑算法对比,验证了本发明的计算优势和计算效能。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1是本发明的流程图。
图2是本发明在机器人SLAM仿真环境设置场景中的应用实例,描述了本发明平滑计算获得的载体运动平滑曲线与实际运动轨迹的对比图,图中曲线如图注所示。
图3是本发明对机器人SLAM系统的处理结果,上图为处理后获得的机器人位置误差,下图为机器人x和y两个方向的标准差数据曲线。
图4是利用扩展Kalman平滑算法计算机器人SLAM系统的平滑运动轨迹曲线和真实运动轨迹曲线对比图。
图5是利用扩展Kalman平滑算法计算机器人SLAM系统的位置误差和机器人x和y方向的标准差数据曲线。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有付出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
本发明提出的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,经由系统非线性函数的Chebyshev正交多项式扩展逼近计算获得非线性系统函数的最优平滑,以SLAM问题系统离散化的随机状态空间模型为:
其中,xk∈Rn是第k步的系统状态变量,yk∈Rm是第k步的系统观测变量,qk-1~N(0,Qk-1)和rk~N(0,Rk)分别表示系统高斯过程噪声和观测噪声,Qk-1表示系统状态变量第k-1步的过程噪声方差,Rk表示系统观测变量第k步的过程噪声方差,f(·)和h(·)表示系统模型的动态模型函数和观测模型函数。
平滑算法一般都是用于改善滤波结果的再平滑计算过程,平滑计算在数学上就是在给定所有观测向量{yk}(T>k)基础上的系统状态变量的条件概率分布p(xk|y1:T)(k=1,2,…,T)的集合。本发明中利用高斯分布假设推理过程,也就是系统状态变量满足高斯正态分布:
p(xk|y1:T)=N(xk|mk|T,Pk|T) (2)
其中,mk|T是系统状态变量的均值向量,Pk|T是第k步平滑计算获得的逼近高斯分布的系统状态变量的方差矩阵。但是在滤波和平滑操作中都需要计算系统非线性函数的一阶和二阶矩,如针对满足高斯分布的泛指非线性向量函数g(x),
其中,表示非线性函数g(x)的期望,。但是对于大多数的非线性函数g(x)而言,式(3)、(4)和(5)这三类矩的精确计算几乎是不可能的,因此,很多文献提出了多种逼近计算方法,如Taylor级数扩展法、Sigma点法、Gauss-Hermite数值积分逼近法、中心差分逼近法和容积点逼近法以及粒子群逼近计算法等等。本发明给出Bayesian平滑的理论计算步骤,假设已经获得了第k步系统状态变量的滤波解的高斯逼近,
p(xk|y1:k)=N(xk|mk|k,Pk|k) (6)
那么泛指非线性函数的高斯一阶矩(期望值)可表示为,
平滑算法从时刻k=T的滤波结果开始平滑计算,直至k=0为止,从而高斯分布的Bayesian平滑计算框架为:
其中,Pk|k为估计方差矩阵,Pk+1|k为预测方差矩阵,Pk+1|T表示第k+1步Chebyshev多项式拟合逼近状态变量的平滑方差矩阵,mk+1|k为计算平滑算法的预测均值,mk|k为估计均值、mk+1|T表示第k+1步Chebyshev多项式拟合逼近的状态变量的平滑均值,Ck+1|k为计算平滑算法的协方差矩阵,Gk为Chebyshev多项式拟合逼近计算的平滑增益,表示平滑增益矩阵Gk的转置。
最终获得系统状态变量逼近计算的平滑分布:
p(xk|y1:T)≈N(xk|mk|T,Pk|T) (9)
本发明利用Chebeshev多项式逼近计算实现非线性系统函数扩展RTS平滑计算,有效减小了基于Taylor级数的扩展表达式计算的复杂性和计算量;假设非线性SLAM问题系统方程为f(x)是非线性非多项式函数,若存在着一个N阶Chebyshev多项式函数g(x),其自变量x满足x∈[-1,1],Chebyshev多项式定义为
TN(x)=cos(Narccos(x)),N=1,2,… (10)
相邻的三个切比雪夫多项式具有递推关系,可以表达为
切比雪夫多项式具有一个重要性质就是正交性,切比雪夫多项式是采用N次的切比雪夫多项式加权代数式来逼近任意的非线性函数,这些加权多项式满足正交特性,在Chebyshev多项式零点上,其正交性可表达为
切比雪夫多项式T(x)具有奇偶性,满足
TN(-x)=(-1)NTN(x), (13)
切比雪夫多项式满足T(x)∈[-1,1]取值区间,并且在此区间中T(x)具有N+1个不同的实数零点,可由这些零点实施切比雪夫多项式插值逼近计算操作。
根据Chebyshev多项式的奇偶性质及取值特性,Chebyshev多项式还可以写为
其中,αN,i表示N次Chebyshev多项式的第i阶单项式的系数,αN,N-2i也是同样的意义,表示N次Chebyshev多项式的第N-2i阶单项式的系数,式中[N/2]表示取整数,从而也可以获得两项Chebyshev多项式的乘积表达式为:
同时切比雪夫多项式函数在区间[-1,1]上交替出现N+1个极值点组,其最大值为1,最小值为-1。切比雪夫多项式的最高次项系数为2N-1,(N=1,2,···)。从而切比雪夫多项式存在着与零偏差最小的特性,并且其偏差为依此特性可以在切比雪夫多项式逼近非线性函数过程中获得多项式插值余项的极小值,这有助于有效改善扩展RTS平滑算法的计算精度。
在实际应用过程中,一般取有限N项Chebyshev多项式插值逼近非线性系统函数,有限N项的Chebyshev多项式逼近为
其中,在自变量x的取值区间[-1,1]上,N+1个Chebyshev多项式系数cj可写为
那么
另外一般情况下系统状态变量的取值区间不在[-1,1]区间,这时需要做系统状态变量的变量替换,若系统状态变量取值区间为[a,b],一般可采用变量替换表达式
那么相应的插值零点变换为:
由此可以开展Chebyshev多项式逼近非线性函数的一阶矩期望和二阶矩方差矩阵的计算工作。
假设利用Chebyshev多项式g(x)拟合非线性泛指函数y=f(x),那么f(x)可简单表达为,
这里c0:N=[c0,c1,…,cN]是由式(17)计算获得的Chebyshev级数多项式的系数项,而AN矩阵是由Chebyshev系数组成,则表示矩阵AN的第n列向量,其定义为
其中,列向量是由n阶Chebyshev多项式中的包括第N个单项式的所有系数组成。
假设自变量x满足高斯分布那么自变量x可以表示为其中Δx是满足分布的零均值高斯随机变量,自变量的幂次项表示为
其中是两项式的系数项,从而利用式(22)可以整理式(20)为
那么非线性函数y的一阶矩期望均值可写为
其中,mi表示随机变量Δx的第i阶矩,那么由式(21)可以整理式(24)为
应该说明的是, 是一个下三角矩阵,为了减小计算量,将其整理为
其中,
其中,⊙表示Hadamard乘积符号。应该说明的是,利用这种分解计算,是一个常值矩阵,对于矩阵,仅通过构建其第N+1行向量将RN的行向量截断来获取矩阵的其他行向量。
利用Chebyshev多项式的性质实施函数y=f(x)的方差矩阵的逼近计算。若Chebyshev多项式的系数可表示为:
那么函数y=f(x)的n阶逼近中心矩(n>1),如可以表达为:
其中,系数向量V0:nN可由第一算法来实现。第一算法可简单描述:首先利用Chebyshev多项式的系数赋值V0:nN←b0:N;从i=1到n-1,计算卷积V0:nN=Conv(V0:nN,b0:N);经过n-1步迭代计算,最后获得系数向量V0:nN。其中,符号Conv(·)表示卷积计算。
其实在滤波和平滑计算中使用的是一阶矩和二阶矩,在这里考虑n=2时的方差矩阵计算,此时那么经由第一算法解算整理式(28)可以获得二阶矩方差的计算,
那么协方差也可以计算为,
如图1所示,一种基于Chebyshev正交多项式扩展RTS Kalman平滑方法,首先建立非线性SLAM问题系统的状态方程和观测方程;实施非线性SLAM问题系统状态模型的状态变量参数的Kalman滤波计算,经由k=T步迭代计算获得第T步的最优滤波结果估计均值mT|T和估计方差矩阵PT|T,并存储各个时刻的估计数据mk|k和Pk|k;然后从k=T-1开始,基于Chebyshev多项式拟合逼近SLAM问题系统非线性的状态方程和观测方程,实施Chebyshev多项式逼近计算操作,计算平滑算法的预测均值mk+1|k、预测方差矩阵Pk+1|k和协方差矩阵Ck+1|k;接着计算Chebyshev多项式拟合逼近计算的平滑增益Gk,获取非线性系统的状态方程和观测方程的Chebyshev多项式拟合逼近计算的平滑均值mk|T以及平滑方差矩阵Pk|T,从而完成SLAM问题系统状态参数变量的估计计算任务。具体步骤为:
步骤一:建立SLAM系统的非线性的状态模型,包括状态方程和非线性观测方程。
构建SLAM系统的模型为式(1),其参数如前所述。
步骤二:对非线性SLAM系统状态模型的状态变量参数的Kalman滤波计算,经由k=T步迭代计算获得第T步的最优滤波结果,并存储各个时刻的估计数据。
经由k=T步的Kalman滤波迭代计算获得系统状态变量的最优滤波结果为估计均值mT|T和估计方差矩阵PT|T,并存储各个时刻的估计数据mk|k和Pk|k
步骤三:从k=T开始,对步骤二获取的各个时刻滤波的估计数据开展逆向平滑操作,基于Chebyshev多项式拟合逼近SLAM系统的非线性状态方程,实施Chebyshev多项式逼近计算操作,根据第T次滤波数据,计算平滑算法的预测均值、预测方差矩阵和协方差矩阵。
从k=T-1开始,基于Chebyshev多项式拟合逼近SLAM问题系统非线性状态方程和观测方程,实施Chebyshev多项式逼近计算操作,计算平滑算法的预测均值mk+1|k,预测方差矩阵Pk+1|k和协方差矩阵Ck+1|k
从k=T-1开始,利用第k步的估计数据mk|k和Pk|k,利用Chebyshev多项式拟合逼近计算第k+1步组合导航系统状态变量的预测均值为:
mk+1|k=E[f(xk)|mk|k,Pk|k] (31)
利用式(16)整理可得到
利用式(23)-(26),进一步整理可以获得预测均值表达式为,
其中,表示第k步已知均值mk|k表达的矩阵,类似于式(26)可表示为,
利用Chebyshev多项式拟合逼近组合导航系统状态变量的平滑预测方差矩阵,利用式(29),经由第一算法的卷积计算获得系数矩阵V0:nN,系统状态变量的预测方差矩阵表达为:
其中,Qk表示系统状态变量的过程噪声方差;V0:2N=Conv(b0:N,b0:N),b0:N如式(27)所述,在第k步具体表达为
b0:N=[a0-mk+1|k,a1,a2,…,aN] (36)
相应的协方差矩阵Ck+1|k可由式(30)计算为,
步骤四:计算Chebyshev多项式拟合逼近计算的平滑增益,获取非线性系统方程的Chebyshev多项式拟合逼近计算的平滑均值及平滑方差矩阵。
接着计算Chebyshev多项式拟合逼近计算的平滑增益Gk,获取非线性系统方程和观测方程的Chebyshev多项式拟合逼近计算的平滑均值mk|T,以及平滑方差矩阵Pk|T
考虑平滑增益Gk计算:
计算第k+1步的平滑均值为,
mk+1|k+1=mk|k+Gk(xk+1-mk+1|k) (39)
相应的平滑方差矩阵为,
最终经由T步迭代平滑计算获得系统状态变量的平滑均值和方差矩阵为,
mk|T=mk|k+Gk(xk+1|T-mk+1|k) (41)
从而完成RTS固定区间平滑算法的迭代计算过程。
步骤五:令k=T-1,根据T-1次的估计数据再次开展Chebyshev多项式RTS平滑计算,获得T-1次的Chebyshev多项式拟合逼近的平滑均值和方差;进而计算直至k=0时的Chebyshev多项式拟合逼近平滑数据,从而完成SLAM问题系统状态参数变量的平滑计算任务。
具体实例:考虑机器人运动载体的即时定位与地图构建(SimultaneousLocalization And Mapping,SLAM)SLAM问题,在笛卡尔坐标系下可以给出载体运动方程为,
这里SLAM系统状态向量为xk=[xk,ykk]T分别表示k时刻的载体位置坐标和方位;V是载体速度,G表示载体转向角,参数WB表示载体轮距,噪声向量vk是高斯过程噪声,vk~N(0,Qk),其中,Qk表示噪声方差。
机器人运动载体配备了距离和方位传感器,能够在方位角±30°范围内感知距离30m范围内的目标物体,由此可以获得机器人SLAM系统的观测方程为,
这里(ri,x,ri,y)是传感器感知的路标位置坐标,wk是观测白噪声,满足分布wk~N(0,Rk),其中Rk表示观测噪声方差。那么SLAM系统初始参数设置为:初始速度V0=3m/s,G=±30°,WB=4m,速度标准差σV=0.3m/s,转向角标准差σG=3°,距离标准差σr=0.1m,方位标准差σB=1°。初始状态向量x0=0,初始方差P0=diag{10-10,10-10,10-10}。由此展开仿真验证工作,并且把Chebyshev多项式扩展RTS平滑算法和扩展Kalman平滑算法计算效能对比,如图2、3、4、5所示。对比图2和图4,很明显两种算法中Chebyshev RTS平滑器算法与机器人真实运行轨迹数据拟合程度较好,而扩展Kalman平滑器算法的拟合程度比较差些;从图3和图5可以对比看出,Chebyshev RTS平滑器算法的计算标准差比较小,误差数据曲线平滑稳定,而扩展Kalman平滑器算法获得的位置误差数据变化剧烈,甚至出现数据发散现象,明显误差数据比较大,相应的标准偏差数据较大,通过利用这两种平滑器算法开展机器人SLAM系统仿真实验,获得的实验数据说明,Chebyshev RTS平滑器算法计算效能优于常规的扩展Kalman平滑算法,表现出了Chebyshev RTS平滑器算法的计算优势。
以上所述仅为本发明的较佳实施例而已,并不用以限制本发明,凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。

Claims (7)

1.一种基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,其步骤如下:
步骤一:建立SLAM系统的非线性的状态模型,包括状态方程和非线性观测方程;
步骤二:对非线性SLAM系统状态模型的状态变量参数的Kalman滤波计算,经由k=T步迭代计算获得第T步的最优滤波结果,并存储各个时刻的估计数据;
步骤三:从k=T开始,对步骤二获取的各个时刻滤波的估计数据开展逆向平滑操作,基于Chebyshev多项式拟合逼近SLAM系统的非线性状态方程,实施Chebyshev多项式逼近计算操作,根据第T步滤波数据,计算平滑算法的预测均值、预测方差矩阵和协方差矩阵;
步骤四:计算Chebyshev多项式拟合逼近计算的平滑增益,获取非线性系统方程的Chebyshev多项式拟合逼近计算的平滑均值及平滑方差矩阵;
步骤五:令k=T-1,根据T-1步的估计数据再次开展Chebyshev多项式RTS平滑计算,获得T-1步的Chebyshev多项式拟合逼近的平滑均值和方差;进而计算直至k=0时的Chebyshev多项式拟合逼近平滑数据,从而完成SLAM问题系统状态参数变量的平滑计算任务。
2.根据权利要求1所述的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,所述步骤一中SLAM系统的非线性的状态模型为
其中,xk∈Rn是第k步的系统状态变量,yk∈Rm是第k步的系统观测变量,Rn和Rm分别表示n和m维的实数空间,qk-1~N(0,Qk-1)和rk~N(0,Rk)分别表示高斯过程噪声和观测噪声,Qk-1表示系统状态变量第k-1步的过程噪声方差,Rk表示系统观测变量第k步的过程噪声方差,f(·)和h(·)分别表示系统模型的动态模型函数和观测模型函数。
3.根据权利要求2所述的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,所述步骤三中从k=T-1开始,利用第k步的估计数据mk|k和Pk|k,利用Chebyshev多项式拟合逼近计算第k+1步组合导航系统状态变量的预测均值为:
mk+1|k=E[f(xk)|mk|k,Pk|k]
利用有限N项的Chebyshev多项式逼近整理可得到:
进一步整理可以获得预测均值表达式为,
其中,表示第k步已知均值mk|k表达的矩阵矩阵表示为:
利用Chebyshev多项式拟合逼近组合导航系统状态变量的平滑预测方差矩阵,由第一算法的卷积计算获得系数矩阵V0:nN,系统状态变量的预测方差矩阵表达为:
其中,Qk表示系统状态变量的过程噪声方差;V0:2N=Conv(b0:N,b0:N),b0:N为Chebyshev多项式第k步的系数:b0:N=[a0-mk+1|k,a1,a2,…,aN];
协方差矩阵Ck+1|k计算为:
4.根据权利要求3所述的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,所述第一算法为:
首先利用Chebyshev多项式的系数赋值V0:nN←b0:N
从i=1到n-1,计算卷积V0:nN=Conv(V0:nN,b0:N);其中,Conv(·)表示卷积计算;
经过n-1步迭代计算获得系数向量V0:nN
5.根据权利要求3所述的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,所述步骤三中利用Chebyshev多项式拟合逼近SLAM系统非线性状态方程开展多项式拟合逼近计算操作:利用有限N项Chebyshev多项式开展SLAM系统状态方程拟合逼近计算为:
其中,y=f(x)表示SLAM系统状态方程的非线性项,g(x)表示有限N项Chebyshev多项式函数;且在自变量x的取值区间[-1,1]上,N+1个Chebyshev多项式系数cj为:
那么,
所述有限N项Chebyshev多项式整理获得:
其中,c0:N=[c0,c1,…,cN]是Chebyshev级数多项式的系数项,AN矩阵由Chebyshev系数组成,表示矩阵AN的第n列向量:
其中,列向量是由N阶Chebyshev多项式中的包括第N个单项式的所有系数组成。
6.根据权利要求4所述的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,所述步骤三中计算SLAM系统状态方程拟合逼近的预测均值:非线性函数y的一阶矩期望均值可写为
其中,mi表示随机变量Δx的第i阶矩,上式可整理为:
系数 是一个下三角矩阵,将其整理为
其中,
其中,⊙表示Hadamard乘积符号,是一个常值矩阵,通过构建第N+1行向量将RN行向量截断来获取矩阵的其他行向量。
7.根据权利要求6所述的基于Chebyshev正交多项式扩展RTS Kalman平滑方法,其特征在于,所述步骤四中平滑增益Gk为:
计算第k+1步的平滑均值为:
mk+1|k+1=mk|k+Gk(xk+1-mk+1|k);
相应的平滑方差矩阵为:
最终经由T步迭代平滑计算获得系统状态变量的平滑均值和方差矩阵分别为:
mk|T=mk|k+Gk(xk+1|T-mk+1|k)
CN201810309568.XA 2018-04-09 2018-04-09 基于Chebyshev正交多项式扩展RTS Kalman平滑方法 Active CN108681621B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810309568.XA CN108681621B (zh) 2018-04-09 2018-04-09 基于Chebyshev正交多项式扩展RTS Kalman平滑方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810309568.XA CN108681621B (zh) 2018-04-09 2018-04-09 基于Chebyshev正交多项式扩展RTS Kalman平滑方法

Publications (2)

Publication Number Publication Date
CN108681621A true CN108681621A (zh) 2018-10-19
CN108681621B CN108681621B (zh) 2021-11-19

Family

ID=63800803

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810309568.XA Active CN108681621B (zh) 2018-04-09 2018-04-09 基于Chebyshev正交多项式扩展RTS Kalman平滑方法

Country Status (1)

Country Link
CN (1) CN108681621B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112464149A (zh) * 2020-12-15 2021-03-09 北京百奥智汇科技有限公司 数据的概率密度分布的确定方法、装置、设备和介质
CN113324547A (zh) * 2021-05-25 2021-08-31 江苏科技大学 基于迭代扩展rts平滑滤波算法的多auv协同定位方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20030172728A1 (en) * 2000-04-12 2003-09-18 Fredrik Gustafsson Tire pressure estimation
CN102540228A (zh) * 2012-03-02 2012-07-04 重庆九洲星熠导航设备有限公司 一种单频gps高精度单点定位系统及方法
CN102968552A (zh) * 2012-10-26 2013-03-13 郑州威科姆科技股份有限公司 一种卫星轨道数据预估与修正方法
CN104038181A (zh) * 2014-06-05 2014-09-10 北京航空航天大学 一种基于nlms算法的自适应滤波器的构建方法
CN105222780A (zh) * 2015-09-07 2016-01-06 郑州轻工业学院 一种基于Stirling插值多项式逼近的椭球集员滤波方法
CN106197424A (zh) * 2016-06-28 2016-12-07 哈尔滨工业大学 遥测数据驱动的无人机飞行状态识别方法
CN106767780A (zh) * 2016-11-28 2017-05-31 郑州轻工业学院 基于Chebyshev多项式插值逼近的扩展椭球集员滤波方法
CN107181474A (zh) * 2017-07-14 2017-09-19 西安交通大学 一种基于函数展开的核自适应滤波器算法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20030172728A1 (en) * 2000-04-12 2003-09-18 Fredrik Gustafsson Tire pressure estimation
CN102540228A (zh) * 2012-03-02 2012-07-04 重庆九洲星熠导航设备有限公司 一种单频gps高精度单点定位系统及方法
CN102968552A (zh) * 2012-10-26 2013-03-13 郑州威科姆科技股份有限公司 一种卫星轨道数据预估与修正方法
CN104038181A (zh) * 2014-06-05 2014-09-10 北京航空航天大学 一种基于nlms算法的自适应滤波器的构建方法
CN105222780A (zh) * 2015-09-07 2016-01-06 郑州轻工业学院 一种基于Stirling插值多项式逼近的椭球集员滤波方法
CN106197424A (zh) * 2016-06-28 2016-12-07 哈尔滨工业大学 遥测数据驱动的无人机飞行状态识别方法
CN106767780A (zh) * 2016-11-28 2017-05-31 郑州轻工业学院 基于Chebyshev多项式插值逼近的扩展椭球集员滤波方法
CN107181474A (zh) * 2017-07-14 2017-09-19 西安交通大学 一种基于函数展开的核自适应滤波器算法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
MOUSSA YAHIA 等: "Use of Chebyshev Polynomial Kalman Filter for pseudo-blind demodulation of CD3S signals", 《INTERNATIONAL JOURNAL OF CONTROL, AUTOMATION AND SYSTEMS》 *
宫轶松: "粒子滤波算法研究及其在GPS/DR组合导航中的应用", 《中国博士学位论文全文数据库 基础科学辑》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112464149A (zh) * 2020-12-15 2021-03-09 北京百奥智汇科技有限公司 数据的概率密度分布的确定方法、装置、设备和介质
CN113324547A (zh) * 2021-05-25 2021-08-31 江苏科技大学 基于迭代扩展rts平滑滤波算法的多auv协同定位方法

Also Published As

Publication number Publication date
CN108681621B (zh) 2021-11-19

Similar Documents

Publication Publication Date Title
CN109597864B (zh) 椭球边界卡尔曼滤波的即时定位与地图构建方法及系统
CN108267731B (zh) 无人机目标跟踪系统的构建方法及应用
CN110208740A (zh) Tdoa-imu数据自适应融合定位装置及方法
CN109724597B (zh) 一种基于函数迭代积分的惯性导航解算方法及系统
CN101871782A (zh) 基于set2fnn的gps/mems-ins组合导航系统定位误差预测方法
CN105222780A (zh) 一种基于Stirling插值多项式逼近的椭球集员滤波方法
CN108520233A (zh) 一种扩展全对称多胞形集员Kalman混合滤波方法
CN104462015A (zh) 处理非高斯Lévy噪声的分数阶线性离散系统状态更新方法
Pelka et al. Introduction, discussion and evaluation of recursive Bayesian filters for linear and nonlinear filtering problems in indoor localization
CN108681621A (zh) 基于Chebyshev正交多项式扩展RTS Kalman平滑方法
CN108566178A (zh) 一种非稳态随机机会网络特征值滤波方法
CN109764876B (zh) 无人平台的多模态融合定位方法
CN111798494A (zh) 广义相关熵准则下的机动目标鲁棒跟踪方法
CN113709662B (zh) 一种基于超宽带的自主式三维反演定位方法
CN110968961A (zh) 一种连续回转电液伺服马达参数辨识方法
CN114370878A (zh) 一种基于stackf的多auv协同定位方法
CN109115228A (zh) 一种基于加权最小二乘容积卡尔曼滤波的目标定位方法
Kang et al. Finite-memory-structured online training algorithm for system identification of unmanned aerial vehicles with neural networks
Eren et al. Implementation of the spline method for mobile robot path control
CN107886058B (zh) 噪声相关的两阶段容积Kalman滤波估计方法及系统
Rathore et al. Robust steering control of autonomous underwater vehicle: Based on PID tuning evolutionary optimization technique
He et al. A SLAM algorithm of fused EKF and particle filter
CN116047495A (zh) 一种用于三坐标雷达的状态变换融合滤波跟踪方法
CN109858191A (zh) 一种广义混沌同步系统构建与电路设计方法
CN110610513A (zh) 自主移动机器人视觉slam的不变性中心差分滤波器方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant