CN101246011B - 一种基于凸优化算法的多目标多传感器信息融合方法 - Google Patents

Abstract

一种基于凸优化算法的多目标多传感器信息融合方法，涉及一种捷联惯性/GPS/大气数据计算机/三轴磁强计组合导航系统信息融合方法，该方法首先使用大气数据计算机、三轴磁强计的速度、姿态数据将惯性导航系统镇定，构成高可靠、高自主性的捷联惯性/大气数据计算机/三轴磁强计组合导航内系统；其次根据现代控制理论最优控制思想，使用凸优化算法设计鲁棒多目标滤波器融合高精度的GPS位置、速度量测信息，构成外系统对内系统进行混合H₂/H_∞/保性能多目标滤波来提高整体系统的导航精度。充分发挥传统导航手段的高可靠性、可用性和卫星导航的高精度，最终实现基于多传感器的高可靠、高精度组合导航技术。

Description

一种基于凸优化算法的多目标多传感器信息融合方法

技术领域

本发明涉及一种多目标多传感器信息融合方法，特别是一种基于凸优化算法的多目标多传感器信息融合方法，适用于航空导航技术领域，可应用于航空器多传感器导航，也可适用于空基、海基导航定位。

背景技术

随着国民经济的发展，导航定位技术已日益成为交通管制、救灾救护、物流管理、土地普查、勘测等领域不可或缺的技术手段；军事领域武器发射平台效率的实现、火力武器的精确打击都依赖于高精度、高可靠的导航手段。目前现有的导航手段较多，但任何一种导航手段都无法单独达到要求，结合多种技术手段充分发挥其优势，避免劣势的组合导航技术是当前导航技术发展的主流方向。如何将多传感器信息有效融合既实现系统的高精度又要保证高可靠性是当前研究的难点。

目前航空导航技术领域用于组合的导航手段比较丰富，主要有惯性导航、卫星导航、天文导航、地磁导航、大气数据系统等。惯性导航是一种全球全天候的导航手段，由于其实现原理使得他们隐蔽性强不易被干扰，在军事应用领域尤其重要，但是它的误差会随时间累积。大气数据系统和磁传感姿态系统可以提供的精度不发散的速度、姿态观测量，它们与惯性导航组合可以构成高自主、高可靠、高隐蔽、精度稳定的导航方式。但这种组合定位、测速精度相对较低，而卫星导航系统也是一种全球全天候的导航系统，它具有较高的精度但自主性较差容易受人为的影响，如果在GPS信号可用时引入GPS信息进行滤波必然可以大幅提高组合导航系统的精度。将多种导航手段结合构成组合导航系统是解决现有单一导航手段不足的有效方法。现有的组合导航研究大都集中在惯性/GPS组合导航kalman滤波设计中，通过Kalman滤波器对导航误差进行估计，但是由于模型和噪声的不确定性使得Kalman滤波效果较差，虽然有EKF、UKF等改进方法但仍然受制于假设条件太强使得实际应用情况不理想。GPS信号不稳定，容易受到干扰而且在大气层某些高度上GPS信号是接收不到的，在GPS信号出现坏点的情况下Kalman滤波表现不佳。粒子滤波可以解决非线性、非高斯情况下的滤波问题且对于观测值中的坏点容忍性较强，但它仍然需要噪声的先验知识，这在实际工程应用中也是比较难以获得的，而且其计算量太大，实时性不能保证。

发明内容

本发明的技术解决问题是：克服现有技术的不足，提出一种基于凸优化算法的多目标多传感器信息融合方法，解决导航系统冗余信息的融合，在模型不确定、噪声非高斯且统计特性未知的情况下改善导航精度的同时提高系统的可靠性，降低工程实施难度。

本发明的技术解决方案为：一种基于凸优化算法的多目标多传感器信息融合方法，步骤如下：

(1)将地球坐标系捷联惯导误差状态方程中参数的变化作为系统的不确定性，并使用不确定系统的多胞模型来描述地球坐标系捷联惯导误差状态方程；

(2)以大气数据计算机、三轴磁强计输出的速度、姿态数据和捷联惯导输出速度、姿态数据的差值作为观测量建立内系统的量测方程；以GPS系统输出的位置、速度数据和捷联惯导位置、速度数据的差值为观测量建立外系统的量测方程；

(3)对内系统量测，采用内系统的鲁棒静态输出反馈镇定器将不稳定的捷联惯导误差状态方程镇定；

(4)对外系统量测，采用鲁棒H₂/H_∞/保性能滤波器，使得从干扰到导航误差满足所要求的H₂/H_∞/保性能要求；

(5)将地球坐标系捷联惯导误差状态方程、内外系统量测方程、内系统镇定器、外系统鲁棒多目标滤波器联立构成闭环系统；同时为平衡对不同噪声输入通道进入的噪声的抑制效果，对噪声输入矩阵做加权处理；

(6)利用凸优化算法求解内系统鲁棒静态输出反馈镇定器和外系统鲁棒混合H₂/H_∞/保性能滤波器。

本发明的原理是：首先，根据捷联惯导误差状态方程时变性主要来源于载体比力变化且是线性依赖的特征，给定载体比力变化域顶点将捷联惯导误差状态方程转化为不确定系统的多胞型描述；其次，将捷联惯导误差状态方程、内系统镇定器和外系统滤波器联立构成闭环系统，由精度要求提出H₂/H_∞/保性能等各项性能指标的可调输出，同时对噪声输入通道做加权处理，根据鲁棒控制理论H₂/H_∞/保性能优化控制和混合多目标控制方法，基于线性矩阵不等式将组合导航系统的多目标信息融合问题转化为凸优化问题；最后，利用Matlab线性矩阵不等式工具箱LMI Toolbox求解该凸优化问题，并通过一系列代数变换从凸优化问题可行解中解出组合导航系统内系统镇定器和外系统多目标滤波器。

