CN103576137B - 一种基于成像策略的多传感器多目标定位方法 - Google Patents

Abstract

本发明公开了一种基于成像策略的多传感器多目标的定位方法，它是通过将传感器网络建模为一个二维稀疏阵，把基于传感器网络的多目标定位问题转变为成像问题，采用像BP算法投影策略，用来解决繁杂的数据关联问题；通过贪婪算法搜寻三维图像的峰值，并利用剔除算法清除虚假目标定位多个目标。与传统多目标定位方法相比，本发明从成像的角度对目标的多个回波进行相参积累，从而能够更加提高目标信噪比，方便、准确的解算多目标；在成像空间方便的考虑多目标定位，从而提高多目标定位的概率。

Description

一种基于成像策略的多传感器多目标定位方法

技术领域

本发明属于雷达组网系统中的多目标定位的技术领域，它特别涉及到了传感器网络下的多目标定位的技术领域。

背景技术

相比于传统的双站雷达，传感器网络能够在一个较大的空间范围获得目标的多个双站延时，这在军事侦察领域具有非常大的优势。利用这些延时，通过解算有关目标位置的非线性距离方程(原理上类似于GPS系统)，我们能够在一个脉冲周期中定位目标，非常有利于跟踪任务。但是由于多目标，噪声以及杂波的存在，传感器网络必须在定位前分配包括噪声和杂波在内的不同目标的混合的延时，即完成目标和传感器之间的数据关联问题，从而创建正确的距离方程，这在实际处理中是一个极其复杂的任务。围绕上述多目标的数据关联的问题，一些大量精密的工作如联合概率数据关联(详见“XU Can,Li Zhi.Study on Joint Probability Denstiy Algorithm InMulti-Sensor Data Fusion,Computer Science and Automation Engineering(CSAE),2012IEEE International Conference onvol.3,May.2012,pp.:32-37”)等，不断的被探索。然而当目标的监视环境变得复杂，目标个数增多时，通过传统的数据关联算法来解决问题在计算量方面是非常棘手的。

发明内容：

为了克服传感器网络中多目标的数据关联问题，本发明提出了一种基于成像策略的多传感器多目标的定位方法，其特点是通过将传感器网络建模为一个二维稀疏阵，从而把基于传感器网络的多目标定位问题转变为成像问题，其中基于成像角度的多目标定位的示意图如图1所示；采用像BP算法(即后向投影算法)投影策略，用来解决繁杂的数据关联问题，其中BP算法的示意图如图2所示；通过贪婪算法搜寻三维图像的峰值，并利用剔除算法清除虚假目标定位多个目标；相对传统多目标定位方法，本发明能够在成像空间实现多目标定位，对目标的回波进行相参积累以提高信噪比，从而提高多目标定位的概率。

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

定义1、双站延时

由发射机发射信号开始计时，该信号经目标散射回接收机时终止计时，即发射机-目标-接收机，这段从发射到回波接收的时间即为该目标对此接收机的双站延时。双站延时可由安装在发射机和接收机上的GPS或者北斗定位系统同步得到。

定义2、双站距离

双站延时乘以光速所得到的距离差即为由发射机-目标-接收机的双站距离。

定义3、双站距离方程

假设附图1的原理图中目标位置为x，发射机位置为T,接收机的数目为N,接收机的位置为r_i(i＝1，2，…，N),则信号从发射机经目标到达第i个接收机的距离为

||T-x||₂+||x-r_i||₂＝v*τ_i

其中，i＝1，2，…，N，v表示光速，*表示乘号，τ_i表示信号从发射机经目标到达第i个接收机的双站延时，||·||₂表示向量的2-范数，上式即为关于目标位置的双站距离方程。

定义4、距离单元编号

在雷达系统中，将雷达测距的范围划分成若干小的区域并将其编号，雷达根据目标回波信号落入的区域编号计算目标与雷达之间的距离。在此定位系统中，根据接收机所测得的各个目标的双站距离，结合距离向分辨率可得到各目标回波的量化编号。

定义5、合成孔径雷达反向投影算法

反向投影算法是基于匹配滤波原理的合成孔径雷达成像算法，其主要通过相干累加实现合成孔径雷达数据的聚焦成像。详细内容可参考文献：“Research on A novel fast back projection algorithmfor strip map bistatic SAR imaging”，Huang Yulin等。

定义6、接收波门

接收的数据从起始时刻到终止时刻的最短时间长度，即发射信号从发射到接收的时间大于这一最短时间长度时才对数据进行接收，否则不接收数据，这一最短时间长度即为接收时间波门。与接收的时间波门对应的则有接收距离波门，即从发射到接收的路程大于某一最小的距离波门时接收数据，否则不接收数据，这一最小的距离波门即为接收距离波门。

定义7、地理空间

地理空间指目标及杂波、噪声等所在的实际的三维场景空间。

定义8、成像空间

雷达成像空间是指将地理空间中的散射点投影到距离向－方位向的二维空间坐标系，该空间由合成孔径雷达成像空间中的两个相互正交的坐标基确定。目前典型合成孔径雷达的成像空间包括距离向－方位向投影空间。

定义9、笛卡尔坐标系

在数学里，笛卡儿坐标系(Cartesian坐标系)，也称直角坐标系，是一种正交坐标系。二维的直角坐标系是由两条相互垂直、0点重合的数轴构成的。在平面内，任何一点的坐标是根据数轴上对应的点的坐标设定的。在平面内，任何一点与坐标的对应关系，类似于数轴上点与坐标的对应关系。三维笛卡儿坐标系是在二维笛卡儿坐标系的基础上根据右手定则增加第三维坐标(即Z轴)而形成的。

定义10、数据关联

数据关联指针对多目标多雷达的定位，需要分清各个被动站的测量数据中那些是来源于同一目标的，并把属于同一目标的测量数据组合起来，即完成测量数据关联。详见“王成，李少洪，黄槐.多站被动雷达系统的多目标测量数据关联.电子学报.vol.30.No.12.Dec.2002”。

定义11、信噪比

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

定义12、单载频脉冲信号

