CN102608590B - 一种基于动态规划和后像投影算法的相参积累方法 - Google Patents

Abstract

本发明公开了一种基于动态规划和后像投影算法的相参检测前跟踪积累方法，它是通过动态规划得到目标的预测航迹，然后通过后像投影算法进行能量积累来提高输出信噪比，实现对低可观测目标的检测。与非相参TBD算法相比，本发明在帧间积累时利用了回波信号的相位信息，对多帧回波进行相参积累来提高输出信噪比进而提高了检测效率；在时域上直接进行距离走动补偿，所以不需要对目标距离进行近似展开和考虑目标的多普勒走动；减少了积累误差，提高了检测概率和检测精度，降低了虚警概率。

Description

一种基于动态规划和后像投影算法的相参积累方法

技术领域：

本发明属于雷达系统中目标检测和跟踪的技术领域，它特别涉及到了低信噪比条件下雷达检测和跟踪低可观测目标(亦称弱目标)的技术领域。

背景技术：

随着目标隐身技术的发展、采用新型吸波材料和改变物体几何外形等隐身技术的不断发展和完善，雷达探测目标的雷达目标反射面积(RCS)降低了几个数量级。这对于雷达检测目标和提高自身的生存能力提出了严峻的挑战。因此对低可观测目标的检测和跟踪成为雷达的一个重要研究方向。从以前的研究中可知，可以通过增加脉冲积累个数的方法来提高信噪比从而提高雷达对低可观测目标的检测能力。现有预警雷达体系一般选用较低的脉冲重复频率以防止距离模糊的产生。在低重频扫描雷达体制下，雷达在一个方位向上发射的脉冲数很少，这导致了在一帧雷达回波数据中目标的回波脉冲数很少，这使得回波数据进行长时间的脉冲积累变得很难实现。

检测前跟踪算法是通过先存储多帧未经过门限处理的雷达回波原始数据，然后进行能量积累来对低可观测目标进行检测和跟踪的技术。由于在低信噪比条件下，单帧数据无法得出检测结果，检测前跟踪方法通过对多帧回波数据的处理，利用目标的运动特性，沿目标轨迹进行能量积累，提高了信噪比，达到了检测低可观测目标的目的。目前用于低可观测目标检测的非相参检测前跟踪算法主要有基于Hough变换的检测前跟踪算法、基于改进Hough变换的检测前跟踪算法、基于粒子滤波的检测前跟踪算法和基于动态规划的检测前跟踪算法等；用于低可观测目标检测的相参检测前跟踪算法主要有基于径向速度估计的相参积累算法(详见“王瑞军，张晓玲，樊玲，多帧相参积累的检测前跟踪方法，计算机工程与应用，2011，47(33).”)和基于Keystone变换的相参积累方法(详见“Wang Kun，Zhang Xiaoling，A TBD method using multi-frame coherent integration.Synthetic Aperture Radar(APSAR)，2011 3rd International Asia-Pacific Conference.26-30 Sept.2011，pages：1-4.”)。在雷达信号处理中脉冲的能量积累包括非相参积累和相参积累，前者指的是仅仅对数据的幅度(也可能是幅度的平方或者幅度对数)进行积累，而后者是指对复数据(即包含幅度和相位的数据)进行积累。由于非相参积累算法在帧间积累时都没有有效利用目标回波的相位信息，因此在低信噪比下积累效率比帧间相参积累时低。传统的非相参检测前跟踪算法为达到好的检测性能需要比相参检测前跟踪算法更多帧的回波数据进行积累，因此数据处理量较大。同时，现有的相参积累算法存在速度模糊问题和距离走动校正的问题，积累结果存在一定的误差。

发明内容：

为了能在低信噪比条件下，使用较少帧数的雷达回波数据就能得到较高的对微弱目标的检测概率，本发明提出了一种基于动态规划和后像投影算法的相参检测前跟踪积累方法，由于它克服了速度模糊问题，并且在时域上直接进行距离走动补偿，所以不需要对目标距离进行近似展开，从而不需要考虑目标的多普勒走动；与现有相参TBD方法相比，减少了积累误差，提高了检测概率和检测精度，减少了虚假航迹，降低了虚警概率。

为了方便描述本发明的内容，首先作以下术语定义：

定义1、检测前跟踪

在雷达系统中，检测前跟踪是指雷达在获得一个扫描周期的回波数据后，先不进行处理，不设检测门限和不宣布检测结果，而是将接收到的每一个扫描时刻的回波数据数字化后存储起来，等达到设定的数据量时在各扫描时刻之间对假设路径包含的点作几乎没有信息损失的相关处理，从而估计出目标的运动轨迹，最后检测结果和目标轨迹同时宣布。详见文献“黄红平.强杂波下微弱目标检测算法研究.电子科技大学硕士学位论文.2010”。

定义2、距离单元

在雷达系统中，将雷达测距的范围划分成若干小的区域并将其编号，雷达根据目标回波信号落入的区域编号计算目标与雷达之间的距离。

定义3、方位向

将雷达扫描空间均匀划分为若干等分，每一等分为一个方位向。

定义4、虚警门限

雷达系统中的参数，当统计值超过虚警门限时雷达报告发现目标，当统计值未超过虚警门限时雷达不报告发现目标。虚警门限值一般由虚警概率计算得到。详见文献“丁鹭飞，耿富录.雷达原理(第三版).西安电子科技大学出版社.2009.8”。

定义5、一帧回波数据

在本发明中，一帧回波数据是指在一个雷达扫描周期内，雷达接收机所接收、采样并存储的在这一个雷达扫描周期内所有发射脉冲的回波数据。

定义6、距离-慢时间二维数据矩阵

本发明中距离-慢时间二维数据矩阵的行代表距离向，其行的数目为雷达对每个回波采样的点数；矩阵的列代表方位向，其列的数目为雷达在每个方位向发射脉冲的序号。假设雷达扫描空间被分为N个方位向，每个方位向发射1个脉冲，雷达对每个发射脉冲的回波采样M次，则在一个雷达扫描周期内雷达连续发射N个脉冲并按方位向编号将采样数据存储为N行M列的二维矩阵S_NxM，如图1所示。

定义7、动态规划算法

动态规划是运筹学的一个分支，是求解决策过程最优化的数学方法。20世纪50年代初美国数学家R.E.Bellman等人在研究多阶段决策过程的优化问题时，提出了著名的最优化原理，把多阶段过程转化为一系列单阶段问题，利用各阶段之间的关系，逐个求解，创立了解决这类过程优化问题的新方法——动态规划。1957年出版了他的名著《DynamicProgramming》，这是该领域的第一本著作。

定义8、后像投影算法

BP算法起源于计算机断层扫描(Computer-Aided Tomography，CAT)技术，是一种精确的时域成像方法。详见文献“LARS M.H.ULANDER，HANS HELLSTEN.GUNNAR STENSTRO¨M.Synthetic-Aperture RadarProcessing Using Fast Factorized Back-Projection.IEEE TRANSACTIONS ON AEROSPACE AND ELECTRONIC SYSTEMSVOL.39.NO.3 JULY 2003”。

定义9、信噪比

信噪比是指信号的功率与环境噪声功率的比值。详见文献“丁鹭飞，耿富录.雷达原理(第三版).西安电子科技大学出版社.2009.8”。

定义10、脉冲压缩

脉冲压缩是一种现代雷达信号处理技术，简单来说就是雷达发射宽脉冲，然后再接收端“压缩”为窄脉冲，从而改善雷达的两种性能：作用距离和距离分辨率。详见“皮亦鸣，杨建宇，付毓生，杨晓波.合成孔径雷达成像原理.第一版.电子科技大学出版社.2007.3”

定义11、包络检波

包络检波是从调幅波包络中提取调制信号的过程：先对调幅波进行整流，得到波包络变化的脉动电流，再以低通滤波器滤除去高频分量，便得到调制信号。详见文献“丁鹭飞，耿富录.雷达原理(第三版).西安电子科技大学出版社.2009.8”。

