CN108519587B - 一种实时的空中目标运动模式识别及参数估计方法 - Google Patents

Abstract

本发明公开一种实时的空中目标运动模式识别及参数估计方法，通过对上报的空中目标的位置点序列进行插值和滤波平滑处理得到平滑的运动航迹，利用计算得到的航向角度信息结合高度信息，判断空中目标在空间上的运动模式。利用计算得到的速度数据，通过滤波和平滑处理判断速度的变化模式。采用最小二乘拟合的方法，拟合计算得到空中目标运动的角速度和加速度等参数估计值。此外，在圆弧运动模式下，根据指定时间段内的航迹点拟合其所在的平面方程和球面方程，进而可以计算得到圆弧运动的圆心坐标和圆弧半径。本发明可以能够很好地解决原始位置数据存在的数据丢失和噪声干扰的影响。

Description

一种实时的空中目标运动模式识别及参数估计方法

技术领域

本发明属于空中目标运动模式识别领域，特别涉及一种实时的空中目标运动模式识别及参数估计方法。

背景技术

随着国内航空运输业高速发展，对空中运输的安全性有了更高的要求，对空中目标的运动模式识别和参数估计的研究也因此变得十分重要。借助ADS-B技术，监控系统可以自动地从相关机载设备获取目标的位置、高度、速度、航向、识别号等信息，通过研究空中目标的运动轨迹和参数，工作人员可以建立一套完备的飞行数据管理平台，可以实现数据收集、数据分析、数据存储、数据应用及数据查询的综合管理功能，极大地提高飞行及航行教学的质量和管理水平。

对空中目标的运动模式的识别及运动参数的估计以实时获取的ADS-B位置数据为基础，接收到的初始的位置数据掺杂了错误数据以及大量的高频噪声干扰，在数据采集以及数据传输的过程当中还会存在数据丢失的情况。

根据空中目标的实时位置信息判断其运动模式主要问题在于实时地对原始数据以及计算得到的数据结果进行插值、去噪和平滑滤波处理。目前的低通滤波等方法的处理结果会产生滞后，难以满足实时处理的要求，采用低通滤波和零相移滤波结合的方法，可以很好消除低通滤波的滞后，并达到系统的实时性要求。在数据读取的过程中进行插值运算，对实时更新的滤波结果按照一定窗口大小做最小二乘拟合，同样也保证了整个运动模式识别和参数计算过程的实时性。

发明内容

发明目的：针对现有技术中接收到的初始的位置数据掺杂了错误数据以及大量的高频噪声干扰的问题，提供一种法结构简单，易于实现，并且可以达到很好的模式识别和参数估算效果的实时的空中目标运动模式识别及参数估计方法。

技术方案：为解决上述技术问题，本发明提供一种实时的空中目标运动模式识别及参数估计方法，包括如下步骤：

(1)读取原始数据，计算当前数据与对应的航迹的末数据的时间间隔；

(2)判断时间间隔是否超过截断时间间隔阀值，如果不超过则进入步骤(3)，如果超过则对航迹数据进行截断处理并进入步骤(13)；

(3)依据相邻航迹点时间间隔条件进行相应的插值处理，将处理后的航迹点存入待插入航迹点缓存区；

(4)依次将待插入航迹点缓存区的航迹点插入到航迹点容器；

(5)航迹点容器中的航迹点每满m个航迹点，即对航迹点进行零相移的巴特沃斯低通滤波处理，并取最后m个航迹点待处理；

(6)计算航迹点之间的高度差，判断竖直方向的运动模式；

(7)计算航迹点之间的航向角和速度；

(8)对最新计算出的m个航向角和速度进行滤波处理；

(9)计算滤波后航向角的标准差，判断是否是圆弧运动，如果是则进入步骤(10)，如果不是则进入步骤(12)；

(10)拟合航向角关于时间的一次曲线，计算角速度；

(11)根据最新处理过的m个航迹点，计算圆心坐标并进入步骤(12)；

(12)拟合速度关于时间的一次曲线，计算加速度；

(13)数据导入至数据库。

