[发明专利]一种惯导系统误差仿真及修理辅助分析方法

权利要求书

1.一种惯导系统误差仿真及修理辅助分析方法，其特征在于：包括以下步骤：

步骤1、建立惯性导航系统误差方程；

步骤2、根据步骤1所建立的惯性导航系统误差方程，设置初始条件，对惯导系统误差进行仿真及修理辅助分析，并根据仿真结果对影响惯导系统精度的决定性因素进行分析。

2.根据权利要求1所述的一种惯导系统误差仿真及修理辅助分析方法，其特征在于：所述步骤1的具体步骤包括：

(1)建立静基座误差方程：

(2)写成状态方程矩阵形式并进行拉氏变换得到：

3.根据权利要求1所述的一种惯导系统误差仿真及修理辅助分析方法，其特征在于：所述步骤2的具体步骤包括：

(1)设置初始条件

初始条件：设陀螺的常值漂移为0.1°/h，加速度计的常值零偏误差为0.0001g，当地纬度39°；起始条件误差为：速度误差0.1m/s，位置误差为0.0005°，水平姿态误差为20角秒，方位姿态误差为5角分。运行时间t＝24小时；

(2)对惯导系统误差进行仿真及修理辅助分析，推导出陀螺漂移引起的系统误差，并画出误差传播曲线；推导出加速度计零偏引起的系统误差，并画出误差传播曲线；推导出初始条件误差引起的系统误差，并画出误差传播曲线；

(3)根据仿真结果对影响惯导系统精度的决定性因素进行分析。

4.根据权利要求3所述的一种惯导系统误差仿真及修理辅助分析方法，其特征在于：所述步骤2第(2)步的具体步骤包括：

①建立系数阵为：

用Matlab求特征矩阵C＝(sI-F)^-1有：

C₁₁＝(R*s)/(R*s^2+g)

C₁₂＝0