本发明提供了用于检测和跟踪低可观测目标的一种基于动态规划和后像投影算法的相参积累方法，该方法的步骤如下：

步骤1、基于动态规划和后向投影算法的积累方法相关参数的初始化

初始化的参数均为已知，且初始化的参数如下：所有的位置坐标信息都是以极坐标形式给出；雷达发射线性调频信号，其发射脉冲的载频为F_c；雷达扫描周期为T；雷达发射脉冲的带宽B；雷达发射脉冲的调频斜率b；雷达发射脉冲的持续时间T_p；雷达的距离分辨率δ_r；相参处理的雷达回波数据帧数为K，(K为正整数)；K帧回波数据中第i帧回波数据为

(t＝1，2，3…K)；雷达扫描空间被划分成的方位向个数N，记方位向序号为n，(n＝1、2、3…N)；雷达在距离向上的采样频率F_s；雷达在距离向上的采样点数M，记距离向序号为m，m＝1、2、3…M；雷达第一虚警概率P_fa1；雷达第二虚警概率P_fa2；根据目标的运动特性可知目标运动速度的上下限分别为V_min和V_max，目标最大加速的限制为a_max。斜率允许误差ε，速度允许误差ε_v和加速度允许误差ε_a。

步骤2、计算目标信息

假设目标直线运动，速度为V，加速度为a。那么目标的速度应该满足V_min≤V≤V_max，目标的加速度应该满足-a_max≤a≤a_max，即雷达只对满足V_min≤V≤V_max与-a_max≤a≤a_max的目标进行检测和跟踪。

由速度限制条件V_min和V_max、加速度限制条件a_max与雷达扫描周期为T计算得到目标在第i帧与第i+1帧之间运动的理论最长距离和理论最短距离

$\mathrm{\Δ}{R}_{\min }=V_{\min }T-\frac{1}{2}{a}_{\max }{T}^2.$

假设背景噪声为复高斯白噪声，其均值为零，方差为σ²。由雷达第一虚警概率P_fa1、雷达第二虚警概率P_fa2以及噪声方差σ²计算得到虚警门限：雷达第一虚警门限值 ${\mathrm{\η}}_1=\sqrt{-2{\mathrm{\σ}}^2\ln \left(P_{\mathrm{fa}1}\right)};$ 雷达第二虚警门限值 ${\mathrm{\η}}_2=\sqrt{-2K{\mathrm{\σ}}^2\ln \left(P_{\mathrm{fa}2}\right)}.$

步骤3、雷达回波数据在距离向上进行脉冲压缩

取出步骤1中初始的K帧回波数据

对K帧回波数据

的N个方位向采用传统的脉冲压缩方法分别进行脉冲压缩，得到脉冲压缩后的K帧回波数据矩阵

然后将脉冲压缩后的K帧回波数据矩阵

存储到计算机中。

步骤4、将回波数据进行初步处理

取出步骤3中存储在计算机中的K帧回波数据矩阵并对K帧回波数据矩阵

中的每个数据取绝对值，得到K帧数据矩阵(i＝1、2、3…K)；然后将K帧数据矩阵

(i＝1、2、3…K)中每个数据与第一虚警门限η₁作比较，将小于第一虚警门限η₁的数据置零，最后得到经过初步处理后的K帧数据矩阵

$Z_{N\×M}^1,Z_{N\×M}^2,Z_{N\×M}^3\·\·\·Z_{N\×M}^K.$

步骤5、用动态规划进行数据关联

动态规划分为如下21个小步骤：

步骤5.1取出步骤4中经过初步处理后的K帧数据矩阵

则经过初步处理后的K帧数据矩阵

中任意一帧即第i帧数据矩阵为

(i＝1、2、3…K)。

定义：第i帧数据矩阵

中大于零的数据个数为P_i，(i＝1、2、3…K，P_i为非负整数)。

创建目标轨迹动态存储矩阵St，具体方法是：

定义：目标轨迹动态存储矩阵St为K+3行，tl列的矩阵，(tl为非负整数)。目标轨迹动态存储矩阵St的每1列表示1个目标轨迹，每1列第1行的数据表示目标的速度和加速度，每1列第2行的数据表示目标运动轨迹的斜率，每1列第3行至每1列第K+2行的数据表示目标从第1帧到第K帧的位置数据，每1列第K+3行的数据表示相参积累后的数值。目标轨迹动态存储矩阵St的格式如图3所示。

定义目标轨迹动态存储矩阵St的列tl＝1。将目标轨迹动态存储矩阵St置为K+3行tl列的零矩阵。

初始化i＝1，即从第1帧数据矩阵

开始进行动态规划。

统计第i帧数据矩阵

中大于零的数据个数P_i。

如果P_i＞0，那么定义第i帧数据矩阵

中第q个大于零的数据位置为(R_i，θ_i)_q，(q＝1、2、3…P_i)。令q＝1，即从第i帧数据矩阵

中第1个大于零的数据位置(R_i，θ_i)₁开始进行动态规划，转到步骤5.2；

如果P_i＝0，那么令目标轨迹动态存储矩阵St为空矩阵，转到步骤7。

步骤5.2取出第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q。

如果目标轨迹动态存储矩阵St是零矩阵，那么将第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St第tl列的第3行，存储格式为R_i+j*θ_i(j表示虚数单位)，转到步骤5.4。

如果目标轨迹动态存储矩阵St不是零矩阵，那么在目标轨迹动态存储矩阵St最后一列的后面增加一列，并将新增加这一列的所有行置为零。将第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St新增加这一列的第3行，存储格式为R_i+j*θ_i；更新tl＝tl+1，转到步骤5.4。

步骤5.3找出目标轨迹动态存储矩阵St第t列第i+1行表示第i-1帧的数据R_i-1+j*θ_i-1与目标轨迹动态存储矩阵St第t列第i+2行表示第i帧的数据R_i+j*θ_i，t表示表示目标轨迹动态存储矩阵St的列序号，(t＝1、2、3…L)。从第i-1帧的数据R_i-1+j*θ_i-1中提取出需要的信息R_i-1和θ_i-1，从第i帧的数据R_i+j*θ_i中提取出需要的信息R_i和θ_i，转到步骤5.4。

步骤5.4由步骤5.2或步骤5.3得到的R_i和θ_i计算出理论上目标在第i+1帧的最小角度

与理论上目标在第i+1帧的最大角度

那么理论上目标在第i+1帧的方位向范围为[θ_min，θ_max]，记该范围为

由于目标在第i+1帧与雷达的距离比目标在第i帧与雷达的距离近，所以得到理论上目标在第i+1帧距离向上的范围为[0，R_i)，记该范围为

转到步骤5.5。

步骤5.5取步骤5.1中得到的第i+1帧数据矩阵统计第i+1帧数据矩阵

中大于零的数据个数P_i+1。

如果P_i+1＞0，定义第i+1帧数据矩阵

中第d个大于零的数据位置为(R_i+1，θ_i+1)_d，(d＝1、2、3…P_i+1)。令d＝1，即取出第i+1帧数据矩阵

中第1个大于零的数据位置(R_i+1，θ_i+1)₁，转到步骤5.6；

如果P_i+1＝0，那么令目标轨迹动态存储矩阵St为空矩阵，转到步骤7。

步骤5.6取步骤5.1中得到的第i+1帧数据矩阵从第i+1帧数据

第d个大于零的数据位置(R_i+1，θ_i+1)_d中提取出R_i+1和θ_i+1。如果同时满足

与

两个条件，转到步骤5.7；如果不满足

与任意一个条件，转到步骤5.16。

步骤5.7由步骤5.2或步骤5.3得到的R_i和θ_i，定义第i帧数据点的位置为(R_i，θ_i)，由步骤5.6中得到的R_i+1和θ_i+1，定义第i+1帧数据点的位置为(R_i+1，θ_i+1)。计算出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的距离ΔR：

