CN112379366B - 一种基于粒子群滤波的雷达目标跟踪算法 - Google Patents

Abstract

本发明提出了一种新的基于粒子群滤波的雷达目标跟踪算法，将标准粒子滤波算法融入到雷达信号处理与数据处理算法的整体流程中，并结合典型的雷达信号处理硬件架构进行算法优化，解决了标准粒子滤波算法计算量过大和收敛速度慢的问题，提高雷达的目标跟踪精度和稳定度。本发明突破传统卡尔曼理论框架下的滤波算法，基于粒子群滤波设计一种搜索最优值的目标跟踪算法，可以在强闪烁噪声或者目标快速机动的情况下有效跳出局部最优点，搜索到理想的粒子最优值，提高雷达机动目标的跟踪精度。

Description

一种基于粒子群滤波的雷达目标跟踪算法

技术领域

本发明属雷达接收技术领域，涉及一种基于粒子群滤波的雷达目标跟踪计算方法。

背景技术

目标跟踪雷达是指能够连续跟踪特定目标、不断对目标坐标进行精确测量并输出目标坐标位置(方位、仰角、斜距、径向速度等)信息的雷达，主要用于精密测量、火控、侦察校射、空间探测和导弹制导等领域。由于雷达系统内部或者外部噪声及干扰信号的影响，跟踪雷达对目标进行探测的过程中检测到真实目标的同时可能产生一些虚假目标。因此，需要对雷达信号处理产生的目标数据进行滤波，以排除虚假目标的干扰并锁定真实目标位置。常规的滤波算法大多是在卡尔曼理论框架下的方法，在闪烁噪声较大或者目标机动性较强时，容易丢失目标或者跟踪精度较差。

发明内容

为解决上述技术问题，本发明提出了一种基于粒子群滤波的雷达目标跟踪计算方法，所述方法包括以下步骤：

步骤1：通过引导信息中的历史点迹，预测目标跟踪初始位置向量，计算预测目标位置所在的距离门和多普勒门；

步骤2：以二维数据中行列为中心，截取距离门和多普勒门的矩形搜索区域；

步骤3：随机生成并初始化N个种群粒子；

步骤4：分别计算N个种群粒子当前所处位置的打分值；

步骤5：初始化全局最优位置，所有粒子的最高打分值作为全局最高打分值，初始化个体最优位置；

步骤6：按照粒子状态转移函数更新粒子群各粒子新的搜索位置；

步骤7：重新评估粒子群各粒子在新的搜索位置上的打分值，如果打分值超过历史最高分，则更新个体最优位置；

步骤8：更新全局最高打分值，并更新全局最优位置

步骤9：重复步骤6到步骤8，直到搜索次数i≥M，则搜索结束，M根据实际算法仿真的收敛速度来确定；

步骤10：全局最优位置为目标跟踪最优位置时，换算目标距离和径向速度。

进一步的，步骤1中利用P_i(R_i，φ_i，θ_i)，预测目标初始位置T₀(r₀，φ₀，θ₀，v₀)，的计算方法为：

v₀＝(R₃-R₁)/(2·t₀)

r₀＝R₃+v₀·τ

θ₀＝θ₃+((θ₃-θ₁)/(2·t₀))·τ

上式中t₀为警戒(搜索)雷达数据刷新周期，τ为引导信息传输系统延时，r₀表示初始目标斜距，

表示初始目标方位角，θ₀表示初始目标俯仰角，v₀表示初始目标径向速度，R_i表示目标斜距，/>

表示目标方位角，θ_i表示目标俯仰角。

进一步的，步骤1中由预测目标初始位置T₀(r₀，φ₀，θ₀，v₀)计算预测目标位置所在的距离门r_cell₀和多普勒门d_cell₀的计算公式为：

r_cell₀＝ceil(2×(r₀-r_bz)/c×t_s)

d_cell₀＝ceil(v₀/(c·λ·n/PRT))

