首页> 中国专利> 一种结合目标角度的外辐射源雷达多目标跟踪方法

一种结合目标角度的外辐射源雷达多目标跟踪方法

摘要

本发明属于外辐射源雷达多目标跟踪技术领域,特别涉及一种结合目标角度的外辐射源雷达多目标跟踪方法。该结合目标角度的外辐射源雷达多目标跟踪方法包括以下步骤:对外辐射源雷达接收的数字电视信号依次进行解调、重编码,重构出各个发射站的直达波信号;利用每个发射站的直达波信号对外辐射源雷达接收的数字电视信号进行自适应杂波相消处理,得到杂波相消后信号;根据所述杂波相消后信号、以及重构出的各个发射站的直达波信号,得到各个目标的观测量;根据每个目标相对于外辐射源雷达的角度,采用概率假设密度粒子滤波方法对外辐射源雷达矩形观测区域内的每个目标实施定位跟踪,每个目标的定位跟踪方法为到达时间定位方法。

著录项

  • 公开/公告号CN104077498A

    专利类型发明专利

  • 公开/公告日2014-10-01

    原文格式PDF

  • 申请/专利权人 西安电子科技大学;

    申请/专利号CN201410351145.6

  • 发明设计人 王俊;王海环;

    申请日2014-07-22

  • 分类号

  • 代理机构西安睿通知识产权代理事务所(特殊普通合伙);

  • 代理人惠文轩

  • 地址 710071 陕西省西安市太白南路2号

  • 入库时间 2023-12-17 01:49:17

法律信息

  • 法律状态公告日

    法律状态信息

    法律状态

  • 2017-06-20

    授权

    授权

  • 2014-10-29

    实质审查的生效 IPC(主分类):G06F19/00 申请日:20140722

    实质审查的生效

  • 2014-10-01

    公开

    公开

说明书

技术领域

本发明属于外辐射源雷达多目标跟踪技术领域,特别涉及一种结合目标 角度的外辐射源雷达多目标跟踪方法。本发明利用外辐射源雷达系统中采用 多发一收结构,将目标双基地距离、多普勒速度与目标角度信息相结合,同 时对多个目标进行定位跟踪,本发明可在杂波背景下,实现对多个目标的准 确跟踪,运算量小,跟踪精度高。

背景技术

随着现代电子技术的发展,传统的有源雷达在电子干扰、反辐射导弹、 超低空突防及隐身技术等威胁中面临严重的生存危机。为了应对现代电子战 争中的四大威胁,外辐射源雷达成为雷达领域的研究热点。外辐射源雷达本 身不发射电磁波,而是利用已存在的民用信号(如FM调频信号、GMS(Global  System for Mobile communication,全球移动通讯系统)、电视信号等)作为第 三方照射源,对目标进行探测和定位,具有生存能力强、反隐身、抗低空突 防及成本低等优势。但是由于民用信号的发射功率较低,目标回波通常被直 达波和多径覆盖,因此实现目标探测和跟踪前需要先对回波信号进行处理, 提取所需的目标参数作为接收站观测量。

利用外辐射源雷达实现多目标跟踪是雷达跟踪领域中的研究热点和难点。 一方面多目标跟踪在民用和军事上都有重要应用;另一方面外辐射源雷达多 目标跟踪要求能从包含杂波且数目时变的观测量中提取出目标个数不定的目 标状态,为实现这一目的,需要配对目标与观测量、发射站与观测量。传统 的处理方法基于数据关联,存在计算量随时间增长的问题。为避免数据关联, Mahler提出概率假设密度(PHD)滤波器,随后概率假设密度粒子滤波器 (PF‐PHD)的提出将PHD算法由理论层面带入工程实践中。

虽然PF‐PHD可以避免数据关联,但需要大量粒子以保证跟踪性能,计算 量太大。多种观测量的结合会提高跟踪精度,在外辐射源雷达跟踪中最常用 的观测量是目标的双基地距离和多普勒速度,但只利用这两种观测量概率假 设密度粒子滤波器(PF‐PHD)仍然需要大量粒子来保证跟踪性能,运算量大, 实时性差。

发明内容

本发明的目的在于提出一种结合目标角度的外辐射源雷达多目标跟踪方 法。本发明以数字电视信号作为外辐射源雷达的第三方照射源,从数字电视 信号的信道估计参数中直接获取目标时延和多普勒频移,简化目标跟踪的前 期信号处理过程;在PF‐PHD滤波器中,通过增加目标角度信息,优化新生粒 子的分布,实现用少量粒子对多个目标进行高精度跟踪。

为实现上述技术目的,本发明采用如下技术方案予以实现。

一种结合目标角度的外辐射源雷达多目标跟踪方法包括以下步骤:

步骤1,利用外辐射源雷达接收数字电视信号,对外辐射源雷达接收的数 字电视信号依次进行解调、重编码,重构出各个发射站的直达波信号;

步骤2,利用每个发射站的直达波信号对外辐射源雷达接收的数字电视信 号进行自适应杂波相消处理,得到杂波相消后信号;根据所述杂波相消后信 号、以及重构出的各个发射站的直达波信号,得到各个目标的观测量;

步骤3,根据每个目标相对于外辐射源雷达的角度,采用概率假设密度 粒子滤波方法对外辐射源雷达矩形观测区域内的每个目标实施定位跟踪,每 个目标的定位跟踪方法为到达时间定位方法。

本发明的特点和进一步改进在于:

所述步骤2的具体子步骤为:

(2.1)利用每个发射站的直达波信号对外辐射源雷达接收的数字电视信 号进行自适应杂波相消处理,得到杂波相消后信号;

(2.2)对所述杂波相消后信号、以及重构出的各个发射站的直达波信号 进行距离—多普勒二维相关运算,得到各个目标的观测量;每个目标的观测 量包括:对应目标的双基地距离和、对应目标的多普勒速度、以及对应目标 相对于外辐射源雷达的角度。

所述步骤3的具体子步骤为:

(3.1)设置参数k,k=0,1,2...;当k=0时,在k时刻用L0个粒子表示对 应目标状态的先验分布,所述L0个粒子表示为表示0时刻所述L0个粒子中第i0个粒子的权值,表示所述L0个粒子中第i0个粒子在0时刻的粒子状态;

(3.2)当k≥1时,将k-1时刻的粒子集为Lk-1为k-1时刻 的粒子数,ik-1取1至Lk-1,为k-1时刻第ik-1个粒子的权值,为k-1时 刻第ik-1个粒子的粒子状态;在k-1时刻每个粒子处产生一组西格玛点,k-1时 刻每个粒子处产生的一组西格玛点包括至以及至n为粒子状态的维数;

将k-1时刻每个粒子处产生的一组西格玛点带入对应目标状态转移方程 中,得到k时刻的一步预测状态值j取0至2n;对至进行加权求和,得到k时刻粒子状态的初步一步预测值计算出k时刻 的初步一步预测状态协方差矩阵

将代入目标观测方程中,得到k时刻的一步预测观测值所述目标观测方程为zk=h(xkk),其中,νk为观测噪声,h(xkk)表示对应 目标的观测函数;对至进行加权求和,得到k时刻观测量的初步 一步预测值

在k时刻对应目标的观测集Zk中,选取出与具有最近欧式距离的对 应目标观测量zk,k时刻对应目标的观测集Zk为:由k时刻对应目标的观测 量组成的集合;用zk对粒子状态的初步一步预测值和的初步一步预 测状态协方差矩阵做进一步修正,得到初步一步预测均值及的 初步一步预测状态协方差矩阵在均值为协方差为的高斯分 布上采样,得到k时刻第ik-1个粒子的粒子状态的一步预测值

根据的概率密度qk,并按照下式计算出所对应的权值

ωk|k-1(ik-1)=ωk-1(ik-1)ps,kfk|k-1(xk|k-1(ik-1)|xk|k-1(ik-1))/qk,ik-1=1,2,···,Lk-1

其中,ps,k为k-1时刻的对应目标在k时刻仍然存在的概率,为对应目标状态转移概率密度函数;

(3.3)当k≥1时,得出k时刻每个发射站对应的每个目标的双基地双基 地椭圆,在k时刻每个发射站对应的每个目标的双基地椭圆曲线上,每个点对 应的双基地距离和为一固定值,每个点对应的双基地距离和指:每个点与对 应发射站的距离、以及对应点与外辐射源雷达的距离之和;

根据k-1时刻每个目标相对于外辐射源雷达的角度、以及k时刻每个目标 相对于外辐射源雷达的角度,在k时刻每个发射站对应的对应目标的双基地椭 圆曲线与外辐射源雷达的矩形观测区域边缘的交点处,分布k时刻新生粒子, 并计算k时刻每个新生粒子的权值,k时刻新生粒子的个数表示为Jk,k时刻 第inew个新生粒子的粒子状态表示为k时刻第inew个新生粒子的权值表 示为inew取Lk-1+1至Lk-1+Jk

(3.4)当k≥1时,根据子步骤(3.2)和子步骤(3.3),得到预测的k时 刻的粒子集其中,表示预测的k时刻的粒子集中第i'k个粒子的粒子状态,表示预测的k时刻的粒子集中第i'k个粒子的权值,i'k取1至Lk-1+Jk;利用k时刻对应目标的观测集Zk,对预测的k时刻的粒子集 中第i'k个粒子的权值进行修正,得出预测的k时刻的粒子集中第i'k个粒 子的修正后权值

(3.5)按照预测的k时刻的粒子集中第i'k个粒子的修正后权值对 预测的k时刻的粒子集进行重采样,得出k时刻的粒子集;对k时刻的粒子集 进行聚类分析,得到每个目标k时刻的状态。

在子步骤(3.3)中,首先,根据k-1时刻目标的观测集Zk-1,获取k-1时刻 每个目标相对于外辐射源雷达的角度,将k-1时每个目标相对于外辐射源雷达 的角度组成k-1时刻目标观测量角度集合θk-1,k-1时刻目标观测量角度集合 θk-1中元素的个数表示为<θk-1>,k-1时刻目标观测量角度集合θk-1中第αk-1个元 素表示为αk-1取1至<θk-1>;根据k时刻目标的观测集Zk,获取k时刻 每个目标相对于外辐射源雷达的角度,将k时刻每个目标相对于外辐射源雷达 的角度组成k时刻目标观测量角度集合θk,k时刻目标观测量角度集合θk中元 素的个数表示为<θk>,k时刻目标观测量角度集合θk中第αk个元素表示为 αk取1至<θk>;

判断k时刻是否产生新生粒子,如果k时刻不产生新生粒子,则k时刻新 生粒子的个数Jk为零;如果k时刻产生新生粒子,则得出k时刻产生的每个新 生粒子的分布区域中心点方位角,k时刻产生的每个新生粒子的分布区域中心 点方位角指:k时刻产生的对应新生粒子的位置的纵坐标与横坐标的比值的反 正切值;

判断k时刻是否产生新生粒子的过程为:将αk-1从1至<θk-1>进行遍历,αk-1每取一个数值时,在k时刻目标观测量角度集合θk中找出满足第一新生粒子存 在条件的元素,所述第一新生粒子存在条件为:k时刻目标观测量角度集合θk中的元素与的差值的绝对值不小于γθ,γθ为设定的角度门限;在将αk-1从1至<θk-1>进行遍历的过程中,若k时刻目标观测量角度集合θk中不存在满 足第一新生粒子存在条件的元素,则k时刻不产生新生粒子;反之,k时刻产 生新生粒子;

