本发明涉及一种近场信号源定位技术，尤其是涉及一种基于累积量域逐级解耦的空时近场源三维参数估计方法。
背景技术：
1、信号源的doa(direction of arrival，到达方向)定位问题是阵列信号处理领域的一个重要课题，该定位问题在无线通信、雷达、声纳等诸多领域均有广泛的应用。doa估计根据信号源与阵列之间的距离，可分为近场和远场测向两种情况。远场测向情况中，信号源离参考阵元的距离满足r＞＞2d2/λ，其中，r表示信号源距离参考阵元的距离，d表示阵列孔径，λ为信号波长，此时信号传播可视为平面波，定位问题与入射角度有关。然而，随着阵列孔径的增大以及信号频率的提升，瑞利距离2d2/λ增大，在信号源离参考阵元的距离满足的情况下，信号源将处于阵列的近场区域，若继续采用基于远场的测向方法，则参数估计的性能会降低。因此，需要对近场区域的阵列信号模型进行重新建模并设计新的算法。此时，信号传播不可以再近似为平面波，而是球面波，定位问题与角度和距离均有关。
2、已有的大部分近场测向算法并未考虑信号幅度衰减的影响，当信号源距离阵列中的传感器较近时，不同传感器接收到的信号幅度往往差异比较大，若忽略阵元接收信号时信号幅度衰减差异的影响，则可能会导致信号模型的失配，从而降低参数估计的性能。已经有许多近场信号源定位的算法被提出，如极大似然法、二阶统计量方法，但是上述方法大多只考虑了近场源的二维参数，即方位角和距离，然而，对于三维场景下的近场信号源的doa定位问题，需联合估计三维参数，即方位角、俯仰角和距离。并且，已有的大部分近场信号源的doa估计算法使用的是信号的二阶统计特性，然而使用高阶累积量可以获得接收信号的更多信息，并且能够很好地抑制高斯噪声，获得更好的参数估计性能。因此，考虑信号幅度衰减的影响，研究累积量域近场信号源三维参数估计的方法势必关键。
技术实现思路
1、本发明所要解决的技术问题是提供一种基于累积量域逐级解耦的空时近场源三维参数估计方法，其考虑信号幅度衰减的影响，并使用高阶累积量，实现三维参数估计，且估计精度更高。
2、本发明解决上述技术问题所采用的技术方案为：一种基于累积量域逐级解耦的空时近场源三维参数估计方法，其特征在于包括如下步骤：
3、步骤1：在三维坐标系上建立一个十字交叉阵列：该十字交叉阵列作为信号接收阵列，其包括位于三维坐标系的x轴上的均匀线阵和位于三维坐标系的y轴上的均匀线阵，位于x轴上的均匀线阵由mx个阵元均匀分布组成，位于y轴上的均匀线阵由my个阵元均匀分布组成，两个均匀线阵共享原点处的中心阵元，并将该中心阵元作为十字交叉阵列的参考阵元，每个均匀线阵以参考阵元为中心对称，每个均匀线阵中相邻两个阵元之间的间距为d；其中，mx表示位于x轴上的均匀线阵的总阵元数，mx＝2nx+1，nx表示位于x轴上的均匀线阵的单边阵元数，my表示位于y轴上的均匀线阵的总阵元数，my＝2ny+1，ny表示位于y轴上的均匀线阵的单边阵元数，十字交叉阵列的总阵元数为m，m＝mx+my-1，d＝λ/8，λ表示信号波长，阵元为传感器；
4、步骤2：设定存在k个空时不相关的窄带近场源同时入射信号到十字交叉阵列中的各个阵元上，每个阵元接收到达其上的所有信号；将k个空时不相关的窄带近场源各自的三维参数分别表示为(θ1，ψ1，r0，1)，…，(θk，ψk，r0，k)，…，(θk，ψk，r0，k)，将位于x轴上的均匀线阵中位置为(m，0)的阵元接收到的信号记为xm(t)，
5、
6、将位于y轴上的均匀线阵中位置为(0，n)的阵元接收到的信号记为yn(t)，
7、
8、其中，nx≥k且ny≥k，k≥1，1≤k≤k，θk表示由第k个窄带近场源入射的信号的俯仰角，ψk表示由第k个窄带近场源入射的信号的方位角，r0，k表示第k个窄带近场源与参考阵元的中心之间的距离，m∈{-nx，-nx+1，…，nx-1，nx}，n∈{-ny，-ny+1，…，ny-1，ny}，t为时间变量，rm，k表示第k个窄带近场源与位置为(m，0)的阵元的中心之间的距离，αk表示由第k个窄带近场源入射的信号与x轴的夹角，αk与θk、ψk的关系为cos(αk)＝sin(θk)cos(ψk)，rn，k表示第k个窄带近场源与位置为(0，n)的阵元的中心之间的距离，βk表示由第k个窄带近场源入射的信号与y轴的夹角，βk与θk、ψk的关系为sin(βk)＝sin(θk)sin(ψk)，p表示先验已知的信号幅度衰减指数值，e表示自然常数，j表示虚数单位，sk(t)表示由第k个窄带近场源入射的信号被任一个阵元接收并解调的基带信号，wm(t)表示位置为(m，0)的阵元上的加性高斯白噪声，wn(t)表示位置为(0，n)的阵元上的加性高斯白噪声，γm，k和γn，k均为引入的中间量，
9、然后对公式(1)中的相移rm，k-r0，k采取菲涅尔近似，将公式(1)改写为：
10、
11、同样，对公式(2)中的相移rn，k-r0，k采取菲涅尔近似，将公式(2)改写为：
12、
13、其中，
14、
15、再在公式(3)的基础上，将位于x轴上的均匀线阵中的所有阵元的接收信号向量记为x(t)，x(t)＝axs(t)+wx(t)；同样，在公式(4)的基础上，将位于y轴上的均匀线阵中的所有阵元的接收信号向量记为y(t)，y(t)＝ays(t)+wy(t)；其中，ax表示位于x轴上的均匀线阵的阵列流形矩阵，ax＝[ax(α1，r0，1)，...，ax(αk，r0，k)，...，ax(αk，r0，k)]，ax(αk，r0，k)表示由第k个窄带近场源入射的信号与x轴所形成的导向矢量，[]t表示矩阵或向量的转置，ay表示位于y轴上的均匀线阵的阵列流形矩阵，ay＝[ay(β1，r0，1)，...，ay(βk，r0，k)，...，ay(βk，r0，k)]，ay(βk，r0，k)表示由第k个窄带近场源入射的信号与y轴所形成的导向矢量，s(t)表示基带信号向量，s(t)＝[s1(t)，...，sk(t)，...，sk(t)]t，wx(t)表示位于x轴上的均匀线阵上的噪声向量，wy(t)表示位于y轴上的均匀线阵上的噪声向量，
16、步骤3：利用均匀线阵的对称性，构造与m，n有关且位置为(-m，0)的阵元和位置为(0，-n)的阵元接收到的信号延时τ后的四阶累积量，记为c(m，n，τ)，
17、
18、进一步获得延时τ后的接收信号四阶累积量矩阵，记为c(τ)，
19、
20、并将c(τ)表示为然后在x(t)和y(t)的长度均为j＝n+l的前提下，对x(t)和y(t)采用最大似然估计，得到延时τ＝lts后的虚拟延时接收数据，记为r(τ)，
21、
22、r(τ)是对c(τ)的最大似然近似；再通过l次延时为τ＝ts，2ts，...，lts的采样，得到基于延时的虚拟伪快拍接收数据矩阵，记为r，r＝[vec(r(ts))，vec(r(2ts))，...，vec(r(lts))]；其中，τ为延时变量，表示求xm(t)的累积量，上标“*”表示取共轭，x-m(t+τ)表示位于x轴上的均匀线阵中位置为(-m，0)的阵元接收到的信号x-m(t)延时τ后得到的信号，y-n(t+τ)表示位于y轴上的均匀线阵中位置为(0，-n)的阵元接收到的信号y-n(t)延时τ后得到的信号，e[]表示求数学期望，表示sk(t)延时τ后的四阶累积量，表示求sk(t)，的累积量，sk(t+τ)表示sk(t)延时τ后得到的信号，r-m，k表示第k个窄带近场源与位置为(-m，0)的阵元的中心之间的距离，r-n，k表示第k个窄带近场源与位置为(0，-n)的阵元的中心之间的距离，和均为引入的中间量，cs(τ)表示延时τ后的基带信号四阶累积量矩阵，“⊙”表示khatri-rao积，n表示快拍数，l表示伪快拍数，l＝1，2，...，l，ts表示抽样时间，表示mx阶的单位矩阵，表示my阶的单位矩阵，上标“h”表示共轭转置，x(t+lts)表示x(t)延时lts后得到的信号向量，y(t+lts)表示y(t)延时lts后得到的信号向量，vec(·)表示对矩阵进行向量化操作；
23、步骤4：计算r的协方差矩阵z，z＝e(rrh)＝rrh/l；然后对z进行特征分解，得到其中，usignal表示z的信号子空间，对z进行特征分解得到的所有特征值按从大到小进行排序，usignal由前k个最大特征值对应的特征向量组成，∑signal为信号对角阵，unoise表示z的噪声子空间，unoise由除前k个最大特征值外的其余所有特征值对应的特征向量组成，∑noise为噪声对角阵；
24、步骤5：根据γ-m，k＝γm，k的对称性，将由第k个窄带近场源入射的信号与x轴所形成的导向矢量ax(αk，r0，k)表示为：
25、
26、
27、同样，根据γ-n，k＝γn，k的对称性，将由第k个窄带近场源入射的信号与y轴所形成的导向矢量ay(βk，r0，k)表示为：
28、
29、步骤6：令α，β，r对应表示由k个窄带近场源入射的信号与x轴的夹角组成的向量变量、由k个窄带近场源入射的信号与y轴的夹角组成的向量变量、k个窄带近场源与参考阵元的中心之间的距离组成的向量变量；根据降秩music算法，依次构造β，α，r各自对应的谱峰搜索器，并对谱峰搜索器进行一维搜索，依次获得β，α，r各自的估计值，具体过程如下：
30、步骤6.1：获取β的估计值：令然后对进行分解，得到：
31、接着根据降秩music算法，构造β对应的谱峰搜索器，表示为：再对β对应的谱峰搜索器进行一维搜索，选择最大的k个谱峰值对应的角度作为β的估计值，记为其中，表示kronecker积，ay(βk，r0，k)中的参数对应替换为β，r得到ay(β，r)，ax(αk，r0，k)中的参数对应替换为α，r得到ax(α，r)，ξy(βk)中的参数替换为β得到ξy(β)，υy(βk，r0，k)中的参数对应替换为β，r得到υy(β，r)，ξx(αk)中的参数替换为α得到ξx(α)，υx(αk，r0，k)中的参数对应替换为α，r得到υx(α，r)，p(β)表示β对应的谱峰搜索器，det()表示取行列式；
32、步骤6.2：获取α的估计值：对进一步分解，得到：；然后根据降秩music算法，并结合步骤6.1获取的构造αk，k＝1，2，…，k对应的谱峰搜索器，表示为：再对αk，k＝1，2，…，k对应的谱峰搜索器进行一维搜索，选择最大的谱峰值对应的角度作为αk，k＝1，2，…，k的估计值，记为即得到α的估计值，记为其中，p(αk)表示αk对应的谱峰搜索器，g(α)中的参数替换为αk得到g(αk)，b(β)中的参数替换为得到
33、步骤6.3：获取r的估计值：根据降秩music算法，并结合步骤6.1获取的和步骤6.2获取的构造r0,k，k＝1，2，…，k对应的谱峰搜索器，表示为：；再对r0,k，k＝1，2，…，k对应的谱峰搜索器进行一维搜索，选择最大的谱峰值对应的距离值作为r0,k，k＝1，2，…，k的估计值，记为即得到r的估计值，记为其中，p(r0,k)表示r0，k对应的谱峰搜索器，中的参数对应替换为r0，k得到
34、步骤7：根据cos(αk)＝sin(θk)cos(ψk)和sin(βk)＝sin(θk)sin(ψk)，以及和计算得到θ1，...，θk，...，θk各自的估计值，记为同时计算得到ψ1，...，ψk，...，ψk各自的估计值，记为进而得到k个窄带近场源各自的三维参数的估计值，分别为
35、与现有技术相比，本发明的优点在于：
36、1)本发明方法利用十字交叉阵列接收信号，并通过菲涅尔近似处理相移；考虑信号幅度衰减，联合十字交叉阵列接收到的信号，构造四阶累积量矩阵，并通过逐级解耦将三维参数估计问题分解为一维参数估计问题，首先，将入射的信号的俯仰角和方位角分解为与x轴和y轴的夹角，再分别采用降秩music算法构造的角度谱峰搜索器进行估计；随后，构造距离谱峰搜索器估计距离；最终，根据入射的信号与x轴和y轴的夹角估计值计算俯仰角和方位角，实现三维参数估计。本发明方法仅需多次一维搜索，计算复杂度低，并且具有良好的参数估计性能。
37、2)通过大量蒙特卡洛仿真，验证了本发明方法能够有效估计近场源入射的信号的俯仰角、方位角，以及近场源与参考阵元的中心之间的距离，并且二维角度和一维距离参数能够自动匹配。
38、3)相比于已有的近场源参数估计方法，本发明方法采用的四阶累积量构造能够有效抑制高斯白噪声的影响，估计性能更优。
1.一种基于累积量域逐级解耦的空时近场源三维参数估计方法，其特征在于包括如下步骤：
