CN113962057B - 基于时序交会的远程导弹主动段运动参数修正方法 - Google Patents

基于时序交会的远程导弹主动段运动参数修正方法 Download PDF

Info

Publication number
CN113962057B
CN113962057B CN202110727620.5A CN202110727620A CN113962057B CN 113962057 B CN113962057 B CN 113962057B CN 202110727620 A CN202110727620 A CN 202110727620A CN 113962057 B CN113962057 B CN 113962057B
Authority
CN
China
Prior art keywords
missile
ballistic
point
satellite
observation
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.)
Active
Application number
CN202110727620.5A
Other languages
English (en)
Other versions
CN113962057A (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.)
Nanjing University of Aeronautics and Astronautics
Original Assignee
Nanjing University of Aeronautics and Astronautics
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 Nanjing University of Aeronautics and Astronautics filed Critical Nanjing University of Aeronautics and Astronautics
Priority to CN202110727620.5A priority Critical patent/CN113962057B/zh
Publication of CN113962057A publication Critical patent/CN113962057A/zh
Application granted granted Critical
Publication of CN113962057B publication Critical patent/CN113962057B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/20—Analysis of motion
    • G06T7/277—Analysis of motion involving stochastic approaches, e.g. using Kalman filters
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30241—Trajectory

Abstract

本发明提出了一种基于时序交会的远程导弹主动段运动参数修正方法。该方法包括：根据弹道导弹目标点坐标到预警卫星红外传感器像点坐标的转换关系，提取在地心固定坐标系下卫星中心到目标弹道导弹的方向单位观测视线矢量；对各个弹道平面间存在的几何关系进行分析，并用这些弹道平面对上述提取的观测视线矢量进行切割，建立弹道面切割模型；构建基于动量守恒的主动段动力学模型，并对备选弹道平面进行筛选；根据图像信息设计转换矩阵，将所有像平面进行转换，得到多个虚拟像平面坐标；对多个虚拟像平面上的像点进行卡尔曼预测；根据预测结果所形成的观测视线对约束进行修正。本发明能有效实现在单星观测背景下弹道导弹主动段参数的估计。

Description