进一步的，所述步骤(2)对于判断时间间隔是否超过截断时间间隔阀值的具体步骤如下：

(2.1)设定时间间隔的合理取值范围为(t_min,t_max]，需截断处理的时间间隔阈值为t_T；

(2.2)依据不同的相邻航迹点时间间隔，对航迹点处理方式归纳为4类：

1)时间间隔Δt满足Δt≤t_min，不进行插值运算，并忽略当前读取的航迹点；

2)时间间隔Δt满足Δt∈(t_min,t_max]，保留当前读取的航迹点，但不需要进行插值处理；

3)时间间隔Δt满足Δt∈(t_max,t_T]，需要对航迹点进行插值处理；

4)时间间隔Δt满足Δt＞t_T，不进行插值处理，并对航迹段进行截断。

进一步的，所述步骤(5)中对航迹点进行巴特沃斯低通滤波处理的具体步骤如下：利用巴特沃斯低通滤波器与零相位滤波算法结合的方法实现对迟滞现象的消除；巴特沃斯滤波器的系数线性差分方程的形式如式(2.1)：

其中，x(n)序列为滤波前的信号序列，a_k,β_j为H(z)系统函数分母与分子的系统数组，求出的y_btw(n)即为滤波后的信号序列；

首先选择合适的截止频率和通频带计算确定巴特沃斯低通滤波器系统函数的分子分母系数，确定零相移巴特沃斯数字低通滤波器的形式，由输入信号得到滤波结果；每当经插值处理后的生成的航迹点满m个时，即对已缓存的航迹数据进行一次零相移低通滤波运算，保留计算结果末尾的m个航迹点，作为初次滤波结果实时缓存。

进一步的，所述步骤(6)中判断竖直方向的运动模式的具体步骤如下：空中目标在竖直方向上的运动模式有爬升运动、下降运动和水平运动三种模式，实时读取经过滤波平滑后的航迹点的高度数值，并计算两两航迹点之间的高度差Δh,设定爬升时高度变化的阈值h_up和下降时高度变化阈值h_down，其中h_up为正值，h_down为负值；当Δh＞h_up时，目标为爬升运动；当Δh<h_down时，目标为下降运动；当Δh介于h_down和h_up时，目标为水平运动。

进一步的，所述步骤(7)中计算航迹点之间的航向角的具体步骤如下：

由两两位置点可以唯一确认空中目标运动的航向角度，在计算角度之前首先需要将位置点的经纬高度坐标转换至以地心为原点的ECEF坐标系下；ECEF坐标系是固定在WGS84参考椭球模型上的随地球旋转的地心地固坐标系，坐标转换过程如下:

其中，lon为纬度，lat为经度，h为距离椭球表面的高度，R为地球的曲率半径；任意两个航迹点p1,p2经过坐标转换得到其在ECEF系下的坐标分别为(x_p1,y_p1,z_p1)和(x_p2,y_p2,z_p2)，那么由航迹点i和j确定的方向矢量的水平投影角

和竖直方向的角度θ计算公式为：

其中Δx＝x_p2-x_p1,Δy＝y_p2-y_p1,Δz＝z_p2-z_p1；

在实时获得的经过滤波平滑后的航迹点数据的基础上计算两两航迹点间的航向角度值，每当有航迹点数据完成滤波平滑处理，即实时读取航迹点位置，并计算其相对前一个航迹点的航向角，计算结果更新至航向角数据的缓存区。

进一步的，所述步骤(9)中判断是否是圆弧运动的具体步骤如下：当滤波后航向角数据的缓存区发生更新时，即从缓存区读取角度数据，每当读取的角度数据个数满m时，计算m个角度值的标准差

设定标准差的变化阈值s_T,比较计算出的标准差与阈值的大小，若

则目标为直线运动模式，若

则目标为圆弧运动模式。

进一步的，所述步骤(12)中计算加速度的具体步骤如下：计算得到平滑的速度数据之后，按m个速度数据为单位，采用最小二乘拟合的方法拟合出速度关于时间的线性函数，目标实时运动的速度v＝v₀+at，那么线性方程的斜率a_v即可近似为目标运动的加速度；计算由所选取的速度数据的首尾端点值确定的速度关于时间的变化率a_v′，以a_v和a_v′的平均值

