本发明公开一种基于红外预警图像的弹道导弹三维轨迹估计方法，对预先获取的连续时间序列的导弹红外遥感图像进行点目标检测，并提取各幅图像点目标质心像素坐标；依据上述提取的目标点的时间及空间信息，跟踪目标并得到运动轨迹；对弹道导弹运动特点分析进行弹道导弹轨迹筛选；将导弹目标的姿态数据与轨道数据统一在一个坐标系下；依据上述的像平面二维轨迹，通过共线方程以及z方向上的约束条件预测目标三维轨迹；根据导弹的运动学模型，及导弹主动段的约束条件，修正误差轨迹。本发明由二维估计三维，为导弹的跟踪预警提供了新的方案。
1.一种基于红外预警图像的弹道导弹三维轨迹估计方法，其特征在于，包括以下步骤：(1)对预先获取的连续时间序列的导弹红外遥感图像进行点目标检测，并提取各幅图像点目标质心像素坐标；
(2)依据上述提取的目标质心的时间及空间信息，跟踪目标并得到运动轨迹；
(3)对弹道导弹运动特点分析进行弹道导弹轨迹筛选，得到像平面二维轨迹；
(4)将导弹目标的姿态数据与轨道数据统一在一个坐标系下；
(5)依据上述的像平面二维轨迹，通过共线方程以及z方向上的约束条件预测目标三维轨迹；
(6)根据导弹的运动学模型，及导弹主动段的约束条件，修正误差轨迹；
所述步骤(5)包括如下步骤：(51)将J2000坐标系转换为WGS84坐标系，利用导弹的质心运动方程和共线条件方程来对z方向进行约束，质心运动方程：其中，V是导弹速度，θ是弹道倾角，σ是弹道偏角；
(52)选取两个点 两个点通过共线方程求出各自X、Y的坐标，这些坐标都是关于h0的函数：(53)通过上述方程Vx、Vy求出弹道倾角和导弹运动速度：则Vz：
(54)通过二维像平面轨迹求出弹道偏角，并推导出z方向上的约束条件，将二维轨迹上的点代入公式中进行迭代，则可求出各个点的三维坐标。
2.根据权利要求1所述的一种基于红外预警图像的弹道导弹三维轨迹估计方法，其特征在于，所述步骤(1)包括以下步骤：(11)对原始图像进行Top‑Hat滤波抑制背景，定义为：TH(f)＝‑open(f,e)，其中，f为输入图像，e为结构元素；
(12)选取与目标特征相同的结构元素对抑制后图像进行卷积处理，并设置阈值分割出可能的目标点，则目标分割后的图像：其中，f(i,j,k)为原始图像，fH(i,j,k)为背景抑制后图像，E(i,j,k)为结构元素卷积后图像，α选取与虚警率密切相关；
(13)使用sobel算子计算图像梯度，聚类分析得目标点个数及质心坐标。
3.根据权利要求1所述的一种基于红外预警图像的弹道导弹三维轨迹估计方法，其特征在于，所述步骤(2)包括如下步骤：(21)在步骤(1)提取后图象中选取n帧图像进行帧积分，并使用meanshift算法收敛到密度较大的区域即目标点运动的区域；
(22)对上述区域内目标进行管道滤波滤除管道外噪声干扰，根据目标运动的连续性，如果第k帧中在像素(x,y)处有目标，则该目标在第k+1帧中必然会出现在像素(x,y)的一个小邻域内；
(23)管道滤波后图像进入卡尔曼滤波器跟踪，得到目标运动轨迹：卡尔曼滤波器时间更新方程：
4.根据权利要求1所述的一种基于红外预警图像的弹道导弹三维轨迹估计方法，其特征在于，所述步骤(3)包括以下步骤：(31)引入导弹质量M，火箭发动机秒耗量a，记k＝a/M为相对秒耗量；并引入高温燃气相对导弹的恒定有效喷射速度uc，建立弹道导弹主动段飞行动量模型：其中，p(t)表示导弹在t时刻的位置，即弹道模型，p0,v0分别表示位置与速度的初始条件；g为重力加速度，下标c表示假定的恒定值；
(32)在步骤(31)的基础上，使用观测的位置数据，即导弹观测的角度信息，建立弹道平面切割模型；单星探测得到的弹道视线观测是二维量，可用方位角a和俯仰角e表示，其中x，y，z是目标在卫星测量坐标系下的坐标值：a＝arctan(y/x)
(33)构建单星观测视线对目标导弹参数的估计算法，通过迭代计算得到一组切割弹道，即至此完成对输入的二维轨迹的参数进行估计，对每条轨迹可计算出一组切割弹道；
(34)利用模型合理性约束体系对结果进行筛选，通过约束筛，排除不符合导弹特征的轨迹对象，即完成筛选出符合导弹运动特征的目标的目的。
5.根据权利要求1所述的一种基于红外预警图像的弹道导弹三维轨迹估计方法，其特征在于，所述步骤(4)包括以下步骤：(41)对J2000坐标系下的三个姿态角进行从J2000坐标系到WGS84坐标系的转换，采用基于CIO的无旋转原点转换，坐标转换公式可写为：其中， 为常值偏差矩阵， 为岁差章动矩阵，为地球自转矩阵， 为极移矩阵；
(42)岁差、章动模型选取iau1976/1980岁差章动模型；
iau1976岁差模型可表示为：其中，T为离参考时刻J2000.0的儒略世纪数， 以世纪为单位；
iau1980章动模型可表示为：式中，幅角 Ai为黄经章动中的振幅项，A′i为黄经章动中振幅项的变化率，Bi为交角章动中的振幅项，B′i为交角章动中振幅项的变化率。
6.根据权利要求1所述的一种基于红外预警图像的弹道导弹三维轨迹估计方法，其特征在于，所述步骤(6)包括如下步骤：(61)建立三维弹道主动段标准模型，根据导弹主动段运动规律，建立主动段弹道三维运动学模型，得到以时间为自变量，导弹一系列的三维坐标：Sx,n+1＝Sx,n+Vx,n+1·dtSy,n+1＝Sy,n+Vy,n+1·dtSz,n+1＝Sz,n+Vz,n+1·dt；
(62)根据约束条件对误差轨迹进行修正，导弹在主动段的运动服从一定规律，因此根据这些约束条件找到轨迹中的误差点，并对其进行修正；约束条件包括：弹道曲线应在地球外侧，即弹道曲线的高度应大于零：弹道曲线的高度变化应满足递增特性：弹道曲线的最大线速度不能超过第一宇宙速度：弹道曲线的各点速度应满足递增特性：弹道曲线的最大加速度不能超过15g：(63)判断修正后轨迹与标准模型拟合效果，当拟合效果符合标准时，输出修正后轨迹，得到最终轨迹结果。
一种基于红外预警图像的弹道导弹三维轨迹估计方法
技术领域
[0001] 本发明属于红外图像处理技术、摄影测量技术及三维轨迹建模技术领域，具体涉及一种基于红外预警图像的弹道导弹三维轨迹估计方法。
背景技术
[0002] 红外预警系统以卫星搭载的红外相机为载体获得红外预警图像，探测的目标为低空飞机和掠海导弹等。远距离时，这些目标在焦平面上成像的面积很小，目标在图像中呈现
为点状没有任何纹理信息可以利用，且受到大气辐射、云层背景等影响信噪比极低，目标很
容易被噪声所淹没。基于以上图像特点，近年来有学者提出基于形态学的滤波算法，通过选
取合适的结构元素有效的提高图像信噪比，检测目标。
[0003] 目标跟踪算法主要是根据目标建模或者对目标特征进行跟踪，主要方法有：其一，基于目标模型建模的方法；通过对目标外观模型进行建模，然后在之后的帧中找到目标。其
二，基于搜索的方法；随着研究的深入，人们将预测算法加入跟踪中，在预值附近进行目标
搜索，减少搜索范围。常见一类预测算法包括Kalman滤波、粒子滤波。另一种减小搜索范围
的方法是内核方法,运用最速下降法的原理向梯度下降方向对目标模板逐步迭代，直到最
优位置，诸如meanshift算法。
[0004] 弹道导弹的三维建模有两种方式，一种是利用先验知识，基于弹道模版库的弹道建模，称为Profile‑dependent Model方法；由于弹道导弹类型有限，因此可以事先将不同
类型目标主动段弹道模版存储在数据库中，称为标称弹道数据库(Nominal Ballistic
Profile)。另一种是不考虑先验弹道信息，称为Profile‑free Model方法。该方法不需要先
验知识，只通过单星观测结果进行弹道主动段参数估计。
发明内容
[0005] 发明目的：本发明提供一种基于红外预警图像的弹道导弹三维轨迹估计方法，能够通过加入z坐标估计约束的方法由二维图像解算三维轨迹，最后建立弹道导弹三维动力
学模型进行轨迹修正。
[0006] 技术方案：本发明所述的一种基于红外预警图像的弹道导弹三维轨迹估计方法，包括以下步骤：
[0007] (1)对预先获取的连续时间序列的导弹红外遥感图像进行点目标检测，并提取各幅图像点目标质心像素坐标；
[0008] (2)依据上述提取的目标质心的时间及空间信息，跟踪目标并得到运动轨迹；
[0009] (3)对弹道导弹运动特点分析进行弹道导弹轨迹筛选；
[0010] (4)将导弹目标的姿态数据与轨道数据统一在一个坐标系下；
[0011] (5)依据上述的像平面二维轨迹，通过共线方程以及z方向上的约束条件预测目标三维轨迹；
[0012] (6)根据导弹的运动学模型，及导弹主动段的约束条件，修正误差轨迹。
[0013] 进一步地，所述步骤(1)包括以下步骤：
[0014] (11)对原始图像进行Top‑Hat滤波抑制背景，定义为：TH(f)＝‑open(f,e)，其中，f为输入图像，e为结构元素；
[0015] (12)选取与目标特征相同的结构元素对抑制后图像进行卷积处理，并设置阈值分割出可能的目标点，则目标分割后的图像：
[0016]
[0017] 其中，f(i,j,k)为原始图像，fH(i,j,k)为背景抑制后图像，E(i,j,k)为结构元素卷积后图像，α选取与虚警率密切相关；
[0018] (13)使用sobel算子计算图像梯度，聚类分析得目标点个数及质心坐标。
[0019] 进一步地，所述步骤(2)包括如下步骤：
[0020] (21)在步骤(1)提取后图象中选取n帧图像进行帧积分，并使用meanshift算法收敛到密度较大的区域即目标点运动的区域；
[0021] (22)对上述区域内目标进行管道滤波滤除管道外噪声干扰，根据目标运动的连续性，如果第k帧中在像素(x,y)处有目标，则该目标在第k+1帧中必然会出现在像素(x,y)的
一个小邻域内；
[0022] (23)管道滤波后图像进入卡尔曼滤波器跟踪，得到目标运动轨迹：
[0023] 卡尔曼滤波器时间更新方程：
[0024] 卡尔曼滤波器状态更新方程：
[0025] 进一步地，所述步骤(3)包括以下步骤：
[0026] (31)引入导弹质量M，火箭发动机秒耗量a，记k＝a/M；为相对秒耗量；并引入高温燃气相对导弹的恒定有效喷射速度uc，建立弹道导弹主动段飞行动量模型：
[0027]
[0028] 其中，p(t)表示导弹在t时刻的位置，即弹道模型，p0,v0分别表示位置与速度的初始条件；g为重力加速度，下标c表示假定的恒定值；
[0029] (32)在步骤(31)的基础上，使用观测的位置数据，即导弹观测的角度信息，建立弹道平面切割模型；单星探测得到的弹道视线观测是二维量，可用方位角a和俯仰角e表示，其
中x，y，z是目标在卫星测量坐标系下的坐标值：
[0030] a＝arctan(y/x)
[0031]
[0032] (33)构建单星观测视线对目标导弹参数的估计算法，通过迭代计算得到一组切割弹道，即至此完成对输入的二维轨迹的参数进行估计，对每条轨迹可计算出一组切割弹道；
[0033] (34)利用模型合理性约束体系对结果进行筛选，通过约束筛，排除不符合导弹特征的轨迹对象，即完成筛选出符合导弹运动特征的目标的目的。
[0034] 进一步地，所述步骤(4)包括以下步骤：
[0035] (41)对J2000坐标系下的三个姿态角进行从J2000坐标系到WGS84坐标系的转换，采用基于CIO的无旋转原点转换，坐标转换公式可写为：
[0036]
[0037] 其中， 为常值偏差矩阵， 为岁差章动矩阵， 为地球自转矩阵，为极移矩阵；
[0038] (42)岁差、章动模型选取iau1976/1980岁差章动模型；
[0039] iau1976岁差模型可表示为：
[0040] 其中，T为离参考时刻J2000.0的儒略世纪数， 以世纪为单位；
[0041] iau1980章动模型可表示为：
[0042]
[0043]
[0044] 式中，幅角 Ai为黄经章动中的振幅项，A′i为黄经章动中振幅项的变化率，Bi为交角章动中的振幅项，B′i为交角章动中振幅项的变
化率。
[0045] 进一步地，所述步骤(5)包括如下步骤
[0046] (51)将J2000坐标系转换为WGS84坐标系，利用导弹的质心运动方程和共线条件方程来对z方向进行约束，质心运动方程：
[0047]
[0048]
[0049]
[0050] 其中，V是导弹速度，θ是弹道倾角，σ是弹道偏角；
[0051] (52)选取两个点 两个点通过共线方程可以求出各自X、Y的坐标，这些坐标都是关于h0的函数：
[0052]
[0053]
[0054]
[0055]
[0056] (53)通过上述方程Vx、Vy可以求出弹道倾角和导弹运动速度：
[0057]
[0058]
[0059] 则Vz就可以求出：
[0060]
[0061] (54)通过二维像平面轨迹求出弹道偏角，并推导出z方向上的约束条件，将二维轨迹上的点代入公式中进行迭代，则可求出各个点的三维坐标。
[0062] 进一步地，所述步骤(6)包括如下步骤：
[0063] (61)建立三维弹道主动段标准模型，根据导弹主动段运动规律，建立主动段弹道三维运动学模型，得到以时间为自变量，导弹一系列的三维坐标：
[0064]
[0065]
[0066]
[0067] Sx,n+1＝Sx,n+Vx,n+1·dt
[0068] Sy,n+1＝Sy,n+Vy,n+1·dt
[0069] Sz,n+1＝Sz,n+Vz,n+1·dt；
[0070] (62)根据约束条件对误差轨迹进行修正，导弹在主动段的运动服从一定规律，因此根据这些约束条件可以找到轨迹中的误差点，并对其进行修正；约束条件主要包括：弹道
曲线应在地球外侧，即弹道曲线的高度应大于零：
[0071]
[0072] 弹道曲线的高度变化应满足递增特性：
[0073]
[0074] 弹道曲线的最大线速度不能超过第一宇宙速度：
[0075]
[0076]
[0077] 弹道曲线的各点速度应满足递增特性：
[0078]
[0079]
[0080] 弹道曲线的最大加速度不能超过15g：
[0081]
[0082]
[0083] (63)判断修正后轨迹与标准模型拟合效果，当拟合效果符合标准时，输出修正后轨迹，得到最终轨迹结果。
[0084] 有益效果：与现有技术相比，本发明的有益效果：1、本发明与传统的基于天基预警系统的导弹弹道估计方法，主要利用星载红外传感器测得的角测量来估计主动段运动状态
不同；利用导弹质心运动方程与共线方程来对缺失的z方向进行约束从而解算三维坐标，有
效避免了由于状态方程及观测方程的非线性所导致的误差不收敛；2、仿真弹道导弹三维轨
迹模型进行轨迹修正，能有效的缩小误差提高轨迹精度。
附图说明
[0085] 图1为本发明的流程图。
具体实施方式
[0086] 下面结合附图对本发明作进一步详细描述。本发明提供一种基于红外预警图像的弹道导弹三维轨迹估计方法，如图1所示，包括如下步骤：
[0087] 步骤一：对仿真所得的一组连续时间序列的导弹红外遥感图像进行点目标检测，并提取各幅图像点目标质心像素坐标。具体包括以下步骤：
[0088] (1)对原始图像进行Top‑Hat滤波抑制背景，定义为：TH(f)＝‑open(f,e)。其中，f为输入图像，e为结构元素。
[0089] (2)选取与目标特征相同的结构元素对抑制后图像进行卷积处理，并设置阈值分割出可能的目标点。则目标分割后的图像：
[0090]
[0091] 其中，f(i,j,k)为原始图像，fH(i,j,k)为背景抑制后图像，E(i,j,k)为结构元素卷积后图像，α选取与虚警率密切相关。
[0092] (3)使用sobel算子计算图像梯度，聚类分析得目标点个数及质心坐标。
[0093] 步骤二：依据上述提取的目标点的时间及空间信息，跟踪目标并得到运动轨迹。具体包括以下步骤：
[0094] (1)在所述提取后图形中选取n帧图像进行帧积分，并使用meanshift算法收敛到密度较大的区域即目标点运动的区域。
[0095] (2)对上述区域内目标进行管道滤波滤除管道外噪声干扰。根据目标运动的连续性，如果第k帧中在像素(x,y)处有目标，则该目标在第k+1帧中必然会出现在像素(x,y)的
一个小邻域内。
[0096] (a)对于先进入“与”管道的第k帧二值图像进行数学形态学的膨胀运算，得到一个掩膜帧，使原二值图像为“1”的像素的小邻域都为“1”；
[0097] (b)将第k+1帧二值图像和掩膜帧进行逻辑“与”运算，得到连续性滤波后的第k+1帧。
[0098] 如此，把连续性滤波管道同检测管道连接起来就构成一个完整的管道滤波器结构。
[0099] (3)管道滤波后图像进入卡尔曼滤波器跟踪，得到目标运动轨迹。
[0100] 卡尔曼滤波器时间更新方程：
[0101] 卡尔曼滤波器状态更新方程：
[0102] 步骤三：对弹道导弹运动特点分析进行弹道导弹轨迹筛选。
[0103] 单星观测使用红外传感器，无法获得直接获得导弹的位置信息，对于导弹定位来说属于不完备观测，也就是说要课题实现最终目的，需要对所得数据进行参数估计，从而建
立尽可能准确的弹道模型。
[0104] 在步骤1和步骤2的基础上，输入为前一过程跟踪拟合所得的像平面二维轨迹数据，即导弹以及干扰物体在卫星观测二维像平面上，随时间t变化的位置坐标。输出则是符
合导弹运动特征的目标，经过修正后的像平面轨迹。具体包括以下步骤：
[0105] (1)引入导弹质量M，火箭发动机秒耗量a，记k＝a/M；为相对秒耗量；并引入高温燃气相对导弹的恒定有效喷射速度uc，建立弹道导弹主动段飞行动量模型：
[0106]
[0107] 其中，p(t)表示导弹在t时刻的位置，即弹道模型，p0,v0分别表示位置与速度的初始条件；g为重力加速度，下标c表示假定的恒定值。
[0108] (22)使用观测的位置数据，即导弹观测的角度信息，建立弹道平面切割模型。单星探测得到的弹道视线观测是二维量，可用方位角a和俯仰角e表示，其中x，y，z是目标在卫星
测量坐标系下的坐标值。
[0109] a＝arctan(y/x)
[0110]
[0111] (3)上述两种模型可构建单星观测视线对目标导弹参数的估计算法，通过迭代计算得到一组切割弹道，即至此完成对输入的二维轨迹的参数进行估计，对每条轨迹可计算
出一组切割弹道。
[0112] (4)最后利用模型合理性约束体系对结果进行筛选，例如导弹目标的射程要求对应于切割弹道的射程下界约束；导弹关机速度受限于第一宇宙速度等。通过约束筛，排除不
符合导弹特征的轨迹对象，即完成筛选出符合导弹运动特征的目标的目的。
[0113] 步骤四：将导弹目标的姿态数据与轨道数据统一在一个坐标系下。具体包括以下步骤：
[0114] (1)对J2000坐标系下的三个姿态角进行从J2000坐标系到WGS84坐标系的转换，采用基于CIO的无旋转原点转换，坐标转换公式可写为：
[0115]
[0116] 其中， 为常值偏差矩阵， 为岁差章动矩阵， 为地球自转矩阵，为极移矩阵。
[0117] (2)岁差、章动模型选取iau1976/1980岁差章动模型；
[0118] iau1976岁差模型可表示为：
[0119] 其中，T为离参考时刻J2000.0的儒略世纪数， 以世纪为单位。
[0120] iau1980章动模型可表示为：
[0121]
[0122]
[0123] 式中，幅角 Ai为黄经章动中的振幅项，A′i为黄经章动中振幅项的变化率，Bi为交角章动中的振幅项，B′i为交角章动中振幅项的变
化率。
[0124] F1＝l＝月球的平近点角
[0125] F2＝l′＝太阳的平近点角
[0126] F3＝F＝L‑Ω(L为地球的平黄经)
[0127] F4＝D＝日月间的平角距
[0128] F5＝Ω＝月球升交点的平黄经
[0129] 步骤五：依据上述的像平面二维轨迹，通过共线方程以及z方向上的约束条件预测目标三维轨迹。具体包括以下步骤：
[0130] (1)将J2000坐标系转换为WGS84坐标系，弹道导弹三维轨迹估计方法的核心思想是利用导弹的质心运动方程和共线条件方程来对z方向进行约束。
[0131] 质心运动方程：
[0132]
[0133]
[0134]
[0135] 其中，V是导弹速度，θ是弹道倾角，σ是弹道偏角。
[0136] (2)选取两个点 两个点通过共线方程可以求出各自X、Y的坐标，这些坐标都是关于h0的函数。
[0137]
[0138]
[0139]
[0140]
[0141] (3)通过上述方程Vx、Vy可以求出弹道倾角和导弹运动速度：
[0142]
[0143]
[0144] 则Vz就可以求出：
[0145]
[0146] (4)弹道偏角可以通过二维像平面轨迹求出，则z方向上的约束条件则可以推导出来。将二维轨迹上的点代入公式中进行迭代，则可求出各个点的三维坐标。
[0147] 步骤六：根据导弹的运动学模型，还有导弹主动段的约束条件，修正误差轨迹。具体包括以下步骤：
[0148] (1)建立三维弹道主动段标准模型，根据导弹主动段运动规律，建立主动段弹道三维运动学模型，得到以时间为自变量，导弹一系列的三维坐标。
[0149]
[0150]
[0151]
[0152] Sx,n+1＝Sx,n+Vx,n+1·dt
[0153] Sy,n+1＝Sy,n+Vy,n+1·dt
[0154] Sz,n+1＝Sz,n+Vz,n+1·dt
[0155] (2)根据约束条件对误差轨迹进行修正。导弹在主动段的运动服从一定规律，因此根据这些约束条件可以找到轨迹中的误差点，并对其进行修正。约束条件主要包括：弹道曲
线应在地球外侧，即弹道曲线的高度应大于零；弹道曲线的高度变化应满足递增特性；弹道
曲线的最大线速度不能超过第一宇宙速度；弹道曲线的各点速度应满足递增特性；弹道曲
线的最大加速度不能超过15g等。
[0156] (3)判断修正后轨迹与标准模型拟合效果，主要从以下方面进行判断：
[0157] 均方根误差RMSE：
[0158]
[0159] 模型效率EF：
[0160]
[0161] 整群剩余系数CRM：
[0162]
[0163] 平均差MD：
[0164]
[0165] 以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下还可以作出若干改进，这些改进也应视为本发明的
保护范围。