$\mathrm{\ΔR}=\sqrt{{\left(R_i\right)}^2+{\left(R_{i+1}\right)}^2-2*R_i*R_{i+1}*\cos ({\mathrm{\θ}}_{i+1}-{\mathrm{\θ}}_i)}.$

比较ΔR与步骤2中得到的ΔR_min、ΔR_max的大小。如果ΔR_min≤ΔR≤ΔR_max，转到步骤5.8；如果不满足ΔR_min≤ΔR≤ΔR_max，转到步骤5.16。

步骤5.8由步骤5.7中得到的第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)计算出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的斜率

如果i＝1，转到步骤5.9；如果i＞1，转到步骤5.10。

步骤5.9如果目标轨迹动态存储矩阵St第tl列的第i+3行为零，则直接将步骤5.8中算出的l_i存入目标轨迹动态存储矩阵St第tl列的第2行；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第tl列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St第tl列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零。将步骤5.8中算出的l_i存入目标轨迹动态存储矩阵St新增加这一列的第2行；将第i帧数据第q个非零数据的位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St新增加这一列的第i+2行，存储格式为R_i+j*θ_i；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增加这一列的第i+3行，存储格式为(R_i+1+j*θ_i+1，转到步骤5.16。

步骤5.10取出目标轨迹动态存储矩阵St中第t列第2行的数据，并记为l_i-1。比较l_i-1与步骤5.8中计算得到的斜率l_i大小。

如果满足-ε≤l_i-1-l_i≤ε(ε为斜率允许误差)，那么由步骤5.3中得到的R_i-1和θ_i-1，定义第i-1帧数据点的位置为(R_i-1，θ_i-1)；从步骤5.7中取出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)。通过第i-1帧数据点的位置(R_i-1，θ_i-1)、第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)来计算目标的运动信息： $\mathrm{\Δ}{R}_{i-1,i}=\sqrt{{(R_{i-1}\cos {\mathrm{\θ}}_{i-1}-R_i\cos {\mathrm{\θ}}_i)}^2+{(R_{i-1}\sin {\mathrm{\θ}}_{i-1}-R_i\sin {\mathrm{\θ}}_i)}^2},$ 其中ΔR_i-1，i表示目标第i-1帧数据点的位置(R_i-1，θ_i-1)与第i帧数据点的位置(R_i，θ_i)之间的距离间隔； $\mathrm{\Δ}{R}_{i,i+1}=\sqrt{{(R_i\cos {\mathrm{\θ}}_i-R_{i+1}\cos {\mathrm{\θ}}_{i+1})}^2+{(R_i\sin {\mathrm{\θ}}_i-R_{i+1}\sin {\mathrm{\θ}}_{i+1})}^2},$ ΔR_i，i+1表示目标第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的距离间隔。计算出ΔR_i-1，i与ΔR_i，i+1后，可以通过ΔR_i-1，i与ΔR_i，i+1计算出目标的速度：

和加速度 $a=\frac{\mathrm{\Δ}{R}_{i,i+1}-\mathrm{\Δ}{R}_{i-1,i}}{T^2},$ 转到步骤5.11；

如果不满足-ε≤l_i-1-l_i≤ε，转到步骤5.16。

步骤5.11由步骤5.10算出速度V和加速度a后，分别将它们与步骤1中已知的速度限制和加速度限制作比较。如果同时满足V_min≤V≤V_max与-a_max≤a≤a_max两个条件，转到步骤5.12；如果不满足V_min≤V≤V_max与-a_max≤a≤a_max任意一个条件，转到步骤5.16。

步骤5.12如果i＝2，转到步骤5.13；如果i＞2，转到步骤5.14。

步骤5.13如果目标轨迹动态存储矩阵St第t列的第i+3行为零，则直接将步骤5.10中算出的速度V和加速度a存入目标轨迹动态存储矩阵St第t列的第1行，存储格式为V+j*a；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第t列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St中第t列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零。将目标轨迹动态存储矩阵St第t列的第2行至目标轨迹动态存储矩阵St第t列的第i+2行的数据复制到目标轨迹动态存储矩阵St新增这一列的第2行至目标轨迹动态存储矩阵St新增这一列的第i+2行；将步骤5.10中算出的速度V和加速度a存入目标轨迹动态存储矩阵St新增这一列的第1行，存储格式为V+j*a；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增这一列的第i+3行，存储格式为R_i+1+j *θ_i+1，转到步骤5.16。

步骤5.14取出目标轨迹动态存储矩阵St第t列第1行的数据，提取出需要的速度信息并令其为V₀，提取出需要的加速度信息并令其为a₀。将步骤5.10中算出的速度V和V₀作比较；将步骤5.10中算出的加速度a和a₀作比较。

如果同时满足-ε_v≤V-V₀≤ε_v(ε_v为速度允许误差)和-ε_a≤a-a₀≤ε_a(ε_a为加速度允许误差)两个条件，转到步骤5.15；

如果不满足-ε_v≤V-V₀≤ε_v和-ε_a≤a-a₀≤ε_a中任意一个条件，转到步骤5.16。

步骤5.15如果目标轨迹动态存储矩阵St第t列的第i+3行为零，则直接将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第t列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St第t列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零。将目标轨迹动态存储矩阵St第t列的第1行至目标轨迹动态存储矩阵St第t列的第i+2行数据复制到目标轨迹动态存储矩阵St新增这一列的第1行至目标轨迹动态存储矩阵St新增这一列的第i+2行；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增这一列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16。

步骤5.16更新d＝d+1；如果d≤P_i+1，转到步骤5.6；如果d＞P_i+1，转到步骤5.17。

步骤5.17如果i＝1，转到步骤5.18；如果i＞1，转到步骤5.19。

步骤5.18更新q＝q+1。剔除目标轨迹动态存储矩阵St中第i+3行为零的列，统计目标轨迹动态存储矩阵St的列数，并且更新目标轨迹动态存储矩阵St的列数tl的值。如果q≤P_i，转到步骤5.2；如果q＞P_i，转到步骤5.20。

步骤5.19更新t＝t+1。如果t≤L，L表示目标轨迹动态存储矩阵St的列数，(L为非负整数)，转到步骤5.3；如果t＞L，转到步骤5.20。

步骤5.20更新i＝i+1。剔除目标轨迹动态存储矩阵St中第i+2行为零的列，然后统计目标轨迹动态存储矩阵St的列数，并令目标轨迹动态存储矩阵St的列数为L，(L为非负整数)。若i≤K-1时，转到步骤5.21；若i＞K-1，转到步骤6。

步骤5.21判断由步骤5.20得到的目标轨迹动态存储矩阵St列数L的大小。

如果L＞0，那么用t来表示目标轨迹动态存储矩阵St的列序号，(t＝1、2、3…L)。令t＝1，转到步骤5.3；

如果L＝0，转到步骤7。

步骤6、用后像投影算法进行相参积累

后像投影算法分为如下4个小步骤：

步骤6.1统计步骤5得到的目标轨迹动态存储矩阵St的列数，并定义目标轨迹动态存储矩阵St的列数为LL(LL为非负整数)。由步骤5.1对目标轨迹动态存储矩阵St的定义可知：目标轨迹动态存储矩阵St的每1列表示1个目标轨迹。于是得到目标轨迹的总个数为LL，转到步骤6.2。

步骤6.2目标轨迹动态存储矩阵St里面的第w个目标轨迹为目标轨迹动态存储矩阵St第w列第3行至目标轨迹动态存储矩阵St第w列第K+2行存储的数据，(w＝1、2、3…LL)。这K个数据为R_i+j*θ_i，(i＝1、2、3…K)，提取出R_i与θ_i，即为(R₁，θ₁)，(R₂，θ₂)，(R₃，θ₃)…，(R_K，θ_K)这K个数据点位置。令w＝1，转到步骤6.3。

步骤6.3从目标轨迹动态存储矩阵St第w列取出步骤6.2中的K个数据点位置(R₁，θ₁)，(R₂，θ₂)，(R₃，θ₃)…，(R_K，θ_K)，这样就得到了这K个数据离雷达的距离R₁、R₂、R₃…R_K。然后通过R₁、R₂、R₃…R_K利用后像投影算法算出相位补偿因子