单一载频脉冲信号可以表示为矩形信号与余弦信号的乘积,即其中f_c是脉冲的载频，T_P是脉冲宽度。详见“皮亦鸣,杨建宇,付毓生,杨晓波.合成孔径雷达成像原理.第一版.电子科技大学出版社.2007.3”。

定义13、相参积累

相参积累指在雷达系统中将来自同一目标的回波进行相位的补偿并将回波数据进行累加以提高信噪比，增加信息量的过程。详见“王瑞军.基于非起伏目标的相参积累TBD方法研究.电子科技大学硕士学位论文.2011”。

定义14：动态存储矩阵

动态存储矩阵指能够根据实际的需要进行动态的更新(存储和清除)矩阵中各元素的值的矩阵。

本发明提出了一种基于成像策略的多传感器多目标定位方法，它包括如下步骤：

步骤1、相关参数的初始化

初始化的参数均已知，如下所示：雷达的空间距离分辨率为δ_r；雷达在距离向上的最大距离单元编号为L，L取值为整数；电磁波传播速度记为v；以发射机位置为原点，在目标，发射机，接收机所在的地理空间建立三维笛卡尔坐标系，则所有位置坐标都是以直角坐标表示；发射机位置记为P_t；空中目标个数记为Ntar，目标的散射截面积记为σ，第i个空中目标的位置记为Ptar(i)(i＝1，2，…，Ntar)，i为目标序号；地面接收机个数记为Nrec，第j个地面接收机位置记为Prec(j)(j＝1，2，…，Nrec)，j为接收机序号；每个接收机的接收距离波门均为Rrec0；每个接收机通道的信噪比记为SNR；雷达发射信号为单载频脉冲信号，其脉冲的载频为f_c，脉冲宽度为T_p；单载频信号由发射机发射经目标i散射，被第j个接收机接收，经过去载频和离散采样后的回波信号记为Echo(i，j)，此回波信号表示为Echo(i，j)＝σ_i，j*exp(-pkR_i，j)，其中σ_i，j表示第i个目标对第j个接收机的雷达散射截面积(RCS)，R_i，j＝||P_t-Ptar(i)||₂+||Ptar(i)-Prec(j)||₂，表示发射信号经目标i散射到达接收机j的双站距离，exp(·)为自然指数为底的指数函数，||·||₂表示向量的2-范数，*表示乘号，π为圆周率，v为光速,f_c为载频，p为虚数单位，即p²＝-1。

步骤2、创建回波存储矩阵