本发明与现有技术相比的优点在于：

(1)本发明利用完全自主的大气数据系统和磁传感姿态系统输出的速度、姿态信息阻尼惯性导航系统构成内系统，获得一定精度的高可靠、高隐蔽性、精度稳定的导航方案；同时利用高精度的卫星导航定位信息对内系统滤波构成外系统，以获取更高精度的导航数据。内系统的鲁棒静态输出反馈镇定器和外系统鲁棒混合H₂/H_∞/保性能滤波器通过凸优化技术对综合性能指标寻优求取最优解，这样在卫星导航可用时通过滤波可以获得高精度的导航数据，在卫星导航不可用时由于大气数据系统和磁传感姿态系统对捷联惯导系统的阻尼，内系统仍然可以稳定可靠工作。

(2)基于凸优化算法得到了对捷联惯导误差状态方程鲁棒的内系统镇定器和外系统H₂/H_∞/保性能多目标滤波器，克服了Kalman滤波和粒子滤波方法递推过程运算量大的缺点。

(3)通过基于凸优化的多目标设计方法保证了系统的混合性能指标要求，克服了Kalman滤波、粒子滤波等单一性能指标滤波方法的不足，提高了系统的可靠性。

(4)克服了Kalman滤波和粒子滤波对噪声统计特性要求较高的问题，避免了在模型不确定、噪声非高斯、噪声统计特性难以获得情况下Kalman滤波和粒子滤波出现的精度恶化，提高系统的可靠性。

(5)通过内系统和外系统的实现，利用了冗余的机载导航设备，达到了整体系统高可靠性、高精度的要求，避免了在GPS信号不可用时惯性/GPS组合导航系统无法正常工作的问题。

附图说明

图1为本发明的一种基于凸优化算法的多目标多传感器组合导航方法的设计流程图。

具体实施方式

如图1所示，本发明具体实现步骤如下：

1.将地球坐标系捷联惯导误差状态方程中参数的变化作为系统的不确定性，并使用不确定系统的多胞模型来描述地球坐标系捷联惯导误差状态方程。

其中建立捷联惯性导航系统误差状态方程的多胞型描述如下：

(1)离散化地球坐标系捷联惯导误差状态方程：

x(k+1)＝A·x(k)+B₁·w₁(k)＝A·x(k)+B·w(k)

$A=\left[\begin{array}{ccc}I+\frac{T^2{N}^e}{2}& \mathrm{TI}+\frac{T^3{N}^e}{6}-T^2{\mathrm{\Ω}}_{\mathrm{ie}}^e& -\frac{T^2{F}^e}{2}+\frac{T^3{\mathrm{\Ω}}_{\mathrm{ie}}^e{F}^e}{3}+\frac{T^3{F}^e{\mathrm{\Ω}}_{\mathrm{ie}}^e}{6}-\frac{T^4{N}^e{F}^e}{24}\\ {\mathrm{TN}}^e& I+\frac{T^2{N}^e}{2}-2T{\mathrm{\Ω}}_{\mathrm{ie}}^e& -{\mathrm{TF}}^e+T^3{\mathrm{\Ω}}_{\mathrm{ie}}^e{F}^e+\frac{T^2{F}^e{\mathrm{\Ω}}_{\mathrm{ie}}^e}{2}-\frac{T^3{N}^e{F}^e}{6}\\ 0& 0& I-T{\mathrm{\Ω}}_{\mathrm{ie}}^e\end{array}\right]$

$B_1=\left[\begin{array}{cc}\frac{T^{2}I}{2}+\frac{T^4{N}^e}{24}-\frac{T^3{\mathrm{\Ω}}_{\mathrm{ie}}^e}{3}& -\frac{T^3{F}^e}{6}+\frac{T^4{\mathrm{\Ω}}_{\mathrm{ie}}^e{F}^e}{12}+\frac{T^4{F}^e{\mathrm{\Ω}}_{\mathrm{ie}}^e}{24}-\frac{T^5{N}^e{F}^e}{120}\\ \mathrm{TI}+\frac{T^3{N}^e}{6}-T^2{\mathrm{\Ω}}_{\mathrm{ie}}^e& -\frac{T^2{F}^e}{2}+\frac{T^3{\mathrm{\Ω}}_{\mathrm{ie}}^e{F}^e}{3}+\frac{T^3{F}^e{\mathrm{\Ω}}_{\mathrm{ie}}^e}{6}-\frac{T^4{N}^e{F}^e}{24}\\ 0& \mathrm{TI}-\frac{T^2{\mathrm{\Ω}}_{\mathrm{ie}}^e}{2}\end{array}\right]$

$N^e=\frac{\mathrm{kM}}{r^3}\left[\begin{array}{ccc}-1+\frac{3x^{e^2}}{r^2}& \frac{3x^e{y}^e}{r^2}& \frac{3x^e{z}^e}{r^2}\\ \frac{3x^e{y}^e}{r^2}& -1+\frac{3y^{e^2}}{r^2}& \frac{{3y}^e{z}^e}{r^2}\\ \frac{3x^e{z}^e}{r^2}& \frac{3y^e{z}^e}{r^2}& -1+\frac{3z^{e^2}}{r^2}\end{array}\right]-{\mathrm{\Ω}}_{\mathrm{ie}}^e\·{\mathrm{\Ω}}_{\mathrm{ie}}^e$