作为空中目标运动的加速度估计值。

进一步的，所述步骤(10)中计算角速度的具体步骤如下：与加速度的计算过程类似，角度在短时间内的变化规律也是线性变化的过程，θ＝θ₀+ωt；当目标的运动模式判定为圆弧运动时，以m个角度数据为单位，拟合角度关于时间的变化斜率ω，计算首尾端点处角度关于时间的变化率ω′，以ω和ω′的平均值

作为目标运动的角速度估计值。

进一步的，所述步骤(11)中计算圆心坐标的具体步骤如下：球面方程为：

Ax+By+Cz-D＝x²+y²+z²

m个数据点所在的空间平面为：

A′x+B′y+C′z+D′＝0

由m组数据经过最小二乘法的计算，得到球面方程和平面方程的各项系数，并计算球面与平面的交线，从而确定圆弧运动所在空间圆圆心和其半径。

与现有技术相比，本发明的优点在于：

本发明提供的实时的空中目标运动模式识别及参数估计方法能够很好地解决原始位置数据存在的数据丢失和噪声干扰的影响，可以实时地判别空中目标的运动模式并估计其运动参数，本发明提供的方法结构简单，易于实现，并且可以达到很好的模式识别和参数估算效果。

在对目标运动模式识别和参数计算过程中，需要对原始数据进行插值处理，由于原始数据十分粗糙，在模式识别和参数估计的运算过程中采用多次滤波处理的方法，对航向角数据以及速度数据进行滤波平滑，并保持其原本的变化趋势。对原始数据分段采用巴特沃斯低通滤波结合零相移滤波算法，可以很好消除低通滤波的滞后。初次滤波的结果再经过卡尔曼滤波和移动平均平滑，可以达到很好的滤波效果。经滤波处理后的数据具有很好的平滑特性并与原始数据的变化趋势保持一致，在滤波结果的基础上，采用最小二乘拟合的方法可以很好地拟合出速度等数据关于时间的一次函数，从而估算出各项运动参数。

附图说明

图1为本发明的方法流程图；

图2为实时插值过程的示意图；

图3为零相移滤波器结构示意图；

图4为移动平均平滑计算过程的示意图；

图5为空间平面与空间球面交线示意图。

具体实施方式

下面结合附图和具体实施方式，进一步阐明本发明。本发明描述的实施例仅仅是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所得到的其他实施例，都属于本发明所保护的范围。

本发明所述的一种实时的空中目标运动模式识别及参数估计方法，包括如下步骤：

1)在时间间隔较大的相邻位置点间进行插值处理，插值过程采用三次样条插值方法，生成的插值后的位置点序列实时更新缓存，以备进一步的处理。

2)对插值后的位置数据进行零相移的巴特沃斯低通滤波处理，滤波结果保留末尾部分的数据，作为初步的滤波结果，在初步滤波结果的基础上再移动平均平滑处理，得到最终的空中目标位置点数据。计算相邻位置点之间的高度差，设定阈值，判断目标在竖直方向上的运动模式。

3)根据空中目标位置点的经纬高度信息，将坐标转换为ECEF坐标系下的空间直角坐标，并计算两两位置点间的飞行航向角度。对航向角度数据进行滤波处理，根据航向角度的偏转特征判别目标运动的轨迹类型。当判断结果为圆弧运动模式时，利用最小二乘拟合的方法，计算角速度的值以及圆弧运动的圆心坐标及圆弧半径。

根据空中目标两两位置点间的空间距离和时间间隔，计算得到飞行速度数据，对速度数据进行滤波平滑处理，得到各位置点处的速度参数，滤波处理过程与对位置点数据的处理过程类似，由于速度的数据波动较大，在经过初次的低通处理过后再进行一次卡尔曼滤波处理，以保证滤波结果足够平滑。在速度数据的基础上，利用最小二乘拟合的方法，估算其加速度参数，判断运动目标的速度变化模式。

