一种优化的地球椭球模型下极区最优阻尼方法

本发明属于惯性导航，涉及惯性导航系统阻尼方法，特别涉及一种优化的地球椭球模型下极区最优阻尼方法，适用于舰船的长航时极区航行。

背景技术：

1、极地地区尤其是北极地区在资源、科研、航道等方面具有重要战略价值，以美国、俄罗斯为代表的各个大国加快极地地区的战略布局及考察研究。由于惯性导航拥有极好的自主性，并且不受极地地区恶劣环境的影响，因此惯性导航已经成为极区重要的导航手段。

2、舰船极区航行导航通常要求导航设备具备长航时导航能力，然而惯性导航系统的误差随时间积累，长航时导航必然导致定位精度下降。此外，纯惯性导航系统的水平通道是无阻尼的，系统产生了三种周期性振荡误差，即舒勒振荡、傅科振荡、地球自转周期振荡，在长航时导航中振荡误差严重影响导航定位精度。根据现代控制理论分析，纯惯性导航的系统误差方程极点存在于复平面的虚轴上，为无阻尼振荡系统，需要进行极点配置实现系统阻尼，提高导航精度。

3、传统机械编排在极区会出现计算溢出、失去航向参考等难题，横向坐标系导航方案常被用于解决这极区导航中的难题，然而横向坐标系导航方案无法解决纯惯性导航系统的振荡误差，极地纯惯性导航同样需要系统阻尼。目前，大多横向坐标系阻尼方案的研究采用phi角误差模型并且预设固定的阻尼网络，然而在地球椭球模型的基础上构建系统误差方程使得该方案所涉及的参数误差求解十分繁琐，固定的阻尼网络设计使得该阻尼系统在面临导航系统需求变化时出现系统阻尼超调。

4、本发明针对目前存在的问题，面向长航时极区航海导航，必须重点解决：1.优化地球椭球模型下横向坐标系机械编排方案；2.变阻尼网络的设计，实现系统自动配置合适极点，使惯性导航系统工作稳定；因此本发明提出一种优化的地球椭球模型下极区最优阻尼方法：通过引入测速仪为外速度参考，构建了基于psi角误差模型的导航方案，在地球椭球模型下使得误差模型更加简洁，使得极点配置更加简单；通过卡尔曼最优估计设计变阻尼网络，实现荡误差的阻尼，抑制误差积累，提高阻尼系统稳定性，具有十分重要的工程意义。

技术实现思路

1、在极区长航时航海导航领域，更高的自主导航精度需要更精准的地球模型和惯性阻尼系统。地球椭球模型致使以phi角误差模型为基础的系统误差参数计算繁琐，系统误差分析难度增大。此外传统阻尼系统常常基于定阻尼网络设计，因此该系统往往不能适应导航需求的变化进而产生超调振荡误差。本发明要解决的技术问题就在于：针对现有技术的不足，本发明提出一种优化的地球椭球模型下极区最优阻尼方法，利用psi角误差模型定义在计算坐标系的特点，简化椭球模型下的系统误差方程；采用现代控制理论对系统状态空间进行分析，配置合适的极点使系统工作在稳定状态，采用卡尔曼最优估计方法实现变阻尼系数，提高系统稳定性。

2、为解决上述技术问题，本发明提出的解决方案为：

3、一种优化的地球椭球模型下极区最优阻尼方法，所述方法包括以下步骤：

4、(1)定义横地球坐标系，定义横向极点，定义横向经度和横向纬度，确定横向位置表示方式：所述横地球坐标系e′的原点位于地心，x轴沿着地球自转轴指向北极，y轴指向本初子午线与赤道的交点，z轴穿过东经90°子午线与赤道的交点；定义(0°,90°e)为横向北极点、(0°,90°w)为横向南极点；定义0°经线和180°经线组成的大椭圆为横向赤道；定义90°e和90°w北半球部分组成的半个大椭圆为0°横经线，且横向本初子午线为地理经度90°e所在的子午圈的北半球部分，横向子午线为过横向极点的平面与地球表面相交的轮廓线；定义地球表面上一点的地理法线与横向赤道面交角为该点的横向纬度；定义该点所在的横向子午面与横向本初子午面的交角为横向经度；根据构建的横经纬网络，将舰船在横地球坐标系中位置表示为(lt，λt，h)，其中，lt表示横纬度，λt表示横经度，h表示高度；

5、(2)定义横地理坐标系：横地理坐标系t的原点位于载体中心，y轴沿横向经线的切线指向横向北极点，z轴垂直于当地水平面指向天向，x轴与y轴和z轴构成右手坐标系，且为“横东-横北-天向”定义；

6、(3)确定坐标系之间的转换关系，步骤如下：

7、根据所述步骤(1)中横地球坐标系定义，确定地球坐标系e到横地球坐标系e′的方向余弦矩阵为：

8、

9、确定地球坐标系e到地理坐标系g的方向余弦矩阵为：