上式中，c为光速，λ为载波波长，PRT为脉冲重复周期，n为脉冲积累个数，r_bz为脉冲雷达收发切换导致的探测盲区距离，t_s为信号处理机ADC的采样周期，ceil表示向上取整。

进一步的，步骤4中种群粒子所在位置的距离信杂比r_snr_xy和多普勒信杂比d_snr_xy的计算方法为：

r_snr_xy＝sum_A/sum_B

d_snr_xy＝sum_C/sum_D

即将信号区幅度数据求和除以噪声区幅度数据求和，S’{s}_m×n是目标回波数据的二维矩阵，x为粒子在二维矩阵的行号，y为粒子在二维矩阵的列号，w0为信号区半宽度，p0为保护区半宽度，q0为噪声区半宽度，s为正交检波采集点的IQ数据，m为单个脉冲采集的点数，n为积累的脉冲数。

进一步的，步骤4中评价函数f(r_snr_xy，d_snr_xy,amp_ab_xy)的定义为：

f(r_snr_xy,d_snr_xy,amp_ab_xy)＝c₁×r_snr_xy+c₂×d_snr_xy+c₃×amp_ab_xy

上式中c₁、c₂和c₃代表调节系数，分别取值为0.25、0.25和0.5，amp_abxy表示种群粒子所在位置的最大绝对信号幅度比。

进一步的，步骤6中的粒子状态转移函数p_tran(p_(i-1)j,pg_(i-1),sg_j)的定义为：

p_ij＝k₀×p_(i-1)j+k₁×rand()×(sg_j-p_(i-1)j)+k₂×rand()×(pg_(i-1)-p_(i-1)j)

上式中k₀为惯性权重，k₁和k₂为加速常数或者学习因子，k₀、k₁和k₂分别取值为0.7、1.5和1.5，第j个种群粒子初始化个体最优位置sg_j＝(x_0j,y_0j)，更新粒子群各粒子新的搜索位置p_ij(x_ij,y_ij)，x_ij,y_ij分别表示位置坐标，其中i为第几轮搜索，j为种群粒子编号，rand()为随机数生成函数。

进一步的，步骤10中由由pg_m(x_m,y_m)换算目标距离r₁和径向速度v_r的计算公式为：

r₁＝x_m×c×t_s/2+r_bz

v_r＝y_m×(c·λ·n/PRT)

上式中，c为光速，λ为载波波长，PRT为脉冲重复周期，n为脉冲积累个数，r_bz为脉冲雷达收发切换导致的探测盲区距离，t_s为信号处理机ADC的采样周期，x_m表示第m个粒子在二维回波矩阵中的行号，y_m表示第m个粒子在二维回波矩阵中的列号。

本发明突破传统卡尔曼理论框架下的滤波算法，基于粒子群滤波设计一种利用“群体智慧”来搜索最优值的目标跟踪算法，可以在强闪烁噪声或者目标快速机动的情况下有效跳出局部最优点，搜索到理想的粒子最优值，提高雷达机动目标的跟踪精度。

附图说明

图1算法流程框架；

图2粒子群迁移影响因素示意图；

图3信杂比能量统计范围示意图。

具体实施方式

本发明提出了一种新的基于粒子群滤波的雷达目标跟踪算法，将标准粒子滤波算法融入到雷达信号处理与数据处理算法的整体流程中，并结合典型的雷达信号处理硬件架构进行算法优化，解决了标准粒子滤波算法计算量过大和收敛速度慢的问题，提高雷达的目标跟踪精度和稳定度。

以下结合附图对本发明的具体实施方式作出详细说明。

正常情况下，跟踪雷达需要与警戒(搜索)雷达配合使用。警戒(搜索)雷达提供待跟踪目标的三个历史点迹作为跟踪雷达的目标引导信息，历史点迹表示为

i＝1,2,3；其中R_i表示目标斜距，/>