计算k时刻每个发射站对应的每个目标的双基地椭圆曲线与外辐射源雷 达的矩形观测区域边缘的交点的位置;然后,根据k时刻每个发射站对应的每 个目标的双基地椭圆曲线与外辐射源雷达的矩形观测区域边缘的交点的位置, 在外辐射源雷达的矩形观测区域边缘中的上边缘、下边缘、左边缘或右边缘 中,找出至少一个满足以下条件的边缘:与k时刻每个发射站对应的每个目标 的双基地椭圆曲线存在交点;将找出的满足以上条件的边缘记为:k时刻目标 的候选新生粒子分布边缘;

得出k时刻目标的新生粒子分布边缘;得出k时刻目标的新生粒子分布边 缘的过程为:当k时刻目标的候选新生粒子分布边缘的个数等于1时,k时刻 目标的新生粒子分布边缘为k时刻目标的候选新生粒子分布边缘;当k时刻目 标的候选新生粒子分布边缘的个数大于1时,在k时刻目标的所有候选新生粒 子分布边缘中找出k时刻目标的新生粒子分布边缘;

当k时刻目标的候选新生粒子分布边缘的个数等于1时,在得出k时刻目 标的新生粒子分布边缘之后,在k时刻目标的新生粒子分布边缘、以及k时刻 每个发射站对应的每个目标的双基地椭圆曲线的所有交点中,找出相距最近 的l个交点,这l个交点与l个发射站一一对应;得出所述相距最近的l个交点 的位置的均值,将所述相距最近的l个交点的位置的均值为:k时刻共同交点; 根据目标观测方程,针对k时刻共同交点,得出k时刻共同交点对应的目标速 度;

当k时刻目标的候选新生粒子分布边缘的个数大于1时,在k时刻目标的 每个候选新生粒子分布边缘、以及k时刻每个发射站对应的每个目标的双基地 椭圆曲线的所有交点中,找出相距最近的l个交点,这l个交点与l个发射站一 一对应;得出所述相距最近的l个交点的位置的均值,将所述相距最近的l个 交点的位置的均值为:对应候选新生粒子分布边缘的k时刻共同交点;根据目 标观测方程,针对每个候选新生粒子分布边缘的k时刻共同交点,得出对应候 选新生粒子分布边缘的k时刻共同交点对应的目标速度;

当k时刻目标的候选新生粒子分布边缘的个数大于1时,找出k时刻目标 的新生粒子分布边缘的过程为:判断k时刻是否产生新生粒子,如果k时刻不 产生新生粒子,则没有k时刻目标的新生粒子分布边缘;如果k时刻产生新生 粒子,针对k时刻的每个新生粒子,判断每个候选新生粒子分布边缘的k时刻 共同交点的方位角是否满足第二新生粒子存在条件,所述第二新生粒子存在 条件为:每个候选新生粒子分布边缘的k时刻共同交点的方位角与k时刻目标 观测量角度集合θk中满足第一新生粒子存在条件的元素的差值的绝对值不大 于设定的方位角门限值;如果每个候选新生粒子分布边缘的k时刻共同交点的 方位角不满足第二新生粒子存在条件,则令k时刻新生粒子的个数Jk取0;如 果至少一个候选新生粒子分布边缘满足第二新生粒子存在条件,则在满足第 二新生粒子存在条件的候选新生粒子分布边缘中,选择其k时刻共同交点的方 位角与k时刻的对应新生粒子的分布区域中心点方位角最接近的候选新生粒 子分布边缘,则k时刻目标的新生粒子分布边缘为:在满足第二新生粒子存在 条件的候选新生粒子分布边缘中选择的候选新生粒子分布边缘,此时,k时刻 共同交点为:k时刻目标的新生粒子分布边缘的k时刻共同交点,k时刻共同 交点对应的目标速度为:k时刻目标的新生粒子分布边缘的k时刻共同交点对 应的目标速度;

在k时刻共同交点处采样Jk个点,得出k时刻Jk个新生粒子;k时刻每 个新生粒子的粒子状态包括:k时刻对应新生粒子的方位状态、k时刻对应新 生粒子的速度状态;

按下式计算k时刻第inew个新生粒子的权值

ωk|k-1(inew)=Qx·qx(xk(inew)|zk)·Qy·qy(yk(inew)|zk)·Qx··qx·(x·k(inew)|zk)·Qy··qy·(y·k(inew)|zk)

其中,inew取Lk-1+1至Lk-1+Jk,为k时刻第inew个新生粒子的位置 在x方向上的概率密度函数,为k时刻第inew个新生粒子的位置在 y方向上的概率密度函数,为k时刻第inew个新生粒子的速度在x 方向上的概率密度函数,为k时刻第inew个新生粒子的速度在y方 向上的概率密度函数;当k时刻目标的新生粒子分布边缘为外辐射源雷达的矩 形观测区域边缘中的上边缘或下边缘时,当k时 刻目标的新生粒子分布边缘为外辐射源雷达的矩形观测区域边缘中的左边缘 或右边缘时,Qx=Qx·=2,Qy=Qy·=1.

在子步骤(3.4)中,当k≥1时,当k≥1时,对预测的k时刻的粒子集中 第i'k个粒子的权值进行修正的过程包括以下子步骤:

a)设置迭代参数j,j=1,2,...当j=1时,执行子步骤b);

b)根据以下公式计算第j个发射站在处的检测概率

pd,k(xk|k-1(ik),j)=Q[2SNR(xk|k-1(ik),j),2ln(1/pFA)]

其中,表示已知的第j个发射站对应的处的信噪比,pFA为 设定的虚警概率,Q[a,b]=bxexp(-x2+a22)I0(ax)dx,I0表示零阶贝塞尔函 数;

c)按照下式计算k时刻杂波平均数λk和杂波概率密度ck

λk=L1l1·L2l2·pFA,ck=1/(L1L2)

其中,L1为外辐射源雷达接收信号的距离范围,L2为外辐射源雷达接收信号 的多普勒频率范围,l1外辐射源雷达接收信号的距离分辨率,l2为外辐射源雷 达接收信号的多普勒分辨率;

d)在k时刻目标的观测集Zk中,利用与第j个发射站对应的k时刻每个 目标的观测量,组成k时刻第j个发射站的目标的观测集然后利用中 的观测量计算似然函数