基于时序交会的远程导弹主动段运动参数修正方法
技术领域
本发明涉及摄影测量技术、红外图像处理技术和轨迹预测技术领域，具体涉及一种基于时序交会的远程导弹主动段运动参数修正方法。
背景技术
天基红外预警系统主要采用星载红外传感器，实现对弹道导弹目标发射的监视、跟踪和参数估计等功能。红外传感器通过被动探测，获得弹道导弹的角观测量，无法进行距离与速度测量，对导弹定位来说属于不完备观测。且在高轨天基预警系统的发展过程中，由于预警星座部署的阶段性再加上系统建设前期的验证与试验，大部分情况下预警探测处于单星探测。
主动段是弹道导弹的空间运动参数估计的关键，对于弹道导弹主动参数估计模型，主要的研究方法是利用先验知识，基于弹道模版库的弹道建模。由于弹道导弹类型有限，因此可以事先将不同类型目标主动段弹道模版存储在数据库中，称为标准弹道数据库(Nominal Ballistic Profile)。现有方法中备选弹道的参数通常根据经验或者遍历来设置，存在难以确定备选弹道平面以及计算效率低的问题。
发明内容
发明目的：本发明提供一种基于时序交会的远程导弹主动段运动参数修正方法，运用卡尔曼滤波和时序交会原理，解决现有技术中心难以确定备选弹道平面的问题。
技术方案：本发明所述的一种基于时序交会的远程导弹主动段运动参数修正方法，包括以下步骤：
(1)根据弹道导弹目标点坐标到预警卫星红外传感器预警图像像点坐标的转换关系，提取在地心固定坐标系下卫星中心到目标弹道导弹的方向单位观测视线矢量；
(2)对各个弹道平面间存在的几何关系进行分析，并用这些弹道平面对上述提取的观测视线矢量进行切割，建立弹道面切割模型；
(3)根据弹道导弹在主动飞行时的受力情况、运动特性以及其共性特征参数构建基于动量守恒的主动段动力学模型，并对备选弹道平面进行筛选；
(4)根据时序转换像平面坐标系，得到多个虚拟像平面坐标系；
(5)利用卡尔曼滤波对像平面坐标系的目标信息进行时序二维轨迹预测；
(6)根据预测二维轨迹的时序观测视线修正轨道。
进一步地，所述步骤(1)包括以下步骤：
(11)基于共线条件方程，对预警卫星红外传感器预警图像像点坐标原始图像建立严密几何成像模型：
Figure BDA0003138093620000021
其中，(x y)为成像模型得到的像平面坐标，(x0 y0)为预警卫星红外传感器预警图像像主点坐标，f为相机焦距，(Xs Ys Zs)为预警卫星在地心固定坐标系下的坐标，(X Y Z)为弹道导弹在地心固定坐标系下的坐标，ak,bk,ck(k＝1,2,3)为预警卫星系统参数的函数，具体表示如下：
Figure BDA0003138093620000022
式中，RJ2000-WGS84为J2000坐标系到WGS84坐标系的旋转矩阵；Rorbit-J2000为预警卫星的轨道测量坐标系到J2000坐标系的旋转矩阵；Rbody-orbit为预警卫星本体坐标系到轨道测量坐标系的转换矩阵；Rcamera-body为传感器坐标系到本体坐标系的变换矩阵；
(12)建立观测卫星中心到目标弹道导弹的方向单位观测视线矢量(u v w)T与导弹的像平面坐标(x,y)之间的数学函数关系：
Figure BDA0003138093620000023
其中，m、λ为比例系数；(x,y)T为弹道导弹成像的像平面坐标；(x0,y0,-f)T为像点对应成像时的内方位元素；根据该关系式得到观测视线矢量(u v w)T，u,v,w分别为X方向、Y方向、Z方向上的观测视线矢量。
进一步地，所述步骤(2)包括如下步骤：
(21)引入预警卫星首发点和终点在地心固定坐标系下的坐标，分别为(X1 Y1 Z1)T和(XN YN ZN)T，所述首发点为预警卫星首次探测发现弹道导弹的位置，所述终点为最后探测发现弹道导弹的位置，则有：
Figure BDA0003138093620000031
其中，(u1 v1 w1)T为首发点对应的观测视线矢量，(uN vN wN)T为终点对应的观测视线矢量，(XS YS ZS)T为预警卫星在地心固定坐标系下的坐标，re为平均地球半径，h1和hN分别为首发点和终点高程，假定观测卫星的星载红外传感器监视时，探测到主动段弹道导弹发动机的尾焰，共有N组观测值，分别对应于观测时刻t1,t2,……,tN；
(22)根据弹道基本约束中的通视约束解算出(X1 Y1 Z1)T和(XN YN ZN)T，所述通视约束指的是弹道平面切割时观测视线矢量不能先于弹道导弹穿过地球，进而确定一个过地心的弹道切割平面。
进一步地，所述步骤(3)包括以下步骤：
(31)忽略空气作用以及柯氏惯性力和牵连惯性力的影响，根据弹道导弹在主动段动量的改变为所受合外力的冲量，构建弹道导弹的加速度模型：
Figure BDA0003138093620000032
其中，p″(t)为导弹在t时刻的加速度，M(t)为t时刻的导弹质量，M′(t)为t时刻导弹质量的变化率，g0为重力加速度，u为导弹高温燃气的相对有效喷射速度；
(32)在步骤(31)的基础上，引入火箭发动机秒耗量a，记s＝a/M，建立弹道导弹主动段动力学模型：
Figure BDA0003138093620000033
其中，p(t)表示导弹在t时刻的位置，p0,v0分别表示位置与速度的初始条件；假定高温燃气相对于导弹以恒定的有效喷射速度为uc，gc为重力加速度常量，下标c表示假定的恒定值；
(33)利用弹道导弹的主动段动力学模型对备选弹道平面进行筛选，将不符合该动力学模型的弹道进行排除，得到符合导弹特征的轨迹对象。
进一步地，所述步骤(4)中利用如下转换矩阵进行坐标系转换：
Figure BDA0003138093620000041
其中各参数分别为：
Figure BDA0003138093620000042
Figure BDA0003138093620000043
Figure BDA0003138093620000044
b1＝cosωsinκ
b2＝cosωcosκ
b3＝-sinω
Figure BDA0003138093620000045
Figure BDA0003138093620000046
Figure BDA0003138093620000047
其中
Figure BDA0003138093620000048
ω、κ分别为卫星的三个姿态角：俯仰角、滚动角和偏航角。
进一步地，所述步骤(5)包括以下步骤：
(51)在步骤(4)所得的虚拟像平面上进行卡尔曼预测，卡尔曼滤波递推公式如下：
状态预测方程
Figure BDA0003138093620000049
状态测量方程
Figure BDA00031380936200000410
均方差预测方程
Figure BDA00031380936200000411
卡尔曼增益方程
Figure BDA00031380936200000412
状态修改方程
Figure BDA00031380936200000413
均方差修改方程
Figure BDA00031380936200000414
其中k代表时间点，矩阵
Figure BDA0003138093620000051
表示k时间点不考虑估计误差的状态，矩阵
Figure BDA0003138093620000052
表示k时间点的状态，矩阵zk表示k时间点测量得到的数据，矩阵Ak表示状态由k-1时间点到k时间点的状态转移矩阵，矩阵
Figure BDA0003138093620000053
表示k时间点的均方差矩阵，矩阵Hk表示卡尔曼增益矩阵；
(52)定义卡尔曼滤波跟踪器系统状态
Figure BDA0003138093620000054
为一个四维向量(Px,Vx,Py,Vy)T，Px，Vx，Py，Vy分别是目标在像平面坐标系下X轴和Y轴方向上的位置与速度；定义观测向量zk＝(Px,Py)T；
定义状态转移矩阵Ak：
Figure BDA0003138093620000055
其中ΔT为卫星相机拍摄时间间隔；
由系统状态和观测状态的关系可知，观测矩阵Ck为：
Figure BDA0003138093620000056
令递推公式中：
Figure BDA0003138093620000057
其中delta＝0.01，r＝100；
卡尔曼滤波根据某一个虚拟像平面内所有像平面坐标进行多次计算后，依据最终Xk得到一个预测的像平面坐标点。
进一步地，所述步骤(6)包括如下步骤：
(61)记卡尔曼滤波计算出的预测点为(x0,y0,z0)，观测到的实际点为(x′0,y′0,z′0)，将观测矢量设为参数方程：
Figure BDA0003138093620000061
其中i1、j1、k1、i2、j2、k2为卫星-导弹向量参数，由卫星与导弹的坐标求出；
(62)根据前方交会原理，预测的时序观测向量与真实的时序观测向量的交点即为导弹位置，为了准确计算两空间直线交点，设两条直线的公垂线与两直线的交点为A(i1m+x0,j1m+y0,k1m+z0)，B(i2n+x′0,j2n+y′0,k2n+z′0)，其中m、n为待求参数，根据公垂线性质，得到：
Figure BDA0003138093620000062
(63)将i1、j1、k1、i2、j2、k2、预测点(x0,y0,z0)和实际点(x′0,y′0,z′0)代入步骤(62)中的计算公式，求出参数m、n的值；将m、n代入步骤(61)中的计算公式，式求出A、B两点三维信息；
(64)计算交点A、B的中点C，所有时间点的C为一个坐标集Ci，约束后轨道与Ci进行比对后修改参数，实现修正轨迹。
有益效果：
1、本发明根据切割弹道平面与观测视线之间几何关系以及动量定理建立导弹运动物理模型，构建强约束的弹道导弹主动段参数估计模型，解析出预警卫星观测空间目标过程中的几何变化规律，可以摆脱参数估计模型对先验信息的依赖，只依靠像平面坐标和卫星信息即可计算得出许多较为准确的备选弹道平面。
2、本发明运用卡尔曼滤波和时序交会原理，利用卡尔曼滤波空间轨迹修正模型，对弹道轨迹进行修正，突破缺少先验信息支持造成参数估计不准确的难题。
附图说明
图1为本发明的基于时序交会的远程导弹主动段运动参数修正方法流程图。
具体实施方式
下面结合附图对本发明作进一步详细描述。本发明提供一种基于时序交会的远程导弹主动段运动参数修正方法，如图1所示，包括如下步骤：
步骤一：根据弹道导弹目标点坐标到预警卫星红外传感器预警图像像点坐标的转换关系，提取在地心固定坐标系下卫星中心到目标弹道导弹的方向单位观测视线矢量。具体包括以下步骤：
(11)基于共线条件方程，对预警卫星红外传感器预警图像像点坐标原始图像建立严密几何成像模型：
Figure BDA0003138093620000071
其中，(x y)为成像模型得到的像平面坐标，(x0 y0)为影像像主点坐标，f为相机焦距，(Xs Ys Zs)为预警卫星中心在地心固定坐标系下的坐标，该地固坐标由测轨设备测量获取，这个对应于平台的轨道模型，涉及轨道动力学和卫星摄动等方面；(X Y Z)为弹道导弹在地心固定坐标系下的坐标，ak,bk,ck(k＝1,2,3)为预警卫星系统参数的函数，具体表示如下：
Figure BDA0003138093620000072
式中，RJ2000-WGS84为J2000坐标系到WGS84坐标系的旋转矩阵；Rorbit-J2000为预警卫星的轨道测量坐标系到J2000坐标系的旋转矩阵，由姿态测量设备获取；Rbody-orbit为预警卫星本体坐标系到轨道测量坐标系的转换矩阵；Rcamera-body为传感器坐标系到本体坐标系的变换矩阵，即相机安装矩阵；
(12)建立观测卫星中心到目标弹道导弹的方向单位观测视线矢量(u v w)T(即测向信息)与导弹的像平面坐标(x,y)之间的数学函数关系：
Figure BDA0003138093620000081
其中，m、λ为比例系数，相当于放缩系数；(x,y)T为弹道导弹成像的像平面坐标；(x0,y0,-f)T为像点对应成像时的内方位元素。即可得到观测视线矢量(u v w)T，其中，u,v,w分别为X方向、Y方向、Z方向上的观测视线矢量。
步骤二：对各个弹道平面间存在的几何关系进行分析，并用这些弹道平面对上述提取的观测视线矢量进行切割，建立弹道面切割模型。具体包括以下步骤：
(21)引入GEO预警卫星首次探测发现弹道导弹(即首发点)和最后探测发现弹道导弹(即终点)在地心固定坐标系下的坐标，分别为(X1 Y1 Z1)T和(XN YN ZN)T，则有：
Figure BDA0003138093620000082
其中，(u1 v1 w1)T为首发点对应的观测视线矢量，(uN vN wN)T为终点对应的观测视线矢量，(XS YS ZS)T为GEO预警卫星在地心固定坐标系下的坐标，re为平均地球半径，h1和hN分别为首发点和终点高程，假定GEO观测卫星的星载红外传感器监视时，探测到主动段弹道导弹发动机的尾焰，共有N组观测值，分别对应于观测时刻t1,t2,……,tN；
(22)根据弹道基本约束中的通视约束，即弹道平面切割时，观测视线矢量不能先于弹道导弹穿过地球，可解算出(X1 Y1 Z1)T和(XN YN ZN)T，进而可以确定一个过地心的弹道切割平面。
基于弹道平面遍历切割的快速切割模型方法假设所观测的弹道导弹没有机动，则弹道导弹主动段一直在一个近似过地球球心的弹道平面内运动，即满足平面假设。用这些弹道平面切割观测视线矢量(ui vi wi)T,i＝1,...N，所有切割出的曲线均有可能为弹道曲线，即为备选弹道。
(a)确定首发点高程和终点高程的区间范围。首发点高程区间[h1min，h1max]可根据首发点预警图像的数据以及当地气象数据得到；终点高程区间的下限边界可选取与首发高程区间的上限边界相同，即hNmin＝h1max，其上界可由hNmax＝-g0·T2/2+T·umax得到。式中，T＝tN-t1为预警卫星从首发点到终点的总观测时长，g0为重力加速度，umax为最大燃气喷速值。
(b)设置h1和hN的遍历步长，遍历所有可能的候选弹道切割平面，并求解出每一个切割平面内的备选弹道轨迹数据。
步骤三：对弹道导弹运动特点进行分析来筛选弹道导弹轨迹。
(31)忽略空气作用以及柯氏惯性力和牵连惯性力的影响，根据弹道导弹在主动段动量的改变为所受合外力的冲量，构建弹道导弹的加速度模型：
Figure BDA0003138093620000091
其中，p″(t)为导弹在t时刻的加速度，M(t)为t时刻导弹质量，M′(t)为t时刻导弹质量的变化率，变化率指的是在t时刻，经过时间微元后，导弹质量的变化量与时间微元的比值。g0为重力加速度，u为导弹高温燃气的相对有效喷射速度。
(32)在步骤(31)的基础上，引入火箭发动机秒耗量a，记s＝a/M，建立弹道导弹主动段动力学模型：
Figure BDA0003138093620000092
其中，p(t)表示导弹在t时刻的位置，即弹道模型，p0,v0分别表示位置与速度的初始条件；假定高温燃气相对于导弹以恒定的有效喷射速度为uc，gc为重力加速度常量，下标c表示假定的恒定值；在弹道导弹主动段，飞行时间段，可忽略因弹道导弹空间位置变化所引起的重力加速度可忽略不计，可将重力加速度做常量处理，记作gc。
(3)利用弹道导弹的动力学模型对备选弹道平面进行严格约束，将不符合弹道导弹动力学模型的弹道进行排除，得到符合导弹特征的轨迹对象。
步骤四：根据图像信息设计转换矩阵，将所有像平面进行转换，得到多个虚拟像平面坐标。
转换多时间点像平面坐标系，得到多个虚拟像平面坐标系。具体转换矩阵为：
Figure BDA0003138093620000101
其中各参数分别为：
Figure BDA0003138093620000102
Figure BDA0003138093620000103
Figure BDA0003138093620000104
b1＝cosωsinκ
b2＝cosωcosκ
b3＝-sinω
Figure BDA0003138093620000105
Figure BDA0003138093620000106
Figure BDA0003138093620000107
其中
Figure BDA0003138093620000108
ω、κ分别为卫星的三个姿态角：俯仰角、滚动角和偏航角，三个旋转角可以根据实际需要修改。
根据实际虚拟像平面需求修改
Figure BDA0003138093620000109
ω、κ的值，即可使步骤一中像平面坐标[x0 y0 -f]T转为虚拟像平面坐标[x′ y′ -f]T：
Figure BDA00031380936200001010
每一个虚拟像平面都有相同数量的坐标点，数值大约为10-20之间，具体数量视实际情况而定。
步骤五：利用卡尔曼滤波对像平面坐标系的目标信息进行二维轨迹预测。具体包括以下步骤：
(51)在步骤四所得到虚拟像平面上进行卡尔曼预测。卡尔曼滤波是一种利用线性系统状态方程，通过系统输入输出观测数据，对系统状态进行最优估计的算法，是目前应用最为广泛的滤波方法，在通信、导航、制导与控制等多领域得到了较好的应用。焦平面上的成像可近似为一个点，点的位置即为卡尔曼系统的系统输入。通过卡尔曼滤波，可以抑制系统噪声，并且得到最优位置估计。
卡尔曼滤波不会保存过去的数据，每次获得新的数据后，根据新的数据和上一时刻的参数估计值，接着借助于系统本身的状态转移方程，按照一套推导出来的递推公式，就能够算出新的参数估计值。卡尔曼滤波装置的存储量和计算量因为这个使用算法特性而大大减少，并且突破了平稳随机过程的限制。
一个离散的动态系统的n维状态方程和m维测量方程分别为：
xk＝Akxk-1+wk-1 (1.1)
yk＝Ckxk+vk (1.2)
式(1.1)和式(1.2)中，矩阵Ak表示状态由k-1时间点到k时间点的状态转移矩阵，Ck称作观测矩阵，yk是测量得到的数据，wk-1和vk均为噪声误差矩阵。Ak，Ck，yk是已知的，利用xk-1的估计值
Figure BDA0003138093620000111
求得k时间点的估计值
Figure BDA0003138093620000112
假如没有wk-1与vk，则从式(1.1)和式(1.2)可以立即求得xk，也就代表着不存在估计问题。估计问题的出现正是因为信号和噪声叠加在一起，所以要估计其中信号的真值。暂时不考虑wk-1和vx，此时按照式(1.1)和(1.2)得到的xk和yk分别用
Figure BDA0003138093620000113
和zk表示，则有：
Figure BDA0003138093620000114
Figure BDA0003138093620000115
将zk与yk的实际观测值作比较，它们的差用
Figure BDA0003138093620000116
表示，有:
Figure BDA0003138093620000117
zk偏离yk而产生
Figure BDA0003138093620000118
明显是由于忽略wk-1和vk引起的，也就是说
Figure BDA0003138093620000119
隐含了wk-1和vk的信息，换句话说
Figure BDA00031380936200001110
隐含了当前观察值yk的信息。
如果将
Figure BDA00031380936200001111
乘以增益Hk来修正原来的
Figure BDA00031380936200001112
值，会得到更好的估计：
Figure BDA00031380936200001113
因此可以说
Figure BDA00031380936200001114
和真值xk的均方误差是一个误差方阵。如果能求得这个误差在最小条件下的Hk，然后将此Hk代入式(1.6)，则所得的
Figure BDA00031380936200001115
就是对xk的线性最优估计。
则卡尔曼滤波递推公式如下：
状态预测方程
Figure BDA0003138093620000121
状态测量方程
Figure BDA0003138093620000122
均方差预测方程
Figure BDA0003138093620000123
卡尔曼增益方程
Figure BDA0003138093620000124
状态修改方程
Figure BDA0003138093620000125
均方差修改方程
Figure BDA0003138093620000126
此处k代表时间点，矩阵
Figure BDA0003138093620000127
表示k时间点不考虑估计误差的状态，矩阵
Figure BDA0003138093620000128
表示k时间点的状态，矩阵zk表示k时间点测量得到的数据，矩阵Ak表示状态由k-1时间点到k时间点的状态转移矩阵，矩阵
Figure BDA0003138093620000129
表示k时间点的均方差矩阵，矩阵Hk表示卡尔曼增益矩阵。
(52)定义卡尔曼滤波跟踪器系统状态
Figure BDA00031380936200001210
是一个四维向量(Px,Vx,Py,Vy)T。Px，Vx，Py，Vy分别是目标在像平面坐标系下X轴和Y轴方向上的位置与速度。由于在图像上只能确定目标的二维位置，所以在此定义观测向量zk＝(Px,Py)T。
由于把感兴趣目标在单位时间间隔内的运动视作匀速运动，所以状态转移矩阵Ak定义为：
Figure BDA00031380936200001211
其中ΔT为卫星相机拍摄时间间隔。
由系统状态和观测状态的关系可知，观测矩阵Ck为：
Figure BDA00031380936200001212
此外，假设递推公式中：
Figure BDA0003138093620000131
其中delta＝0.01，r＝100。
卡尔曼滤波根据某一个虚拟像平面内所有像平面坐标进行多次计算后，即可依据最终Xk的得到一个预测的像平面坐标点。
步骤六：根据预测二维轨迹的观测视线修正轨道。具体包括以下步骤：
(61)假定根据时序的卡尔曼滤波计算出的预测点为(x0,y0,z0)，观测到的实际点为(x′0,y′0,z′0)。则可把观测矢量设为参数方程：
Figure BDA0003138093620000132
其中i1、j1、k1、i2、j2、k2为卫星-导弹向量参数，可由卫星与导弹的坐标求出。
(62)根据前方交会原理，预测的时序观测向量与真实的时序观测向量的交点即为导弹位置。为了准确计算两空间直线交点，设两条直线的公垂线与两直线的交点为A(i1m+x0,j1m+y0,k1m+z0)，B(i2n+x′0,j2n+y′0,k2n+z′0)。其中m、n为待求参数。根据公垂线性质，即可得到：
Figure BDA0003138093620000133
(63)将i1、j1、k1、i2、j2、k2、预测点(x0,y0,z0)和实际点(x′0,y′0,z′0)代入式(11)，即可求出参数m、n的值。将m、n代入式(10)即可求出A、B两点三维信息。
(64)计算交点A、B的中点C，C即为(A+B)/2。所有虚拟像平面计算所得C的所有集合为Ci。根据三点确定一个平面的方法，可以从Ci集合中任选三个坐标定义一个平面D，所有不同平面D的集合为Di。
平面D定义为：
CS＝OX+PY+LZ (12)
其中O、P、L、CS为平面参数。
如果步骤三备选弹道平面缺少平面集合Di的平面，则把缺少的平面加入备选弹道平面。
本发明研究单星视轴下的弹道平面切割预估模型，解析弹道主动段运动特性，利用卡尔曼滤波构造约束平面，确定弹道空间轨迹，可以解决现有技术中心难以确定备选弹道平面的问题。
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下还可以作出若干改进，这些改进也应视为本发明的保护范围。