i＝1、2、3…K。取出步骤3存储在计算机中的K帧回波数据矩阵

于是第i帧回波数据矩阵为

(i＝1、2、3…K)。在第i帧回波数据矩阵

的第(R_i，θ_i)位置处提取出该位置的数据并定义为τ_i。这样就得到了K个数据τ₁、τ₂、τ₃…τ_K。将这K个数据乘以相位补偿因子后直接相加得到相参积累值

将积累值μ_w取绝对值后得到最终积累值ρ_w，并将最终积累值ρ_w存储到该轨迹在目标轨迹动态存储矩阵St中第w列的第K+3行，转到步骤6.4。

步骤6.4更新w＝w+1。如果w≤LL，转到步骤6.3；如果w＞LL，转到步骤7。

步骤7、输出结果

步骤7.1判断目标轨迹动态存储矩阵St是否为空矩阵。如果目标轨迹动态存储矩阵St为空矩阵，转到步骤7.3；如果目标轨迹动态存储矩阵St为非空矩阵，转到步骤7.2。

步骤7.2将目标轨迹动态存储矩阵St中第K+3行的值与雷达第二虚警门限值η₂作比较，删除目标轨迹动态存储矩阵St中第K+3行的值小于η₂的列，转到步骤7.3。

步骤7.3判断目标轨迹动态存储矩阵St是否为空矩阵。

如果目标轨迹动态存储矩阵St为空矩阵，那么认为没有目标存在；

如果目标轨迹动态存储矩阵St为非空矩阵，那么就认为有目标存在。计算目标轨迹动态存储矩阵St的列数，并定义列数为LLL(LLL为正整数)。由步骤5.1对目标轨迹动态存储矩阵St的定义可知：目标轨迹动态存储矩阵St的每1列表示1个目标轨迹。于是得到目标轨迹的总个数为LLL，同时输出LLL个目标的轨迹、速度和加速度。

本发明的创新点在于：它是通过动态规划得到目标的预测航迹，然后通过后像投影算法进行能量积累来提高输出信噪比，实现对低可观测目标的检测。针对传统非相参TBD方法在帧间积累时没有利用目标回波的相位信息，针对现有相参TBD方法存在速度模糊和距离走动问题，本发明提供一种基于动态规划和后向投影算法的相参积累方法。该方法在帧间积累时利用了回波信号的相位信息，对多帧回波进行相参积累来提高输出信噪比进而提高了检测效率。

本发明的优点：本发明在多帧积累时利用了回波信号的相位信息，与非相参TBD算法相比，大大提高检测概率，减少了虚假航迹。由于该算法克服了速度模糊问题，，并且由于该算法在时域上直接进行距离走动补偿，所以不需要对目标距离进行近似展开，从而不需要考虑目标的多普勒走动；同时不仅仅局限于对径向运动目标的检测，与现有相参TBD算法相比，减少了积累误差，提高了检测概率和检测精度，降低了虚警概率。

附图说明

图1为一帧回波数据的存储格式

其中，纵坐标代表距离-慢时间二维矩阵的慢时间向即方位向，慢时间单元代表雷达发射脉冲的时刻。横坐标代表距离-慢时间二维矩阵的距离向，距离单元代表目标与雷达之间的距离。其中N表示雷达每个扫描周期将扫描空间划分的方位向个数；M为雷达距离单元个数。

图2为目标的运动模型

其中，R_i为目标在第i帧的位置与雷达之间的距离，θ_i为目标在第i帧的位置与雷达之间的距离，R_i+1目在第i+1帧的位置与雷达之间的距离，θ_i+1为目标在第i+1帧的位置与雷达之间的距离，ΔR为目标在两帧之间实际运动的距离，ΔR_min为目标在两帧之间理论运动的最大距离，ΔR_max为目标在两帧之间理论运动的最大距离，V为目标运动速度，β为运动方向。

图3为目标轨迹动态存储矩阵St

其中，每一列代表一个完整轨迹，如果有L条航迹，那么该矩阵就有L列。目标轨迹动态存储矩阵St每一列有K+3行：第1行数据表示目标的速度和加速度，第2行数据表示目标运动轨迹的斜率，第3行至第K+2行数据表示目标从第1帧到第K帧的位置数据，第K+3行存储相参积累后的数值。

图4为本发明的流程图。

具体实施方式

本发明主要采用计算机仿真的方法进行验证，所有步骤、结论都在MATLAB-R2010b上验证正确。具体实施步骤如下：

步骤1、基于动态规划和后像投影算法的多帧相参积累目标检测前跟踪方法相关参数的初始化

该方法需要进行初始化的参数如下：光速C＝300000Km/s；雷达发射脉冲的载频F_c＝100MHZ；雷达扫描周期为T＝2s；雷达发射脉冲的带宽B＝10MHZ；雷达发射脉冲的调频斜率b＝10THZ/s；雷达发射脉冲的持续时间T_p＝1μs；雷达的距离分辨率δ_r＝20m；相参处理的雷达数据帧数为K＝3；雷达扫描空间被划分成的方位向个数N＝80；雷达在距离向上的采样频率F_s＝20MHZ；雷达在距离向上的采样点数M＝200。目标速度的上下限分别为V_min＝900m/s、V_max＝1100m/s与最大加速度a_max＝10m/s²(假设目标的速度在3马赫左右，加速度在10m/s²以内)。第一门限虚警概率为P_fa1＝0.6，第二门限虚警概率为P_fa2＝0.0001。斜率允许误差ε＝0.001，速度允许误差ε_v＝10m/s和加速度允许误差ε_a＝0.5m/s²。仿真使用的背景噪声为复高斯白噪声，其均值为零，方差为σ²＝1。

步骤2、计算相关信息

计算出目标在第i帧与第i+1帧之间运动的理论最短距离ΔR_min＝880m，目标在第i帧与第i+1帧之间运动的理论最长距离ΔR_max＝1080m，雷达第一虚警门限值η₁＝1.0108与雷达第二虚警门限值η₂＝7.4338。

步骤3、雷达回波数据在距离向上进行脉冲压缩

取出步骤1中初始的3帧回波数据

对3帧回波数据

的80个方位向采用传统的脉冲压缩方法分别进行脉冲压缩，得到脉冲压缩后的3帧回波数据矩阵

然后将脉冲压缩后的3帧回波数据矩阵

存储到计算机中。

步骤4、将回波数据进行初步处理

取出步骤3中存储在计算机中的3帧回波数据矩阵

并对3帧回波数据矩阵

中的每个数据取绝对值，得到3帧数据矩阵(i＝1、2、3)；然后将3帧数据矩阵(i＝1、2，3)中每个数据与第一虚警门限η₁＝1.0108作比较，将小于第一虚警门限η₁＝1.0108的数据置零，最后得到经过初步处理后的3帧数据矩阵

步骤5、用动态规划进行数据关联

动态规划分为如下21个小步骤：

步骤5.1取出步骤4中经过初步处理后的3帧数据矩阵

则经过初步处理后的3帧数据矩阵中任意一帧即第i帧数据矩阵为

(i＝1、2、3)。

定义：第i帧数据矩阵

中大于零的数据个数为P_i，(i＝1、2、3，P_i为非负整数)。

创建目标轨迹动态存储矩阵St，具体方法是：

定义：目标轨迹动态存储矩阵St为6行，tl列的矩阵，(tl为非负整数)。目标轨迹动态存储矩阵St的每1列表示1个目标轨迹，每1列第1行的数据表示目标的速度和加速度，每1列第2行的数据表示目标运动轨迹的斜率，每1列第3行至每1列第5行的数据表示目标从第1帧到第3帧的位置数据，每1列第6行的数据表示相参积累后的数值。目标轨迹动态存储矩阵St的格式如图3所示。

