CN102621543A - 基于粒子滤波和航迹管理的微弱目标检测前跟踪方法 - Google Patents

Abstract

本发明公开了一种基于粒子滤波和航迹管理的微弱目标检测前跟踪方法，属于雷达数据处理领域。基于粒子滤波的微弱目标检测前跟踪方法不能提供目标航迹信息，且在目标信噪比较低时，容易出现漏检的现象。本发明提出的基于粒子滤波和航迹管理的微弱目标检测前跟踪方法立足于解决此类问题。本发明具有结构简单，计算快速，易于硬件实现等优点，同时克服了基于粒子滤波的微弱目标检测前跟踪方法的局限性，通过有效解决漏检问题，进而保证目标航迹的连续性，因此具有较强的工程应用价值和推广前景。

Description

基于粒子滤波和航迹管理的微弱目标检测前跟踪方法

技术领域

本发明涉及一种雷达数据处理方法，特别是涉及一种雷达对微弱目标的检测与跟踪方法。

背景技术

随着目标隐身技术的日益成熟和目标飞行高度的降低，各种隐身飞机、掠海导弹以及巡航导弹等微弱目标大量出现，严重的影响着雷达目标检测跟踪性能的发挥，同时给雷达的生存带来严峻的威胁。因此，实现雷达对微弱目标的有效检测与跟踪，对提高其作战效能和生存能力具有重要意义。目前，针对微弱目标的检测与跟踪问题，已有很多的专家和学者进行了大量的研究，其中，基于粒子滤波的微弱目标检测前跟踪方法由于具有处理方法简单、易于硬件实现、可适应于非线性非高斯系统等优点，在微弱目标检测与跟踪领域受到广泛重视。该方法主要通过以下步骤实现：

(1)将雷达接收机输出的目标回波信号进行A/D变换，送雷达数据处理计算机执行以下步骤；

(2)利用目标状态转移方程对粒子集进行预测；

(3)利用新量测对各粒子权重进行更新；

(4)对权重更新后的粒子集进行重采样，并得到目标状态估计。

基于粒子滤波的微弱目标检测前跟踪方法具有以下两个缺陷：(1)不能提供目标的航迹信息；(2)在目标信噪比较低的情况时容易出现漏检。

发明内容

本发明的目的是提出一种基于粒子滤波和航迹管理的微弱目标检测前跟踪方法，解决基于粒子滤波的微弱目标检测前跟踪方法无法提供目标航迹信息和容易出现漏检的问题。

本发明提出的基于粒子滤波和航迹管理的微弱目标检测前跟踪方法的技术方案包括以下步骤：

步骤1：变量初始化

K是总仿真时间；

T是雷达扫描周期；

N是滤波器采用的粒子数；

SNR是目标信噪比；

I_min和I_max是目标可能的最小和最大强度；

v_min和v_max是目标可能的最小和最大速度；

γ是目标存在判决门限；

η₀是目标初始存在概率；

L×H表示雷达监测区域的大小；

∏_k是目标存在状态概率转移矩阵；

F_k是目标状态转移矩阵；

Q_k是过程噪声协方差；

步骤2：令k＝0，进行粒子初始化，得到初始粒子集其中，

包含目标的位置

速度

和强度信息，

是目标存在状态变量，i＝1，2，…，N

(1)令i＝1；

(2)目标位置初始化

$x_k^l=L\×\mathrm{rand}\left(1\right)---\left(1\right)$

$y_k^l=H\×\mathrm{rand}\left(1\right)---\left(2\right)$

rand(1)表示按照均匀分布产生一个[0，1]区间上的随机数；

(3)目标速度初始化

${\stackrel{\·}{x}}_k^l=v_{\min }+(v_{\max }-v_{\min })\×\mathrm{rand}\left(1\right)---\left(3\right)$

${\stackrel{\·}{y}}_k^l=v_{\min }+(v_{\max }-v_{\min })\×\mathrm{rand}\left(1\right)---\left(4\right)$

(4)目标强度初始化

${\stackrel{\·}{I}}_k^l=v_{\min }+(v_{\max }-v_{\min })\×\mathrm{rand}\left(1\right)---\left(5\right)$

(5)根据目标初始存在概率η₀初始化目标存在状态变量

(6)粒子权重初始化

$w_k^l=\frac{1}{N}---\left(6\right)$

(7)令i＝i+1，若i≤N，转(2)；

步骤3：令k＝k+1，获得k时刻的雷达量测