表示目标方位角，θ_i表示目标俯仰角。跟踪雷达根据该目标引导信息，控制伺服转台转向预测的当前目标方位，并向目标预测位置方向发射雷达波束，然后雷达采集回波数据，进行信号处理和数据处理，得出目标最新跟踪位置并根据位置偏移控制伺服转台锁定跟随目标，重复执行上述雷达波束发射、回波采集、信号处理和数据处理以及转台控制的闭环过程实现对目标的跟踪探测。目标回波数据表示为S{s}_m×n，回波数据为m行n列的二维矩阵，s为正交检波采集点的IQ数据，m为单个脉冲采集的点数，n为积累的脉冲数，即m为距离门个数、n为多普勒门个数。本发明为通过采集的目标回波数据得出待跟踪目标坐标位置信息的算法过程。算法总体流程框架如附图1所示，具体步骤为：

步骤1：通过引导信息中的

预测目标跟踪初始位置向量

其中r₀表示初始目标斜距，/>

表示初始目标方位角，θ₀表示初始目标俯仰角，v₀表示初始目标径向速度，计算预测目标位置所在的距离门r_cell₀和多普勒门d_cell₀；

步骤2：以二维数据S{s}_m×n中第d_cell₀行第r_cell₀列为中心，截取距离门±w个单元和多普勒门±h个单元的(2w+1)×(2h+1)矩形搜索区域

i＝1,2,3；w和h根据系统对目标输出刷新率要求以及信号处理硬件的总体资源量来评估选取适当值；

步骤3：随机生成并初始化N个种群粒子P_0j(x_0j,y_0j)，x_0j和y_0j分别为粒子当前位置的横纵坐标，x_0j∈(0,2w+1),y_0j∈(0,2h+1),j为种群粒子编号；为提高信号处理硬件资源的利用率，N可设置为信号处理运算单元的核心数量，如采用TI公司6678DSP芯片作为信号处理的核心处理器，N取为8(单片6778DSP有8个核)；

步骤4：分别计算N个种群粒子当前所处位置的打分值score_0j，0表示第0轮搜索，为种群粒子编号，评价函数为f(r_snr_xy，d_snr_xy,amp_ab_xy),r_snr_xy为种群粒子所在位置的距离信杂比，d_snr_xy为种群粒子所在位置的多普勒信杂比，amp_ab_xy种群粒子所在位置的最大绝对信号幅度比；

步骤5：初始化全局最优位置pg_j＝argmax{score_0j(x_0j,y_0j)}，score_m₀＝max{score_0j}，即将当前所有粒子所处位置打分值最高的粒子位置作为全局最优位置，所有粒子的最高打分值作为全局最高打分值，初始化个体最优位置sg_j＝(x_0j,y_0j)；

步骤6：按照粒子状态转移函数p_tran(p_(i-1)j,pg_(i-1),sg_j)更新粒子群各粒子新的搜索位置p_ij(x_ij,y_ij)，i为第几轮搜索，j为种群粒子编号，粒子群迁移影响因素示意图如附图2所示；

步骤7：重新评估粒子群各粒子在新的搜索位置上的打分值score_ij，i为第几轮搜索，j为种群粒子编号，如果打分值超过历史最高分，则更新个体最优位置sg_j＝(x_ij,y_ij)；

步骤8：如果max{score_ij(x_ij,y_ij)}>score_m_i-1，则更新全局最高打分值score_m_i＝max{score_ij(x_ij,yi_j)}，并更新全局最优位置pg_i＝argmax{score_ij(x_ij,y_ij)}，i为第几轮搜索，j为种群粒子编号，否则score_m_i＝score_m_i-1且pg_i＝pg_(i-1)；

步骤9：重复步骤6—步骤8，直到搜索次数i≥M，则搜索结束，M根据实际算法仿真的收敛速度来确定；

步骤10：如果score_m_m-score_m_m-1≤S_T，则全局最优位置pg_m为目标跟踪最优位置，根据pg_m(x_m,y_m)换算目标距离r₁和径向速度v_r；否则表示目标跟踪算法流程未收敛，目标丢失。

下面对上述步骤中部分参数的具体计算方法进行说明。