定义目标轨迹动态存储矩阵St的列tl＝1。将目标轨迹动态存储矩阵St置为6行tl列的零矩阵。

初始化i＝1，即从第1帧数据矩阵

开始进行动态规划。

统计第i帧数据矩阵

中大于零的数据个数P_i。

如果P_i＞0，那么定义第i帧数据矩阵

中第q个大于零的数据位置为(R_i，θ_i)_q，(q＝1、2、3…P_i)。令q＝1，即从第i帧数据矩阵

中第1个大于零的数据位置(R_i，θ_i)₁开始进行动态规划，转到步骤5.2；

如果P_i＝0，那么令目标轨迹动态存储矩阵St为空矩阵，转到步骤7。

步骤5.2取出第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q。

如果目标轨迹动态存储矩阵St是零矩阵，那么将第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St第tl列的第3行，存储格式为R_i+j*θ_i(j表示虚数单位)，转到步骤5.4。

如果目标轨迹动态存储矩阵St不是零矩阵，那么在目标轨迹动态存储矩阵St最后一列的后面增加一列，并将新增加这一列的所有行置为零。将第i帧数据矩阵中第q个大于零的数据位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St新增加这一列的第3行，存储格式为R_i+j*θ_i；更新tl＝tl+1，转到步骤5.4。

步骤5.3找出目标轨迹动态存储矩阵St第t列第i+1行表示第i-1帧的数据R_i-1+j *θ_i-1与目标轨迹动态存储矩阵St第t列第i+2行表示第i帧的数据R_i+j*θ_i，t表示目标轨迹动态存储矩阵St的列序号(t＝1、2、3…L)。从第i-1帧的数据R_i-1+j*θ_i-1中提取出需要的信息R_i-1和θ_i-1，从第i帧的数据R_i+j*θ_i中提取出需要的信息R_i和θ_i，转到步骤5.4。

步骤5.4由步骤5.2或步骤5.3得到的R_i和θ_i计算出理论上目标在第i+1帧的最小角度

与理论上目标在第i+1帧的最大角度

那么理论上目标在第i+1帧的方位向范围为[θ_min，θ_max]，记该范围为

由于目标在第i+1帧离雷达的距离比目标在第i帧离雷达的距离近，所以得到理论上目标在第i+1帧距离向上的范围为[0，R_i)，记该范围为

转到步骤5.5。

步骤5.5取步骤5.1中得到的第i+1帧数据矩阵

统计第i+1帧数据矩阵

中大于零的数据个数P_i+1。

如果P_i+1＞0，定义第i+1帧数据矩阵中第d个大于零的数据位置为(R_i+1，θ_i+1)_d，(d＝1、2、3…P_i+1)。令d＝1，即取出第i+1帧数据矩阵

中第1个大于零的数据位置(R_i+1，θ_i+1)₁，转到步骤5.6；

如果P_i+1＝0，那么令目标轨迹动态存储矩阵St为空矩阵，转到步骤7。

步骤5.6取步骤5.1中得到的第i+1帧数据矩阵

从第i+1帧数据

第d个大于零的数据位置(R_i+1，θ_i+1)_d中提取出R_i+1和θ_i+1。如果同时满足

与

两个条件，转到步骤5.7；如果不满足

与

任意一个条件，转到步骤5.16。

步骤5.7由步骤5.2或步骤5.3得到的R_i和θ_i，定义第i帧数据点的位置为(R_i，θ_i)，由步骤5.6中得到的R_i+1和θ_i+1，定义第i+1帧数据点的位置为(R_i+1，θ_i+1)。计算出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的距离ΔR：

$\mathrm{\ΔR}=\sqrt{{\left(R_i\right)}^2+{\left(R_{i+1}\right)}^2-2*R_i*R_{i+1}*\cos ({\mathrm{\θ}}_{i+1}-{\mathrm{\θ}}_i)}.$

比较ΔR与步骤2中得到的ΔR_min、ΔR_max的大小。如果ΔR_min≤ΔR≤ΔR_max，转到步骤5.8；如果不满足ΔR_min≤ΔR≤ΔR_max，转到步骤5.16。

步骤5.8由步骤5.7中得到的第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)计算出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的斜率如果i＝1，转到步骤5.9；如果i＞1，转到步骤5.10。

步骤5.9如果目标轨迹动态存储矩阵St第tl列的第i+3行为零，则直接将步骤5.8中算出的l_i存入目标轨迹动态存储矩阵St第tl列的第2行；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第tl列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St第tl列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零。将步骤5.8中算出的l_i存入目标轨迹动态存储矩阵St新增加这一列的第2行；将第i帧数据

第q个非零数据的位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St新增加这一列的第i+2行，存储格式为R_i+j*θ_i；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增加这一列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16。

步骤5.10取出目标轨迹动态存储矩阵St中第t列第2行的数据，并记为l_i-1。比较l_i-1与步骤5.8中计算得到的斜率l_i大小。

如果满足-0.001≤l_i-1-l_i≤0.001，那么由步骤5.3中得到的R_i-1和θ_i-1，定义第i-1帧数据点的位置为(R_i-1，θ_i-1)；从步骤5.7中取出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)。通过第i-1帧数据点的位置(R_i-1，θ_i-1)、第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)来计算目标的运动信息：

$\mathrm{\Δ}{R}_{i-1,i}=\sqrt{{(R_{i-1}\cos {\mathrm{\θ}}_{i-1}-R_i\cos {\mathrm{\θ}}_i)}^2+{(R_{i-1}\sin {\mathrm{\θ}}_{i-1}-R_i\sin {\mathrm{\θ}}_i)}^2},$ 其中ΔR_i-1，i表示目标第i-1帧数据点的位置(R_i-1，θ_i-1)与第i帧数据点的位置(R_i，θ_i)之间的距离间隔； $\mathrm{\Δ}{R}_{i,i+1}=\sqrt{{(R_i\cos {\mathrm{\θ}}_i-R_{i+1}\cos {\mathrm{\θ}}_{i+1})}^2+{(R_i\sin {\mathrm{\θ}}_i-R_{i+1}\sin {\mathrm{\θ}}_{i+1})}^2},$ ΔR_i，i+1表示目标第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的距离间隔。计算出ΔR_i-1，i与ΔR_i，i+1后，可以通过ΔR_i-1，i与ΔR_i，i+1计算出目标的速度：

和加速度 $a=\frac{\mathrm{\Δ}{R}_{i,i+1}-\mathrm{\Δ}{R}_{i-1,i}}{4},$ 转到步骤5.11；

如果不满足-0.001≤L_i-1-l_i≤0.001，转到步骤5.16。

步骤5.11由步骤5.10算出速度V和加速度a后，分别将它们与步骤1中已知的速度限制和加速度限制作比较。如果同时满足900≤V≤1100与-10≤a≤10两个条件，转到步骤5.12；如果不满足900≤V≤1100与-10≤a≤10任意一个条件，转到步骤5.16。

步骤5.12如果i＝2，转到步骤5.13；如果i＞2，转到步骤5.14。

步骤5.13如果目标轨迹动态存储矩阵St第t列的第i+3行为零，则直接将步骤5.10中算出的速度V和加速度a存入目标轨迹动态存储矩阵St第t列的第1行，存储格式为V+j*a；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第t列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St中第t列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零。将目标轨迹动态存储矩阵St第t列的第2行至目标轨迹动态存储矩阵St第t列的第i+2行的数据复制到目标轨迹动态存储矩阵St新增这一列的第2行至目标轨迹动态存储矩阵St新增这一列的第i+2行；将步骤5.10中算出的速度V和加速度a存入目标轨迹动态存储矩阵St新增这一列的第1行，存储格式为V+j*a；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增这一列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16。

步骤5.14取出目标轨迹动态存储矩阵St第t列第1行的数据，提取出需要的速度信息并令其为V₀，提取出需要的加速度信息并令其为a₀。将步骤5.10中算出的速度V和V₀作比较；将步骤5.10中算出的加速度a和a₀作比较。