C₁₃＝(R*g*s*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₁₄＝(R*g*s*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₁₅＝-(R*g*(s^2-wie^2*sin(L)^2+wie^2))/((s^2+wie^2)*(R*s^2+g))

C₁₆＝-(R*g*wie^2*sin(2*L))/(2*(s^2+wie^2)*(R*s^2+g))

C₂₁＝0

C₂₂＝(R*s)/(R*s^2+g)

C₂₃＝-(R*g*wie^2)/((s^2+wie^2)*(R*s^2+g))

C₂₄＝(R*g*s^2)/((s^2+wie^2)*(R*s^2+g))

C₂₅＝(R*g*s*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₂₆＝-(R*g*s*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

C₃₁＝0

C₃₂＝1/(R*s^2+g)

C₃₃＝(R*s^3+(R*wie^2+g)*s)/((s^2+wie^2)*(R*s^2+g))

C₃₄＝(g*s)/((s^2+wie^2)*(R*s^2+g))

C₃₅＝(g*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₃₆＝-(g*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

C₄₁＝0

C₄₂＝-1/(R*s^2+g)

C₄₃＝-(R*s*wie^2)/((s^2+wie^2)*(R*s^2+g))

C₄₄＝(R*s^3)/((s^2+wie^2)*(R*s^2+g))

C₄₅＝(R*s^2*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₄₆＝-(R*s^2*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

C₅₁＝1/(R*s^2+g)

C₅₂＝0

C₅₃＝-(R*s^2*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₅₄＝-(R*s^2*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

C₅₅＝(R*s*(s^2-wie^2*sin(L)^2+wie^2))/((s^2+wie^2)*(R*s^2+g))

C₅₆＝(R*s*wie^2*sin(2*L))/(2*(s^2+wie^2)*(R*s^2+g))

C₆₁＝tan(L)/(R*s^2+g)

C₆₂＝0

C₆₃＝(wie*(R*s^2*cos(L)^2+g))/(cos(L)*(s^2+wie^2)*(R*s^2+g))

C₆₄＝(wie*(R*s^2*cos(L)^2+g))/(cos(L)*(s^2+wie^2)*(R*s^2+g))

C₆₅＝-(s*sin(L)*(g-R*wie^2*cos(L)^2))/(cos(L)*(s^2+wie^2)*(R*s^2+g))

C₆₆＝(R*s^3+(R*wie^2*sin(L)^2+g)*s)/((s^2+wie^2)*(R*s^2+g))

各式中，常数值g表示重力加速度，ws表示舒勒角频率，wie表示地球自转角速度，R表示地球半径，L表示当地纬度；

Dx表示东向加速度计常值零偏，Dy表示北向加速度计常值零偏，Ex表示东向陀螺漂移，Ey表示北向陀螺漂移，Ez表示方位陀螺漂移；

Dx0表示初始东向速度误差，Dy0表示初始北向速度误差，L0表示初始纬度位置误差，N0表示初始经度位置误差，Qx0表示初始东向姿态误差，Qy0表示初始北向姿态误差，Qz0表示初始方位姿态误差；

dVx表示系统东向速度误差，dVy表示系统北向速度误差；dL表示纬度误差，dN表示经度误差；Qx、Qy、Qz表示平台误差；

②计算东向陀螺漂移εx引起的系统误差及误差传播曲线；

dVxs＝(Ex*R*g*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

dVxt＝(Ex*R*g*sin(t*wie)*sin(L))/(g-R*wie^2)-(Ex*R^(3/2)*g^(1/2)*wie*sin((g^(1/2)*t)/R^(1/2))*sin(L))/(g-R*wie^2)

dVys＝(Ex*R*g*s)/((s^2+wie^2)*(R*s^2+g))

dVyt＝(Ex*R*g*cos(t*wie))/(g-R*wie^2)-(Ex*R*g*cos((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2)

dLs＝(Ex*g)/((s^2+wie^2)*(R*s^2+g))

dLt＝(Ex*g*sin(t*wie))/(wie*(g-R*wie^2))-(Ex*R^(1/2)*g^(1/2)*sin((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2)

dNs＝(Ex*g*wie*sin(L))/(cos(L)*(s^2+wie^2)*(R*s^2+g))

dNt＝((Ex*g*sin(t*wie)*sin(L))/(g-R*wie^2)-(Ex*R^(1/2)*g^(1/2)*wie*sin((g^(1/2)*t)/R^(1/2))*sin(L))/(g-R*wie^2))/cos(L)

Qxs＝(Ex*R*s^2)/((s^2+wie^2)*(R*s^2+g))

Qxt＝(Ex*R^(1/2)*g^(1/2)*sin((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2)-(Ex*R*wie*sin(t*wie))/(g-R*wie^2)

Qys＝-(Ex*R*s*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

Qyt＝(Ex*R*wie*cos((g^(1/2)*t)/R^(1/2))*sin(L))/(g-R*wie^2)-(Ex*R*wie*cos(t*wie)*sin(L))/(g-R*wie^2)

Qzs＝(Ex*wie*(R*s^2*cos(L)^2+g))/(s*cos(L)*(s^2+wie^2)*(R*s^2+g))

Qzt＝-((Ex*cos(t*wie)*(g-R*wie^2*cos(L)^2))/(wie*(g-R*wie^2))-Ex/wie+(Ex*R*wie*cos((g^(1/2)*t)/R^(1/2))*(cos(L)-1)*(cos(L)+1))/(g-R*wie^2))/cos(L)

③计算北向陀螺漂移εy引起的系统误差及误差传播曲线；

dVxs＝-(Ey*R*g*(s^2-wie^2*sin(L)^2+wie^2))/(s*(s^2+wie^2)*(R*s^2+g))

dVxt＝Ey*R*(sin(L)-1)*(sin(L)+1)+(Ey*R*cos((g^(1/2)*t)/R^(1/2))*(R*wie^2*sin(L)^2-R*wie^2+g))/(g-R*wie^2)-(Ey*R*g*cos(t*wie)*sin(L)^2)/(g-R*wie^2)

dVys＝(Ey*R*g*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

dVyt＝(Ey*R*g*sin(t*wie)*sin(L))/(g-R*wie^2)-(Ey*R^(3/2)*g^(1/2)*wie*sin((g^(1/2)*t)/R^(1/2))*sin(L))/(g-R*wie^2)

dLs＝(Ey*g*wie*sin(L))/(s*(s^2+wie^2)*(R*s^2+g))

dLt＝(Ey*sin(L))/wie+(Ey*R*wie*cos((g^(1/2)*t)/R^(1/2))*sin(L))/(g-R*wie^2)-(Ey*g*cos(t*wie)*sin(L))/(wie*(g-R*wie^2))

dNs＝-(Ey*g*(s^2-wie^2*sin(L)^2+wie^2))/(s*cos(L)*(s^2+wie^2)*(R*s^2+g))

dNt＝(Ey*(sin(L)-1)*(sin(L)+1)+(Ey*cos((g^(1/2)*t)/R^(1/2))*(R*wie^2*sin(L)^2-R*wie^2+g))/(g-R*wie^2)-(Ey*g*cos(t*wie)*sin(L)^2)/(g-R*wie^2))/cos(L)

Qxs＝(Ey*R*s*wie*sin(L))/((s^2+wie^2)*(R*s^2+g))

Qxt＝(Ey*R*wie*cos(t*wie)*sin(L))/(g-R*wie^2)-(Ey*R*wie*cos((g^(1/2)*t)/R^(1/2))*sin(L))/(g-R*wie^2)

Qys＝(Ey*R*(s^2-wie^2*sin(L)^2+wie^2))/((s^2+wie^2)*(R*s^2+g))

Qyt＝(Ey*R^(1/2)*sin((g^(1/2)*t)/R^(1/2))*(R*wie^2*sin(L)^2-R*wie^2+g))/(g^(1/2)*(g-R*wie^2))-(Ey*R*wie*sin(t*wie)*sin(L)^2)/(g-R*wie^2)

Qzs＝-(Ey*sin(L)*(g-R*wie^2*cos(L)^2))/(cos(L)*(s^2+wie^2)*(R*s^2+g))

Qzt＝-((Ey*sin(t*wie)*sin(L)*(g-R*wie^2*cos(L)^2))/(wie*(g-R*wie^2))-(Ey*R^(1/2)*sin((g^(1/2)*t)/R^(1/2))*sin(L)*(g-R*wie^2*cos(L)^2))/(g^(1/2)*(g-R*wie^2)))/cos(L)

④计算方位陀螺漂移εz引起的系统误差及误差传播曲线；

dVxs＝-(Ez*R*g*wie^2*sin(2*L))/(s*(R*s^2+g)*(2*s^2+2*wie^2))

dVxt＝(Ez*R*g*sin(2*L)*cos(t*wie))/(2*g-2*R*wie^2)-(Ez*R^2*wie^2*sin(2*L)*cos((g^(1/2)*t)/R^(1/2)))/(2*g-2*R*wie^2)-(Ez*R*sin(2*L))/2

dVys＝-(Ez*R*g*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

dVyt＝(Ez*R^(3/2)*g^(1/2)*wie*cos(L)*sin((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2)-(Ez*R*g*cos(L)*sin(t*wie))/(g-R*wie^2)

dLs＝-(Ez*g*wie*cos(L))/(s*(s^2+wie^2)*(R*s^2+g))

dLt＝(Ez*g*cos(t*wie)*cos(L))/(wie*(g-R*wie^2))-(Ez*R*wie*cos((g^(1/2)*t)/R^(1/2))*cos(L))/(g-R*wie^2)-(Ez*cos(L))/wie

dNs＝-(Ez*g*wie^2*sin(2*L))/(s*cos(L)*(R*s^2+g)*(2*s^2+2*wie^2))

dNt＝-(Ez*sin(2*L)-(Ez*g*sin(2*L)*cos(t*wie))/(g-R*wie^2)+(Ez*R*wie^2*sin(2*L)*cos((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2))/(2*cos(L))

Qxs＝-(Ez*R*s*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

Qxt＝(Ez*R*wie*cos((g^(1/2)*t)/R^(1/2))*cos(L))/(g-R*wie^2)-(Ez*R*wie*cos(t*wie)*cos(L))/(g-R*wie^2)

Qys＝(Ez*R*wie^2*sin(2*L))/((R*s^2+g)*(2*s^2+2*wie^2))

Qyt＝(Ez*R*wie*sin(2*L)*sin(t*wie))/(2*g-2*R*wie^2)-(Ez*R^(3/2)*wie^2*sin(2*L)*sin((g^(1/2)*t)/R^(1/2)))/(g^(1/2)*(2*g-2*R*wie^2))

Qzs＝(Ez*(R*s^3+(R*wie^2*sin(L)^2+g)*s))/(s*(s^2+wie^2)*(R*s^2+g))

Qzt＝(Ez*sin(t*wie)*(R*wie^2*sin(L)^2-R*wie^2+g))/(wie*(g-R*wie^2))-(Ez*R^(3/2)*wie^2*sin((g^(1/2)*t)/R^(1/2))*sin(L)^2)/(g^(1/2)*(g-R*wie^2))

⑤计算东向加速度计零偏引起的系统误差及误差传播曲线；

dVxs＝(Dx*R)/(R*s^2+g)

dVxt＝(Dx*R^(1/2)*sin((g^(1/2)*t)/R^(1/2)))/g^(1/2)

dVys＝0

dVyt＝0

dLs＝0

dLt＝0

dNs＝Dx/(cos(L)*(R*s^2+g))

dNt＝(Dx*sin((g^(1/2)*t)/R^(1/2)))/(R^(1/2)*g^(1/2)*cos(L))

Qxs＝0

Qxt＝0

Qys＝Dx/(s*(R*s^2+g))

Qyt＝Dx/g-(Dx*cos((g^(1/2)*t)/R^(1/2)))/g

Qzs＝(Dx*tan(L))/(s*(R*s^2+g))

Qzt＝(Dx*tan(L))/g-(Dx*cos((g^(1/2)*t)/R^(1/2))*tan(L))/g

⑥计算初始方位姿态误差Φ_z0引起的系统误差及误差传播曲线；

dVxs＝-(Qz0*R*g*wie^2*sin(2*L))/((R*s^2+g)*(2*s^2+2*wie^2))

dVxt＝(Qz0*R^(3/2)*g^(1/2)*wie^2*sin(2*L)*sin((g^(1/2)*t)/R^(1/2)))/(2*g-2*R*wie^2)-(Qz0*R*g*wie*sin(2*L)*sin(t*wie))/(2*g-2*R*wie^2)

dVys＝-(Qz0*R*g*s*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

dVyt＝(Qz0*R*g*wie*cos((g^(1/2)*t)/R^(1/2))*cos(L))/(g-R*wie^2)-(Qz0*R*g*wie*cos(t*wie)*cos(L))/(g-R*wie^2)

dLs＝-(Qz0*g*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

dLt＝(Qz0*R^(1/2)*g^(1/2)*wie*cos(L)*sin((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2)-(Qz0*g*cos(L)*sin(t*wie))/(g-R*wie^2)

dNs＝-(Qz0*g*wie^2*sin(2*L))/(cos(L)*(R*s^2+g)*(2*s^2+2*wie^2))

dNt＝-((Qz0*g*wie*sin(2*L)*sin(t*wie))/(g-R*wie^2)-(Qz0*R^(1/2)*g^(1/2)*wie^2*sin(2*L)*sin((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2))/(2*cos(L))

Qxs＝-(Qz0*R*s^2*wie*cos(L))/((s^2+wie^2)*(R*s^2+g))

Qxt＝(Qz0*R*wie^2*cos(L)*sin(t*wie))/(g-R*wie^2)-(Qz0*R^(1/2)*g^(1/2)*wie*cos(L)*sin((g^(1/2)*t)/R^(1/2)))/(g-R*wie^2)

Qys＝(Qz0*R*s*wie^2*sin(2*L))/((R*s^2+g)*(2*s^2+2*wie^2))

Qyt＝(Qz0*R*wie^2*sin(2*L)*cos(t*wie))/(2*g-2*R*wie^2)-(Qz0*R*wie^2*sin(2*L)*cos((g^(1/2)*t)/R^(1/2)))/(2*g-2*R*wie^2)

Qzs＝(Qz0*(R*s^3+(R*wie^2*sin(L)^2+g)*s))/((s^2+wie^2)*(R*s^2+g))

Qzt＝(Qz0*cos(t*wie)*(R*wie^2*sin(L)^2-R*wie^2+g))/(g-R*wie^2)-(Qz0*R*wie^2*cos((g^(1/2)*t)/R^(1/2))*sin(L)^2)/(g-R*wie^2)。