10、

11、其中l表示舰船所处的纬度，λ表示舰船所处的经度；

12、确定横地球坐标系e′到横地理坐标系t的方向余弦矩阵

13、

14、根据链式法则，确定地理坐标系g到横地理坐标系t的方向余弦矩阵

15、

16、式中表示为方向余弦矩阵的转置；σ表示横地理坐标系t与地理坐标系g之间的夹角，具体表示为：

17、

18、

19、(4)利用惯性导航获得载体姿态、速度、位置相关信息，确定横地理坐标系下的姿态更新方程、速度更新方程、位置更新方程，具体步骤如下：

20、(4.1)确定横地理坐标系下的姿态更新方程：

21、

22、式中，表示从载体坐标系b到横地理坐标系t的方向余弦矩阵；表示载体坐标系b相对于惯性坐标系i的旋转角速度在载体坐标系b下的投影；表示横地理坐标系t相对于惯性坐标系i的旋转角速度在横地理坐标系t下的投影；

23、其中：表示地球坐标系e相对于惯性坐标系i的旋转角速度在横地理坐标系t下的投影，表示横地理坐标系t相对于地球坐标系e的旋转角速度在横地理坐标系t下的投影，具体表示为：

24、

25、

26、其中，表示地球自转角速度在地球坐标系e下的投影，ωie表示地球自转角速度的大小；表示从地球坐标系e到横地理坐标系t的方向余弦矩阵；表示地理坐标系g相对于地球坐标系e的旋转角速度在地理坐标系g下的投影；表示横地理坐标系t相对于地理坐标系g的旋转角速度在横地理坐标系t下的投影；

27、和具体表示为：

28、

29、

30、其中，和分别表示地理坐标系下载体的北向速度和东向速度；rm表示载体处的子午圈半径，rn表示载体处的卯酉圈半径，具体表示为：

31、

32、其中，re表示地球长半轴半径，ρ表示地球的偏心率；

33、确定

34、

35、式中，和分别为横地理坐标系t下载体东向和北向的速度；表示横地理坐标系下载体处的扭曲率，和分别为横地理东向和横地理北向的曲率，具体表示为：

36、

37、

38、

39、(4.2)确定横地理坐标系下的速度vt的更新方程：

40、

41、式中，vt表示横地理坐标系t下的载体速度；fb表示载体坐标系b下表示的比力；gt表示横地理坐标系t下表示的重力矢量；

42、(4.3)确定横地理坐标系下的位置更新方程：

43、步骤(4.1)所述的参数其变化由横向经度、横向纬度的变化所引起，具体表示为：

44、

45、对比步骤(4.1)所确定的参数确定横向经度、横向纬度微分方程：

46、

47、

48、高度变化由天向速度所引起，确定高度微分方程：

49、

50、式中，表示横地理坐标系t下载体的天向速度；

51、(5)确定计算坐标系、平台坐标系与真实导航坐标系之间的转换关系：

52、确定计算坐标系c到平台坐标系p的方向余弦矩阵为：

53、

54、确定真实导航坐标系tg到平台坐标系p的方向余弦矩阵为：

55、

56、确定真实导航坐标系tg到计算地理坐标系c的方向余弦矩阵为：

57、

58、式中i3×3表示3×3的单位矩阵；ψ为漂移误差角，φ为姿态误差角，δθ为位置误差角；确定漂移误差角、姿态误差角、位置误差角之间的关系为：

59、φ＝ψ+δθ

60、(6)确定计算横地理坐标系下阻尼系统的误差模型：

61、(6.1)确定计算横地理坐标系下阻尼系统的状态方程：

62、

63、其中xc'为系统状态向量；f为状态转移矩阵；g为系统噪声分配矩阵；w为系统噪声向量；k为阻尼系数矩阵；u为阻尼系统反馈矩阵；

64、将系统状态向量xc'表示为：

65、

66、其中，表示三维漂移误差角矢量在计算横地理坐标系c'下的投影，各分量分别为计算横地理坐标系c'下东向、北向、天向的漂移误差角；表示三维速度误差矢量在计算横地理坐标系c'下的投影，各分量分别为计算横地理坐标系c'下东向、北向、天向的速度误差；表示三维位置误差矢量在计算横地理坐标系c'下的投影，各分量分别为计算横地理坐标系c'下东向、北向、天向的位置误差；表示陀螺的零偏矢量，各分量分别为x、y、z轴陀螺的零偏；表示加速度计的零偏矢量，各分量分别为x、y、z轴加速度计的零偏；δk表示测速仪标度因数误差；δη、δγ表示测速仪的俯仰角安装误差、方位角安装误差；

67、(6.2)确定计算横地理坐标系下阻尼系统的误差方程：

68、(6.2.1)确定计算横坐标系下惯性导航系统的姿态、速度和位置误差方程：

69、

70、

71、