具体表现为：

1原始数据的插值处理

从ADS-B系统获取到的原始位置点航迹数据的时间间隔没有固定的采样周期，其时间间隔可能会在较大范围内波动。为了保证航迹点之间的时间间隔尽量保持一致，需要对时间间隔较大的航迹段进行插值处理。

1.1插值操作的航迹点时间间隔要求

由于插值处理生成的航迹点会存在误差，在本发明提供的方法中，尽可能的使用原始的航迹数据，设定一个时间间隔的合理取值范围(t_min,t_max]以及需截断处理的时间间隔阈值t_T，在取值范围以内的航迹点均予以保留。对于时间间隔小于t_min的相邻航迹点，在其之间进行插值会导致数据过于密集，此时将当前读取到的航迹点视为错误点并忽略，继续读取下一个航迹点再次判断。若航迹点的时间间隔在(t_max,t_T]区间内，则需要航迹点进行插值运算。若相邻航迹点时间间隔过大，即时间间隔大于t_T，进行插值操作的实际意义不大，此时应对航迹进行截断处理，将整条航迹进行分割，将先前运算的结果保留，并从下一个航迹点开始重新进行处理。

综上，依据不同的相邻航迹点时间间隔，对航迹点处理方式归纳为4类：

1)时间间隔Δt满足Δt≤t_min，不进行插值运算，并忽略当前读取的航迹点；

2)时间间隔Δt满足Δt∈(t_min,t_max]，保留当前读取的航迹点，但不需要进行插值处理；

3)时间间隔Δt满足Δt∈(t_max,t_T]，需要对航迹点进行插值处理；

4)时间间隔Δt满足Δt＞t_T，不进行插值处理，并对航迹段进行截断；

1.2航迹数据读取过程

在对航迹点的插值处理过程当中，为了保证插值的航迹点能够与原航迹具有相同的变化趋势，需要利用到待插值航迹段端点前后的航迹点数据。本发明处理方法是为每一个航班创建一个大小为n的航迹点缓存区，n为偶数，缓存区结构如图2所示。

航迹点处理程序从原始数据源中依次读取航迹点数据，在读取过程中判断当前读取的航迹点与上一个航迹点的时间间隔是否小于截断时间阈值t_T，如果时间间隔满足条件，则将数据存入到航迹点缓存区。当航迹点缓存区存储的航迹点个数达到n时对缓冲内的数据进行进一步的插值处理。每一次的插值处理会删除最早存入缓存区的航迹点数据，此时继续从数据源中读取数据存入至航迹点缓存区，保证插值处理过程中航迹点缓存区的大小为n。

1.3航迹点的插值处理过程

对航迹点的插值处理过程如图2所示。每当航迹点缓存区存入新的航迹点数据时，即对缓存区内下标为

及

的航迹点的时间间隔进行判断，若时间满足插值条件，利用缓冲区的n个航迹点计算三次样条插值曲线，在第

和第

的航迹点之间按照固定时间间隔计算产生插值数据并存入到待插入航迹点缓存区以备之后依次插入到航迹点容器进行进一步的处理，插值处理结束后清除航迹点缓存区第0航迹点数据，并继续读入新的数据。在航迹点的读取过程中，已经筛除时间间隔过大的情况，若第

和第

的航迹点间隔小于t_min，则忽略第

个航迹点，即从航迹点缓存区删除该航迹点，并继续读取下一个航迹点。若第

和第

的航迹点间隔在t_min和t_max之间，无需进行插件处理，清除航迹点缓存区第0航迹点数据，并继续读入新的数据。

2航迹点的滤波处理

插值后的航迹点序列是以大量的原始数据为主体的航迹点序列，由于受到采样设备以及外界干扰因素等的影响，原始数据包含大量的噪声信息，对运动模式的判别以及参数的估算造成很大的干扰，因此在还需对插值后的航迹点进行滤波平滑。

2.1零相移的巴特沃斯低通滤波处理