由安装在发射机和接收机上的全球定位系统(GPS或者北斗定位系统，收集到空中目标对接收机的双站延时τ_i，j及相应延时的回波信号Echo(i，j)，τ_i，j表示第i个目标对第j个接收机的双站延时，Echo(i，j)表示第i个目标对第j个接收机相应延时的回波信号。由收集的延时τ_i，j和对应的回波Echo(i，j)建立一个存储目标回波的动态存储矩阵，回波的动态存储矩阵建立过程如下：

步骤2.1初始化矩阵

根据步骤1中的空中目标i、地面接收机j，定义一个L行Nrec列的动态存储矩阵，记为矩阵的每一单元记为 ${\mathrm{Echo}}_L^{\mathrm{Nrec}}(l,c),l=\mathrm{1,2},...,L,c=\mathrm{1,2},...,\mathrm{Nrec},$ 其中l表示矩阵的行号，c表示矩阵的列号；下标L为距离向上的最大距离单元编号，由所有接收机的最大探测距离而定；上标Nrec表示Nrec个接收机通道，即用于存储接收到的所有目标的回波，转到步骤2.2；

步骤2.2存储目标回波

针对步骤2.1中定义的动态存储矩阵进行动态的回波存储。由各个接收机所获得的所有目标的双站延时τ_i，j，计算所有目标的双站延时τ_i，j的距离单元编号，将所有目标的延时τ_i，j所对应的目标回波Echo(i，j)存储在动态存储矩阵的相应单元中；具体存储步骤见以下步骤2.2.1和步骤2.2.2；

步骤2.2.1根据延时信息计算相应回波的双站距离

根据接收机j所获得的目标i的延时τ_i，j及对应的目标回波信号Echo(i，j)，由双站距离公式R_i，j＝v*τ_i，j计算目标i对第j个接收机的双站距离，转到步骤2.2.2；

步骤2.2.2根据双站距离计算相应回波的距离单元编号

根据步骤2.2.1中所获得的双站距离R_i，j计算得到相应回波的距离单元编号，记为ID_i，j，计算公式：其中ID_i，j表示目标i的回波在接收机j中的距离单元编号，1≤ID_i，j≤L，且ID_i，j取值为正整数；其中L为最大的距离单元编号，转到步骤2.2.3；

步骤2.2.3根据距离单元编号将相应的回波进行存储

根据步骤2.2.2中得到的目标i对接收机j的距离单元编号ID_i，j，将步骤1提供的去载频和离散采样后的回波信号Echo(i，j)，存储于动态存储矩阵中的单元 ${\mathrm{Echo}}_L^{\mathrm{Nrec}}(l,c)={\mathrm{ID}}_{i,j},c=j,$ 即 ${\mathrm{Echo}}_L^{\mathrm{Nrec}}(l,c)=\mathrm{Echo}(i,j)(l={\mathrm{ID}}_{i,j},c=j),$ 采用传统的遍历方法遍历所有的目标序号i和接收机序号j，得到Ntar个目标对Nrec个接收机的回波信号Echo(i，j)的距离单元编号ID_i，j，将去载频和离散采样后的回波信号Echo(i，j)根据编号ID_i，j存储于动态存储矩阵的相应单元中，得到完整的Nrec个接收机存储Ntar个目标回波的动态存储矩阵回波的动态存储矩阵示意图如附图3所示；

步骤3、栅格划分

在步骤1中所建立笛卡尔坐标系中，设定一个包含所有目标的三维长方体搜索区域Ω_Sear,长方体的长、宽、高分别平行于坐标系的x轴、y轴、z轴。长方体的长为10×δ_r×Nx,长方体的宽为10×δ_r×Ny，长方体的高为10×δ_r×Nz,其中Nx，Ny，Nz分别为长方体沿x轴，y轴，z轴的采样点数，长方体搜索区域Ω_Sear的中心记为P0，根据步骤1中初始的空间距离分辨率δ_r将长方体搜索区域划分为一组微小的栅格，每一栅格大小为10δ_r×10δ_r×10δ_r，并为每一个栅格单元分配一个代表点，代表点坐标记为P_xyz＝P0+[x-Nx/2，y-Ny/2，z-Nz/2]*10*δ_r，x＝1，2…Nx，y＝1，2…Ny，z＝1，2…Nz，x，y，z分别是长方体沿x轴，y轴，z轴的采样点序号，即P_xyz表示采样点序号为x，y，z的代表点的三维坐标，相邻代表点之间相距10δ_r。

步骤4、信息投影

采用传统的遍历法在步骤3建立的整个搜索区域Ω_Sear中，寻找每一栅格代表点所对应各个接收机的回波，采用传统的相参积累方法对每一栅格代表点所对应各个接收机的回波进行相参积累；具体步骤如下：

步骤4.1

初始化定义一个动态存储目标成像信息的三维矩阵，记为Ima_Nx×Ny×Nz，矩阵的元素单元记为Ima_Nx×Ny×Nz(k1，k2，k3)，k1、k2和k3分别表示矩阵Ima_Nx×Ny×Nz行号、列号和层号，Ima_Nx×Ny×Nz(k1，k2，k3)表示矩阵的任一元素，k1＝1，2，…Nx，k2＝1，2，…Ny，k3＝1，2，…Nz，Nx表示矩阵的行数，Ny表示矩阵的列数，Nz表示矩阵的层数，Nx×Ny×Nz表示矩阵的大小。初始化接收机序号j＝1，长方体沿x轴，y轴，z轴的采样点序号x＝1，y＝1，z＝1，转到步骤4.2；

步骤4.2

针对步骤3中提供的栅格代表点P_xyz＝P0+[x-Nx/2，y-Ny/2，z-Nz/2]*10*δ_r和接收机j，计算从发射机P_t经栅格代表点P_xyz到达接收机j的双站距离，记为R_xyz，j，则R_xyz，j＝||P_t-P_xyz||₂+||P_xyz-Prec(j)||₂，即R_xyz，j表示发射信号经采样点序号为x，y，ｚ的代表点P_xyz散射回接收机j的双站距离，并由R_xyz，j计算代表点P_xyz对接收机j的距离单元编号，记为ID_xyz，j，则即ID_xyz，j表示采样点序号为x，y，z的代表点P_xyz对接收机j的距离单元编号，Rrec0和δ_r分别为步骤1中初始化已知的接收机的接收距离波门和雷达的空间距离分辨率，转到步骤4.3；

步骤4.3

取出步骤2所建立的动态存储矩阵的第j列第ID_xyz，j行中的数据作为代表点P_xyz所对应的回波，记为Echo(xyz，j)，Echo(xyz，j)表示采样点序号为x，y，z的代表点对接收机j的回波，并将回波Echo(xyz，j)在成像空间相参积累，即Ima_Nx×Ny×Nz(k1，k2，k3)＝Ima_Nx×Ny×Nz(k1，k2，k3)+Echo(xyz，j)*exp(pkR_xyz，j)，k1＝x，k2＝y，k3＝z，接收机序号j增加1，并转到步骤4.4；

步骤4.4

如果j≤Nrec，则转到步骤4.2；

如果j＞Nrec，则j置1，x增加1，转到步骤4.5；

步骤4.5

如果x≤Nx，则转到步骤4.2；

如果x＞Nx，则x置1，y增加1，转到步骤4.6；

步骤4.6

如果y≤Ny，则转到步骤4.2；

如果y＞Ny，则y置1，z增加1，转到步骤4.7；

步骤4.7

如果z≤Nz，则转到步骤4.2；

如果z＞Nz，则遍历结束，得到成像空间中一个完整的三维图像矩阵Ima_Nx×Ny×Nz，转到步骤5；

步骤5、提取目标位置

首先，在步骤4得到的三维图像Ima_Nx×Ny×Nz中找出三维图像Ima_Nx×Ny×Nz的Ntar个最大值，记录找到的Ntar个最大值的序号为n(n＝1，2…，Ntar)；然后从三维图像Ima_Nx×Ny×Nz中剔除Ntar个最大值对应的全部信息；具体步骤如下：

步骤5.1

找出步骤4中得到的三维图像矩阵Ima_Nx×Ny×Nz的最大值，记为Max_n，转至步骤5.2；

步骤5.2

记录下最大值Max_n对应的三维矩阵中的位置，记为(k1_{Max_n}，k2_{Max_n}，k3_{Max_n})，即Ima_Nx×Ny×Nz(k1，k2，k3)＝Max_n，k1＝k1_{Max_n}，k2＝k2_{Max_n}，k3＝k3_{Max_n}，k1_{Max_n}，k2_{Max_n}，k3_{Max_n}分别为三维图像矩阵中最大值Max_n所对应的行号，列号，层号，则该三维矩阵位置(k1_{Max_n}，k2_{Max_n}，k3_{Max_n})所对应的地理空间位置记为P_n＝P0+[k1_{Max_n}-Nx/2，k2_{Max_n}-Ny/2，k3_{Max_n}-Nz/2]*10*δ_r，即P_n表示获取的第n个目标的地理空间位置，初始化接收机序号j＝1，转到步骤5.3；

步骤5.3

选择第n个位置P_n＝P0+[k1_{Max_n}-Nx/2，k2_{Max_n}-Ny/2，k3_{Max_n}-Nz/2]*10*δ_r和接收机j，由双站距离公式计算从发射机P_t经位置P_n到达接收机j的双站距离，记为R_n，j，则R_n，j＝||P_t-P_n||₂+||P_n-Prec(j)||₂，即R_n，j表示第n个目标对接收机j的双站距离，并由R_n，j计算目标n的回波距离单元编号ID_n，j表示目标n对接收机j的回波的距离单元编号，转到步骤5.4；

步骤5.4

在步骤2建立的动态存储矩阵的第j列第ID_n，j个距离单元中找出目标n对应的回波，记为Echo(n，j)，即 $\mathrm{Echo}(n,j)={\mathrm{Echo}}_L^{\mathrm{Nrec}}(l,c)(l={\mathrm{ID}}_{n,j},c=j),$ Echo(n，j)表示目标n对接收机j的回波，并将此回波从图像矩阵中剔除，即

Ima_Nx×Ny×Nz(k1，k2，k3)＝Ima_Nx×Ny×Nz(k1，k2，k3)-Echo(n，j)*exp(pkR_n，j)

其中k1＝k1_{Max_n}，k2＝k2_{Max_n}，k3＝k3_{Max_n}，j增加1，并转到步骤5.5；

步骤5.5

如果j≤Nrec，则转到步骤5.3；

如果j≥Nrec，则得到消除第n个最大值信息后的三维矩阵Ima_Nx×Ny×Nz，n增加1，并转到步骤5.6；

步骤5.6

如果n≤Ntar，则转到步骤5.1；

如果n＞Ntar，则转到步骤6；

步骤6

多目标定位结束，得到Ntar个目标在三维地理空间的位置P_n(n＝1，2，…，Ntar)。

本发明的创新点在于传统的基于传感器网络的多目标定位中，需要进行复杂的数据关联，而且又需要创建并解算复杂的多目标定位的非线性距离方程组，本发明提供一种基于成像策略的多传感器多目标的定位方法。该方法从成像的角度对目标的多个回波进行相参积累，从而能够更加提高目标信噪比，并方便、准确的解算多目标。

本发明的实质是通过将传感器网络建模为一个二维稀疏阵，从而把基于传感器网络的多目标定位问题看作成像问题，其中基于成像角度的多目标定位的示意图如图1所示；我们采用像BP算法，(详见“Jun,S.；Long,M.；Xiaoling,Z.,Streaming BP for Non-Linear Motion Compensation SAR ImagingBased on GPU,IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing,vol.PP,Issue:99,2013,pp.1–16”)，即后向投影算法这样的投影策略，能够很容易的被用来解决繁杂的数据关联问题，其中BP算法的示意图如图2所示；另外，投影方法的优势是三维图像的值包含了一些额外的目标信息，譬如目标的雷达散射截面积(RCS)或者存在概率等等，而这些信息在传统的基于数据关联的有关目标位置的双站距离方程组的创建中被舍弃了；通过贪婪算法搜寻三维图像的峰值，并利用剔除算法清除虚假目标，我们能够定位多个目标。

本发明的优点：本发明从成像的角度考虑多目标定位问题，并利用了目标多个方向的回波信号的相位信息，进行相参积累以提高信噪比，相对传统多目标定位方法，本发明能够在成像空间方便的考虑多目标定位，并能够利用成像的知识对目标的回波进行相参积累以提高信噪比，从而提高多目标定位的概率。

附图说明

图1为基于成像角度的多目标定位的模型

其中，多个传感器可看作接收目标散射信号的稀疏阵列。

图2为多传感器定位中的BP算法示意图

其中，A和B表示不同的目标的回波经过距离压缩后脉冲，F表示虚假目标的回波脉冲。N表示传感器接收机的个数，t表示时间轴。针对同一目标在不同传感器中的回波在成像空间进行相参累加，从而进行能量的积累。

图3为多目标的回波存储矩阵示意图

其中，L表示最大的距离单元编号，ID表示距离单元序号，回波存储矩阵为L行Nrec列，ID表示距离单元的序号，n表示接收机通道的序号。每一列对应一个接收机通道，用于存储各个目标的回波，。同一目标对不同传感器的散射回波放在不同的列，不同的目标对同一传感器的散射信号放在同一列的不同距离单元中。如目标1对接收机机1的回波放在回波矩阵的第一列对应的距离单元编号中。

图4为本发明方法的流程图

具体实施方式

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

步骤1、基于成像策略的传感器网络多目标定位方法相关参数的初始化

初始化目标定位的参数如下：雷达的空间距离分辨率为δ_r＝10m；雷达在距离向上的距离单元个数为L＝10000；电磁波传播速度v＝3×10⁸m/s；以发射机位置为原点，在目标，发射机，接收机所在的地理空间建立笛卡尔坐标系；发射机位置记为P_t＝[0，0，0]；空中目标个数Ntar＝10，目标的散射截面积σ＝(30×rand(1)+1)m²，目标随机分布在以[100，200，10]×1000为中心，边长两公里的立方体区域中，得第i个目标位置记为Ptar(i)＝[100，200，10]×1000+[round((rand-0.5)/0.1)，round((rand-0.5)/0.1)，round((rand-0.5)/0.1)]，i＝1，2，…，10；地面接收机个数Nrec＝40，接收机随机分布在以[100，50，0]×1000为中心的，边长20km的正方形区域中，得到第个接收机位置记为Prec(j)＝[100，50，0]×1000+[rand(1)-0.5，rand(1)-0.5，0]×20×1000，j＝1，2，…，40；每个接收机的接收距离波门为Rrec0＝300km；每个接收机通道的信噪比为SNR＝-20db；雷达发射信号为单载频信号，脉冲的载频f_c＝2GHz；脉冲宽度T_p＝20μs；单载频信号由发射机发射经目标i散射，被第j个接收机接收，经过去载频和离散采样后的回波信号记为Echo(i，j)， $\mathrm{Echo}(i,j)=(30\×\mathrm{rand}\left(1\right)+1)\×\exp (-p\frac{2\mathrm{\π}\×2\×{10}^9}{3\×{10}^8}{R}_{i,j}),$ R_i，j＝||P_t-Ptar(i)||₂+||Ptar(i)-Prec(j)||_２，表示发射信号经目标i散射到达接收机j的双站距离，i＝1，2，…，10，j＝1，2，…，40，exp(·)为自然指数为底的指数函数，||·||₂表示向量的2-范数，π为圆周率，v为光速，p为虚数单位，即p²＝-1。

步骤2、创建回波存储矩阵

由安装在发射机和接收机上的全球定位系统(GPS)或者北斗定位系统，收集到目标对接收机的双站延时τ_i，j及相应延时的回波信号Echo(i，j)，i＝1，2，…，10，j＝1，2，…，40，i为目标序号，j为接收机序号，τ_i，j表示第i个目标对第j个接收机的双站延时，Echo(i，j)表示第i个目标对第j个接收机散射的回波信号。由收集到的延时τ_i，j和对应的回波Echo(i，j)建立一个存储目标回波的动态存储矩阵。回波的动态存储矩阵建立过程如下：

步骤2.1初始化矩阵

根据步骤1中的空中目标i(i＝1，2，…10)，地面接收机j(j＝1，2，…40)以及距离向上的最大距离单元编号L＝10000，定义一个10000行40列的动态存储矩阵，记为矩阵的每一单元记为 ${\mathrm{Echo}}_{10000}^{40}(l,c),l=\mathrm{1,2},...,10000,c=\mathrm{1,2},...,40,$ 其中l表示矩阵的行号，c表示矩阵的列号，转到步骤2.2；

步骤2.2存储目标回波

针对步骤2.1中定义的矩阵进行动态的回波存储。数据处理中心根据各个接收机所获得的所有目标的双站延时τ_i，j，计算此延时的距离单元编号，并据此编号在相应接收机通道中将延时τ_i，j所对应的目标回波Echo(i，j)存储在矩阵的相应单元中，i＝1，2，…10；j＝1，2，…40，其中i为目标序号，j为接收机序号，具体步骤见以下步骤2.2.1和步骤2.2.2；

步骤2.2.1根据延时信息计算相应回波的双站距离

数据处理中心根据接收机j所获得的目标i的延时τ_i，j及对应的回波信号Echo(i，j)，计算得到目标i对第j个接收机的双站距离R_i，j＝3×10⁸×τ_i，j，R_i，j则表示第i个目标对第j个接收机的双站距离，i＝1，2，…10；j＝1，2，…40，其中i为目标序号，j为接收机序号，10为目标数目，40为接收机数目，转到步骤2.2.2；

步骤2.2.2根据双站距离计算相应回波的距离单元编号

根据步骤2.2.1中所获得的双站距离R_i，j计算得到相应回波的距离单元编号，记为ID_i，j，则ID_i，j表示目标i的回波在接收机j中的距离单元编号，1≤ID_i，j≤10000，且ID_i，j取值为正整数；i＝1，2，…，10，j＝1，2，…，40，其中i为目标序号，j为接收机序号，10为目标数目40为接收机数目，转到步骤2.2.3；

步骤2.2.3根据距离单元编号将相应的回波进行存储

根据步骤2.2.2中得到的目标i对接收机j的距离单元编号ID_i，j，将数据处理中心所收集得到的相应回波Echo(i，j)，存储于矩阵中的单元即遍历所有的目标序号i(i＝1，2，…，10)和接收机序号j(j＝1，2，…，40)，得到10个目标对40个接收机的回波Echo(i，j)的距离单元编号ID_i，j(i＝1，2，…，10，j＝1，2，…，40)，将回波Echo(i，j)根据编号ID_i，j存储于动态存储矩阵的相应单元中，得到最终完整的10个接收机存储40个目标回波的动态存储矩阵其中1≤ID_i，j≤10000且ID_i，j取值为正整数，i＝1，2，…，10，j＝1，2，…，40，其中i为目标序号，j为接收机序号，目标数目为10个，接收机数目为40个，回波存储矩阵示意图如附图3所示；

步骤3、栅格划分

在步骤1中所建立笛卡尔坐标系中设定一个包含所有目标的三维长方体搜索区域Ω_Sear，长方体的长，宽，高分别平行于坐标系的x轴，y轴，z轴。长方体沿三个坐标轴的采样点数皆为40，长方体的长为10×10×40，长方体的宽为10×10×40，长方体的高为10×10×40，其长方体搜索区域Ω_Sear的中心记为P0＝[100，200，10]×1000，根据步骤1中的空间分辨率δ_r＝10m将搜索区域划分为一组微小的栅格，每一栅格大小为100×100×100，并为每一个栅格单元分配一个代表点，代表点坐标记为P_xyz＝[100，200，10]×1000+[x-40/2，y-40/2，z-40/2]×10×10，x＝1，2…40，y＝1，2…40，z＝1，2…40，x，y，z分别是立方体沿x轴，y轴，z轴的采样点序号，即P_xyz表示采样点为x，y，z的代表点的三维坐标，相邻代表点之间相距100m。

步骤4、信息投影

采用传统的遍历法寻找步骤3整个搜索区域Ω_Sear中的每一栅格代表点针对各个接收机的回波并相参积累；

步骤4.1

初始化定义一个动态存储目标成像信息的三维矩阵，记为Ima_N40×40×40，矩阵的元素单元记为Ima_40×40×40(k1，k2，k3)，k1，k2和k3分别表示矩阵Ima_40×40×40行号，列号和层号，Ima_40×40×40(k1，k2，k3)表示矩阵的任一元素，k1＝1，2，…40，k2＝1，2，…40，k3＝1，2，…40，40×40×40表示矩阵的大小。初始化接收机序号j＝1，长方体沿x轴，y轴，z轴的采样序号为x＝1，y＝1，z＝1，转到步骤4.2；

步骤4.2

选择点P_xyz＝[100，200，10]×1000+[x-40/2，y-40/2，z-40/2]×10×10和接收机j，计算从发射机P_t＝[0，0，0]经代表点P_xyz到达接收机j的双站距离R_xyz，j＝||P_t-P_xyz||₂+||P_xyz-Prec(j)||₂，即R_xyz，j表示发射信号经采样序号为x，y，z的代表点P_xyz散射回接收机j的双站距离，||·||₂表示向量的2-范数，并由此R_xyz，j计算代表点P_xyz对接收机j的距离单元编号，记为ID_xyz，j，则即ID_xyz，j表示样序号为x，y，z的代表点P_xyz对接收机j的距离单元编号，转到步骤4.3；

步骤4.3

如果0＜ID_xyz，j≤10000，则取出步骤2建立的回波矩阵的第j列第ID_xyz，j行的数据作为代表点P_xyz所对应的回波，记为Echo(xyz，j)，Echo(xyz，j)表示采样点序号为x，y，z的代表点对接收机j的回波，并将回波Echo(xyz，j)在成像空间相参积累，即 $\begin{array}{c}{\mathrm{Ima}}_{40\×40\×40}(k1,k2,k3)=\\ {\mathrm{Ima}}_{40\×40\×40}(k1,k2,k3)+\mathrm{Echo}(\mathrm{xyz},j)*\exp \left(p\frac{2\mathrm{\π}\×2\×{10}^9}{3\×{10}^8}{R}_{\mathrm{xyz},j}\right)\end{array},,$ k1＝x，k2＝y，k3＝z，R_xyz，j表示发射信号经采样序号为x，y，z的代表点P_xyz散射回接收机j的双站距离，*表示乘号，p为虚数单位，即p²＝-1，exp(·)为自然指数为底的指数函数，j增加1并转到步骤4.4；

如果ID_xyz，j≤0或者ID_xyz，j＞10000,则接收机序号j增加1，转到步骤4.4；

步骤4.4

如果j≤40，则转到步骤4.2；

如果j＞40，则j置1，x增加1，转到步骤4.5；

步骤4.5

如果x≤40，则转到步骤4.2；

如果x＞40，则x置1，y增加1，转到步骤4.6；

步骤4.6

如果y≤40，则转到步骤4.2；

如果y＞40，则y置1，z增加1，转到步骤4.7；

步骤4.7

如果z≤40，则转到步骤4.2；

如果z＞40，则遍历结束，得到成像空间中一个完整的三维图像Ima_40×40×40，转到步骤5；

步骤5、提取目标位置

逐步找出步骤4中得到的三维图像的10个最大值，并从图像中剔除最大值对应的全部信息，记录找到的10个最大值的序号为n(i＝1，２…，10),并初始化n＝1,转到步骤5.1：

步骤5.1

找出步骤4中得到的最终完整的三维图像矩阵Ima_40×40×40的最大值，记为Max_n，转至步骤5.2；

步骤5.2

记录下最大值Max_n对应的三维矩阵中的位置，记为(k1_{Max_n}，k2_{Max_n}，k3_{Max_n})，即Ima_40×40×40(k1，k2，k3)＝Max_n，k1＝k1_{Max_n}，k2＝k2_{Max_n}，k3＝k3_{Max_n}，k1_{Max_n}，k2_{Max_n}，k3_{Max_n}分别为成像矩阵中最大值Max_n所对应的行号，列号，层号，则该三维矩阵位置(k1_{Max_n}，k2_{Max_n}，k3_{Max_n})所对应的地理空间位置记为P_n＝[100，200，10]×1000+[k1_{Max_n}-40/2，k2_{Max_n}-40/2，k3_{Max_n}-40/2]×100，即P_n表示获取的第n个目标的地理空间位置，初始化接收机序号j＝1，转到步骤5.3；

步骤5.3

选择第n个位置P_n＝[100，200，10]×1000+[k1_{Max_n}-Nx/2，k2_{Max_n}-Ny/2，k3_{Max_n}-Nz/2]×100和接收机j，计算从发射机P_t＝[0，0，0]经位置P_n到达接收机j的双站距离，记为R_n，j，R_n，j＝||P_t-P_Max||₂+||P_Max-Prec(j)||₂，即R_n，j表示第n个目标对接收机j的双站距离，||·||₂表示向量的2-范数，并由R_n，j计算目标n的回波距离单元编号ID_n，j表示目标n对接收机j的回波的距离单元编号，转到步骤5.4；

步骤5.4

如果0＜ID_n，j≤10000，则在步骤2建立的回波矩阵的第j列第ID_n，j个距离单元中找出目标n对应的回波，记为Echo(n，j)，即 $\mathrm{Echo}(n,j)={\mathrm{Echo}}_{10000}^{40}(l,c)(l={\mathrm{ID}}_{n,j},c=j),$ Echo(n，j)表示目标n对接收机j的回波，并将此回波从图像矩阵中剔除，即

$\begin{array}{c}{\mathrm{Ima}}_{40\×40\×40}(k1,k2,k3)\\ ={\mathrm{Ima}}_{40\×40\×40}(k1,k2,k3)-\mathrm{Echo}(n,j)*\exp \left(p\frac{2\mathrm{\π}\×2\×{10}^9}{3\×{10}^8}{R}_{n,j}\right)\end{array}$

其中k1＝k1_{Max_n}，k2＝k2_{Max_n}，k3＝k3_{Max_n}，R_n，j表示第n个目标对接收机j的双站距离，*表示乘号，p为虚数单位，即p²＝-1，exp(·)为自然指数为底的指数函数，j增加1并转到步骤5.5；

如果ID_n，j≤0或者ID_n，j＞10000,则j增加1，转到步骤5.5；

步骤5.5

如果j≤40，则转到步骤5.3；

如果j＞40，则得到消除第n个最大值信息后的三维矩阵Ima_40×40×40，n增加1，并转到步骤5.6；

步骤5.6

如果n≤10，则j置1，求取三维矩阵Ima_40×40×40的最大值赋给Max_n，转到步骤5.2；

如果n＞10，则转到步骤6；

步骤6

多目标定位结束，并得到Ntar个目标在三维地理空间的位置P_n(n＝1，２，…，10)。

通过本发明的具体实施可以看出，本发明通过成像的角度对基于传感器网络多目标回波进行了处理，较传统的基于数据关联建立非线性方程的方法快捷而且方便。该方法在基于上述初始化参数，并在接收机通道具有-20db的信噪比的情况下，经过多次蒙特卡洛实验，多目标定位的统计概率可以达到0.8265，可见该方法在低信噪比的情况下依然可以较为准确的对多目标进行定位。

Claims (1)

1.一种基于成像策略的多传感器多目标定位方法，其特征是该方法包括如下步骤：

步骤1、相关参数的初始化

初始化的参数均已知，如下所示：雷达的空间距离分辨率为δ_r；雷达在距离向上的最大距离单元编号为L，由所有接收机的最大探测距离而定，L取值为整数；电磁波传播速度记为v；以发射机位置为原点，在目标，发射机，接收机所在的地理空间建立三维笛卡尔坐标系，则所有位置坐标都是以直角坐标表示；发射机位置记为P_t；空中目标个数记为Ntar，目标的散射截面积记为σ，第i个空中目标的位置记为Ptar(i)，i＝1，2，…，Ntar，i为目标序号；地面接收机个数记为Nrec，第j个地面接收机位置记为Prec(j)，j＝1，2，…，Nrec，j为接收机序号；每个接收机的接收距离波门均为Rrec0；每个接收机通道的信噪比记为SNR；雷达发射信号为单载频脉冲信号，其脉冲的载频为f_c，脉冲宽度为T_p；单载频信号由发射机发射经目标i散射，被第j个接收机接收，经过去载频和离散采样后的回波信号记为Echo(i，j)，此回波信号表示为Echo(i，j)＝σ_i，j*exp(-pkR_i，j)，其中σ_i，j表示第i个目标对第j个接收机的雷达散射截面积(RCS)，R_i，j＝||P_t-Ptar(i)||₂+||Ptar(i)-Prec(j)||₂，表示发射信号经目标i散射到达接收机j的双站距离，exp(·)为自然指数为底的指数函数，||·||₂表示向量的2-范数，*表示乘号，π为圆周率，v为光速,f_c为载频，p为虚数单位，即p²＝-1；

步骤2、创建回波存储矩阵

由安装在发射机和接收机上的全球定位系统(GPS)或者北斗定位系统，收集到空中目标对接收机的双站延时τ_i，j及相应延时的回波信号Echo(i，j)，τ_i，j表示第i个目标对第j个接收机的双站延时，Echo(i，j)表示第i个目标对第j个接收机相应延时的回波信号；由收集的延时τ_i，j和对应的回波Echo(i，j)建立一个存储目标回波的动态存储矩阵；回波的动态存储矩阵建立过程如下：

步骤2.1初始化矩阵

根据步骤1中的空中目标i，地面接收机j，以及距离向上的最大距离单元编号L，定义一个L行Nrec列的动态存储矩阵，记为矩阵的每一单元记为 l＝1，2，…，L，c＝1，2，…，Nrec，其中l表示矩阵的行号，c表示矩阵的列号；下标L为最大的距离单元编号，由所有接收机的最大探测距离而定；上标Nrec表示Nrec个接收机通道，即用于存储接收到的所有目标的回波，转到步骤2.2；

步骤2.2存储目标回波

针对步骤2.1中定义的动态存储矩阵进行动态的回波存储；由各个接收机所获得的所有目标的双站延时τ_i，j，计算所有目标的双站延时τ_i，j的距离单元编号，将所有目标的延时τ_i，j所对应的目标回波Echo(i，j)存储在动态存储矩阵的相应单元中；具体存储步骤见以下步骤2.2.1和步骤2.2.2；

步骤2.2.1根据延时信息计算相应回波的双站距离

根据接收机j所获得的目标i的延时τ_i，j及对应的目标回波信号Echo(i，j)，由双站距离公式R_i，j＝v*τ_i，j计算目标i对第j个接收机的双站距离，转到步骤2.2.2；

步骤2.2.2根据双站距离计算相应回波的距离单元编号

根据步骤2.2.1中所获得的双站距离R_i，j计算得到相应回波的距离单元编号，记为ID_i，j，计算公式：其中ID_i，j表示目标i的回波在接收机j中的距离单元编号，1≤ID_i，j≤L，且ID_i，j取值为正整数；其中L为步骤1中已知的最大的距离单元编号，转到步骤2.2.3；

步骤2.2.3根据距离单元编号将相应的回波进行存储

根据步骤2.2.2中得到的目标i对接收机j的距离单元编号ID_i，j，将步骤1提供的去载频和离散采样后的回波信号Echo(i，j)，存储于动态存储矩阵中的单元 即采用传统的遍历方法遍历所有的目标序号i和接收机序号j，得到Ntar个目标对Nrec个接收机的回波信号Echo(i，j)的距离单元编号ID_i，j，将去载频和离散采样后的回波信号Echo(i，j)根据编号ID_i，j存储于动态存储矩阵的相应单元中，得到完整的Nrec个接收机存储Ntar个目标回波的动态存储矩阵

步骤3、栅格划分

在步骤1中所建立笛卡尔坐标系中，设定一个包含所有目标的三维长方体搜索区域Ω_Sear,长方体的长、宽、高分别平行于坐标系的x轴、y轴、z轴；长方体的长为10×δ_r×Nx,长方体的宽为10×δ_r×Ny,长方体的高为10×δ_r×Nz,其中Nx，Ny，Nz分别为长方体沿x轴，y轴，z轴的采样点数，长方体搜索区域Ω_Sear的中心记为P0，根据步骤1中初始的空间距离分辨率δ_r将长方体搜索区域划分为一组微小的栅格，每一栅格大小为10δ_r×10δ_r×10δ_r，为每一个栅格单元分配一个代表点，代表点坐标记为P_xyz＝P0+[x-Nx/2，y-Ny/2，z-Nz/2]*10*δ_r，x＝1，2…Nx，y＝1，2…Ny，z＝1，2…Nz，x，y，z分别是长方体沿x轴，y轴，z轴的采样点序号，即P_xyz表示采样点序号为x，y，z的代表点的三维坐标，相邻代表点之间相距10δ_r；

步骤4、信息投影

采用传统的遍历法在步骤3建立的整个搜索区域Ω_Sear中，寻找每一栅格代表点所对应各个接收机的回波，采用传统的相参积累方法对每一栅格代表点所对应各个接收机的回波进行相参积累；具体步骤如下：

步骤4.1

初始化定义一个动态存储目标成像信息的三维矩阵，记为Ima_Nx×Ny×Nz，矩阵的元素单元记为Ima_Nx×Ny×Nz(k1，k2，k3)，k1、k2和k3分别表示矩阵Ima_Nx×Ny×Nz行号、列号和层号，Ima_Nx×Ny×Nz(k1，k2，k3)表示矩阵的任一元素，k1＝1，2，…Nx，k2＝1，2，…Ny，k3＝1，2，…Nz,Nx表示矩阵的行数，Ny表示矩阵的列数，Nz表示矩阵的层数，Nx×Ny×Nz表示矩阵的大小；初始化接收机序号j＝1，长方体沿x轴，y轴，z轴的采样点序号x＝1，y＝1，z＝1，转到步骤4.2；

步骤4.2

针对步骤3中提供的栅格代表点P_xyz＝P0+[x-Nx/2，y-Ny/2，z-Nz/2]*10*δ_r和接收机j，计算从发射机P_t经栅格代表点P_xyz到达接收机j的双站距离，记为R_xyz，j，则R_xyz，j＝||P_t-P_xyz||₂+||P_xyz-Prec(j)||₂，即R_xyz，j表示发射信号经采样点序号为x，y，z的代 表点P_xyz散射回接收机j的双站距离，并由R_xyz，j计算代表点P_xyz对接收机j的距离单元编号，记为ID_xyz，j，则即ID_xyz，j表示采样点序号为x，y，z的代表点p_xyz对接收机j的距离单元编号，Rrec0和δ_r分别为步骤1中初始化已知的接收机的接收距离波门和雷达的空间距离分辨率，转到步骤4.3；

步骤4.3

取出步骤2所建立的动态存储矩阵的第j列第ID_xyz，j行中的数据作为代表点P_xyz所对应的回波，记为Echo(xyz，j)，Echo(xyz，j)表示采样点序号为x，y，z的代表点对接收机j的回波，并将回波Echo(xyz，j)在成像空间相参积累，即 

k1＝x，k2＝y，k3＝z；接收机序号j增加1，并转到步骤4.4；

步骤4.4

如果j≤Nrec，则转到步骤4.2；

如果j＞Nrec，则j置1，x增加1，转到步骤4.5；

步骤4.5

如果x≤Nx，则转到步骤4.2；

如果x＞Nx，则x置1，y增加1，转到步骤4.6；

步骤4.6

如果y≤Ny，则转到步骤4.2；

如果y＞Ny，则y置1，z增加1，转到步骤4.7；

步骤4.7

如果z≤Nz，则转到步骤4.2；

如果z＞Nz，则遍历结束，得到成像空间中一个完整的三维图像矩阵Ima_Nx×Ny×Nz，转到步骤5；

步骤5、提取目标位置

首先，在步骤4得到的三维图像Ima_Nx×Ny×Nz中找出三维图像Ima_Nx×Ny×Nz的Ntar个最大值，记录找到的Ntar个最大值的序号为n，n＝1，2…，Ntar；然后从三维图像Ima_Nx×Ny×Nz中剔除Ntar个最大值对应的全部信息；具体步骤如下：

步骤5.1

找出步骤4中得到的三维图像矩阵Ima_Nx×Ny×Nz的最大值，记为Max_n，转至步骤5.2；

步骤5.2

记录下最大值Max_n对应的三维矩阵中的位置，记为(k1_{Max_n}，k2_{Max_n}，k3_{Max_n})，即Ima_Nx×Ny×Nz(k1，k2，k3)＝Max_n，k1＝k1_{Max_n}，k2_{Max_n}，k3＝k3_{Max_n}，k1_{Max_n}，k2_{Max_n}，k3_{Max_n}分别为三维图像矩阵中最大值Max_n所对应的行号，列号，层号，则该三维矩阵位置(k1_{Max_n}，k2_{Max_n}，k3_{Max_n})所对应的地理空间位置记为P_n＝P0+[k1_{Max_n}-Nx/2，k2_{Max_n}-Ny/2，k3_{Max_n}-Nz/2]*10*δ_r，即P_n表示获取的第n个目标的地理空间位置，初始化接收机序号j＝1，转到步骤5.3；

步骤5.3

选择第n个位置P_n＝P0+[k1_{Max_n}-Nx/2，k2_{Max_n}-Ny/2，k3_{Max_n}-Nz/2]*10*δ_r和接收机j，由双站距离公式计算从发射机P_t经位置P_n到达接收机j的双站距离，记为R_n，j，则R_n，j＝||P_t-P_n||₂+||P_n-Prec(j)||₂，即R_n，j表示第n个目标对接收机j的双站距离，并由R_n，j计算目标n的回波距离单元编号ID_n，j表示目标n对接收机j的回波的距离单元编号，转到步骤5.4；

步骤5.4

在步骤2建立的动态存储矩阵的第j列第ID_n，j个距离单元中找出目标n对应的回波，记为Echo(n，j)，即Echo(n，j)表示目标n对接收机j的回波，并将此回波从图像矩阵中剔除，即

Ima_Nx×Ny×Nz(k1，k2，k3)＝Ima_Nx×Ny×Nz(k1，k2，k3)-Echo(n，j)*exp(pkR_n，j)

其中k1＝k1_{Max_n}，k2＝k2_{Max_n}，k3＝k3_Max__n，j增加1，并转到步骤5.5；

步骤5.5

如果j≤Nrec，则转到步骤5.3；

如果j＞Nrec，则得到消除第n个最大值信息后的三维矩阵Ima_Nx×Ny×Nz，n增加1，并转到步骤5.6；

步骤5.6

如果n≤Ntar，则转到步骤5.1；

如果n＞Ntar，则转到步骤6；

步骤6

多目标定位结束，得到Ntar个目标在三维地理空间的位置P_n，n＝1，2，…，Ntar。