如果同时满足-10≤V-V₀≤10和-0.5≤a-a₀≤0.5两个条件，转到步骤5.15；

如果不满足-10≤V-V₀≤10和-0.5≤a-a₀≤0.5中任意一个条件，转到步骤5.16。

步骤5.15如果目标轨迹动态存储矩阵St第t列的第i+3行为零，则直接将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第t列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St第t列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零。将目标轨迹动态存储矩阵St第t列的第1行至目标轨迹动态存储矩阵St第t列的第i+2行数据复制到目标轨迹动态存储矩阵St新增这一列的第1行至目标轨迹动态存储矩阵St新增这一列的第i+2行；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增这一列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16。

步骤5.16更新d＝d+1；如果d≤P_i+1，转到步骤5.6；如果d＞P_i+1，转到步骤5.17。

步骤5.17如果i＝1，转到步骤5.18；如果i＞1，转到步骤5.19。

步骤5.18更新q＝q+1。剔除目标轨迹动态存储矩阵St中第i+3行为零的列，统计目标轨迹动态存储矩阵St的列数，并且更新目标轨迹动态存储矩阵St的列数tl的值。如果q≤P_i，转到步骤5.2；如果q＞P_i，转到步骤5.20。

步骤5.19更新t＝t+1。如果t≤L，L表示目标轨迹动态存储矩阵St的列数，(L为非负整数)，转到步骤5.3；如果t＞L，转到步骤5.20。

步骤5.20更新i＝i+1。剔除目标轨迹动态存储矩阵St中第i+2行为零的列，然后统计目标轨迹动态存储矩阵St的列数，并令目标轨迹动态存储矩阵St的列数为L，(L为非负整数)。若i≤2时，转到步骤5.21；若i＞2，转到步骤6。

步骤5.21判断由步骤5.20得到的目标轨迹动态存储矩阵St列数L的大小。

如果L＞0，那么用t来表示目标轨迹动态存储矩阵St的列序号，(t＝1、2、3…L)；令t＝1，转到步骤5.3；

如果L＝0，转到步骤7。

步骤6、用后像投影算法进行相参积累

后像投影算法分为如下4个小步骤：

步骤6.1统计步骤5得到的目标轨迹动态存储矩阵St的列数，并定义目标轨迹动态存储矩阵St的列数为LL(LL为非负整数)。由步骤5.1对目标轨迹动态存储矩阵St的定义可知：目标轨迹动态存储矩阵St的每1列表示1个目标轨迹。于是得到目标轨迹的总个数为LL，转到步骤6.2。

步骤6.2目标轨迹动态存储矩阵St里面的第w个目标轨迹为目标轨迹动态存储矩阵St第w列第3行至目标轨迹动态存储矩阵St第w列第5行存储的数据，(w＝1、2、3…LL)。这3个数据为R_i+j*θ_i，(i＝1、2、3)，提取出R_i与θ_i，即为(R₁，θ₁)，(R₂，θ₂)，(R₃，θ₃)这3个数据点位置。令w＝1，转到步骤6.3。

步骤6.3从目标轨迹动态存储矩阵St第w列取出步骤6.2中的3个数据点位置(R₁，θ₁)(R₂，θ₂)，(R₃，θ₃)，这样就得到了这3个数据离雷达的距离R₁、R₂、R₃。然后通过R₁、R₂、R₃利用后像投影算法算出相位补偿因了

(i＝1、2、3)。取出步骤3存储在计算机中的3帧回波数据矩阵于是第i帧回波数据矩阵为

(i＝1、2、3)。在第i帧回波数据矩阵的第(R_i，θ_i)位置处提取出该位置的数据并定义为τ_i。这样就得到了3个数据τ₁、τ₂、τ₃。将这3个数据乘以相位补偿因子后直接相加得到相参积累值

将积累值μ_w取绝对值后得到最终积累值ρ_w，并将最终积累值ρ_w存储到该轨迹在目标轨迹动态存储矩阵St中第w列的第6行，转到步骤6.4。

步骤6.4更新w＝w+1。如果w≤LL，转到步骤6.3；如果w＞LL，转到步骤7。

步骤7、输出结果

步骤7.1判断目标轨迹动态存储矩阵St是否为空矩阵。如果目标轨迹动态存储矩阵St为空矩阵，转到步骤7.3；如果目标轨迹动态存储矩阵St为非空矩阵，转到步骤7.2。

步骤7.2将目标轨迹动态存储矩阵St中第6行的值与雷达第二虚警门限值η₂＝7.4338作比较，删除目标轨迹动态存储矩阵St中第6行的值小于η₂＝7.4338的列，转到步骤7.3。

步骤7.3判断目标轨迹动态存储矩阵St是否为空矩阵。

如果目标轨迹动态存储矩阵St为空矩阵，那么认为没有目标存在；

如果目标轨迹动态存储矩阵St为非空矩阵，那么就认为有目标存在。计算目标轨迹动态存储矩阵St的列数，并定义列数为LLL(LLL为正整数)。由步骤5.1对目标轨迹动态存储矩阵St的定义可知：目标轨迹动态存储矩阵St的每1列表示1个目标轨迹。于是得到目标轨迹的总个数为LLL，同时输出LLL个目标的轨迹、速度和加速度。

由具体实施可以看出，本发明通过动态规划后运用后像投影算法对多帧回波数据进行相参积累处理，通过仿真可以得到在3帧回波数据的条件下，在检测概率为0.9时相参动态规划算法比非相参动态规划算法信噪比有大约1.8dB的改善，十分接近理论上的数值。与现有相参TBD算法相比，本发明克服了现有相参TBD算法存在速度模糊的问题，并且由于该算法在时域上直接进行距离走动补偿，所以不需要对目标距离进行近似展开，从而不需要考虑目标的多普勒走动，减少了积累误差，提高了检测概率和检测精度，降低了虚警概率。

Claims (1)

1.一种基于动态规划和后像投影算法的相参积累方法，其特征是它包括以下步骤：

步骤1、基于动态规划和后像投影算法的积累方法相关参数的初始化；

初始化的参数均为已知，且初始化的参数如下：所有的位置坐标信息都是以极坐标形式给出；雷达发射线性调频信号，其发射脉冲的载频为F_c;雷达扫描周期为T；雷达发射脉冲的带宽B;雷达发射脉冲的调频斜率b；雷达发射脉冲的持续时间T_p;雷达的距离分辨率δ_r;相参处理的雷达回波数据帧数为K，K为正整数；雷达扫描空间被划分成的方位向个数N，记方位向序号为n，n=1、2、3…N;雷达在距离向上的采样频率F_s;雷达在距离向上的采样点数M，记距离向序号为m，m＝1、２、３…M；K帧回波数据中第i帧回波数据为

i=1、2、3…K;雷达第一虚警概率P_fa1；雷达第二虚警概率P_fa2;根据目标的运动特性可知目标运动速度的上下限分别为V_min和V_max，目标最大加速的限制为a_max；斜率允许误差ε，速度允许误差ε_v和加速度允许误差ε_a；

步骤2、计算目标信息；

假设目标直线运动，速度为V，加速度为a；那么目标的速度应该满足V_min≤V≤V_max，目标的加速度应该满足-a_max≤a≤a_max，即雷达只对满足V_min≤V≤V_max与-a_max≤a≤a_max的目标进行检测和跟踪；

由速度限制条件V_min和V_max、加速度限制条件a_max与雷达扫描周期为T计算得到目标在第i帧与第i+1帧之间运动的理论最长距离

和理论最短距离 ${\mathrm{\ΔR}}_{\min }=V_{\min }T-\frac{1}{2}{a}_{\max }{T}^2;$

假设背景噪声为复高斯白噪声，其均值为零，方差为σ²雷达第一虚警概率P_fa1、雷达第二虚警概率P_fa2以及噪声方差σ²计算得到虚警门限：雷达第一虚警门限值 ${\mathrm{\η}}_1=\sqrt{-2{\mathrm{\σ}}^2\ln \left(P_{\mathrm{fa}1}\right)};$ 雷达第二虚警门限值 ${\mathrm{\η}}_2=\sqrt{-2K{\mathrm{\σ}}^2\ln \left(P_{\mathrm{fa}2}\right)};$