目标运动轨迹的变化相对较为缓慢,而噪声的变化则相对较快,即信号的有用成分和希望去除的成分各自占有不同的频带。因此可在信号频域上分析滤波问题,通过经典滤波器将噪声信号有效地去除。为了满足实时性处理的要求，还需要克服传统低通滤波算法的迟滞效应，在本发明采用的方法是利用巴特沃斯低通滤波器与零相位滤波算法结合的方法实现对迟滞现象的消除。

巴特沃斯滤波器的系数线性差分方程的形式如式(2.1)：

其中，x_in(n)序列为滤波前的信号序列，a_k,β_j为H(z)系统函数分母与分子的系统数组，求出的y_btw(n)即为滤波后的信号序列；

零相移滤波器能够使得经滤波过后的信号相位响应为零，在具体的实现当中，零相移滤波器使用到时域的翻转，即将一个信号序列按时间先后顺序翻转成为另一个信号序列。经过滤波处理后，输入信号序列x(n)的频谱只是在幅度上被频率响应函数有所修改，而相位无变化，即实现了零相移滤波。根据时域翻转前后信号序列的Z变换之间的关系，零相移滤波器的构成框架可如图3所示。在具体的实现过程中先利用巴特沃斯滤波算法对原始数据进行一次正向滤波，之后再对正向滤波后的结果进行一次逆向的处理，最终结果即为滤波结果。

在对航迹点进行滤波处理时，首先选择合适的截止频率和通频带计算确定巴特沃斯低通滤波器系统函数的分子分母系数，确定零相移巴特沃斯数字低通滤波器的形式，由输入信号得到滤波结果。由于航迹点数据是实时更新上报的，为满足实时性要求，每当经插值处理后的生成的航迹点满m个时，即对已缓存的航迹数据进行一次零相移低通滤波运算，保留计算结果末尾的m个航迹点，作为初次滤波结果实时缓存。

2.2平滑滤波处理

为了进一步保证航迹数据的平滑性，本发明提供的方法采用长度为p的滑动窗口进行移动平均平滑，其移动平均平滑过程如图4所示。(a₀,a₁,a₂,…,a_l)构成原始的数据序列，从原始序列的首个元素开始即可进行移动平均平滑处理，平滑结果为(b₀,b₁,b₂,…,b_l)，但是数据序列最开始的两个数据a₀和a₁前后的数据不足构成长度为p的滑动窗口，在处理结果上令b₀＝a₀,

末尾的数据也作相同处理。一般情况下，平滑计算结果b_i为对应原始数据段的均值，在实时运算过程当中由于不能确定当前读取到的航迹点数据是否是某条航迹的最末的航迹点，因此对于长度一定的原始数据序列，经过一次平滑滤波后的结果序列长度相对变短，图3中末尾处的数据暂时不计算其数值，若确定a_l为航迹段的末尾数据即可立即更新。

3空中目标运动模式识别

3.1飞行航向角特征计算

由两两位置点可以唯一确认空中目标运动的航向角度，在计算角度之前首先需要将位置点的经纬高度坐标转换至以地心为原点的ECEF坐标系下。ECEF坐标系是固定在WGS84参考椭球模型上的随地球旋转的地心地固坐标系，坐标转换过程如下:

其中，lon为纬度，lat为经度，h为距离椭球表面的高度，R为地球的曲率半径；任意两个航迹点p1,p2经过坐标转换得到其在ECEF系下的坐标分别为(x_p1,y_p1,z_p1)和(x_p2,y_p2,z_p2)，那么由航迹点i和j确定的方向矢量的水平投影角

和竖直方向的角度θ计算公式为：

其中Δx＝x_b-x_a,Δy＝y_b-y_a,Δz＝z_b-z_a。

在实时获得的经过滤波平滑后的航迹点数据的基础上计算两两航迹点间的航向角度值，每当有航迹点数据完成滤波平滑处理，即实时读取航迹点位置，并计算其相对前一个航迹点的航向角，计算结果更新至航向角数据的缓存区。

3.2航向角数据的实时滤波处理

