[发明专利]基于后数据的旋转制导炮弹快速空中对准方法

摘要

本发明提供一种基于后数据的旋转制导炮弹快速空中对准方法，该方法中对准的位置和速度由卫星导航系统给出，并利用卫星导航输出的速度信息解算出对应时刻的航向角和俯仰角，然后根据设定时刻T到对准时刻T1的俯仰角、航向角变化率、俯仰角变化率和INS数据中的陀螺输出角速度，确定横滚角观测方程的系数矩阵，并通过最小二乘法求解横滚角观测方程，从而实现对惯导初始位置、速度和姿态角的精确解算，即实现自旋制导炮弹的快速空中对准，大大提高了自旋制导炮弹的落点精度；该对准算法简单、精度高，而且对准时间短，对准速度快，为提高自旋制导炮弹的落点精度和缩短打击时间做出了重要铺垫。

主权项

基于后数据的旋转制导炮弹快速空中对准方法，其特征在于包括如下步骤：(1)、发射后的制导炮弹接收卫星导航信号进行导航处理，其中实现卫星导航信号捕获跟踪并输出导航结果的时刻为T0；然后从所述时刻T0到设定的对准时刻T1保存卫星导航系统输出的M组卫星导航结果和INS输出的N组INS数据，其中：TGPS为卫星导航结果输出周期，TINS为INS数据输出周期，且TGPS＝Q×TINS，Q为正整数；所述卫星导航结果包括制导炮弹的速度和位置；所述INS数据包括前向陀螺、左向陀螺和上向陀螺输出的角速度，其中，前向陀螺敏感弹体横滚角角速度，左向陀螺敏感俯仰角角速度，上向陀螺敏感航向角角速度；(2)、根据M组卫星导航结果中的制导炮弹的速度，计算得到相应的M组航向角和俯仰角；(3)、对步骤(2)计算得到的M组航向角和俯仰角的计算结果进行拟合计算，得到航向角和俯仰角在时刻T0到时刻T1之间随时间变换的函数；然后对所述函数进行求导运算，得到航向角变化率和俯仰角变化率的时间函数；(4)、将输出N组INS数据的时刻值代入到步骤(3)确定的4个时间函数中，计算得到输出N组INS数据时的制导炮弹的航向角、俯仰角、航向角变化率和俯仰角变化率；(5)、根据设定时刻T到时刻T1的Np′组俯仰角、航向角变化率、俯仰角变化率和INS数据中的陀螺输出角速度，计算观测方程Z＝H×X中的观测矩阵H和测量矩阵Z；其中X为两维的观测向量，X(1)为制导炮弹横滚角的正弦值，X(2)为制导炮弹横滚角的余弦值；其中，设定时刻T的取值范围为T0≤T<T1，正整数(6)、利用最小二乘法对观测方程Z＝H×X进行求解，得到观测向量X＝(HTH)‑1HTZ；(7)、根据观测向量X计算结果中的制导炮弹横滚角的正弦值和余弦值，计算得到制导炮弹的横滚角；(8)、将步骤(7)计算得到的横滚角，以及时刻T1的卫星导航系统输出的速度、位置，以及根据所述速度计算得到的航向角和俯仰角，作为空中对准结果，输出到制导炮弹的导航系统，用于对所述制导炮弹进行导航和控制；在步骤(5)中，根据时刻T～T1的Np′组俯仰角、航向角变化率和俯仰角变化率，以及INS数据中的陀螺输出角速度，计算观测矩阵H和测量矩阵Z，具体计算过程如下：(5a)、对观测矩阵H和测量矩阵Z进行初始化，得到初始观测矩阵H0和测量矩阵Z0：如果初始化H0＝[a(T1‑TINS)b(T1‑TINS)]，则Z0＝z(T1‑TINS)；如果初始化H0＝[‑b(T1‑TINS)a(T1‑TINS)]，则Z0＝z′(T1‑TINS)；其中：a(T1-TINS)=cos[(ωx(T1-TINS)+▿φgz(T1-TINS)sin(φgx(T1-TINS)))×TINS];]]>b(T1-TINS)=sin[(ωx(T1-TINS)+▿φgz(T1-TINS)sin(φgx(T1-TINS)))×TINS];]]>z(T1-TINS)=ωy(T1-TINS)▿φgz(T1-TINS)cos(φgx(T1-TINS))-ωz(T1-TINS)▿φgx(T1-TINS)▿φgx(T1-TINS)2+(▿φgz(T1-TINS)cos(φgx(T1-TINS)))2;]]>z′(T1-TINS)=ωz(T1-TINS)▿φgz(T1-TINS)cos(φgx(T1-TINS))+ωy(T1-TINS)▿φgx(T1-TINS)▿φgx(T1-TINS)2+(▿φgz(T1-TINS)cos(φgx(T1-TINS)))2;]]>其中：ωx(T1‑TINS)、ωy(T1‑TINS)和ωz(T1‑TINS)分别为时刻T1‑TINS前向陀螺、左向陀螺和上向陀螺输出的角速度测量值；φgx(T1‑TINS)为时刻T1‑TINS的俯仰角；和分别为时刻T1‑TINS的航向角变化率和俯仰角变化率；(5b)、在时刻Tn′＝T1‑(n+1)×TINS，n＝1、2、…Np′‑1，按照如下的迭代公式对观测矩阵H和测量矩阵Z进行迭代更新，得到时刻Tn′的观测矩阵Hn和测量矩阵Zn；如果Hn＝[Hn‑1；a(Tn′),b(Tn′)]，则Zn＝[Zn‑1；z(Tn′)]；如果Hn＝[Hn‑1；‑b(Tn′),a(Tn′)]，则Zn＝[Zn‑1；z′(Tn′)]；其中：a(Tn′)=cos[Σm=0n(ωx(Tm′)+▿φgz(Tm′)sin(φgx(Tm′)))×TINS];]]>b(Tn′)=sin[Σm=0n(ωx(Tm′)+▿φgz(Tm′)sin(φgx(Tm′)))×TINS];]]>z(Tn′)=ωy(Tn′)▿φgz(Tn′)cos(φgx(Tn′))-ωz(Tn′)▿φgx(Tn′)▿φgx(Tn′)2+(▿φgz(Tn′)cos(φgx(Tn′)))2;]]>z′(Tn′)=ωz(Tn′)▿φgz(Tn′)cos(φgx(Tn′))+ωy(Tn′)▿φgx(Tn′)▿φgx(Tn′)2+(▿φgz(Tn′)cos(φgx(Tn′)))2;]]>其中：ωx(Tm′)、ωy(Tm′)和ωz(Tm′)分别为时刻Tm′前向陀螺、左向陀螺和上向陀螺输出的角速度测量值，φgx(Tm′)和分别为时刻Tm′的俯仰角和航向角变化率，m＝0～n且T0′＝T1‑TINS；φgx(Tn′)为时刻Tn′的俯仰角；和分别为时刻Tn′的航向角变化率和俯仰角变化率；ωx(Tn′)、ωy(Tn′)和ωz(Tn′)分别为时刻Tn′前向陀螺、左向陀螺和上向陀螺输出的角速度测量值；(5c)、将步骤(5b)迭代得到的时刻T＝T1‑Np′×TINS的观测矩阵和测量矩阵作为观测矩阵H和测量矩阵Z的最终计算结果，用于步骤(6)的最小二乘计算。