g(zk(j)|xk|k-1(ik))=12πRvexp(-12(zk(j)-h(xk|k-1(ik),vk))Rv-1(zk-h(xk|k-1(ik),vk))T)

其中,Rv是观测噪声νk的协方差矩阵,上标‐1表示矩阵的逆,上 标T表示矩阵或向量的转置;

e)按照下式计算出预测的k时刻的粒子集中第i'k个粒子的修正后权值

ωk(ik)=ωk|k-1(ik)·[1-pd,k(xk|k-1(ik))+Σzk(j)Zk(j)ψk,z(xk|k-1(ik))λkck+<ωk|k-1,ψk,z>]

ψk,z(xk|k-1(ik))=pd,k(xk|k-1(ik))·g(zk(j)|xk|k-1(ik))

<ωk|k-1,ψk,z>=Σik=1Lk-1+Jkψk,z(xk|k-1(ik))·ωk|k-1(ik)

其中,i'k取1至Lk-1+Jk

然后,令j值自增1;

f)判断j与发射站的个数为l的大小关系,如果j<l,则令然后返回至子步骤b),如果j=l,则迭代终止,预测的k时刻的粒子集中第i'k个粒子的修正后权值为

在子步骤(3.5)中,首先根据预测的k时刻的粒子集中第i'k个粒子的修 正后权值计算k时刻的估计目标数

针对每个目标固定采样m个粒子,对预测的k时刻的粒子集进行重采样, 得到k时刻的粒子集ik取1至Lk,表示k时刻 的粒子集中第ik个粒子的粒子状态,表示k时刻的粒子集中第ik个粒子权 值。

在步骤3之后,针对每个目标0时刻至k时刻的状态,进行航迹关联, 剔除虚假航迹,得出每个目标的航迹。

本发明的有益效果为:

第一,本发明通过数字电视信号的信道估计参数直接得到所需目标状态 参数(包括时延和多普勒频率),不需要用每个发射站的直达波与其相应的回 波做距离—多普勒二维相关处理,从而简化了外辐射源雷达信号处理流程。

第二,本发明除利用目标的双基地距离和多普勒速度外,还结合目标的 角度信息实现目标跟踪。在分布新生粒子前,先利用角度信息判断当前时刻 是否有新目标出现,而不是在每一时刻都产生新生粒子,从而降低了运算量。

第三,本发明在分布新生粒子时,利用角度信息剔除了观测区域边缘上 虚假的共同交点,缩小了新生粒子分布范围,进一步降低跟踪所需粒子数。

第四,本发明利用相对于每个发射站的观测量依次对预测权值进行更新, 有利于剔除对应单个发射站的观测量中可能存在的杂波,从而提高跟踪精度。

附图说明

图1为本发明的使用场景示意图;

图2为本发明应用的多目标到达时间定位原理示意图;

图3为本发明中对每个目标实施定位跟踪的步骤框图;

图4a为仿真实验中每个发射站对应的对应目标的双基地椭圆曲线与外辐 射源雷达的矩形观测区域边缘的交点示意图;

图4b为图4a的局部放大图;

图5a为增加目标角度信息之前和之后的每时刻估计目标数和真实目标数 的对比图;

图5b为增加目标角度信息之前和之后得出的目标方位估计误差的对比图;

图5c为增加目标角度信息之前和之后得出的目标速度估计误差的对比图。

具体实施方式

下面结合附图对本发明作进一步说明:

步骤1,参照图1,为本发明的使用场景示意图。本发明实施例中,利用 外辐射源雷达(即图1中的接收站)接收数字电视信号,所述数字电视信号 包括各个发射站的直达波信号(图1中,3个发射站直接发送到接收站的数字 电视信号)、各个发射站的多径信号以及各个发射站的目标回波信号。本发明 实施例中,发射站的个数为l。

对外辐射源雷达接收的数字电视信号依次进行解调、重编码,重构出各 个发射站的直达波信号。

步骤2,利用每个发射站的直达波信号对外辐射源雷达接收的数字电视信 号进行自适应杂波相消处理,得到杂波相消后信号;根据所述杂波相消后信 号、以及重构出的各个发射站的直达波信号,得到各个目标的观测量。

其具体子步骤为:

(2.1)利用每个发射站的直达波信号对外辐射源雷达接收的数字电视信 号进行自适应杂波相消处理,消除外辐射源雷达接收的数字电视信号中各个 发射站的直达波信号和各个发射站的多径信号,得到杂波相消后信号。

(2.2)根据所述杂波相消后信号、以及重构出的各个发射站的直达波信 号,得到各个目标的观测量。具体地说,对所述杂波相消后信号、以及重构 出的各个发射站的直达波信号进行距离—多普勒二维相关运算,得到各个目 标的观测量。

本发明实施例中,每个目标的观测量包括:对应目标的双基地距离和、 对应目标的多普勒速度、以及对应目标相对于外辐射源雷达(接收站)的角 度。

对应目标的双基地距离和包括:对应目标相对于每个发射站的双基地距 离和,对应目标相对于第l'个发射站的双基地距离和为:对应目标与第l'个发 射站的距离、以及对应目标与外辐射源雷达的距离的和,l'取1至l。对应目 标相对于外辐射源雷达(接收站)的角度指:对应目标的外辐射源雷达连线 方向与水平面的夹角(是否正确),对应目标的外辐射源雷达连线方向为:对 应目标与外辐射源雷达的连线方向。

下面说明对应目标相对于每个发射站的双基地距离和的计算过程、以及 对应目标的速度的计算过程。首先利用数字电视信号的信道估计参数直接得 到对应目标相对于每个发射站的时延、以及对应目标的多普勒频率。对于其 中一个发射站来说,如果对应目标相对于对应发射站的时延为τ,对应目标的 多普勒频率为fd,每个发射站发射信号的波长为γ,则对应目标相对于发射 站的时延为bτ/2,对应目标的多普勒速度为fdγ/2。