对空中目标运动模式的判断从水平投影和竖直方向进行分析，对圆弧运动的判定主要基于水平投影角度特征。计算出的航向角度结果会存在部分噪声信号，因此在对运动模式判别之前还需对航向角数据进行滤波处理。FIR滤波器是一种非递归型滤波器，因此，FIR滤波器在保证满足幅频响应要求的同时还可获得严格的线性相位特性，从而保持稳定。对因果的FIR系统，其系统函数仅有零点(除z＝0的极点外)，因此系数a_k(k≠0)全为零，差分方程就简化为：

经过滤波处理后的航向角数据具有很好的平滑特性。

从航向角数据缓存区实时读取水平投影角数据，每当读取的水平投影角数据个数满m时，即对当前读取的水平投影角度进行一次FIR滤波运算，计算结果更新至滤波后航向角数据的缓存区。

3.3目标运动模式识别

3.3.1竖直方向运动模式判别

空中目标在竖直方向上的运动模式有爬升运动、下降运动和水平运动三种模式，实时读取经过滤波平滑后的航迹点的高度数值，并计算两两航迹点之间的高度差Δh,设定爬升时高度变化的阈值h_up和下降时高度变化阈值h_down，其中h_up为正值，h_down为负值。当Δh＞h_up时，目标为爬升运动；当Δh<h_down时，目标为下降运动；当Δh介于h_down和h_up时，目标为水平运动。

3.3.2目标运动轨迹模式判别

对航向角变化的判断可以根据每段航迹点数据的标准差作为判据。在航向角的变化阶段，计算的出的标准差较稳定阶段有很大的突变，设置合适的阈值条件即可判断出轨迹模式的改变。

当滤波后航向角数据的缓存区发生更新时，即从缓存区读取角度数据，每当读取的角度数据个数满m时，计算m个角度值的标准差

设定标准差的变化阈值s_T,比较计算出的标准差与阈值的大小，若

则目标为直线运动模式，若

则目标为圆弧运动模式。

4空中目标运动参数估计

4.1空中运动速度的滤波处理

4.1.1目标运动速度计算

在计算相邻两个航迹点之间的航向角的同时可以计算得出航迹段的速度值，速度计算方法为相邻两个航迹点间的空间距离与时间间隔的商。在计算飞行航向角过程中已经得到在ECEF坐标系下，相邻航迹点在XYZ轴上的差值，速度计算值为

实时计算得到的速度数据实时存入至速度数据缓存区。

4.1.2初次低通滤波处理

由于原始航迹数据的粗糙性以及插值过程的可能与实际运动情况不一致，计算出的速度数据掺杂了大量的高频噪声。速度数据的滤波处理流程与航迹点的滤波处理过程类似，采用巴特沃斯低通滤波算法和零相移滤波算法过滤掉速度信号中的高频噪声。从速度数据缓存区实时读取数据，每当读取的数据满m个时，即对已读取的速度数据进行一次零相移低通滤波运算，保留计算结果末尾的m个计算结果，作为初次滤波结果实时缓存。

4.1.3卡尔曼滤波处理

为了保证速度数据的良好的平滑特性，对滤波结果再进行卡尔曼滤波平滑。在空中目标的的运动过程当中，其运动状态变化较为平稳，并且数据采集的周期较短，经过初次滤波处理过后，速度数值的变化相对平滑，因此建立的卡尔曼滤波的系统离散方程和观测方程采用最简单的跟踪形式，系统参数和测量系统参数的值取1。根据卡尔曼滤波算法的迭代规则，输入的原始信号可以实时地输出滤波结果，不会影响数据处理的实时性。

4.2加速度和角速度的估算

4.2.1加速度参数估算

在较短的时间周期内，空中目标运动的速度和角速度关于时间的变化可以近似为线性变化。计算得到平滑的速度数据之后，按m个速度数据为单位，采用最小二乘拟合的方法拟合出速度关于时间的线性函数，目标实时运动的速度v＝v₀+at，那么线性方程的斜率a即可近似为目标运动的加速度。在本发明提供的方法中，还计算了由所选取的速度数据的首尾端点值确定的速度关于时间的变化率a′，以a和a′的平均值