${\mathrm{\Ω}}_{\mathrm{ie}}^e=\left[\begin{array}{ccc}0& -{\mathrm{\ω}}_{\mathrm{ie}}& 0\\ {\mathrm{\ω}}_{\mathrm{ie}}& 0& 0\\ 0& 0& 0\end{array}\right]$ $F^e=\left[\begin{array}{ccc}0& -f_z^e& f_y^e\\ f_z^e& 0& -f_x^e\\ -f_y^e& f_x^e& 0\end{array}\right]$

其中，x＝[δx δy δz δV_x δV_y δV_z ε_x ε_y ε_z]^T为捷联惯导误差状态方程的状态向量， $w=\left[\begin{array}{cccccc}\mathrm{\δ}{f}^e& {\mathrm{\δ\ω}}_{\mathrm{ib}}^e& {\mathrm{\δv}}_{\mathrm{ADC}}& {\mathrm{\δa}}_{\mathrm{mag}}& {\mathrm{\δp}}_{\mathrm{GPS}}& {\mathrm{\δv}}_{\mathrm{GPS}}\end{array}\right]$ 为地球坐标系组合导航系统噪声向量。δr＝[δx δy δz]^T为捷联惯导地球坐标系位置计算误差，单位m；δV＝[δV_x δV_y δV_z]^T为捷联惯导地球坐标系速度计算误差，单位m/s；ε＝[ε_x ε_y ε_z]^T为捷联惯导地球坐标系姿态误差角，单位rad。 $\mathrm{\δ}{f}^e={\left[\begin{array}{ccc}\mathrm{\δ}{f}_x^e& \mathrm{\δ}{f}_y^e& {\mathrm{\δf}}_z^e\end{array}\right]}^T$ 为捷联惯导地球坐标系加速度计加速度量测噪声，单位m/s²； ${\mathrm{\δ\ω}}_{\mathrm{ib}}^e={\left[\begin{array}{ccc}{\mathrm{\δ\ω}}_{\mathrm{ibx}}^e& {\mathrm{\δ\ω}}_{\mathrm{iby}}^e& {\mathrm{\δ\ω}}_{\mathrm{ibz}}^e\end{array}\right]}^T$ 为捷联惯导地球坐标系陀螺角速率量测噪声，单位rad/s。 $\mathrm{\δ}{p}_{\mathrm{GPS}}={\left[\begin{array}{ccc}{\mathrm{\δp}}_{{\mathrm{GPS}}_x}& {\mathrm{\δp}}_{{\mathrm{GPS}}_y}& \mathrm{\δ}{p}_{{\mathrm{GPS}}_z}\end{array}\right]}^T$ 为GPS接收机地球坐标系位置量测噪声，单位m， ${\mathrm{\δv}}_{\mathrm{GPS}}={\left[\begin{array}{ccc}\mathrm{\δ}{v}_{{\mathrm{GPS}}_x}& {\mathrm{\δv}}_{{\mathrm{GPS}}_y}& {\mathrm{\δv}}_{{\mathrm{GPS}}_z}\end{array}\right]}^T$ 为GPS接收机地球坐标系速度量测噪声，单位m/s； ${\mathrm{\δv}}_{\mathrm{ADC}}={\left[\begin{array}{ccc}{\mathrm{\δv}}_{{\mathrm{ADC}}_x}& {\mathrm{\δv}}_{{\mathrm{ADC}}_y}& {\mathrm{\δv}}_{{\mathrm{ADC}}_z}\end{array}\right]}^T$ 为大气数据计算机地球坐标系速度量测噪声，单位m/s。 $\mathrm{\δ}{a}_{\mathrm{mag}}={\left[\begin{array}{ccc}{\mathrm{\δ\α}}_{{\mathrm{mag}}_x}& {\mathrm{\δ\α}}_{{\mathrm{mag}}_y}& {\mathrm{\δ\α}}_{{\mathrm{mag}}_z}\end{array}\right]}^T$ 为磁传感姿态系统地球坐标系姿态量测噪声，单位rad。T为离散化周期，单位s。Ω_ie ^e为地球坐标系地球自转角速度的反对称矩阵，ω_ie为地球自转角速度。F^e为地球坐标系比力向量反对称矩阵，r^e＝[x^e y^e z^e]为捷联惯导地球坐标系位置向量，kM为地球万有引力常数与地球质量乘积， $r=\sqrt{x^{e^2}+y^{e^2}+z^{e^2}}.$ 地球坐标系是指固连在地球上的坐标系，原点在地球中心，它相对惯性坐标系以地球自转角速率旋转，Z轴沿地轴指向北极，X轴在赤道平面内指向格林威治子午线，Y轴在赤道平面内指向东经90度的方向；

(2)地球坐标系捷联惯导误差状态方程的多胞型描述的建立

矩阵A、B中变化项包括N^e和F^e，其中N^e变化十分缓慢，近似可以看作常数，则将误差状态方程表示成多胞型主要就是将F^e即地球坐标系载体比力量测值反对称矩阵用一个凸多边形表示。直接在地球坐标系描述F^e不方便，所以首先在当地水平坐标系将比力f^l的变化域用一个凸多面体代表，取该凸多面体的顶点 $f_i^l(i=\mathrm{1,2},\·\·\·)$ 并利用公式 $f^e=R_l^e\·f^l$ (R_l ^e为地理坐标系到地球坐标系坐标变换矩阵)将 $f_i^l(i=\mathrm{1,2},\·\·\·)$ 变换到地球坐标系中得到地球坐标系载体比力变化域顶点 $f_i^e(i=\mathrm{1,2},\·\·\·),$ 将 $f_i^e(i=\mathrm{1,2},\·\·\·)$ 的反对称矩阵F_i ^e代入惯性导航误差状态方程即得到多胞型描述误差传播方程的顶点Ω_i