将雷达接收到的信号进行A/D变换，得到k时刻的雷达量测 $z_k\left\{z_k^{(p,q)}\right|p=\mathrm{1,2},...,m;q=\mathrm{1,2},...,n\},$ 送雷达数据处理计算机；

步骤4：对初始粒子集进行预测，得到预测粒子集

其中，

(1)令i＝1；

(2)根据前一时刻的目标存在状态变量

和目标存在状态概率转移矩阵∏_k预测当前时刻目标存在状态变量

(3)若

且

利用方程(1)-方程(5)生成粒子状态

(4)若

且利用状态转移方程对粒子i的状态进行一步预测

$x_{k|k-1}^l=F_k{x}_{k-1}^l+v_k---\left(7\right)$

其中，v_k为过程噪声，其噪声协方差为Q_k；

(5)令i＝i+1，若i ≤N，转(2)；

步骤5：利用当前时刻的量测对粒子权重进行更新

(1)令i＝1；

(2)确定似然区域的范围

$C_p\left(x_{k|k-1}^l\right)=\{p_0-s,...,p_0-1,p_0,p_0+1,...,p_0+s\}---\left(8\right)$

$C_q\left(y_{k|k-1}^l\right)=\{q_0-s,...,q_0-1,q_0,q_0+1,...,q_0+s\}---\left(9\right)$

其中，p₀表示最接近状态矢量元素

的整数，q₀表示最接近状态矢量元素

的整数，s为预设参数；

(3)粒子权重更新