步骤3,采用概率假设密度粒子滤波器(PF‐PHD)对外辐射源雷达矩形观测 区域内的每个目标实施定位跟踪,定位方法采用到达时间(TOA)定位方法; 参照图2,为本发明应用的多目标到达时间定位原理示意图。参照图3,为本 发明中对每个目标实施定位跟踪的步骤框图。

步骤3具体包括以下子步骤:

(3.1)初始化概率假设密度粒子滤波器(PF‐PHD)中的粒子状态和权值。 本发明实施例中,粒子状态由粒子的方位状态和速度状态组成。

具体地说,设置时刻参数k,k=0,1,2...;

当k=0时,在k时刻用L0个粒子表示目标状态的先验分布,所述L0个粒子 (即0时刻粒子集)表示为L0为大于1的自然数,i0=1,2,…,L0; 表示0时刻所述L0个粒子中第i0个粒子的权值,与估计目标数和粒子 数目有关,表示0时刻设定的目标的个数;表示所述L0个粒子中第i0个粒子在0时刻的粒子状态。

(3.2)当k≥1时,基于贝叶斯滤波原理,利用k时刻目标的观测集Zk, 通过对k时刻存活粒子在k-1时刻的粒子状态进行无迹卡尔曼滤波,得到k时 刻存活粒子的一步预测状态和权值。k-1时刻的粒子集为Lk-1为k-1时刻的粒子数,ik-1取1至Lk-1,为k-1时刻第ik-1个粒子的权值,为k-1时刻第ik-1个粒子的粒子状态。k时刻目标的观测集Zk为:由k时刻每个 目标的观测量组成的集合。

具体地说,在子步骤(3.2)中,首先在k-1时刻每个粒子处产生一组西 格玛点(Sigma Point)n为粒子状态的维数,例如,n=4。

当j=0时,χk-1(ik-1,j)=xk-1(ik-1),ωk-1(ik-1,j)=λ/(n+λ),λ为尺度因子,

λ=α2(n+ε0)-n,(α<1,ε0=0)

其中,α为小于1的常数,ε0=0。

当j取1至n时,χk-1(ik-1,j)=xk-1(ik-1)+(n+λ)Pk-1,ωk-1(ik-1,j)=λ/2(n+λ),Pk-1为 k-1时刻第ik-1个粒子的状态协方差矩阵。

当j取n+1至2n时,χk-1(ik-1,j)=xk-1(ik-1)-(n+λ)Pk-1,ωk-1(ik-1,j)=λ/2(n+λ).

然后,将k-1时刻每个粒子处产生的一组西格玛点带入目标状态转移方程 中,得到k时刻的一步预测状态值所述目标状态转移方程为: xk=f(xk-1,uk),xk-1为k-1时刻目标的状态,xk为k时刻目标的状态,uk为过 程噪声,f(xk-1,uk)表示对应目标的状态转移函数;对至进行加 权求和,得到k时刻粒子状态的初步一步预测值计算出k时刻的初 步一步预测状态协方差矩阵具体操作过程如下:

χk|k-1(ik-1,j)=f(χk-1(ik-1,j),uk)

xk|k-1(ik-1)=Σj=02nωk-1(ik-1,j)χk|k-1(ik-1,j)

Pk|k-1(ik-1)=Σj=02nωk-1(ik-1,j)[χk-1(ik-1,j)-xk|k-1(ik-1)][χk-1(ik-1,j)-xk|k-1(ik-1)]T

其中,上标T表示向量或矩阵的转置。

然后,将代入目标观测方程中,得到k时刻的一步预测观测值 所述目标观测方程为zk=h(xkk),其中,νk为观测噪声,h(xkk)表 示目标的观测函数。对至进行加权求和,得到k时刻观测量的初 步一步预测值具体操作过程如下:

zk|k-1(ik-1,j)=h(χk|k-1(ik-1,j),vk),ηk|k-1(ik-1)=Σj=02nωk-1(ik-1,j)zk|k-1(ik-1,j)

然后,在k时刻目标的观测集Zk中,选取出与具有最近欧式距离的 对应目标观测量zk,用zk对粒子状态的初步一步预测值和的初步一 步预测状态协方差矩阵做进一步修正,得到初步一步预测均值及 的初步一步预测状态协方差矩阵在均值为协方差为的 高斯分布上采样,得到k时刻第ik-1个粒子的粒子状态的一步预测值具 体操作过程如下:

x~k|k-1(ik-1)=xk|k-1(ik-1)+Kk(ik-1)(zk-ηk|k-1(ik-1)),Kk(ik-1)=Gk(ik-1)[Sk(ik-1)]-1

Pk(ik-1)=Pk|k-1(ik-1)-Gk(ik-1)[Sk(ik-1)]-1(Gk(ik-1))T

Sk(ik-1)=Σj=02nωk-1(ik-1,j)(zk|k-1(ik-1,j)-ηk|k-1(ik-1))(zk|k-1(ik-1,j)-ηk|k-1(ik-1))T

Gk(ik-1)=Σj=02nωk-1(ik-1,j)(χk-1(ik-1,j)-xk|k-1(ik-1))(zk|k-1(ik-1,j)-ηk|k-1(ik-1))T

xk|k-1(ik-1)~N(x~k|k-1(ik-1),Pk(ik-1))

其中,上标-1表示矩阵的逆,上标T表示矩阵或向量的转置,表 示均值为协方差为的高斯分布。

然后,根据均值为协方差为的高斯分布在处的概率密度 qk,按照下式计算出一步预测状态为的粒子的权值

ωk|k-1(ik-1)=ωk-1(ik-1)ps,kfk|k-1(xk|k-1(ik-1)|xk-1(ik-1))/qk,ik-1=1,2,···,Lk-1

其中,ps,k为k-1时刻的目标在k时刻仍然存在的概率,为对 应目标状态转移概率密度函数。