(A_i，B_i)。以Ω_i

(A_i，B_i)为顶点的多胞型Ω(A，B)∈Υ为捷联惯性导航系统误差状态方程的多胞型描述。

建立地球坐标系捷联惯导误差状态方程的多胞型描述。首先，根据大气数据系统、磁传感姿态系统、GPS系统的数据输出频率确定组合周期T，一般为1秒。将初始位置代入N^e，将地球自转角速度ω_ie代入Ω_ie ^e。

其次，确定地球坐标系载体比力变化域顶点。假设地理坐标系水平方向机动加速度变化范围均为±a_EN，高度方向加速度变化范围为g±a_U(g为重力加速度)，则地理坐标系加速度计比力变化域顶点集为：

${\mathrm{\Ω}}_{f^l}=\left\{\begin{array}{cccccccc}{a}_{\mathrm{EN}}& a_{\mathrm{EN}}& a_{\mathrm{EN}}& a_{\mathrm{EN}}& {-a}_{\mathrm{EN}}& -a_{\mathrm{EN}}& -a_{\mathrm{EN}}& -a_{\mathrm{EN}}\\ a_U& a_U& -a_U& -a_U& a_U& a_U& -a_U& -a_U\\ g+a_U& g-a_U& g+a_U& g-a_U& g+a_U& g-a_U& g+a_U& g-a_U\end{array}\right\}$

由初始位置L、λ(L和λ代表地理纬度和经度)求初始点地理坐标系到地球坐标系坐标变换矩阵：

$R_l^e=\left[\begin{array}{ccc}-\sin \mathrm{\λ}& -\sin L\cos \mathrm{\λ}& \cos L\cos \mathrm{\λ}\\ \cos \mathrm{\λ}& -\sin L\sin \mathrm{\λ}& \cos L\sin \mathrm{\λ}\\ 0& \cos L& \sin L\end{array}\right]$

将Ω_fl每列分别左乘R_l ^e，所得8列向量为地球坐标系加速度计比力变化域8个顶点。将地球坐标系比力变化域顶点代入矩阵A、B中得到地球坐标系捷联惯性导航误差状态方程多胞型描速的顶点Ω_i (A_i，B_i)(i＝1，2，…8)，地球坐标系捷联惯性导航误差状态方程多胞型描述的顶点集为 $\stackrel{\‾}{\mathrm{\Ω}}=\{{\mathrm{\Ω}}_1,{\mathrm{\Ω}}_2,\·\·\·,{\mathrm{\Ω}}_8\}.$

2.以大气数据计算机、三轴磁强计输出的速度、姿态数据和捷联惯导输出速度、姿态数据的差值作为观测量建立内系统的量测方程；以GPS系统输出的位置、速度数据和捷联惯导位置、速度数据的差值为观测量建立外系统的量测方程。

以大气数据计算机、三轴磁强计输出的速度、姿态数据和捷联惯导输出速度、姿态数据的差值作为观测量建立内系统的量测方程为：

$y_1=\left[\begin{array}{ccc}0& I& 0\\ 0& 0& I\end{array}\right]\·\left[\begin{array}{c}\mathrm{\δr}\\ \mathrm{\δv}\\ \mathrm{\ϵ}\end{array}\right]+\left[\begin{array}{cc}I& 0\\ 0& I\end{array}\right]\·\left[\begin{array}{c}{\mathrm{\δv}}_{\mathrm{ADC}}\\ {\mathrm{\δv}}_a\end{array}\right]=C_1\·x+D_{11}\·w_1=C_1\·x+D_1\·w;$

以GPS系统输出的位置、速度数据和捷联惯导位置、速度数据的差值为观测量建立外系统的量测方程为：

$y_2=\left[\begin{array}{ccc}I& 0& 0\\ 0& I& 0\end{array}\right]\·\left[\begin{array}{c}\mathrm{\δr}\\ \mathrm{\δv}\\ \mathrm{\ϵ}\end{array}\right]+\left[\begin{array}{cc}I& 0\\ 0& I\end{array}\right]\·\left[\begin{array}{c}{\mathrm{\δp}}_{\mathrm{GPS}}\\ \mathrm{\δ}{v}_{\mathrm{GPS}}\end{array}\right]=C_2\·x+D_{21}\·w_2=C_2\·x+D_2\·w.$

3.对内系统量测，采用内系统的鲁棒静态输出反馈镇定器将不稳定的捷联惯导误差状态方程镇定，为保证系统整体性能，镇定器与下面所述滤波器一起使用凸优化方法同时求取。

4.对外系统量测，采用鲁棒H₂/H_∞/保性能滤波器，使得从干扰到导航误差满足所要求的H₂/H_∞/保性能要求。

鲁棒H₂/H_∞/保性能滤波器为：