Claims (7)

1.一种基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，包括以下步骤：
(1)根据弹道导弹目标点坐标到预警卫星红外传感器预警图像像点坐标的转换关系，提取在地心固定坐标系下卫星中心到目标弹道导弹的方向单位观测视线矢量；
(2)对各个弹道平面间存在的几何关系进行分析，并用这些弹道平面对上述提取的观测视线矢量进行切割，建立弹道面切割模型，包括：
(21)引入预警卫星首发点和终点在地心固定坐标系下的坐标，分别为(X1 Y1 Z1)T和(XNYN ZN)T，所述首发点为预警卫星首次探测发现弹道导弹的位置，所述终点为最后探测发现弹道导弹的位置，则有：
Figure FDA0003605614920000011
其中，(u1 v1 w1)T为首发点对应的观测视线矢量，(uN vN wN)T为终点对应的观测视线矢量，(XS YS ZS)T为预警卫星在地心固定坐标系下的坐标，re为平均地球半径，h1和hN分别为首发点和终点高程，假定观测卫星的星载红外传感器监视时，探测到主动段弹道导弹发动机的尾焰，共有N组观测值，分别对应于观测时刻t1,t2,……,tN；
(22)根据弹道基本约束中的通视约束解算出(X1 Y1 Z1)T和(XN YN ZN)T，所述通视约束指的是弹道平面切割时观测视线矢量不能先于弹道导弹穿过地球，进而确定一个过地心的弹道切割平面；
(3)根据弹道导弹在主动飞行时的受力情况、运动特性以及其共性特征参数构建基于动量守恒的主动段动力学模型，并对备选弹道平面进行筛选；
(4)根据时序转换像平面坐标系，得到多个虚拟像平面坐标系；
(5)利用卡尔曼滤波对像平面坐标系的目标信息进行时序二维轨迹预测；
(6)根据预测二维轨迹的时序观测视线修正轨道。
2.根据权利要求1所述的基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，所述步骤(1)包括以下步骤：
(11)基于共线条件方程，对预警卫星红外传感器预警图像像点坐标原始图像建立严密几何成像模型：
Figure FDA0003605614920000021
其中，(x y)为成像模型得到的像平面坐标，(x0 y0)为预警卫星红外传感器预警图像像主点坐标，f为相机焦距，(Xs Ys Zs)为预警卫星在地心固定坐标系下的坐标，(X Y Z)为弹道导弹在地心固定坐标系下的坐标，ak,bk,ck(k＝1,2,3)为预警卫星系统参数的函数，具体表示如下：
Figure FDA0003605614920000022
式中，RJ2000-WGS84为J2000坐标系到WGS84坐标系的旋转矩阵；Rorbit-J2000为预警卫星的轨道测量坐标系到J2000坐标系的旋转矩阵；Rbody-orbit为预警卫星本体坐标系到轨道测量坐标系的转换矩阵；Rcamera-body为传感器坐标系到本体坐标系的变换矩阵；
(12)建立观测卫星中心到目标弹道导弹的方向单位观测视线矢量(u v w)T与导弹的像平面坐标(x,y)之间的数学函数关系：
Figure FDA0003605614920000023
其中，m、λ为比例系数；(x,y)T为弹道导弹成像的像平面坐标；(x0,y0,-f)T为像点对应成像时的内方位元素；根据该关系式得到观测视线矢量(u v w)T，u,v,w分别为X方向、Y方向、Z方向上的观测视线矢量。
3.根据权利要求1所述的基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，所述步骤(22)包括：
(a)确定首发点高程和终点高程的区间范围，其中首发点高程区间[h1min，h1max]根据首发点预警图像的数据以及当地气象数据得到；终点高程区间的下限边界选取与首发高程区间的上限边界相同，即hNmin＝h1max，其上界由hNmax＝-g0·T2/2+T·umax得到，式中，T＝tN-t1为预警卫星从首发点到终点的总观测时长，g0为重力加速度，umax为最大燃气喷速值；
(b)设置h1和hN的遍历步长，遍历所有可能的候选弹道切割平面，并求解出每一个切割平面内的备选弹道轨迹数据。
4.根据权利要求1所述的基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，所述步骤(3)包括以下步骤：
(31)忽略空气作用以及柯氏惯性力和牵连惯性力的影响，根据弹道导弹在主动段动量的改变为所受合外力的冲量，构建弹道导弹的加速度模型：
Figure FDA0003605614920000031
其中，p″(t)为导弹在t时刻的加速度，M(t)为t时刻的导弹质量，M′(t)为t时刻导弹质量的变化率，g0为重力加速度，u为导弹高温燃气的相对有效喷射速度；
(32)在步骤(31)的基础上，引入火箭发动机秒耗量a，记s＝a/M，建立弹道导弹主动段动力学模型：
Figure FDA0003605614920000032
其中，p(t)表示导弹在t时刻的位置，p0,v0分别表示位置与速度的初始条件；假定高温燃气相对于导弹以恒定的有效喷射速度为uc，gc为重力加速度常量，下标c表示假定的恒定值；
(33)利用弹道导弹的主动段动力学模型对备选弹道平面进行筛选，将不符合该动力学模型的弹道进行排除，得到符合导弹特征的轨迹对象。
5.根据权利要求1所述的基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，所述步骤(4)中利用如下转换矩阵进行坐标系转换：
Figure FDA0003605614920000041
其中各参数分别为：
Figure FDA0003605614920000042
Figure FDA0003605614920000043
Figure FDA0003605614920000044
b1＝cosωsinκ
b2＝cosωcosκ
b3＝-sinω
Figure FDA0003605614920000045
Figure FDA0003605614920000046
Figure FDA0003605614920000047
其中
Figure FDA0003605614920000048
ω、κ分别为卫星的三个姿态角：俯仰角、滚动角和偏航角。
6.根据权利要求1所述的基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，所述步骤(5)包括以下步骤：
(51)在步骤(4)所得的虚拟像平面上进行卡尔曼预测，卡尔曼滤波递推公式如下：
状态预测方程
Figure FDA0003605614920000049
状态测量方程
Figure FDA00036056149200000410
均方差预测方程
Figure FDA00036056149200000411
卡尔曼增益方程
Figure FDA00036056149200000412
状态修改方程
Figure FDA00036056149200000413
均方差修改方程
Figure FDA00036056149200000414
其中k代表时间点，矩阵
Figure FDA00036056149200000415
表示k时间点不考虑估计误差的状态，矩阵
Figure FDA00036056149200000416
表示k时间点的状态，矩阵zk表示k时间点测量得到的数据，矩阵Ak表示状态由k-1时间点到k时间点的状态转移矩阵，矩阵
Figure FDA00036056149200000417
表示k时间点的均方差矩阵，矩阵Hk表示卡尔曼增益矩阵；
(52)定义卡尔曼滤波跟踪器系统状态
Figure FDA0003605614920000051
为一个四维向量(Px,Vx,Py,Vy)T，Px，Vx，Py，Vy分别是目标在像平面坐标系下X轴和Y轴方向上的位置与速度；定义观测向量zk＝(Px,Py)T；
定义状态转移矩阵Ak：
Figure FDA0003605614920000052
其中ΔT为卫星相机拍摄时间间隔；
由系统状态和观测状态的关系可知，观测矩阵Ck为：
Figure FDA0003605614920000053
令递推公式中：
Figure FDA0003605614920000054
其中delta＝0.01，r＝100；
卡尔曼滤波根据某一个虚拟像平面内所有像平面坐标进行多次计算后，依据最终Xk得到一个预测的像平面坐标点。
7.根据权利要求1所述的基于时序交会的远程导弹主动段运动参数修正方法，其特征在于，所述步骤(6)包括如下步骤：
(61)记卡尔曼滤波计算出的预测点为(x0,y0,z0)，观测到的实际点为(x′0,y′0,z′0)，将观测矢量设为参数方程：
Figure FDA0003605614920000055
其中i1、j1、k1、i2、j2、k2为卫星-导弹向量参数，由卫星与导弹的坐标求出；
(62)根据前方交会原理，预测的时序观测向量与真实的时序观测向量的交点即为导弹位置，为了准确计算两空间直线交点，设两条直线的公垂线与两直线的交点为A(i1m+x0,j1m+y0,k1m+z0)，B(i2n+x′0,j2n+y′0,k2n+z′0)，其中m、n为待求参数，根据公垂线性质，得到：
Figure FDA0003605614920000061
(63)将i1、j1、k1、i2、j2、k2、预测点(x0,y0,z0)和实际点(x′0,y′0,z′0)代入步骤(62)中的计算公式，求出参数m、n的值；将m、n代入步骤(61)中的计算公式，式求出A、B两点三维信息；
(64)计算交点A、B的中点C，所有时间点的C为一个坐标集Ci，约束后轨道与Ci进行比对后修改参数，实现修正轨迹。
CN202110727620.5A 2021-06-29 2021-06-29 基于时序交会的远程导弹主动段运动参数修正方法 Active CN113962057B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110727620.5A CN113962057B (zh) 2021-06-29 2021-06-29 基于时序交会的远程导弹主动段运动参数修正方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110727620.5A CN113962057B (zh) 2021-06-29 2021-06-29 基于时序交会的远程导弹主动段运动参数修正方法

