一种未知衰减系数下基于多频率代价函数融合的宽带信号直接定位方法

1.本发明属于无源定位技术领域，具体涉及一种未知衰减系数下基于多频率代价函数融合的宽带信号直接定位方法。


背景技术：

2.定位技术在无线通信、信号处理、导航、目标检测等诸多领域发挥了巨大的作用。无源定位技术作为定位技术的重要组成部分，自20世纪40年代以来渐渐受到了广泛的关注。与有源定位不同，要确定辐射源的位置，无源定位设备不需要发送信号，只需要接收辐射源发出的信号。因此，它具有较高的隐蔽性。在军事方面，无源定位主要用于雷达、通信系统和一些信噪比低的飞行器的定位。在民用方面，被动定位在室内导航中得到了广泛的应用。
3.早期无源定位一般采用传统的两步定位法实现：第一步是利用观测站接收到的信号估计中间参数。最常用的参数有到达方向(doa)、到达时间(toa)、到达频差(fdoa)、到达时差(tdoa)、接收信号强度(rss)等；第二步是根据上述参数和观测站分布建立相应的方程。然后通过某种特定的方法求解方程，得到辐射源的位置。
4.传统的两步定位法虽然被广泛应用，但也有其缺点。很明显，估计的性能受到第一步中间参数估计的限制。此外，在两步定位法中，没有考虑来自同一辐射源的不同观测站接收信号之间的相关性。随着无线定位技术在更加复杂的定位场景中得到广泛应用，两步定位法定位精度和分辨率不足的缺点逐渐显现。
5.为了避免上述问题，提高定位性能，直接定位(dpd)方法被提出。该技术利用多阵列观测站将目标的角度信息融合到阵列流形中。dpd无需进行参数估计，可以直接从观测站的接收数据中估计出目标位置。因此，dpd方法可以显著提高定位性能，尤其是在信噪比较低的情况下。然而，由于要进行多维矩阵运算，dpd的算法复杂度较高，计算量往往较大，需要花费较长时间。此外，基于阵列的dpd方法需要花费的硬件成本较高，而由单天线构成的分布式基站成本较低，可以用来实现目前还较少的基于tdoa的dpd方法。


技术实现要素：