(3.3)当k≥1时,得出k时刻每个发射站对应的每个目标的双基地双基 地椭圆,在k时刻每个发射站对应的每个目标的双基地椭圆曲线上,每个点对 应的双基地距离和为一固定值,每个点对应的双基地距离和指:每个点与对 应发射站的距离、以及对应点与外辐射源雷达的距离之和。

根据k-1时刻每个目标相对于外辐射源雷达的角度、以及k时刻每个目标 相对于外辐射源雷达的角度,在k时刻每个发射站对应的对应目标的双基地椭 圆曲线与外辐射源雷达的矩形观测区域边缘的交点处,分布k时刻新生粒子, 并计算k时刻每个新生粒子的权值。k时刻新生粒子的个数表示为Jk,k时刻 第inew个新生粒子的粒子状态表示为k时刻第inew个新生粒子的权值表 示为inew取Lk-1+1至Lk-1+Jk

具体地说,在子步骤(3.3)中,首先,根据k-1时刻目标的观测集Zk-1, 获取k-1时刻每个目标相对于外辐射源雷达的角度,将k-1时每个目标相对于外 辐射源雷达的角度组成k-1时刻目标观测量角度集合θk-1,k-1时刻目标观测量 角度集合θk-1中元素的个数表示为<θk-1>,k-1时刻目标观测量角度集合θk-1中 第αk-1个元素表示为αk-1取1至<θk-1>。根据k时刻目标的观测集Zk, 获取k时刻每个目标相对于外辐射源雷达的角度,将k时刻每个目标相对于外 辐射源雷达的角度组成k时刻目标观测量角度集合θk,k时刻目标观测量角度 集合θk中元素的个数表示为<θk>,k时刻目标观测量角度集合θk中第αk个元 素表示为αk取1至<θk>。

判断k时刻是否产生新生粒子,如果k时刻不产生新生粒子,则k时刻新 生粒子的个数Jk为零。如果k时刻产生新生粒子,则得出k时刻产生的每个新 生粒子的分布区域中心点方位角,k时刻产生的每个新生粒子的分布区域中心 点方位角指:k时刻产生的对应新生粒子的位置的纵坐标与横坐标的比值的反 正切值。

判断k时刻是否产生新生粒子的过程为:将αk-1从1至<θk-1>进行遍历,αk-1每取一个数值时,在k时刻目标观测量角度集合θk中找出满足第一新生粒子存 在条件的元素,所述第一新生粒子存在条件为:k时刻目标观测量角度集合θk中的元素与的差值的绝对值不小于γθ,γθ为设定的角度门限;在将αk-1从1至<θk-1>进行遍历的过程中,若k时刻目标观测量角度集合θk中不存在满 足第一新生粒子存在条件的元素,则k时刻不产生新生粒子;反之,k时刻产 生新生粒子。

计算k时刻每个发射站对应的每个目标的双基地椭圆曲线与外辐射源雷 达的矩形观测区域边缘的交点的位置。然后,根据k时刻每个发射站对应的每 个目标的双基地椭圆曲线与外辐射源雷达的矩形观测区域边缘的交点的位置, 在外辐射源雷达的矩形观测区域边缘中的上边缘、下边缘、左边缘或右边缘 中,找出至少一个满足以下条件的边缘:与k时刻每个发射站对应的每个目标 的双基地椭圆曲线存在交点;将找出的满足以上条件的边缘记为:k时刻目标 的候选新生粒子分布边缘。

然后得出k时刻目标的新生粒子分布边缘。得出k时刻目标的新生粒子分 布边缘的过程为:当k时刻目标的候选新生粒子分布边缘的个数等于1时,k 时刻目标的新生粒子分布边缘为k时刻目标的候选新生粒子分布边缘。当k时 刻目标的候选新生粒子分布边缘的个数大于1时,在k时刻目标的所有候选新 生粒子分布边缘中找出k时刻目标的新生粒子分布边缘。

当k时刻目标的候选新生粒子分布边缘的个数等于1时,在得出k时刻目 标的新生粒子分布边缘之后,在k时刻目标的新生粒子分布边缘、以及k时刻 每个发射站对应的每个目标的双基地椭圆曲线的所有交点中,找出相距最近 的l个交点,这l个交点与l个发射站一一对应。得出所述相距最近的l个交点 的位置的均值,将所述相距最近的l个交点的位置的均值为:k时刻共同交点; 根据目标观测方程,针对k时刻共同交点,得出k时刻共同交点对应的目标速 度。

当k时刻目标的候选新生粒子分布边缘的个数大于1时,在k时刻目标的 每个候选新生粒子分布边缘、以及k时刻每个发射站对应的每个目标的双基地 椭圆曲线的所有交点中,找出相距最近的l个交点,这l个交点与l个发射站一 一对应。得出所述相距最近的l个交点的位置的均值,将所述相距最近的l个 交点的位置的均值为:对应候选新生粒子分布边缘的k时刻共同交点;根据目 标观测方程,针对每个候选新生粒子分布边缘的k时刻共同交点,得出对应候 选新生粒子分布边缘的k时刻共同交点对应的目标速度。例如,k时刻目标的 候选新生粒子分布边缘的个数表示为Wk,Wk为小于或等于4的自然数。在k 时刻目标的第wk个候选新生粒子分布边缘、以及k时刻每个发射站对应的每 个目标的双基地椭圆曲线的所有交点中,找出相距最近的l个交点,wk取1 至Wk。得出所述相距最近的l个交点的位置的均值,将所述相距最近的l个交 点的位置的均值为:第wk个候选新生粒子分布边缘的k时刻共同交点。