作为空中目标运动的加速度估计值。

4.2.2目标速度变化模式判别

根据加速度的取值，结合阈值条件，可以判断目标速度的变化模式：当超出一定阈值判断为加速运动，低于某阈值判断为减速运动，否则即视为匀速运动。设定速度模式变化的加速度阈值a_acc、a_dec，a_acc为正值而a_dec为负值。当加速度估计值

时，运动模式为加速运动；当加速度估计值

时，运动模式为减速运动；当加速度估计值

介于a_dec和a_acc之间时，运动模式判别为匀速运动。

4.2.3角速度参数估算

与加速度的计算过程类似，角度在短时间内的变化规律也是线性变化的过程，θ＝θ₀+ωt。当目标的运动模式判定为圆弧运动时，以m个角度数据为单位，拟合角度关于时间的变化斜率ω，计算首尾端点处角度关于时间的变化率ω′，以ω和ω′的平均值

作为目标运动的角速度估计值。

4.3圆弧运动模式下圆弧半径及圆心坐标的估算

当运动目标的航向角发生改变时，其与运动轨迹近似的看作空间圆弧，空间圆弧为空间圆与空间平面的交线，如图5。

设球面方程为：

Ax+By+Cz-D＝x²+y²+z² (4.1)

设m个数据点所在的空间平面为：

A′x+B′y+C′z+D′＝0 (4.2)

由m组数据经过最小二乘法的计算，得到球面方程和平面方程的各项系数，并计算球面与平面的交线，从而确定圆弧运动所在空间圆圆心和其半径。

Claims (9)

1.一种实时的空中目标运动模式识别及参数估计方法，其特征在于，包括如下步骤：

(1)读取原始数据，计算当前数据与对应的航迹的末数据的时间间隔；

(2)判断时间间隔是否超过截断时间间隔阀值，如果不超过则进入步骤(3)，如果超过则对航迹数据进行截断处理并进入步骤(13)；

(3)依据相邻航迹点时间间隔条件进行相应的插值处理，将处理后的航迹点存入待插入航迹点缓存区；

(4)依次将待插入航迹点缓存区的航迹点插入到航迹点容器；

(5)航迹点容器中的航迹点每满m个航迹点，即对航迹点进行零相移的巴特沃斯低通滤波处理，并取最后m个航迹点待处理；

(6)计算航迹点之间的高度差，判断竖直方向的运动模式；

(7)计算航迹点之间的航向角和速度；

(8)对最新计算出的m个航向角和速度进行滤波处理；

(9)计算滤波后航向角的标准差，判断是否是圆弧运动，如果是则进入步骤(10)，如果不是则进入步骤(12)；

(10)拟合航向角关于时间的一次曲线，计算角速度；

(11)根据最新处理过的m个航迹点，计算圆心坐标并进入步骤(12)；

(12)拟合速度关于时间的一次曲线，计算加速度；

(13)数据导入至数据库。

2.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(2)对于判断时间间隔是否超过截断时间间隔阀值的具体步骤如下：

(2.1)设定时间间隔的合理取值范围为(t_min,t_max]，需截断处理的时间间隔阈值为t_T；

(2.2)依据不同的相邻航迹点时间间隔，对航迹点处理方式归纳为4类：

1)时间间隔Δt满足Δt≤t_min，不进行插值运算，并忽略当前读取的航迹点；

2)时间间隔Δt满足Δt∈(t_min,t_max]，保留当前读取的航迹点，但不需要进行插值处理；

3)时间间隔Δt满足Δt∈(t_max,t_T]，需要对航迹点进行插值处理；

4)时间间隔Δt满足Δt＞t_T，不进行插值处理，并对航迹段进行截断。

3.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(5)中对航迹点进行巴特沃斯低通滤波处理的具体步骤如下：利用巴特沃斯低通滤波器与零相位滤波算法结合的方法实现对迟滞现象的消除；巴特沃斯滤波器的系数线性差分方程的形式如式(2.1)：