${\mathrm{\Σ}}_f:\left\{\begin{array}{c}v(k+1)=A_f\·v\left(k\right)+B_f\·y_2\left(k\right)\\ y_f\left(k\right)=C_f\·v\left(k\right)+D_f\·y_2\left(k\right)\end{array}\right..$

5.内外系统量测方程、内系统镇定器、外系统鲁棒多目标滤波器联立构成闭环系统；同时为平衡对不同噪声输入通道进入的噪声的抑制效果，对噪声输入矩阵做加权处理。

(1)闭环系统为：

${\mathrm{\Σ}}_o:\left\{\begin{array}{c}\left[\begin{array}{c}x(k+1)\\ v(k+1)\end{array}\right]=\left[\begin{array}{cc}A-H\·C_1& 0\\ B_f\·C_2& A_f\end{array}\right]\·\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]+\left[\begin{array}{c}B-H\·D_1\\ B_f\·D_2\end{array}\right]\·w\left(k\right)\\ z_{\∞}\left(k\right)=L_{\∞}\·\left[\begin{array}{cc}I-D_f\·C_2& -C_f\end{array}\right]\·\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]-L_{\∞}\·D_f\·D_2\·w\left(k\right)\\ z_g\left(k\right)=L_g\·\left[\begin{array}{cc}I-D_f\·C_2& -C_f\end{array}\right]\·\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]\\ z_2\left(k\right)=L_2\·\left[\begin{array}{cc}I-D_f\·C_2& -C_f\end{array}\right]\·\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]-L_2\·D_f\·D_2\·w\left(k\right)\end{array}\right.;$

(2)噪声输入通道加权矩阵的确定

噪声统计特性很难获取但噪声强度却比较容易估计。对捷联惯导进行静态测试，采集一组加速度计和陀螺数据并计算其均方差记为E_acc和E_gyro；对大气数据系统进行静态测试，采集一组速度数据并计算其均方差记为E_vADC；对磁传感姿态系统进行静态测试，采集一组姿态数据并计算其均方差记为E_a；对GPS系统进行静态测试，采集一组位置、速度数据并计算其均方差记为E_pGPS、E_vGPS。

噪声输入通道加权矩阵Pow为15维对角矩阵，主对角线元素分别为E_acc、E_acc、E_acc、E_gyro、E_gyro、E_gyro、E_vADC、E_vADC、E_vADC、E_a、E_a、E_a、E_pGPS、E_pGPS、E_pGPS、E_vGPS、E_vGPS、E_vGPS。

6.利用凸优化算法求解内系统鲁棒静态输出反馈镇定器和外系统鲁棒混合H₂/H_∞/保性能滤波器。

(1)H₂、H_∞和保性能可调输出矩阵的选取

闭环系统的H₂、H_∞和保性能可调输出矩阵L₂、L_∞、L_g均为对角矩阵，为简单起见取L₂＝L_∞＝L_g。主对角线元素分别为E_pGPS ^-1、E_pGPS ^-1、E_pGPS ^-1、E_vGPS ^-1、E_vGPS ^-1、E_vGPS ^-1、E_a ^-1、E_a ^-1、E_a ^-1。

(2)综合性能加权参数α₁、α₂、α₃的选取

α₁为1；α₂取值在1-10之间；给定α₁、α₂，α₃在不显著增加系统的H₂和H_∞性能上界的情况下尽量取大。

(3)多目标多传感器组合导航方法存在的LMI条件

给定导航误差初始值x(0)的协方差矩阵的近似值为COV_int，多胞型顶点集

噪声通道加权矩阵，可调输出矩阵，综合性能加权参数，使用matlab里LMI工具箱中的mincx函数解以下LMIs优化问题

min(α₁·γ+α₂·Trace(∏)+α₃·Trace(COV_int·P₂₁·COV_int))

${\mathrm{\&Phi;}}_1=\left[\begin{array}{cccc}{\mathrm{\&Phi;}}_{111}& {\mathrm{\&Phi;}}_{12}& 0& {\mathrm{\&Phi;}}_{14}\\ {{\mathrm{\&Phi;}}_{12}}^T& {\mathrm{\&Phi;}}_{122}& {\mathrm{\&Phi;}}_{123}^T& 0\\ 0& {\mathrm{\&Phi;}}_{123}& -I& {\mathrm{\&Phi;}}_{134}\\ {{\mathrm{\&Phi;}}_{14}}^T& 0& {{\mathrm{\&Phi;}}_{134}}^T& -{\mathrm{\&gamma;}}_{1}I\end{array}\right]<0,$ $\&ForAll;\mathrm{\&Omega;}\stackrel{\mathrm{\&Delta;}}{=}(A,B)\&Element;\stackrel{\&OverBar;}{\mathrm{\&Omega;}}$

${\mathrm{\&Phi;}}_2=\left[\begin{array}{ccc}{\mathrm{\&Phi;}}_{211}& {\mathrm{\&Phi;}}_{12}& 0\\ {\mathrm{\&Phi;}}_{12}^T& {\mathrm{\&Phi;}}_{222}& {\mathrm{\&Phi;}}_{223}^T\\ 0& {\mathrm{\&Phi;}}_{223}& -I\end{array}\right]<0,$ $\&ForAll;\mathrm{\&Omega;}\stackrel{\mathrm{\&Delta;}}{=}(A,B)\&Element;\stackrel{\&OverBar;}{\mathrm{\&Omega;}}$

${\mathrm{\&Phi;}}_3=\left[\begin{array}{ccc}{\mathrm{\&Phi;}}_{311}& {\mathrm{\&Phi;}}_{12}& 0\\ {\mathrm{\&Phi;}}_{12}^T& {\mathrm{\&Phi;}}_{322}& {\mathrm{\&Phi;}}_{323}^T\\ 0& {\mathrm{\&Phi;}}_{323}& -I\end{array}\right]<0,$ $\&ForAll;\mathrm{\&Omega;}\stackrel{\mathrm{\&Delta;}}{=}(A,B)\&Element;\stackrel{\&OverBar;}{\mathrm{\&Omega;}}$