当k时刻目标的候选新生粒子分布边缘的个数大于1时,找出k时刻目标 的新生粒子分布边缘的过程为:判断k时刻是否产生新生粒子,如果k时刻不 产生新生粒子,则没有k时刻目标的新生粒子分布边缘。如果k时刻产生新生 粒子,针对k时刻的每个新生粒子,判断每个候选新生粒子分布边缘的k时刻 共同交点的方位角是否满足第二新生粒子存在条件,所述第二新生粒子存在 条件为:每个候选新生粒子分布边缘的k时刻共同交点的方位角(每个候选新 生粒子分布边缘的k时刻共同交点的位置的纵坐标与横坐标的比值的反正切 值)与k时刻目标观测量角度集合θk中满足第一新生粒子存在条件的元素的差 值的绝对值不大于设定的方位角门限值;如果每个候选新生粒子分布边缘的k 时刻共同交点的方位角不满足第二新生粒子存在条件,则令k时刻新生粒子的 个数Jk取0。如果至少一个候选新生粒子分布边缘满足第二新生粒子存在条件, 则在满足第二新生粒子存在条件的候选新生粒子分布边缘中,选择其k时刻共 同交点的方位角与k时刻的对应新生粒子的分布区域中心点方位角最接近的 候选新生粒子分布边缘,则k时刻目标的新生粒子分布边缘为:在满足第二新 生粒子存在条件的候选新生粒子分布边缘中选择的候选新生粒子分布边缘, 此时,k时刻共同交点为:k时刻目标的新生粒子分布边缘的k时刻共同交点, k时刻共同交点对应的目标速度为:k时刻目标的新生粒子分布边缘的k时刻 共同交点对应的目标速度。

然后,在k时刻共同交点处采样Jk个点,得出k时刻Jk个新生粒子。k时 刻每个新生粒子的粒子状态包括:k时刻对应新生粒子的方位状态、k时刻对 应新生粒子的速度状态。具体地,在k时刻共同交点的位置高斯分布中采样Jk个点,得出k时刻Jk个新生粒子的方位状态,k时刻共同交点的位置高斯分布 指:均值为k时刻共同交点的位置坐标(包括横坐标和纵坐标),标准差为双 基地距离单元(外辐射源雷达的距离分辨率)的高斯分布。在k时刻共同交点 的速度高斯分布采样Jk个点,得出k时刻Jk个新生粒子的速度状态,k时刻共 同交点的速度高斯分布指:均值为k时刻共同交点的速度(包括速度的横坐标 以及速度的纵坐标),标准差为外辐射源雷达的多普勒分辨率的高斯分布。

然后,按下式计算k时刻第inew个新生粒子的权值

ωk|k-1(inew)=Qx·qx(xk(inew)|zk)·Qy·qy(yk(inew)|zk)·Qx··qx·(x·k(inew)|zk)·Qy··qy·(y·k(inew)|zk)

其中,inew取Lk-1+1至Lk-1+Jk,为k时刻第inew个新生粒子的位置 在x方向上的概率密度函数,为k时刻第inew个新生粒子的位置在 y方向上的概率密度函数,为k时刻第inew个新生粒子的速度在x 方向上的概率密度函数,为k时刻第inew个新生粒子的速度在y方 向上的概率密度函数。当k时刻目标的新生粒子分布边缘为外辐射源雷达的矩 形观测区域边缘中的上边缘或下边缘时,当k时 刻目标的新生粒子分布边缘为外辐射源雷达的矩形观测区域边缘中的左边缘 或右边缘时,Qx=Qx·=2,Qy=Qy·=1.

(3.4)当k≥1时,根据子步骤(3.2)和子步骤(3.3),得到预测的k时 刻的粒子集其中,表示预测的k时刻的粒子集中第i'k个粒子的粒子状态,表示预测的k时刻的粒子集中第i'k个粒子的权值,i'k取1至Lk-1+Jk。利用k时刻目标的观测集Zk,对预测的k时刻的粒子集中第 i'k个粒子的权值进行修正,得出预测的k时刻的粒子集中第i'k个粒子的 修正后权值

当k≥1时,对预测的k时刻的粒子集中第i'k个粒子的权值进行修正 的过程包括以下子步骤:

a)设置迭代参数j,j=1,2,...当j=1时,执行子步骤b)。

b)根据以下公式计算第j个发射站在处的检测概率

pd,k(xk|k-1(ik),j)=Q[2SNR(xk|k-1(ik),j),2ln(1/pFA)]

SNR(xk|k-1(ik),j)=KjRT,j2RR2

Kj=PT,jGT,jGRλf,j2σrcsFT,j2FR2(4π)3k0T0(1/CPI)(NF)

其中,表示第j个发射站对应的处的信噪比(本发明实施例 中,为已知数值),pFA为设定的虚警概率,Q[·]表示marcum Q 函数,Q[a,b]=bxexp(-x2+a22)I0(ax)dx,I0表示零阶贝塞尔函数;RT,j为预 测的k时刻的粒子集中第i'k个粒子对应的目标到第j个发射站的距离,RR为 预测的k时刻的粒子集中第i'k个粒子对应的目标到外辐射源雷达(接收站) 的距离;PT,j为第j个发射站的信号发射功率,GT,j为第j个发射站的发射天 线增益,GR为外辐射源雷达的接收天线增益,λf,j为第j个发射站的发射信号 波长,σrcs为预测的k时刻的粒子集中第i'k个粒子对应的目标的双基雷达(指 第j个发射站与外辐射源雷达)反射面积,FT,j为第j个发射站的发射信道中 的信号传播衰减因子,FR为外辐射源雷达的接收信道中的信号传播衰减因子。 k0为玻尔兹曼常数,T0为环境温度,CPI为外辐射源雷达接收信号时的脉冲积 累时间,NF指环境噪声系数(此处包含干扰)。

c)按照下式计算k时刻杂波平均数λk和杂波概率密度ck

λk=L1l1·L2l2·pFA,ck=1/(L1L2)

其中,L1为外辐射源雷达接收信号的距离范围,L2为外辐射源雷达接收信号 的多普勒频率范围,l1外辐射源雷达接收信号的距离分辨率,l2为外辐射源雷 达接收信号的多普勒分辨率。

d)在k时刻目标的观测集Zk中,利用与第j个发射站对应的k时刻每个 目标的观测量,组成k时刻第j个发射站的目标的观测集然后利用中 的观测量计算似然函数