其中，x_in(n)序列为滤波前的信号序列，a_k,β_j为H(z)系统函数分母与分子的系统数组，求出的y_btw(n)即为滤波后的信号序列；

首先选择合适的截止频率和通频带计算确定巴特沃斯低通滤波器系统函数的分子分母系数，确定零相移巴特沃斯数字低通滤波器的形式，由输入信号得到滤波结果；每当经插值处理后的生成的航迹点满m个时，即对已缓存的航迹数据进行一次零相移低通滤波运算，保留计算结果末尾的m个航迹点，作为初次滤波结果实时缓存。

4.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(6)中判断竖直方向的运动模式的具体步骤如下：空中目标在竖直方向上的运动模式有爬升运动、下降运动和水平运动三种模式，实时读取经过滤波平滑后的航迹点的高度数值，并计算两两航迹点之间的高度差Δh,设定爬升时高度变化的阈值h_up和下降时高度变化阈值h_down，其中h_up为正值，h_down为负值；当Δh＞h_up时，目标为爬升运动；当Δh<h_down时，目标为下降运动；当Δh介于h_down和h_up时，目标为水平运动。

5.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(7)中计算航迹点之间的航向角的具体步骤如下：

由两两位置点可以唯一确认空中目标运动的航向角度，在计算角度之前首先需要将位置点的经纬高度坐标转换至以地心为原点的ECEF坐标系下；ECEF坐标系是固定在WGS84参考椭球模型上的随地球旋转的地心地固坐标系，坐标转换过程如下:

其中，lon为纬度，lat为经度，h为距离椭球表面的高度，R为地球的曲率半径；任意两个航迹点p1,p2经过坐标转换得到其在ECEF系下的坐标分别为(x_p1,y_p1,z_p1)和(x_p2,y_p2,z_p2)，那么由航迹点i和j确定的方向矢量的水平投影角

和竖直方向的角度θ计算公式为：

其中Δx＝x_p2-x_p1,Δy＝y_p2-y_p1,Δz＝z_p2-z_p1；

在实时获得的经过滤波平滑后的航迹点数据的基础上计算两两航迹点间的航向角度值，每当有航迹点数据完成滤波平滑处理，即实时读取航迹点位置，并计算其相对前一个航迹点的航向角，计算结果更新至航向角数据的缓存区。

6.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(9)中判断是否是圆弧运动的具体步骤如下：当滤波后航向角数据的缓存区发生更新时，即从缓存区读取角度数据，每当读取的角度数据个数满m时，计算m个角度值的标准差

设定标准差的变化阈值s_T,比较计算出的标准差与阈值的大小，若

则目标为直线运动模式，若

则目标为圆弧运动模式。

7.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(12)中计算加速度的具体步骤如下：计算得到平滑的速度数据之后，按m个速度数据为单位，采用最小二乘拟合的方法拟合出速度关于时间的线性函数，目标实时运动的速度v＝v₀+at，那么线性方程的斜率a_v即可近似为目标运动的加速度；计算由所选取的速度数据的首尾端点值确定的速度关于时间的变化率a_v′，以a_v和a_v′的平均值

作为空中目标运动的加速度估计值。

8.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(10)中计算角速度的具体步骤如下：与加速度的计算过程类似，角度在短时间内的变化规律也是线性变化的过程，θ＝θ₀+ωt；当目标的运动模式判定为圆弧运动时，以m个角度数据为单位，拟合角度关于时间的变化斜率ω，计算首尾端点处角度关于时间的变化率ω′，以ω和ω′的平均值

作为目标运动的角速度估计值。

9.根据权利要求1所述的一种实时的空中目标运动模式识别及参数估计方法，其特征在于，所述步骤(11)中计算圆心坐标的具体步骤如下：球面方程为：

Ax+By+Cz-D＝x²+y²+z²

m个数据点所在的空间平面为：

A′x+B′y+C′z+D′＝0

由m组数据经过最小二乘法的计算，得到球面方程和平面方程的各项系数，并计算球面与平面的交线，从而确定圆弧运动所在空间圆圆心和其半径。

Images