${\mathrm{\&Phi;}}_4=\left[\begin{array}{ccc}{\mathrm{\&Phi;}}_{211}& {\mathrm{\&Phi;}}_{14}& 0\\ {{\mathrm{\&Phi;}}_{14}}^T& -\mathrm{\&Pi;}& {{\mathrm{\&Phi;}}_{223}}^T\\ 0& {\mathrm{\&Phi;}}_{234}& -I\end{array}\right]<0,$ $\&ForAll;\mathrm{\&Omega;}\stackrel{\mathrm{\&Delta;}}{=}(A,B)\&Element;\stackrel{\&OverBar;}{\mathrm{\&Omega;}}$

${\mathrm{\&Phi;}}_5=\left[\begin{array}{cc}{P}_{22}-\mathrm{\&delta;}& 0\\ 0& W-G-G^T+\mathrm{\&delta;}\end{array}\right]<0$

其中：

${\mathrm{\&Phi;}}_{i11}=\left[\begin{array}{cc}{P}_{i11}-G-G^T& V\\ V^T& -W\end{array}\right],$ ${\mathrm{\&Phi;}}_{i22}=\left[\begin{array}{cc}-P_{i11}& -P_{12}\\ -{P_{12}}^T& -P_{22}\end{array}\right],$ i＝1，2，3

${\mathrm{\&Phi;}}_{12}=\left[\begin{array}{cc}\mathrm{GA}-K{C}_1+{\stackrel{\&OverBar;}{B}}_f{C}_2& \mathrm{GA}-K{C}_1+{\stackrel{\&OverBar;}{B}}_f{C}_2+{\stackrel{\&OverBar;}{A}}_f\\ \mathrm{GA}-K{C}_1& \mathrm{GA}-{\mathrm{KC}}_1\end{array}\right],$ ${\mathrm{\&Phi;}}_{14}=\left[\begin{array}{c}\mathrm{GB}-{\mathrm{KD}}_1+{\stackrel{\&OverBar;}{B}}_f{D}_2\\ \mathrm{GB}-K{D}_1\end{array}\right],$

${\mathrm{\&Phi;}}_{123}=\left[\begin{array}{cc}{L}_{\&infin;}-L_{\&infin;}{D}_f{C}_2& L_{\&infin;}-L_{\&infin;}{\stackrel{\&OverBar;}{C}}_f\end{array}\right],$ ${\mathrm{\&Phi;}}_{223}=\left[\begin{array}{cc}{L}_2-L_2{D}_f{C}_2& L_2-L_2{\stackrel{\&OverBar;}{C}}_f\end{array}\right],$

${\mathrm{\&Phi;}}_{323}=\left[\begin{array}{cc}{L}_g-L_g{D}_f{C}_2& L_g-L_g{\stackrel{\&OverBar;}{C}}_f\end{array}\right],$ Ф₁₃₄＝L_∞D_fD₂，Ф₂₃₄＝L₂D_fD₂

P₁₁₁，P₂₁₁，P₃₁₁，P₂₂，W，∏，δ为一组正定的LMI矩阵变量，P₁₂，K，V，G，

D_f为一组LMI矩阵变量，γ为LMI标量。

(4)鲁棒静态输出反馈和滤波器求解：

a.由Cholesky分解求X，其中XX^T＝P₂₂

b.由Cholesky分解求Q₁₁，[Q₁₁X-(X)^-T][Q₁₁X-(X)^-T]^T＝P₂₂ ^-1-G^-1WG^-T

c.将Q₁₁代入式 ${\stackrel{\&OverBar;}{P}}_{12}{Q_{11}}^{-T}=P_{12},$ 得到

d.由LU分解求G₁₂、Q₁₂，其中 $G_{12}\&CenterDot;Q_{12}^T={\stackrel{\&OverBar;}{P}}_{12}-G\&CenterDot;Q_{11}^T-I-V\&CenterDot;G^{-T}$

e.H＝G^-1K

f. $C_f=({\stackrel{\&OverBar;}{C}}_f-D_f{C}_2){Q_{11}}^T{Q_{12}}^{-T}$

g. $B_f={G_{12}}^{-1}{\stackrel{\&OverBar;}{B}}_f$

h. $A_f={G_{12}}^{-1}{\stackrel{\&OverBar;}{A}}_f{Q_{11}}^T{Q_{12}}^{-T}$

(5)系统性能的改进

为进一步改进性能，需要对Pow和综合性能加权参数α₁、α₂、α₃以及可调输出矩阵进行细调并通过实验确定一组合适的值。

a.对Pow和可调输出矩阵的调整。首先在0.1到10倍范围内调整Pow主对角线前三个元素使得其尽量大的同时又未造成综合性能指标明显升高，其次按相同方法调节Pow主对角线第4、5、6三个元素和第7、8、9三个元素。在经过上述调节的基础上，提高位置估计精度可以通过增大Pow主对角线第13、14、15元素或减小可调输出矩阵主对角线前三个元素实现；提高速度估计精度可以通过增大Pow主对角线后三个元素或减小可调输出矩阵主对角线第4、5、6三个元素实现；提高姿态估计精度可以通过增大Pow主对角线第10、11、12三个元素或减小可调输出矩阵主对角线后三个元素实现。