6.本发明所要解决的技术问题是针对上述现有技术的不足，提供一种未知衰减系数下基于多频率代价函数融合的宽带信号直接定位方法，基于分布式监测节点接收辐射源信号，利用信号分段以及将宽带信号分解为窄带信号集合，有利于减少矩阵计算的计算量，能够有效降低定位所需的时间成本，且不需要参数估计，可以有效提升估计精度；此外，利用分布式的单天线观测节点进行定位相比于阵列天线有助于节约成本。
7.为实现上述技术目的，本发明采取的技术方案为：
8.一种未知衰减系数下基于多频率代价函数融合的宽带信号直接定位方法，包括：
9.步骤1：分布式监测节点接收辐射源宽带信号并进行采样，将采样后的宽带信号传
输至中心站点；
10.步骤2：中心站点将所有接收数据整合，将接收数据分段后进行数据预处理，得到宽带信号中每个频率对应的窄带信号；
11.步骤3：对每个窄带信号的协方差矩阵进行特征分解，得到其噪声子空间，根据多频带中的正交关系构造出一个新的代价函数；
12.步骤4：选取定位区域和搜索网格大小，采用新的代价函数，在衰减系数未知情况下结合拉格朗日乘子法，在选定区域内进行网格化的谱峰搜索，最终得到辐射源坐标估计值。
13.为优化上述技术方案，采取的具体措施还包括：
14.上述的步骤1具体步骤包括：
15.(11)假设共有m个未知位置的辐射源，且坐标分别为：p1(x1,y1,z1)，
…
，pm(xm,ym,zm)，分布式监测节点数量为l，坐标分别为s1(x
s1
,y
s1
,z
s1
)，
…
，s
l
(x
sl
,y
sl
,z
sl
)；
16.设采样率为fs，则各监测节点的接收信号表示为：
[0017][0018]
其中，sm(t)表示第m个辐射源在t时刻的发射信号，α
lm
表示传播路径上的信号衰减；
[0019]nl
(t)为相应节点接收信号时产生的噪声，考虑为零均值高斯噪声；
[0020]
τ
lm
是表示信号从辐射源传输到监测站的时延，其表达式是：
[0021][0022]
其中，c表示电磁波传播的速度；
[0023]
(12)对接收信号进行采样，采样间隔为t＝1/fs，设t＝nst，则接收信号的离散形式表示为：
[0024][0025]
其中，ns表示信号的采样点数。
[0026]
上述的步骤2具体步骤包括：
[0027]
(21)中心站点对所有监测节点的信号进行接收与整合，得到行向量对应监测节点接收信号的矩阵，该矩阵的列数为ns，将其按列平均划分为j个子矩阵，则每个子矩阵列数为n＝ns/j，此时宽带信号被划分为j段信号；
[0028]
(22)对每个子矩阵的行向量，即分段后的接收信号进行频域变换，分离出信号中的时延和频率信息，此时分段前后信号频谱包络相同，且每段信号由n个频率点的信息组成；
[0029]
(23)将同一个监测节点的分段信号相同频率对应的频谱信息按序取出，组成一个该频率下的窄带数据，其中第n个信号为xn(j)，此时原本的宽带信号被划分为n个长度为j的窄带信号的集合，每个频率对应信号的矩阵维数为l
×
j。
[0030]
上述的(21)中，把所有观测站收到的信号整合起来，得到合成信号的
表达式：
[0031][0032]
将矩阵y(ns)平均划分为j个子矩阵，每个子矩阵的列数为n＝ns/j，则分段后的信号表示为
[0033]
y(n)＝[y1(n) y2(n)
ꢀ…ꢀ
yj(n)],1≤n≤n
[0034]
其中，式中，其中，s
mj
(n)＝sm(ns),(j-1)n 1≤ns≤jn。
[0035]
上述的(22)中，数据分段后，对yj(n)进行频域变换，得到
[0036][0037]
其中，y
lj
(k)是y
lj
(n)频域形式，且表示发射信号中对应的分段；
[0038]smj
(k)表示s
mj
(n)的dft，nj(k)是分段后噪声的dft；
[0039]
τ
lm
是辐射源位置的函数，可以被表示为τ
l
(pm)；
[0040]
其中，pm表示第m个辐射源的位置；
[0041]
则重写为：
[0042][0043]
其中，矩阵φ(α)包含所有的衰减系数，a(ωk,p)是包含位置信息的时延因子，其表达式分别为：
[0044][0045][0046]
yj(k)包含n组不同频率的信息。
[0047]
上述的步骤3包含以下几个步骤：
[0048]
(31)求出所有窄带信号对应的协方差矩阵，并对其进行特征分解，得到其噪声子空间；
[0049]
(32)根据每个频带中的正交关系写出对应的单频带代价函数，并将所有频带下的代价函数整合，得到多频率融合的新的代价函数。
[0050]
上述的(31)中，对的协方差矩阵进行特征分解，得到：
[0051][0052]
其中，表示第n个窄带信号的信号子空间，表示噪声子空间。
[0053]
上述的(32)中，由噪声子空间和节点与辐射源方向流形的正交关系，得到
[0054][0055]
其中，表示方向向量；
[0056]
表示目标坐标组成的向量；
[0057]
包含各观测站的所有衰减系数；
[0058]
将所有频率对应的等式合成为下式：
[0059][0060]
由上式定义新的代价函数为：
[0061][0062]
上述的步骤4具体实现为：
[0063]
划分出感兴趣的区域作为定位区域，并定义搜索网格的大小，将每一个网格点坐标代入代价函数，根据拉格朗日乘子法，在选定区域内搜索代价函数第一个元素的极大值，极大值对应的坐标即为辐射源的坐标估计值。
[0064]
本发明具有以下有益效果：
[0065]
本发明适用于由多个单天线观测站组成的分布式监测节点定位场景。
[0066]
本发明实现的具体步骤如下：首先将接收信号分段并进行不同频带下的数据预处理，得到不同频带下的噪声子空间；之后利用接收信号多个频带中存在的正交关系给出各自频率下的代价函数，并基于共同的衰减系数进一步得到多频率融合的代价函数；由于信号的衰减系数未知，因此结合拉格朗日乘子法，进行参数分离估计，通过空间搜索最终可以得到辐射源位置的估计。
[0067]
本发明使用单天线的分布式监测节点取代常规的阵列天线，有效节约了花费的成本；将宽带信号转换成多频带窄带信号的集合，并通过数据分段的方式，大大减少了矩阵运算的复杂度，节约了时间成本，提高了定位效率；相比于传统的两步定位法，本发明使用的直接定位方法避免了参数估计，减少了定位误差；不同于大多数直接定位方法基于窄带信号源的假设，本发明可以应用于单个或多个宽带信号源的定位。
[0068]
本发明针对目前越来越广泛应用的宽带信号，通过信号分段和频带划分相关处理，有效减少算法的计算复杂度，同时相比于传统时差定位和直接定位可获得更好的定位性能，相关方面经过了仿真和实测验证；此外，较于基于阵列天线的直接定位方法，本发明使用单天线的监测站，有利于降低定位成本。
附图说明
[0069]
图1为本发明提供的未知衰减系数下基于多频率代价函数融合的宽带信号直接定位方法流程图；
[0070]
图2为步骤2中的数据频带划分示意图；
[0071]
图3为本发明方法与传统直接定位及tdoa方法针对仿真数据的定位误差比较图；
[0072]
图4为本发明提供的实测数据场景图；
[0073]
图5为本发明方法在实测数据下定位区域的热力图；
[0074]
图6为本发明方法与传统直接定位及tdoa方法针对实测数据的定位误差cdf图。
具体实施方式
[0075]
以下结合附图对本发明的实施例作进一步详细描述。
[0076]
本发明提供的一种未知衰减系数下基于多频率代价函数融合的宽带信号直接定位方法的详细流程如图1所示：
[0077]
具体实现如下：
[0078]
步骤1：分布式监测节点接收辐射源宽带信号并进行采样，将采样后的宽带信号传输至中心站点；
[0079]
(11)假设共有m个未知位置的辐射源，且坐标分别为：p1(x1,y1,z1)，
…
，pm(xm,ym,zm)，分布式监测节点数量为l，坐标分别为s1(x
s1
,y
s1
,z
s1
)，
…
，s
l
(x
sl
,y
sl
,z
sl
)；
[0080]
设采样率为fs，则各监测节点的接收信号表示为：
[0081][0082]
其中，sm(t)表示第m个辐射源在t时刻的发射信号，α
lm
表示传播路径上的信号衰减；
[0083]nl
(t)为相应节点接收信号时产生的噪声，考虑为零均值高斯噪声；
[0084]
τ
lm
是表示信号从辐射源传输到监测站的时延，其表达式是：
[0085][0086]
其中，c表示电磁波传播的速度；
[0087]
(12)对接收信号进行采样，采样间隔为t＝1/fs，设t＝nst，则接收信号的离散形式表示为：
[0088][0089]
其中，ns表示信号的采样点数。
[0090]
步骤2：中心站点将所有接收数据整合，将较长的接收数据分段后进行数据预处理，得到宽带信号中每个频率对应的窄带信号：
[0091]
中心站点将接收信号分段并进行数据预处理，目的是将宽带信号转换为若干个不同频率的窄带信号的集合；
[0092]
(21)中心站点对所有监测节点的信号进行接收与整合，得到行向量对应监测节点接收信号的矩阵，该矩阵的列数为ns，将其按列平均划分为k个子矩阵，则每个子矩阵列数为n＝ns/k，此时较长的宽带信号就被划分为k段长度较短的信号，便于后续计算；
[0093]
把所有观测站收到的信号整合起来，得到合成信号的表达式：
[0094][0095]
将矩阵y(ns)平均划分为j个子矩阵，每个子矩阵的列数为n＝ns/j，则分段后的信号表示为
[0096]
y(n)＝[y1(n) y2(n)
ꢀ…ꢀ
yj(n)],1≤n≤n
[0097]
其中，式中，其中，s
mj
(n)＝sm(ns),(j-1)n 1≤ns≤jn；
[0098]
(22)对每个子矩阵的行向量，即分段后的接收信号进行频域变换，分离出信号中的时延和频率信息，此时分段前后信号频谱包络应大致相同，且每段信号由n个频率点的信息组成；
[0099]
数据分段后，对yj(n)进行频域变换，得到
[0100][0101]
其中，y
lj
(k)是y
lj
(n)频域形式，且表示发射信号中对应的分段；
[0102]smj
(k)表示s
mj
(n)的dft，nj(k)是分段后噪声的dft；
[0103]
τ
lm
是辐射源位置的函数，因此，它可以被表示为τ
l
(pm)；
[0104]
其中，pm表示第m个辐射源的位置；
[0105]
则重写为：
[0106][0107]
其中，矩阵φ(α)包含所有的衰减系数，a(ωk,p)是包含位置信息的时延因子，其表达式分别为：
[0108]
[0109][0110]
yj(k)包含n组不同频率的信息。
[0111]
(23)为了避免直接处理宽带信号，将同一个监测节点的分段信号相同频率对应的频谱信息按序取出，组成一个该频率下的窄带数据，此时原本的宽带信号被划分为n个长度为j的窄带信号的集合，每个频率对应信号的矩阵维数为l
×
j。
[0112]
信号的频带划分方式如图2所示。定义第n个信号为xj(j)。此时，j变成了窄带信号的长度。在现实中，接收到的宽带信号在某些频率上的频谱并不理想，在这种情况下，可以选择一些频谱较高的频率作为处理的信号。
[0113]
步骤3：对每个窄带信号的协方差矩阵进行特征分解，得到其噪声子空间，根据多频带中的正交关系构造出一个新的代价函数；
[0114]
通过组合不同频率下的代价函数，可以得到多频带融合的新的代价函数；
[0115]
(31)求出所有窄带信号对应的协方差矩阵，并对其进行特征分解，得到其噪声子空间；对的协方差矩阵进行特征分解，得到：
[0116][0117]
其中，表示第n个窄带信号的信号子空间，表示噪声子空间。
[0118]
(32)由噪声子空间和节点与辐射源方向流形的正交关系，有
[0119][0120]
其中表示方向向量；
[0121]
表示目标坐标组成的向量。
[0122]
包含各观测站的所有衰减系数。
[0123]
由于同一宽带发射信号不同频率的窄带信号对应的未知衰减系数相同，即它们有相同的因此，可以将所有频率对应的等式合成为下式
[0124][0125]
由上式可以定义新的代价函数
[0126]
[0127]
步骤4：选取定位区域和搜索网格大小，采用代价函数，在衰减系数未知情况下结合拉格朗日乘子法，在选定区域内进行网格化的谱峰搜索，最终得到辐射源坐标估计值。
[0128]
将选定区域网格化，求取选定区域内网格点的代价函数值并进行谱峰搜索，得到辐射源位置坐标的估计。
[0129]
首先根据代价函数和衰减系数的关系，可以得到但是由于衰减系数未知，所以利用拉格朗日乘子法解决这个问题。
[0130]
其次选取定位区域，并将其网格化，将网格中每个点的位置代入搜索的代价函数。
[0131]
假设未知的的第一个元素为1，则由拉格朗日乘子法，可以通过寻找的第一行和第一列元素中绝对值最大的m个来得到辐射源的坐标估计。
[0132]
图3为本发明所述方法与传统直接定位及tdoa方法针对仿真数据的定位误差cdf比较图。其中，仿真信号的采样率为fs＝125mhz，lfm仿真信号的时域宽度为t0＝0.1ms；频带宽度b＝20mhz，调频斜率k＝b/t0，在信噪比为-15db的情况下，进行40次仿真，分别用本发明方法和传统直接定位及两步定位中的tdoa方法对仿真数据进行处理，信号分段的段数为320，dft点数为32，取其中15个幅度值较大的频带，将累计分布函数cdf作为误差评价指标。从图3中可以看出，本发明方法的定位精度优于另外两种定位方法。
[0133]
图4为本发明提供的实测数据场景图，其中图(a)、(b)分别为辐射源和监测节点的实物图，图(c)为4个分布式监测节点和单辐射源的分布图。本发明的实测数据来源即通过这几个监测节点接收的辐射源信号。
[0134]
图5为本发明方法在实测数据下定位区域的热力图。实测数据采样率为fs＝125mhz，信号长度为32508，dft点数为256，则分段数为126，选取频谱中幅度值较大的15个频率对应的信号进行定位，得到定位区域的热力图。图5中，整体区域是选取的400m
×
400m的定位范围，网格边长取为1m。图5中的target区域表示热力值较高的区域，即估计的目标位置。从图中可以看出，对于实测数据而言，本发明方法的定位结果较为准确。
[0135]
图6为本发明方法与传统直接定位及tdoa方法针对实测数据的定位误差cdf图。其中，选取200m
×
200m的定位范围，网格边长为1m，dft点数为128，选取频谱中幅度值较大的9个频率对应的信号进行定位对采集的20组实测数据进行估计和误差计算，图6中，横轴表示定位结果的误差值范围，纵轴为cdf值。从图6中可以看出，本发明的方法产生的误差明显小于另外两种传统算法，体现了本发明方法在估计精度上的优越性，以及在实际情况中的良好应用性。
[0136]
以上仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，应视为本发明的保护范围。