步骤1中利用

预测目标初始位置/>

的计算方法为：

v₀＝(R₃-R₁)/(2·t₀)

r₀＝R₃+v₀·τ

θ₀＝θ₃+((θ₃-θ₁)/(2·t₀))·τ

上式中t₀为警戒(搜索)雷达数据刷新周期，τ为引导信息传输系统延时。

步骤1中由预测目标初始位置

计算预测目标位置所在的距离门r_cell₀和多普勒门d_cell₀的计算公式为：

r_cell₀＝ceil(2×(r₀-r_bz)/c×t_s)

d_cell₀＝ceil(v₀/(c·λ·n/PRT))

上式中，c为光速，λ为载波波长，PRT为脉冲重复周期，n为脉冲积累个数，r_bz为脉冲雷达收发切换导致的探测盲区距离，t_s为信号处理机ADC的采样周期，ceil表示向上取整。

步骤4中种群粒子所在位置的距离信杂比r_snr_xy和多普勒信杂比d_snr_xy的计算方法为：

r_snr_xy＝sum_A/sum_B

d_snr_xy＝sum_C/sum_D

即如附图3所示，将信号区幅度数据求和除以噪声区幅度数据求和，噪声区与信号区中间间隔保护区，防止目标信号泄露的能量计入噪声能量。

步骤4中评价函数f(r_snr_xy，d_snr_xy,amp_ab_xy)的定义为：

f(r_snr_xy,d_snr_xy,amp_ab_xy)＝c₁×r_snr_xy+c₂×d_snr_xy+c₃×amp_ab_xy

上式中c₁、c₂和c₃分别取值为0.25、0.25和0.5；

步骤6中的粒子状态转移函数p_tran(p_(i-1)j,pg_(i-1),sg_j)的定义为：

p_ij＝k₀×p_(i-1)j+k₁×rand()×(sg_j-p_(i-1)j)+k₂×rand()×(pg_(i-1)-p_(i-1)j)

上式中k₀为惯性权重，k₁和k₂为加速常数或者学习因子，k₀、k₁和k₂分别取值为0.7、1.5和1.5。

步骤10中由pg_m(x_m,y_m)换算目标距离r₁和径向速度v_r的计算公式为：

r₁＝x_m×c×t_s/2+r_bz

v_r＝y_m×(c·λ·n/PRT)

上式中，c为光速，λ为载波波长，PRT为脉冲重复周期，n为脉冲积累个数，r_bz为脉冲雷达收发切换导致的探测盲区距离，t_s为信号处理机ADC的采样周期。

最后应说明的是，以上实施方式仅用以说明本发明实施例的技术方案而非限制，尽管参照以上较佳实施方式对本发明实施例进行了详细说明，本领域的普通技术人员应当理解，可以对本发明实施例的技术方案进行修改或等同替换都不应脱离本发明实施例的技术方案的精神和范围。

Claims (7)

1.一种基于粒子群滤波的雷达目标跟踪算法，其特征在于，所述算法包括以下步骤：

步骤1：通过引导信息中的历史点迹

预测目标跟踪初始位置向量，计算预测目标位置所在的距离门和多普勒门；

步骤2：以二维数据中行列为中心，截取距离门和多普勒门的矩形搜索区域；

步骤3：随机生成并初始化N个种群粒子；

步骤4：分别计算N个种群粒子当前所处位置的打分值；

步骤5：初始化全局最优位置，所有粒子的最高打分值作为全局最高打分值，初始化个体最优位置；

步骤6：按照粒子状态转移函数更新粒子群各粒子新的搜索位置；

步骤7：重新评估粒子群各粒子在新的搜索位置上的打分值，如果打分值超过历史最高分，则更新个体最优位置；

步骤8：更新全局最高打分值，并更新全局最优位置

步骤9：重复步骤6到步骤8，直到搜索次数i≥M，则搜索结束，M根据实际算法仿真的收敛速度来确定；