b.对综合性能加权参数α₁、α₂、α₃的调整。原则为固定α₁，改变α₂、α₃对α₁的比例，通过实验确定一组合适的综合性能加权参数。

本发明说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。

Claims (8)

1.一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于包括以下步骤：

(1)将地球坐标系捷联惯导误差状态方程中参数的变化作为系统的不确定性，并使用不确定系统的多胞模型来描述地球坐标系捷联惯导误差状态方程；

(2)以大气数据计算机、三轴磁强计输出的速度、姿态数据和捷联惯导输出速度、姿态数据的差值作为观测量建立内系统的量测方程；以GPS系统输出的位置、速度数据和捷联惯导位置、速度数据的差值为观测量建立外系统的量测方程；

(3)对内系统量测，采用内系统的鲁棒静态输出反馈镇定器将不稳定的捷联惯导误差状态方程镇定；

(4)对外系统量测，采用鲁棒H₂/H_∞/保性能滤波器，使得从干扰到导航误差满足所要求的H₂/H_∞/保性能要求；

(5)将地球坐标系捷联惯导误差状态方程、内外系统量测方程、内系统镇定器、外系统混合多目标滤波器联立构成闭环系统；同时为平衡对不同噪声输入通道进入的噪声的抑制效果，对噪声输入矩阵做加权处理；

(6)利用凸优化算法求解内系统鲁棒静态输出反馈镇定器和外系统鲁棒混合H₂/H_∞/保性能滤波器。

2.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(1)中的捷联惯导地球坐标系误差状态方程的多胞型描述方法为：捷联惯导误差状态方程的多胞型描述主要是将时变的捷联惯导误差状态方程中参数的变化看作系统的不确定性，将时变的捷联惯导误差状态方程看作不确定系统，并使用鲁棒控制中常采用的多胞模型来描述；其中主要是确定载体加速度变化域，取其顶点并将其代入捷联惯导误差状态方程得到多胞型描述的捷联惯导误差状态方程的顶点，进而得到多胞型描述的捷联惯导误差状态方程。

3.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(2)中的内系统量测方程建立如下：

内系统量测y1以大气数据计算机、三轴磁强计的速度、姿态数据和捷联惯导速度、姿态数据的差值作为量测值，量测方程为：

$y_1=\left[\begin{array}{ccc}0& I& 0\\ 0& 0& I\end{array}\right]\&CenterDot;\left[\begin{array}{c}\mathrm{\&delta;r}\\ \mathrm{\&delta;v}\\ \mathrm{\&epsiv;}\end{array}\right]+\left[\begin{array}{cc}I& 0\\ 0& I\end{array}\right]\&CenterDot;\left[\begin{array}{c}\mathrm{\&delta;}{v}_{\mathrm{ADC}}\\ \mathrm{\&delta;}{v}_a\end{array}\right]=C_1\&CenterDot;x+D_{11}\&CenterDot;w_1=C_1\&CenterDot;x+D_1\&CenterDot;w$

其中，状态x＝[δr^Tδv^Tε^T]^T为系统的状态向量，δr，δv，ε分别表示位置、速度、姿态的差值；

为系统的噪声变量，δv_ADC，δv_a分别表示速度误差、姿态误差的量测噪声；D₁₁为噪声变量的增益阵，I为适当维数的单位阵；

$C_1=\left[\begin{array}{ccc}0& I& 0\\ 0& 0& I\end{array}\right],$ $D_{11}=\left[\begin{array}{cc}I& 0\\ 0& I\end{array}\right],$ $w=\left[\begin{array}{c}{w}_1\\ w_2\end{array}\right],$ D₁＝[D₁₁ 0]。

4.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(2)中外系统量测方程的建立如下：

外系统量测y₂以GPS系统输出的位置、速度数据和捷联惯导位置、速度数据的差值为量测值，量测方程为：

$y_2=\left[\begin{array}{ccc}I& 0& 0\\ 0& I& 0\end{array}\right]\&CenterDot;\left[\begin{array}{c}\mathrm{\&delta;r}\\ \mathrm{\&delta;v}\\ \mathrm{\&epsiv;}\end{array}\right]+\left[\begin{array}{cc}I& 0\\ 0& I\end{array}\right]\&CenterDot;\left[\begin{array}{c}\mathrm{\&delta;}{p}_{\mathrm{GPS}}\\ \mathrm{\&delta;}{v}_{\mathrm{GPS}}\end{array}\right]=C_2\&CenterDot;x+D_{21}\&CenterDot;w_2=C_2\&CenterDot;x+D_2\&CenterDot;w$

其中，

为系统的噪声变量，δP_GPS，δv_GPS分别表示位置误差、速度误差的量测噪声；D₂₁为噪声变量的增益阵；

$C_2=\left[\begin{array}{ccc}I& 0& 0\\ 0& I& 0\end{array}\right],D_{21}=\left[\begin{array}{cc}I& 0\\ 0& I\end{array}\right],D_2=\left[\begin{array}{cc}0& D_{21}\end{array}\right].$

5.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(3)中的鲁棒静态输出反馈镇定器设计如下：设计静态输出反馈镇定器∑将原系统镇定：

∑：x(k+1)＝A·x(k)+B·w(k)+H·y₁(k)

其中，A为系统阵，B为噪声增益阵，H为输出反馈镇定器∑的增益阵；在此只确定镇定器形式，为保证系统整体性能，镇定器∑与滤波器∑_f一起使用凸优化方法同时求取。

6.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(4)中的鲁棒H₂/H_∞/保性能滤波器∑_f的设计如下：