72、式中，表示计算横地理坐标系c'下地球自转角速度，表示计算横地理坐标系c'下计算横地理坐标系c'相对地球坐标系e的角速度，表示从载体坐标系b到计算横地理坐标系c'的方向余弦矩阵，fc'表示计算横地理坐标系c'下表示的比力；

73、式中分别表示载体坐标系b下的陀螺误差和加速度计误差，表示为：

74、

75、

76、式中，和分别表示陀螺和加速度计的噪声；

77、(6.2.2)确定陀螺零偏、加速度计零偏、测速仪标度因数误差、测速仪俯仰角安装误差、测速仪方位角安装误差的误差方程：

78、

79、

80、

81、式中，τε和分别表示陀螺和加速度计的一阶马尔可夫相关时间，wε和分别表示陀螺和加速度计的高斯白噪声；

82、将计算横地理坐标系下阻尼系统表示为：

83、

84、

85、其中，表示惯性导航系统的速度估计值与测速仪的速度输出值；03×3表示3×3的零矩阵；

86、(7)确定阻尼系统的阻尼系数矩阵k：

87、步骤(6)所述阻尼系数矩阵k与卡尔曼最优估计的增益系数一致，即建立卡尔曼最优估计的状态方程和观测方程即获得阻尼系数矩阵k；由于卡尔曼最优估计是依据惯性导航系统不同时刻不同状态的特征进行估计，因此其增益系数为时变的，进而实现系统的变阻尼网络；

88、(7.1)确定卡尔曼最优估计系统状态方程：

89、

90、其中系统状态向量xc'、状态转移矩阵f、系统噪声分配矩阵g、系统噪声向量w与步骤(6)所述相同；

91、(7.2)确定卡尔曼最优估计系统观测方程：

92、

93、其中，为测速仪安装误差矩阵；vc'表示计算横地理坐标系c'下的速度矢量；υ为测量噪声向量；h为状态观测矩阵，将其表示为：

94、h＝[-[vc'×] i3×3 03×3 03×3 03×3 -vc' h1 h2]

95、式中，h1和h2分别为矩阵的第一列和第三列向量；

96、(7.3)确定阻尼系数矩阵k：

97、根据卡尔曼最优估计方法，其中预测过程为：

98、

99、更新过程为：

100、

101、式中，为系统第n-1时刻的误差状态，为系统第n时刻预测的误差状态，为系统第n时刻的误差状态；fn/n-1为从n-1时刻至n时刻的状态转移矩阵；为系统第n-1时刻误差状态的协方差矩阵，为系统第n时刻预测的误差状态的协方差矩阵，为系统第n时刻误差状态的协方差矩阵；qn-1为系统第n-1时刻噪声向量的协方差矩阵；bn-1为系统第n-1时刻系统噪声分配矩阵；hn为系统第n时刻系统的观测矩阵；rn为系统第n时刻测量噪声向量的协方差矩阵；zn为系统第n时刻系统测量误差向量；kn为系统第n时刻的卡尔曼最优估计增益矩阵，即系统阻尼系数矩阵；i18×18表示18×18的单位矩阵；

102、(7.4)确定系统误差状态修正方式：卡尔曼最优估计的系统状态向量定义在计算横地理坐标系中，通过对系统状态进行修正使其定义在横地理坐标系t中：

103、φt＝ψc'+δθt

104、δvt＝δvc'-δθt×vc'

105、

106、式中，φt表示姿态误差角φ在横地理坐标系t中的投影；δvt为速度误差δv在横地理坐标系t中的投影；表示计算横地理坐标系c'到横地理坐标系t的方向余弦矩阵；δθt表示位置误差角δθ在横地理坐标系t中的投影，具体表示为：

107、

108、式中，分别为δrt的北向分量、东向分量；δrt为位置误差δr在横地理坐标系t中的投影，其中经纬度误差转换为位置误差δrt表示为：

109、δrt＝[(rn+h)δλtcoslt (rm+h)δlt δh]t

110、位置误差δrt转换为经纬度误差表示为：

111、

112、其中为位置误差δrt的天向分量。

113、进一步的，本发明的卡尔曼最优估计对系统的姿态误差、速度误差、位置误差、陀螺和加速度计零偏采用闭环反馈，测速仪标度因数误差和安装误差采用开环反馈，且每次闭环反馈的系统误差状态校正后置0。

114、进一步的，若载体接收到其他传感器的位置信息，包括但不限于gnss位置信息、重力匹配位置信息、地磁匹配位置信息，基于接收到的位置信息对转换关系或进行修正更新。

115、与现有技术相比，本发明具有以下优点：

116、本发明基于psi角误差模型简化了地球椭球模型下横向坐标系导航系统的误差方程，利用卡尔曼最优估计实现阻尼系数的时变估计与极点自动配置，实现了纯惯性导航的振荡误差抑制，提高系统的稳定性；本发明构建的系统误差方程更加简洁，变阻尼网络设计可以增加阻尼系统的稳定性，有利于提高长航时导航精度。