Publications (2)

Publication Number Publication Date
CN113962057A CN113962057A (zh) 2022-01-21
CN113962057B true CN113962057B (zh) 2022-06-24

Family

ID=79460314

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110727620.5A Active CN113962057B (zh) 2021-06-29 2021-06-29 基于时序交会的远程导弹主动段运动参数修正方法

Country Status (1)

Country Link
CN (1) CN113962057B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117052542B (zh) * 2023-10-13 2023-12-08 太仓点石航空动力有限公司 一种航空发动机的推进控制优化方法及系统

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105022035A (zh) * 2015-07-31 2015-11-04 中国电子科技集团公司第三十八研究所 一种基于模型修正的弹道目标发射点估计装置及其方法
CN106643297A (zh) * 2016-12-23 2017-05-10 南京长峰航天电子科技有限公司 一种运动平台矢量脱靶量参数估计修正方法
CN107607947A (zh) * 2017-08-22 2018-01-19 西安电子科技大学 基于卡尔曼滤波的星载雷达成像参数在线估计方法
CN109933847A (zh) * 2019-01-30 2019-06-25 中国人民解放军战略支援部队信息工程大学 一种改进的主动段弹道估计算法
CN111462182A (zh) * 2020-03-31 2020-07-28 南京航空航天大学 一种基于红外预警图像的弹道导弹三维轨迹估计方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DK1316774T3 (da) * 2001-11-28 2006-10-09 Rheinmetall Waffe Munition Projektiler med höj penetrations- og lateralvirkning med integreret sönderdelingsindretning
US8639476B2 (en) * 2008-04-22 2014-01-28 The United States Of America As Represented By The Secretary Of The Navy Process for estimation of ballistic missile boost state

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105022035A (zh) * 2015-07-31 2015-11-04 中国电子科技集团公司第三十八研究所 一种基于模型修正的弹道目标发射点估计装置及其方法
CN106643297A (zh) * 2016-12-23 2017-05-10 南京长峰航天电子科技有限公司 一种运动平台矢量脱靶量参数估计修正方法
CN107607947A (zh) * 2017-08-22 2018-01-19 西安电子科技大学 基于卡尔曼滤波的星载雷达成像参数在线估计方法
CN109933847A (zh) * 2019-01-30 2019-06-25 中国人民解放军战略支援部队信息工程大学 一种改进的主动段弹道估计算法
CN111462182A (zh) * 2020-03-31 2020-07-28 南京航空航天大学 一种基于红外预警图像的弹道导弹三维轨迹估计方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
"Trajectory and Launch Point Estimation for Ballistic Missiles from Boost Phase LOS Measurements";Yicong Li etc.;《1999 IEEE Aerospace Conference. Proceedings (Cat. No.99TH8403)》;20020806;全文 *
"单星无源探测弹道导弹射向估计新方法";申镇 等;《宇航学报》;20110731;全文 *