${\mathrm{\&Sigma;}}_f:\left\{\begin{array}{c}v(k+1)=A_f\&CenterDot;v\left(k\right)+B_f\&CenterDot;y_2\left(k\right)\\ y_f\left(k\right)=C_f\&CenterDot;v\left(k\right)+D_f\&CenterDot;y_2\left(k\right)\end{array}\right.$

其中，v为滤波器状态变量，y_f为滤波器的输出变量；A_f，B_f，C_f，D_f为滤波器的增益阵；使得y_f尽可能逼近导航误差x。

7.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(5)中的闭环系统的联立以及对噪声输入矩阵的加权处理方法如下：

(1)闭环系统状态方程

将地球坐标系捷联惯性导航系统误差状态方程、内外系统量测方程、内系统镇定器、外系统多目标滤波器联立构成闭环系统状态方程：

$\left[\begin{array}{c}x(k+1)\\ v(k+1)\end{array}\right]=\left[\begin{array}{cc}A-H\&CenterDot;C_1& 0\\ B_f\&CenterDot;C_2& A_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]+\left[\begin{array}{c}B-H\&CenterDot;D_1\\ B_f\&CenterDot;D_2\end{array}\right]\&CenterDot;w\left(k\right)$

其中，x为输出反馈系统∑的状态变量，A为系统的系统阵，B为噪声增益阵，H为输出反馈镇定器的增益阵；v为滤波器∑_f的状态变量，A_f，B_f为滤波器的增益阵；

D₁＝[D₁₁0]，D₂＝[0D₂₁]；

(2)闭环系统H₂/H_∞/保性能可调输出

为提高滤波器对捷联惯导导航误差的估计，将滤波器输出与捷联惯导导航误差之差作为闭环系统输出，并协调H₂/H_∞/保性能各项性能指标对滤波精度的不同要求，确定可调输出矩阵，得闭环系统的H₂/H_∞/保性能可调输出：

$\left\{\begin{array}{c}{z}_{\&infin;}\left(k\right)=L_{\&infin;}\&CenterDot;\left[\begin{array}{cc}I-D_f\&CenterDot;C_2& -C_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]-L_{\&infin;}\&CenterDot;D_f\&CenterDot;D_2\&CenterDot;w\left(k\right)\\ z_g\left(k\right)=L_g\&CenterDot;\left[\begin{array}{cc}I-D_f\&CenterDot;C_2& -C_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]\\ z_2\left(k\right)=L_2\&CenterDot;\begin{array}{}\end{array}\left[\begin{array}{cc}I-D_f\&CenterDot;C_2& -C_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]-L_2\&CenterDot;D_f\&CenterDot;D_2\&CenterDot;w\left(k\right)\end{array}\right.$

其中z_∞、z_g、z₂分别为H_∞、保性能和H₂性能的可调输出，L_∞、L_g、L₂分别为H_∞、保性能和H₂性能可调输出矩阵；C_f，D_f为滤波器∑_f的增益阵；

(3)噪声输入通道加权处理

闭环系统状态方程和闭环系统可调输出构成闭环系统∑_o：

${\mathrm{\&Sigma;}}_o:\left[\begin{array}{c}\left[\begin{array}{c}x(k+1)\\ v(k+1)\end{array}\right]=\left[\begin{array}{cc}A-H\&CenterDot;C_1& 0\\ B_f\&CenterDot;C_2& A_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]+\left[\begin{array}{c}B-H\&CenterDot;D_1\\ B_f\&CenterDot;D_2\end{array}\right]\&CenterDot;w\left(k\right)\\ z_{\&infin;}\left(k\right)=L_{\&infin;}\&CenterDot;\left[\begin{array}{cc}I-D_f\&CenterDot;C_2& -C_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]-L_{\&infin;}\&CenterDot;D_f\&CenterDot;D_2\&CenterDot;w\left(k\right)\\ z_g\left(k\right)=L_g\&CenterDot;\left[\begin{array}{cc}I-D_f\&CenterDot;C_2& -C_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]\\ z_2\left(k\right)=L_2\&CenterDot;\left[\begin{array}{cc}I-D_f\&CenterDot;C_2& -C_f\end{array}\right]\&CenterDot;\left[\begin{array}{c}x\left(k\right)\\ v\left(k\right)\end{array}\right]-L_2\&CenterDot;D_f\&CenterDot;D_2\&CenterDot;w\left(k\right)\end{array}\right]$

为平衡对不同噪声输入通道进入的噪声的抑制效果，需预先对噪声输入矩阵做加权处理，即B＝B·Pow，D₁＝D₁·Pow，D₂＝D₂·Pow，Pow为系统噪声输入通道重要性加权矩阵，根据对各通道噪声的干扰抑制需要确定各输入通道的重要性加权系数；Pow为对角阵，对角线相应元素为相应噪声输入通道的重要性加权系数。

8.根据权利要求1所述的一种基于凸优化算法的多目标多传感器信息融合方法，其特征在于：所述步骤(6)中的利用凸优化算法求解组合导航系统多目标多传感器信息融合方法的内系统鲁棒静态输出反馈镇定器和外系统鲁棒混合H₂/H_∞/保性能滤波器为：将问题用一组线性矩阵不等式LMI表示，并编制相应程序使用Matlab的线性矩阵不等式工具箱LMIToolbox中lmivar、lmiterm、mincx函数求解，最后通过代数变换将所得LMI决策变量转化为内系统鲁棒静态输出反馈镇定器和外系统的鲁棒混合H₂/H_∞/保性能滤波器。

Images