步骤3、雷达回波数据在距离向上进行脉冲压缩；

取出步骤1中初始的K帧回波数据

对K帧回波数据

的N个方位向采用传统的脉冲压缩方法分别进行脉冲压缩，得到脉冲压缩后的K帧回波数据矩阵

然后将脉冲压缩后的K帧回波数据矩阵

存储到计算机中；

步骤4、将回波数据进行初步处理；

取出步骤3中存储在计算机中的K帧回波数据矩阵

并对K帧回波数据矩阵

中的每个数据取绝对值，得到K帧数据矩阵

i=1、2、3…K;然后将K帧数据矩阵

中每个数据与第一虚警门限η₁作比较，将小于第一虚警门限

的数据置零，最后得到经过初步处理后的K帧数据矩阵

步骤5、用动态规划进行数据关联；

动态规划分为如下21个小步骤：

步骤5.1 取出步骤4中经过初步处理后的K帧数据矩阵

则经过初步处理后的K帧数据矩阵

中任意一帧即第i帧数据矩阵为

定义：第i帧数据矩阵

中大于零的数据个数为P_i，i=1、2、3…K，P_i为非负整数；

创建目标轨迹动态存储矩阵St，具体方法是：

定义：目标轨迹动态存储矩阵St为K+3行，tl列的矩阵，tl为非负整数；目标轨迹动态存储矩阵St的每1列表示1个目标轨迹，每1列第1行的数据表示目标的速度和加速度，每1列第2行的数据表示目标运动轨迹的斜率，每1列第3行至每1列第K+2行的数据表示目标从第1帧到第K帧的位置数据，每1列第K+3行的数据表示相参积累后的数值；

定义目标轨迹动态存储矩阵St的列tl=1;将目标轨迹动态存储矩阵St置为K+3行tl列的零矩阵；

初始化i=1，即从第1帧数据矩阵

开始进行动态规划；

统计第i帧数据矩阵

中大于零的数据个数P_i；

如果P_i＞0，那么定义第i帧数据矩阵

中第q个大于零的数据位置为(R_i，θ_i)_q，R_i和θ_i为目标在第i帧中极坐标形式的位置信息，q=1、2、3…P_i;令q=1，即从第i帧数据矩阵中第1个大于零的数据位置(R_i，θ_i)₁开始进行动态规划，转到步骤5.2；

如果P_i=0，那么令目标轨迹动态存储矩阵St为空矩阵，转到步骤7；

步骤5.2 取出第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q;

如果目标轨迹动态存储矩阵St是零矩阵，那么将第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St第tl列的第3行，存储格式为R_i+j*θ_i，j表示虚数单位，转到步骤5.4；

如果目标轨迹动态存储矩阵St不是零矩阵，那么在目标轨迹动态存储矩阵St最后一列的后面增加一列，并将新增加这一列的所有行置为零；将第i帧数据矩阵

中第q个大于零的数据位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St新增加这一列的第3行，存储格式为R_i+j*θ_i;更新tl=tl+1，转到步骤5.4；

步骤5.3 找出目标轨迹动态存储矩阵St第t列第i+1行表示第i-1帧的数据R_i-1+j*θ_i-1与目标轨迹动态存储矩阵St第t列第i+2行表示第i帧的数据R_i+j*θ_i，t表示表示目标轨迹动态存储矩阵St的列序号，t=1、2、3…L;从第i-1帧的数据R_i-1+j*θ_i-1中提取出需要的信息R_i-1和θ_i-1，从第i帧的数据R_i+j*θ_i中提取出需要的信息R_i和θ_i，转到步骤5.4；

步骤5.4 由步骤5.2或步骤5.3得到的R_i和θ_i计算出理论上目标在第i+1帧的最小角度

与理论上目标在第i+1帧的最大角度

那么理论上目标在第i+1帧的方位向范围为[θ_min，θ_max]，记该范围为

由于目标在第i+1帧与雷达的距离比目标在第i帧与雷达的距离近，所以得到理论上目标在第i+1帧距离向上的范围为[0，R_i)，记该范围为

转到步骤5.5；

步骤5.5 取步骤5.1中得到的第i+1帧数据矩阵

统计第i+1帧数据矩阵

中大于零的数据个数P_i+1；

如果P_i+1＞0，定义第i+1帧数据矩阵

中第d个大于零的数据位置为(R_i+1，θ_i+1)_d，d=1、2、3…P_i+1；令d=1，即取出第i+1帧数据矩阵

中第1个大于零的数据位置(R_i+1，θ_i+1)₁，转到步骤5.6；

如果P_i+1=0，那么令目标轨迹动态存储矩阵St为空矩阵，转到步骤7；

步骤5.6 取步骤5.1中得到的第i+1帧数据矩阵

从第i+1帧数据

第d个大于零的数据位置(R_i+1，θ_i+1)_d中提取出R_i+1和θ_i+1；如果同时满足

与

两个条件，转到步骤5.7；如果不满足

与

任意一个条件，转到步骤5.16；

步骤5.7 由步骤5.2或步骤5.3得到的R_i和θ_i，定义第i帧数据点的位置为(R_i，θ_i)，由步骤5.6中得到的R_i+1和θ_i+1，定义第i+1帧数据点的位置为(R_i+1，θ_i+1)；计算出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的距离ΔR：

$\mathrm{\ΔR}=\sqrt{{\left(R_i\right)}^2+{\left(R_{i+1}\right)}^2-2*R_i*R_{i+1}*\cos ({\mathrm{\θ}}_{i+1}-{\mathrm{\θ}}_i)};$

比较ΔR与步骤2中得到的ΔR_min、ΔR_max的大小；如果ΔR_min≤ΔR≤ΔR_max，转到步骤5.8；如果不满足ΔR_min≤ΔR≤ΔR_max，转到步骤5.16；

步骤5.8 由步骤5.7中得到的第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)计算出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的斜率如果i=1，转到步骤5.9；如果i＞1，转到步骤5.10；

步骤5.9 如果目标轨迹动态存储矩阵St第tl列的第i+3行为零，则直接将步骤5.8中算出的l_i存入目标轨迹动态存储矩阵St第tl列的第2行；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第tl列的第i+3行，存储格式为R_i+1+j×θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St第tl列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零；将步骤5.8中算出的l_i存入目标轨迹动态存储矩阵St新增加这一列的第2行；将第i帧数据

第q个非零数据的位置(R_i，θ_i)_q中的R_i和θ_i直接存储到目标轨迹动态存储矩阵St新增加这一列的第i+2行，存储格式为R_i+j*θ_i;将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增加这一列的第i+3行，存储格式为R_i+1+j*θ_i+1,转到步骤5.16；

步骤5.10 取出目标轨迹动态存储矩阵St中第t列第2行的数据，并记为l_i-1；比较l_i-1与步骤5.8中计算得到的斜率l_i大小；

如果满足-ε≤l_i-1-l_i≤ε，ε为斜率允许误差，那么由步骤5.3中得到的R_i-1和θ_i-1，定义第i-1帧数据点的位置为(R_i-1，θ_i-1)；从步骤5.7中取出第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)；通过第i-1帧数据点的位置(R_i-1，θ_i-1)、第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)来计算目标的运动信息：

${\mathrm{\ΔR}}_{i-1,i}=\sqrt{{(R_{i-1}\cos {\mathrm{\θ}}_{i-1}-R_i\cos {\mathrm{\θ}}_i)}^2+{(R_{i-1}\sin {\mathrm{\θ}}_{i-1}-R_i\sin {\mathrm{\θ}}_i)}^2},$ 其中ΔR_i-1，i表示目标第i-1帧数据点的位置