步骤10：全局最优位置为目标跟踪最优位置时，换算目标距离和径向速度。

2.如权利要求1所述的方法，其特征在于：步骤1中利用

预测目标初始位置/>

的计算方法为：

v₀＝(R₃-R₁)/(2·t₀)

r₀＝R₃+v₀·τ

θ₀＝θ₃+((θ₃-θ₁)/(2·t₀))·τ

上式中t₀为警戒(搜索)雷达数据刷新周期，τ为引导信息传输系统延时，r₀表示初始目标斜距，

表示初始目标方位角，θ₀表示初始目标俯仰角，v₀表示初始目标径向速度，R_i表示目标斜距，/>

表示目标方位角，θ_i表示目标俯仰角，i表示目标序号。

3.如权利要求1所述的方法，其特征在于：步骤1中由预测目标初始位置

计算预测目标位置所在的距离门r_cell₀和多普勒门d_cell₀的计算公式为：

r_cell₀＝ceil(2×(r₀-r_bz)/c×t_s)

d_cell₀＝ceil(v₀/(c·λ·n/PRT))

上式中，c为光速，λ为载波波长，PRT为脉冲重复周期，n为脉冲积累个数，r_bz为脉冲雷达收发切换导致的探测盲区距离，t_s为信号处理机ADC的采样周期，r₀表示初始目标斜距，v₀表示初始目标径向速度，ceil表示向上取整。

4.如权利要求1所述的方法，其特征在于：步骤4中种群粒子所在位置的距离信杂比r_snr_xy和多普勒信杂比d_snr_xy的计算方法为：

r_snr_xy＝sum_A/sum_B

d_snr_xy＝sum_C/sum_D

即将信号区幅度数据求和除以噪声区幅度数据求和，S’{s}_m×n是目标回波数据的二维矩阵，x为粒子在二维矩阵的行号，y为粒子在二维矩阵的列号，w0为信号区半宽度，p0为保护区半宽度，q0为噪声区半宽度，s为正交检波采集点的IQ数据，m为单个脉冲采集的点数，n为积累的脉冲数。

5.如权利要求1所述的方法，其特征在于：步骤4中评价函数f(r_snr_xy，d_snr_xy,amp_ab_xy)的定义为：

f(r_snr_xy,d_snr_xy,amp_ab_xy)＝c₁×r_snr_xy+c₂×d_snr_xy+c₃×amp_ab_xy

上式中c₁、c₂和c₃代表调节系数，分别取值为0.25、0.25和0.5，amp_abxy表示种群粒子所在位置的最大绝对信号幅度比。

6.如权利要求1所述的方法，其特征在于：步骤6中的粒子状态转移函数p_tran(p_(i-1)j,pg_(i-1),sg_j)的定义为：

p_ij＝k₀×p_(i-1)j+k₁×rand()×(sg_j-p_(i-1)j)+k₂×rand()×(pg_(i-1)-p_(i-1)j)

上式中k₀为惯性权重，k₁和k₂为加速常数或者学习因子，k₀、k₁和k₂分别取值为0.7、1.5和1.5，第j个种群粒子初始化个体最优位置sg_j＝(x_0j,y_0j)，更新粒子群各粒子新的搜索位置p_ij(x_ij,y_ij)，x_ij,y_ij分别表示位置坐标，其中i为第几轮搜索，j为种群粒子编号，rand()为随机数生成函数。

7.如权利要求1所述的方法，其特征在于：步骤10中由由pg_m(x_m,y_m)换算目标距离r₁和径向速度v_r的计算公式为：

r₁＝x_m×c×t_s/2+r_bz

v_r＝y_m×(c·λ·n/PRT)

上式中，c为光速，λ为载波波长，PRT为脉冲重复周期，n为脉冲积累个数，r_bz为脉冲雷达收发切换导致的探测盲区距离，t_s为信号处理机ADC的采样周期，x_m表示第m个粒子在二维回波矩阵中的行号，y_m表示第m个粒子在二维回波矩阵中的列号。

Images