$w_{k|k-1}^l=\left\{\begin{array}{cc}\underset{p\∈C_p\left(x_{k|k-1}^l\right)}{\mathrm{\Π}}\underset{q\∈C_q\left(x_{k|k-1}^l\right)}{\mathrm{\Π}}l\left(z_k^{(p,q)}\right|x_{k|k-1}^l)& E_{k|k-1}^l=1\\ 1& E_{k|k-1}^l=0\end{array}\right.---\left(10\right)$

其中

$l\left(z_k^{(p,q)}\right|x_{k|k-1}^l)=\exp \{-\frac{h_k^{(p,q)}\left(x_{k|k-1}^l\right)(h_k^{(p,q)}\left(x_{k|k-1}^l\right)-{2z}_k^{(p,q)})}{{2\mathrm{\σ}}^2}\}---\left(11\right)$

$h_k^{(p,q)}\left(x_{k|k-1}^l\right)=\frac{{\mathrm{\Δ}}_x{\mathrm{\Δ}}_y{I}_{k|k-1}^l}{2\mathrm{\π}{\mathrm{\Σ}}^2}\exp \{-\frac{{({\mathrm{p\Δ}}_x-x_{k|k-1}^l)}^2+{({\mathrm{q\Δ}}_y-y_{k|k-1}^l)}^2}{2{\mathrm{\Σ}}^2}\}---\left(12\right)$

其中，σ²为量测噪声协方差，Δ_x和Δ_y表示传感器分辨率，∑表示传感器模糊斑点数；

(4)令i＝i+1，若i≤N，转(2)；

(5)权重归一化

$w_k^l=\frac{w_{k|k-1}^l}{\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{w}_{k|k-1}^j},$ i＝1，2，…，N    (13)

步骤6：重采样

对权重更新后的粒子集

进行系统重采样

${\{y_k^l,(1/N)\}}_{l=1}^N=\mathrm{resample}{\left\{y_{k|k-1}^l{w}_k^l\right\}}_{l=1}^N---\left(14\right)$

得到k时刻的粒子集

其中

步骤7：目标状态和状态协方差估计

(1)计算目标存在概率

${\hat{P}}_k=\frac{1}{N}\underset{l=1}{\overset{N}{\mathrm{\Σ}}}{E}_k^l---\left(15\right)$

若

转步骤9；

(2)目标状态估计

${\hat{x}}_k=\frac{\underset{l=1}{\overset{N}{\mathrm{\Σ}}}{x}_k^l{E}_k^l}{\underset{l=1}{\overset{N}{\mathrm{\Σ}}}{E}_k^l}---\left(16\right)$

(3)状态协方差估计

${\hat{B}}_k=\frac{1}{N}\underset{l=1}{\overset{N}{\mathrm{\Σ}}}({\hat{x}}_k-x_k^l){({\hat{x}}_k-x_k^l)}^T---\left(17\right)$

步骤8：航迹更新

(1)将目标状态估计

等效为新的量测值，状态协方差估计

等效为量测误差；

(2)若当前时刻存在目标航迹τ_k-1或τ′_k-1，利用新量测对其进行更新，得到更新的航迹τ_k，转(5)；

(3)若当前时刻不存在目标航迹，将新量测作为航迹头，利用3/4逻辑法进行航迹起始；

(4)若航迹起始成功，记为τ_k，设置航迹维持标记flag＝0；

(5)转步骤10；

步骤9：航迹维持和终结

(1)若当前时刻不存在目标航迹，转步骤10；

(2)若当前时刻存在目标航迹τ_k-1或τ′_k-1，令flag＝flag+1，若flag＞2，转(4)；

(3)对τ_k-1或τ′_k-1进行一步预测，得到τ′_k，转步骤10；

(4)输出航迹τ_k-3并删除航迹τ′_k-1；

步骤10：重复步骤3～步骤9，直至雷达关机。

和背景技术相比，本发明的有益效果说明：

(1)本发明提出的基于粒子滤波和航迹管理的微弱目标检测前跟踪方法，可以在对目标状态估计的同时给出目标航迹，解决了背景技术中基于粒子滤波的检测前跟踪方法不能提供目标航迹信息的问题；(2)本发明提出的方法通过对航迹进行管理，可以有效的解决低信噪比条件下基于粒子滤波的微弱目标检测前跟踪方法容易出现漏检的问题，进而保证目标航迹的连续性。

附图说明

附图1是本发明的基于粒子滤波和航迹管理的微弱目标检测前跟踪方法的整体流程图，附图中各符号的含义与发明内容部分相应符号的含义相同；

附图2是本发明实施例中基于粒子滤波的微弱目标检测前跟踪方法估计的目标存在概率；

附图3是本发明实施例中目标真实航迹与基于粒子滤波的微弱目标检测前跟踪方法估计的目标状态，附图中X位置和Y位置分别表示目标的x坐标和y坐标；

附图4是本发明实施例中目标真实航迹与基于粒子滤波和航迹管理的微弱目标检测前跟踪方法估计的目标航迹，附图中X位置和Y位置分别表示目标的x坐标和y坐标；

附图5是本发明实施例中基于粒子滤波和航迹管理的微弱目标检测前跟踪方法的目标位置均方根误差；

附图6是本发明实施例中基于粒子滤波和航迹管理的微弱目标检测前跟踪方法的目标速度均方根误差；

具体实施方式

下面结合附图对本发明的基于粒子滤波和航迹管理的微弱目标检测前跟踪方法进行详细描述。

不失一般性，设置一个目标匀速直线运动的二维仿真场景。假设目标初始存在概率η₀＝0.05，目标初始强度为I＝20，其可能的最小强度和最大强度分别为I_min＝10和I_max＝30，目标初始位置和初始速度分别为(4.2，7.2)和(0.45，0.25)，可能的最小速度和最大速度分别为v_min＝-1和v_max＝1，目标状态转移概率P_b＝P_d＝0.05，目标存在概率判决门限γ＝0.6，目标由k＝7时出现，直到k＝22消失，过程噪声强度分别为q₁＝0.001和q₂＝0.01；雷达扫描周期T＝1，分辨率Δx＝Δy＝1，分辨单元的格式m＝n＝20，量测噪声标准差σ＝3.25，传感器模糊斑点数∑＝0.7；滤波器采用的粒子数为N＝40000，总仿真时间为K＝30个扫描周期，目标似然区域选择参数s＝2。其步骤如附图1所示。

(1)根据以上仿真条件进行变量初始化，其余变量初始化如下：

目标信噪比

$\mathrm{SNR}=10\lg {\left[\frac{I{\mathrm{\Δ}}_x{\mathrm{\Δ}}_y/2\mathrm{\π}{\mathrm{\Σ}}^2}{\mathrm{\σ}}\right]}^2=6\mathrm{dB}$

雷达监测区域

L＝Δx×m＝20，H＝Δy×n＝20

目标运动状态转移矩阵F_k和过程噪声协方差Q_k分别为

$F_k=\left[\begin{array}{ccccc}1& T& 0& 0& 0\\ 0& 1& 0& 0& 0\\ 0& 0& 1& T& 0\\ 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 1\end{array}\right],$ $Q_k=\left[\begin{array}{ccc}{Q}_0& 0& 0\\ 0& Q_0& 0\\ 0& 0& 0.01T\end{array}\right],$ $Q_0=0.001\left[\begin{array}{cc}\frac{1}{3}{T}^3& \frac{1}{2}{T}^2\\ \frac{1}{2}{T}^2& T^3\end{array}\right]$

目标存在状态转移矩阵

${\mathrm{\Π}}_k=\left[\begin{array}{cc}1-P_b& P_b\\ P_d& 1-P_d\end{array}\right]=\left[\begin{array}{cc}0.95& 0.05\\ 0.05& 0.95\end{array}\right]$

(2)按发明内容部分步骤2所述的方法进行粒子集初始化，得到初始粒子集；

(3)按发明内容部分步骤3所述的方法获取当前时刻的量测；

(4)按发明内容部分步骤4所述的方法对粒子集进行一步预测，得到预测的粒子集；

(5)按发明内容部分步骤5所述的方法对预测的粒子集的权重进行更新；

(6)按发明内容部分步骤6所述的方法对权重更新后的粒子集进行重采样，得到当前时刻的粒子集；

(7)按发明内容部分步骤7所述的方法得到目标的状态和状态协方差估计；

(8)按发明内容部分步骤8和步骤9所述的方法对目标航迹进行管理；

(9)循环执行发明内容部分步骤3～步骤8，直至雷达关机。

本实施例中目标存在概率判决门限γ＝0.6，由附图2可以看出，在k＝17，21，23等处，目标存在概率低于0.6，出现漏检的问题，此外，由附图3可以看出，基于粒子滤波的微弱目标检测前跟踪方法只能给出各时刻的目标状态估计，不能给出目标航迹；而由附图4可以看出，本发明的方法可以有效消除漏检的现象，并给出了目标的航迹信息。附图5和附图6分别给出了本实施例中基于粒子滤波和航迹管理的微弱目标检测前跟踪方法的目标位置和速度均方根误差，是30次Monte Carlo仿真的平均结果，表明算法具有较好的收敛性。

Claims (1)

1.基于粒子滤波和航迹管理的微弱目标检测前跟踪方法，其特征包括以下步骤：步骤1：变量初始化

K是总仿真时间；

T是雷达扫描周期；

N是滤波器采用的粒子数；

SNR是目标信噪比；

I_min和I_max是目标可能的最小和最大强度；

v_min和v_max是目标可能的最小和最大速度；

γ是目标存在判决门限；

η₀是目标初始存在概率；

L×H表示雷达监测区域的大小；

∏_k是目标存在状态概率转移矩阵；

F_k是目标状态转移矩阵；

Q_k是过程噪声协方差；

步骤2：令k＝0，进行粒子初始化，得到初始粒子集

其中，

包含目标的位置速度和强度

信息，是目标存在状态变量，i＝1，2，…，N

(1)令i＝1；

(2)目标位置初始化

$x_k^l=L\×\mathrm{rand}\left(1\right)---\left(1\right)$

$y_k^l=H\×\mathrm{rand}\left(1\right)---\left(2\right)$

rand(1)表示按照均匀分布产生一个[0，1]区间上的随机数；

(3)目标速度初始化

${\stackrel{\·}{x}}_k^l=v_{\min }+(v_{\max }-v_{\min })\×\mathrm{rand}\left(1\right)---\left(3\right)$

${\stackrel{\·}{y}}_k^l=v_{\min }+(v_{\max }-v_{\min })\×\mathrm{rand}\left(1\right)---\left(4\right)$

(4)目标强度初始化

${\stackrel{\·}{I}}_k^l=v_{\min }+(v_{\max }-v_{\min })\×\mathrm{rand}\left(1\right)---\left(5\right)$

(5)根据目标初始存在概率η₀初始化目标存在状态变量

；

(6)粒子权重初始化

$w_k^l=\frac{1}{N}---\left(6\right)$

(7)令i＝i+1，若i≤N，转(2)；

步骤3：令k＝k+1，获得k时刻的雷达量测

将雷达接收到的信号进行A/D变换，得到k时刻的雷达量测 $z_k\left\{z_k^{(p,q)}\right|p=\mathrm{1,2},...,m;q=\mathrm{1,2},...,n\},$ 送雷达数据处理计算机；

步骤4：对初始粒子集进行预测，得到预测粒子集

其中，

(1)令i＝1；

(2)根据前一时刻的目标存在状态变量

和目标存在状态概率转移矩阵∏k预测当前时刻目标存在状态变量

(3)若

且

利用方程(1)-方程(5)生成粒子状态

(4)若

且

利用状态转移方程对粒子i的状态进行一步预测

$x_{k|k-1}^l=F_k{x}_{k-1}^l+v_k---\left(7\right)$

其中，v_k为过程噪声，其噪声协方差为Q_k；

(5)令i＝i+1，若i≤N，转(2)；

步骤5：利用当前时刻的量测对粒子权重进行更新

(1)令i＝1；

(2)确定似然区域的范围

$C_p\left(x_{k|k-1}^l\right)=\{p_0-s,...,p_0-1,p_0,p_0+1,...,p_0+s\}---\left(8\right)$

$C_q\left(y_{k|k-1}^l\right)=\{q_0-s,...,q_0-1,q_0,q_0+1,...,q_0+s\}---\left(9\right)$

其中，p₀表示最接近状态矢量元素的整数，q₀表示最接近状态矢量元素的整数，s为预设参数；

(3)粒子权重更新

$w_{k|k-1}^l=\left\{\begin{array}{cc}\underset{p\∈C_p\left(x_{k|k-1}^l\right)}{\mathrm{\Π}}\underset{q\∈C_q\left(x_{k|k-1}^l\right)}{\mathrm{\Π}}l\left(z_k^{(p,q)}\right|x_{k|k-1}^l)& E_{k|k-1}^l=1\\ 1& E_{k|k-1}^l=0\end{array}\right.---\left(10\right)$

其中

$l\left(z_k^{(p,q)}\right|x_{k|k-1}^l)=\exp \{-\frac{h_k^{(p,q)}\left(x_{k|k-1}^l\right)(h_k^{(p,q)}\left(x_{k|k-1}^l\right)-{2z}_k^{(p,q)})}{{2\mathrm{\σ}}^2}\}---\left(11\right)$

$h_k^{(p,q)}\left(x_{k|k-1}^l\right)=\frac{{\mathrm{\Δ}}_x{\mathrm{\Δ}}_y{I}_{k|k-1}^l}{2\mathrm{\π}{\mathrm{\Σ}}^2}\exp \{-\frac{{({\mathrm{p\Δ}}_x-x_{k|k-1}^l)}^2+{({\mathrm{q\Δ}}_y-y_{k|k-1}^l)}^2}{2{\mathrm{\Σ}}^2}\}---\left(12\right)$

其中，σ²为量测噪声协方差，Δ_x和Δ_y表示传感器分辨率，∑表示传感器模糊斑点数；

(4)令i＝i+1，若i ≤N，转(2)；

(5)权重归一化

$w_k^l=\frac{w_{k|k-1}^l}{\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{w}_{k|k-1}^j},$ i＝1，2，…，N    (13)

步骤6：重采样

对权重更新后的粒子集

进行系统重采样

${\{y_k^l,(1/N)\}}_{l=1}^N=\mathrm{resample}{\{y_{k|k-1}^l,w_k^l\}}_{l=1}^N---\left(14\right)$

得到k时刻的粒子集

其中

步骤7：目标状态和状态协方差估计

(1)计算目标存在概率

${\hat{P}}_k=\frac{1}{N}\underset{l=1}{\overset{N}{\mathrm{\Σ}}}{E}_k^l---\left(15\right)$

若

转步骤9；

(2)目标状态估计

${\hat{x}}_k=\frac{\underset{l=1}{\overset{N}{\mathrm{\Σ}}}{x}_k^l{E}_k^l}{\underset{l=1}{\overset{N}{\mathrm{\Σ}}}{E}_k^l}---\left(16\right)$

(3)状态协方差估计

${\hat{B}}_k=\frac{1}{N}\underset{l=1}{\overset{N}{\mathrm{\Σ}}}({\hat{x}}_k-x_k^l){({\hat{x}}_k-x_k^l)}^T---\left(17\right)$

步骤8：航迹更新

(1)将目标状态估计

等效为新的量测值，状态协方差估计

等效为量测误差；

(2)若当前时刻存在目标航迹τ_k-1或τ′_k-1，利用新量测对其进行更新，得到更新的航迹τ_k，转(5)；

(3)若当前时刻不存在目标航迹，将新量测作为航迹头，利用3/4逻辑法进行航迹起始；

(4)若航迹起始成功，记为τ_k，设置航迹维持标记flag＝0；

(5)转步骤10；

步骤9：航迹维持和终结

(1)若当前时刻不存在目标航迹，转步骤10；

(2)若当前时刻存在目标航迹τ_k-1或τ′_k-1，令flag＝flag+1，若flag＞2，转(4)；

(3)对τ_k-1或τ′_k-1进行一步预测，得到τ′_k，转步骤10；

(4)输出航迹τ_k-3并删除航迹τ′_k-1；

步骤10：重复步骤3～步骤9，直至雷达关机。

Images