Also Published As

Publication number Publication date
CN113962057A (zh) 2022-01-21

Similar Documents

Publication Publication Date Title
CN110081881B (zh) 一种基于无人机多传感器信息融合技术的着舰引导方法
US9073648B2 (en) Star tracker rate estimation with kalman filter enhancement
CN107741229B (zh) 一种光电/雷达/惯性组合的舰载机着舰导引方法
CN104316060B (zh) 空间非合作目标的交会对接方法与装置
Mammarella et al. Machine vision/GPS integration using EKF for the UAV aerial refueling problem
CN109059907B (zh) 轨迹数据处理方法、装置、计算机设备和存储介质
US20080195316A1 (en) System and method for motion estimation using vision sensors
Zhao et al. Vision-aided estimation of attitude, velocity, and inertial measurement bias for UAV stabilization
CN114216454B (zh) 一种gps拒止环境下基于异源图像匹配的无人机自主导航定位方法
US7957899B2 (en) Method for determining the attitude, position, and velocity of a mobile device
Johnson et al. Combining stereo vision and inertial navigation for automated aerial refueling
CN111912430B (zh) 高轨光学卫星的在轨几何定标方法、装置、设备及介质
US9816786B2 (en) Method for automatically generating a three-dimensional reference model as terrain information for an imaging device
CN113962057B (zh) 基于时序交会的远程导弹主动段运动参数修正方法
CN108225276B (zh) 一种单星成像目标运动特性反演方法及系统
Helgesen et al. Tracking of ocean surface objects from unmanned aerial vehicles with a pan/tilt unit using a thermal camera
CN113963035A (zh) 一种单星条件下预警图像弹道导弹射向估计方法
CN116642482A (zh) 基于固态激光雷达和惯性导航的定位方法、设备和介质
Campbell et al. Vision-based geolocation tracking system for uninhabited aerial vehicles
Deschênes et al. Lidar scan registration robust to extreme motions
Moore et al. Vision-only estimation of wind field strength and direction from an aerial platform
WO2012052738A1 (en) Sensor positioning for target tracking
KR102050995B1 (ko) 공간좌표의 신뢰성 평가 장치 및 방법
Veth et al. Tightly-coupled ins, gps, and imaging sensors for precision geolocation
Bashir et al. Kalman Filter Based Sensor Fusion for Altitude Estimation of Aerial Vehicle

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