g(zk(j)|xk|k-1(ik))=12πRvexp(-12(zk(j)-h(xk|k-1(ik),vk))Rv-1(zk-h(xk|k-1(ik),vk))T)

其中,Rv是观测噪声νk的协方差矩阵,上标‐1表示矩阵的逆,上 标T表示矩阵或向量的转置。

e)按照下式计算出预测的k时刻的粒子集中第i'k个粒子的修正后权值

ωk(ik)=ωk|k-1(ik)·[1-pd,k(xk|k-1(ik))+Σzk(j)Zk(j)ψk,z(xk|k-1(ik))λkck+<ωk|k-1,ψk,z>]

ψk,z(xk|k-1(ik))=pd,k(xk|k-1(ik))·g(zk(j)|xk|k-1(ik))

<ωk|k-1,ψk,z>=Σik=1Lk-1+Jkψk,z(xk|k-1(ik))·ωk|k-1(ik)

其中,i'k取1至Lk-1+Jk

然后,令j值自增1。

f)判断j与发射站的个数为l的大小关系,如果j<l,则令然后返回至子步骤b),如果j=l,则迭代终止,预测的k时刻的粒子集中第 i'k个粒子的修正后权值为

(3.5)按照预测的k时刻的粒子集中第i'k个粒子的修正后权值对 预测的k时刻的粒子集进行重采样,得出k时刻的粒子集;对k时刻的粒子集 进行聚类分析,得到每个目标k时刻的状态。

具体地,在子步骤(3.5)中,首先根据预测的k时刻的粒子集中第i'k个 粒子的修正后权值计算k时刻的估计目标数

针对每个目标固定采样m个粒子,对预测的k时刻的粒子集进行重采样, 得到k时刻的粒子集ik取1至Lk,表示k时刻 的粒子集中第ik个粒子的粒子状态,表示k时刻的粒子集中第ik个粒子权 值,ωk(ik)=Nk/Lk.

在步骤3之后,针对每个目标0时刻至k时刻的状态,进行航迹关联,剔 除虚假航迹,得出每个目标的航迹。

本发明的效果可以通过以下仿真实验进行验证。

1)仿真条件

本发明的仿真实验中,设置三个发射站与一个接收站(外辐射源雷达), 本发明的仿真实验软件平台为MATLAB,操作系统为Win XP系统。每个发 射站的信号发射功率为1kW,每个发射站的信号发射频率为780MHz,每个 发射站的信号发射带宽为8MHz,每个发射站的发射天线增益为8dB,接收站 的接收天线增益为13dB,信号在传播中无衰减。设定的虚警概率pFA为10-4, 环境噪声系数(包含干扰)为30dB。目标在观测区域内匀速运动,且每一帧 最多有一个新生目标出现,每个目标对每个发射站的RCS(双基雷达反射面 积)都相同,均为10dBm2。目标存活概率为0.98,每个目标固定分配30个 粒子,每个新生目标分配20个粒子。

2)仿真结果分析

参照图4a,为仿真实验中k时刻每个发射站对应的对应目标的双基地椭 圆曲线与外辐射源雷达的矩形观测区域边缘的交点示意图。参照图4b,为图 4a的局部放大图。图4b中,新生粒子分布在A点和B点的附近区域内。如图 4a和图4b所示,A、B两点都为k时刻共同交点,若只利用目标的距离和信 息和多普勒信息,无法进一步选择。本发明方法利用当前时刻观测量中角度 信息判断新目标更可能出现的位置,从而剔除虚假公共交点。

参照图5a,为增加目标角度信息之前和之后的每时刻估计目标数和真实 目标数的对比图,图5a中,横轴表示时间,纵轴表示目标数,无标示直线代 表真实目标数,以三角形标示的直线代表增加目标角度信息之前的估计目标 数,以菱形标示的增加目标角度信息之后(对应本发明)的估计目标数。

对增加目标角度信息之前和本发明(增加目标角度信息之后)得出的方 位估计误差和速度估计误差进行对比,采用最优子模式分配(Optimal  Sub‐pattern Assignment,OSPA)作为多目标跟踪性能评估标准,其中距离敏感性 参数选择2,关联敏感性参数选择5。参照图5b,为增加目标角度信息之前和 之后得出的目标方位估计误差的对比图。图5b中,横轴表示时间,纵轴表示 目标方位估计误差,单位为m,图5b中,以三角形标示的直线代表增加目标 角度信息之前的目标方位估计误差,以菱形标示的直线代表增加目标角度信 息之后(对应本发明)的目标方位估计误差。参照图5c,为增加目标角度信 息之前和之后得出的目标速度估计误差的对比图。图5c中,横轴表示时间, 纵轴表示目标速度估计误差,单位为m,图5c中,以三角形标示的直线代表 增加目标角度信息之前的目标速度估计误差,以菱形标示的直线代表增加目 标角度信息之后(对应本发明)的目标速度估计误差。从图5a至图5c可以 看出,采用本发明之后,目标数目的估计精度虽然有所下降,但与增加目标 角度信息之前相比,性能相差不大;但是增加目标角度信息后,目标的方位 和速度估计精度得到了显著提高。此外,由于并非每个时刻都产生新生粒子, 且在分布新生粒子时有效剔除了虚假点,增加目标角度信息后,降低了运算 耗时。

显然,本领域的技术人员可以对本发明进行各种改动和变型而不脱离本 发明的精神和范围。这样,倘若本发明的这些修改和变型属于本发明权利要 求及其等同技术的范围之内,则本发明也意图包含这些改动和变型在内。

去获取专利,查看全文>

相似文献

  • 专利
  • 中文文献
  • 外文文献
获取专利

客服邮箱:kefu@zhangqiaokeyan.com

京公网安备:11010802029741号 ICP备案号:京ICP备15016152号-6 六维联合信息科技 (北京) 有限公司©版权所有
  • 客服微信

  • 服务号