与第i帧数据点的位置(R_i，θ_i)之间的距离间隔； ${\mathrm{\ΔR}}_{i,i+1}=\sqrt{{(R_i\cos {\mathrm{\θ}}_i-R_{i+1}\cos {\mathrm{\θ}}_{i+1})}^2+{(R_i\sin {\mathrm{\θ}}_i-R_{i+1}\sin {\mathrm{\θ}}_{i+1})}^2},$ ΔR_i，i+1表示目标第i帧数据点的位置(R_i，θ_i)与第i+1帧数据点的位置(R_i+1，θ_i+1)之间的距离间隔；计算出ΔR_i-1，i与ΔR_i，i+1后，可以通过ΔR_i-1，i与ΔR_i，i+1计算出目标的速度：和加速度

转到步骤5.11；

如果不满足-ε≤l_i-1-l_i≤ε，转到步骤5.16；

步骤5.11 由步骤5.10算出速度V和加速度a后，分别将它们与步骤1中已知的速度限制和加速度限制作比较；如果同时满足V_min≤V≤V_max与-a_max≤a≤a_max两个条件，转到步骤5.12；如果不满足

与-a_max≤a≤a_max任意一个条件，转到步骤5.16；

步骤5.12 如果i=2，转到步骤5.13；如果i＞2，转到步骤5.14；

步骤5.13 如果目标轨迹动态存储矩阵St第t列的第i+3行为零，则直接将步骤5.10中算出的速度V和加速度a存入目标轨迹动态存储矩阵St第t列的第1行，存储格式为V+j*a;将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第t列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St中第t列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零；将目标轨迹动态存储矩阵St第t列的第2行至目标轨迹动态存储矩阵St第t列的第i+2行的数据复制到目标轨迹动态存储矩阵St新增这一列的第2行至目标轨迹动态存储矩阵St新增这一列的第i+2行；将步骤5.10中算出的速度V和加速度a存入目标轨迹动态存储矩阵St新增这一列的第1行，存储格式为V+j*a;将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增这一列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

步骤5.14 取出目标轨迹动态存储矩阵St第t列第1行的数据，提取出需要的速度信息并令其为V₀，提取出需要的加速度信息并令其为a₀；将步骤5.10中算出的速度V和V₀作比较；将步骤5.10中算出的加速度a和a₀作比较；

如果同时满足-ε_v≤V-V₀≤ε_v和-ε_a≤a-a₀≤ε_a两个条件，ε_v为速度允许误差，ε_a为加速度允许误差，转到步骤5.15；

如果不满足-ε_v≤V-V₀≤ε_v和-ε_a≤a-a₀≤ε_a中任意一个条件，转到步骤5.16；

步骤5.15 如果目标轨迹动态存储矩阵St第t列的第i+3行为零，则直接将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St第t列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

如果目标轨迹动态存储矩阵St第t列的第i+3行不为零，则在目标轨迹动态存储矩阵St最后一列的后面新增加一列，并将新增加这一列的所有行置为零；将目标轨迹动态存储矩阵St第t列的第1行至目标轨迹动态存储矩阵St第t列的第i+2行数据复制到目标轨迹动态存储矩阵St新增这一列的第1行至目标轨迹动态存储矩阵St新增这一列的第i+2行；将步骤5.6中得到的R_i+1和θ_i+1存入目标轨迹动态存储矩阵St新增这一列的第i+3行，存储格式为R_i+1+j*θ_i+1，转到步骤5.16；

步骤5.16 更新d=d+1；如果d≤P_i+1，转到步骤5.6；如果d＞P_i+1，转到步骤5.17；

步骤5.17 如果i=1，转到步骤5.18；如果i＞1，转到步骤5.19；

步骤5.18 更新q=q+1；剔除目标轨迹动态存储矩阵St中第i+3行为零的列，统计目标轨迹动态存储矩阵St的列数，并且更新目标轨迹动态存储矩阵St的列数tl的值；如果q≤P_i，转到步骤5.2；如果q＞P_i，转到步骤5.20；

步骤5.19 更新t=t+1；如果t≤L，L表示目标轨迹动态存储矩阵St的列数，L为非负整数，转到步骤5.3；如果t＞L，转到步骤5.20；

步骤5.20 更新i=i+1；剔除目标轨迹动态存储矩阵St中第i+2行为零的列，然后统计目标轨迹动态存储矩阵St的列数，并令目标轨迹动态存储矩阵St的列数为L，L为非负整数；若i≤K-1时，转到步骤5.21；若i＞K-1，转到步骤6；

步骤5.21 判断由步骤5.20得到的目标轨迹动态存储矩阵St列数L的大小；

如果L＞0，那么用t来表示目标轨迹动态存储矩阵St的列序号，t=1、２、３…Ｌ;令t=1，转到步骤5.3；

如果L=0，转到步骤7；

步骤6、用后像投影算法进行相参积累

后像投影算法分为如下4个小步骤：

步骤6.1 统计步骤5得到的目标轨迹动态存储矩阵St的列数，并定义目标轨迹动态存储矩阵St的列数为LL，LL为非负整数；由步骤5.1对目标轨迹动态存储矩阵St的定义可知：目标轨迹动态存储矩阵St的每1列表示1个目标轨迹；于是得到目标轨迹的总个数为LL，转到步骤6.2；

步骤6.2 目标轨迹动态存储矩阵St里面的第w个目标轨迹为目标轨迹动态存储矩阵St第w列第3行至目标轨迹动态存储矩阵St第w列第K+2行存储的数据，w=1、２、３…LＬ;这Ｋ个数据为R_i+j*θ_i，i=1、２、３…Ｋ；提取出R_i与θ_i，即为(R₁，θ₁)，(R_2，θ₂)，(R₃，θ₃)…，(R_K，θ_K)这K个数据点位置；令w=1，转到步骤6.3；

步骤6.3 从目标轨迹动态存储矩阵St第w列取出步骤6.2中的K个数据点位置(R₁，θ₁)，(R₂，θ₂)，(R₃，θ₃)…，(R_K，θ_K)，这样就得到了这K个数据离雷达的距离R₁、R₂、R₃…R_K;然后通过R₁、R₂、R₃…R_K利用后像投影算法算出相位补偿因子

i=1、2、3…K;取出步骤3存储在计算机中的K帧回波数据矩阵

于是第i帧回波数据矩阵为

i=1、2、3…K;在第i帧回波数据矩阵

的第(R_i，θ_i)位置处提取出该位置的数据并定义为τ_i；这样就得到了K个数据

将这K个数据乘以相位补偿因子后直接相加得到相参积累值

将积累值μ_w取绝对值后得到最终积累值ρ_w，并将最终积累值ρ_w存储到该轨迹在目标轨迹动态存储矩阵St中第w列的第K+3行，转到步骤6.4；

步骤6.4 更新w=w+1；如果w≤LL，转到步骤6.3；如果w＞LL，转到步骤7；

步骤7、输出结果；

步骤7.1 判断目标轨迹动态存储矩阵St是否为空矩阵；如果目标轨迹动态存储矩阵St为空矩阵，转到步骤7.3；如果目标轨迹动态存储矩阵St为非空矩阵，转到步骤7.2；

步骤7.2 将目标轨迹动态存储矩阵St中第K+3行的值与雷达第二虚警门限值η₂作比较，删除目标轨迹动态存储矩阵St中第K+3行的值小于η₂的列，转到步骤7.3；

步骤7.3 判断目标轨迹动态存储矩阵St是否为空矩阵；

如果目标轨迹动态存储矩阵St为空矩阵，那么认为没有目标存在；

如果目标轨迹动态存储矩阵St为非空矩阵，那么就认为有目标存在；计算目标轨迹动态存储矩阵St的列数，并定义列数为LLL，LLL为正整数；由步骤5的第一个小步骤对目标轨迹动态存储矩阵St的定义可知：目标轨迹动态存储矩阵St的每1列表示1个目标轨迹；得到目标轨迹的总个数为LLL，同时输出LLL个目标的轨迹、速度和加速度。

Images