一种基于空时信息联合的近场源三维参数欠定估计方法

本发明涉及一种近场信号源定位技术，尤其是涉及一种基于空时信息联合的近场源三维参数欠定估计方法。

背景技术：

1、近场信号源定位在雷达、通信、无线电检测等领域中应用场景不断扩大。传统的高分辨测向算法通常基于远场信号源平面波模型，当信号源位于阵列的近场区域时，阵元接收信号的相位差异是关于阵元间距的非线性函数，远场信号源定位算法无法直接应用，需要进行菲涅尔近似处理，定位精度受到影响。对于近场信号源定位，不同阵元对同一信号的接收数据存在与距离相关的幅度衰减，忽略幅度衰减会导致估计性能的下降甚至无法满足测向精度的需要。现有的大多数近场信号源测向算法，或将时延进行近似处理，或将幅度衰减视作相等处理，算法的计算复杂度有所降低但其测向精度受到限制。将幅度衰减视作相等处理的如经典子空间类算法music算法，其直接对物理阵元接收数据进行协方差计算，通过特征值分解得到噪声子空间，再利用信号子空间与噪声子空间正交的性质对近场信号源参数进行估计，然而测向精度不理想。利用接收数据的时间信息可进行虚拟阵元扩展并提高算法估计性能，基于此，本发明提出一种基于空时信息联合的近场源三维参数欠定估计方法。

技术实现思路

1、本发明所要解决的技术问题是提供一种基于空时信息联合的近场源三维参数欠定估计方法，其结合空时信息计算四阶延时互累积量，通过矢量化操作进行虚拟阵元扩展，计算虚拟接收数据矩阵的协方差并提取噪声子空间特征向量矩阵，利用信号子空间与噪声子空间正交的性质，对伪功率函数进行三维谱峰搜索实现欠定估计近场信号源的方向角、俯仰角、同参考阵元的距离，且估计精度高。

2、本发明解决上述技术问题所采用的技术方案为：一种基于空时信息联合的近场源三维参数欠定估计方法，其特征在于包括如下步骤：

3、步骤1：建立基于l型阵列的三维近场信号源空间传播模型：在三维坐标系的x轴上布置一个由n个阵元构成的阵列作为x阵列，在三维坐标系的y轴上布置一个由m个阵元构成的阵列作为y阵列，x阵列和y阵列共用三维坐标系的原点处阵元构成l型阵列；x阵列为嵌套阵列，x阵列由两个均匀线性阵列构成，其中一个均匀线性阵列的阵元个数为lx且阵元间距为d，另一个均匀线性阵列的阵元个数为n-lx且阵元间距为(lx+1)d；y阵列也为嵌套阵列，y阵列也由两个均匀线性阵列构成，其中一个均匀线性阵列的阵元个数为ly且阵元间距为d，另一个均匀线性阵列的阵元个数为m-ly且阵元间距为(ly+1)d；用dxn表示x阵列的第n个阵元的位置，x阵列的第n∈[1,lx]个阵元的位置为dxn＝(n-1)d，x阵列的第n∈[lx+1,n]个阵元的位置为dxn＝(n-lx)(lx+1)d-d；用dym表示y阵列的第m个阵元的位置，y阵列的第m∈[1,ly]个阵元的位置为dym＝(m-1)d，y阵列的第m∈[ly+1,m]个阵元的位置为dym＝(m-ly)(ly+1)d-d；有k个近场信号源入射到l型阵列中，第k个近场信号源以方向角θk和俯仰角φk、从距离原点rk的位置入射，第k个近场信号源的坐标位置表示为(rk sinφk cosθk,rk sinφk sinθk,rk cosφk)；

4、其中，n＞1，m＞1，l型阵列的阵元个数为n+m-1，n＝1,2,…,n，m＝1,2,…,m，1＜lx＜n，1＜ly＜m，d取λ/2，λ为信号波长，k≥1，k＝1,2,…,k，θk∈(-180°,180°]，φk∈[0°,90°]，

5、步骤2：用rxn,k表示第k个近场信号源到x阵列的第n个阵元的距离，用rym,k表示第k个近场信号源到y阵列的第m个阵元的距离，然后考虑幅度衰减，计算x阵列的第n个阵元接收到第k个近场信号源发出的信号的幅度相位因子，记为axn(θk,φk,rk)，同样，计算y阵列的第m个阵元接收到第k个近场信号源发出的信号的幅度相位因子，记为aym(θk,φk,rk)，再用x(t)表示x阵列的所有阵元接收到的数据，用y(t)表示y阵列的所有阵元接收到的数据，x(t)＝axs(t)+nx(t)，y(t)＝ays(t)+ny(t)；

6、其中，g表示已知的路径损耗衰减指数，e表示自然常数，j为虚部表示，表示x阵列的第n个阵元接收到第k个近场信号源发出的信号相对于参考阵元即原点处阵元接收到第k个近场信号源发出的信号的衰减，表示x阵列的第n个阵元接收到第k个近场信号源发出的信号相对于参考阵元即原点处阵元接收到第k个近场信号源发出的信号的相位差，表示y阵列的第m个阵元接收到第k个近场信号源发出的信号相对于参考阵元即原点处阵元接收到第k个近场信号源发出的信号的衰减，表示y阵列的第m个阵元接收到第k个近场信号源发出的信号相对于参考阵元即原点处阵元接收到第k个近场信号源发出的信号的相位差，x(t)和y(t)均为列向量，t表示时间变量，ax表示x阵列的阵列流型矩阵，ax＝[ax(θ1,φ1,r1),ax(θ2,φ2,r2),…,ax(θk,φk,rk)]，ax的维数为n×k，ax(θk,φk,rk)＝[ax1(θk,φk,rk),ax2(θk,φk,rk),…,axn(θk,φk,rk)]t，ax(θk,φk,rk)的维数为n×1，上标“t”表示向量或矩阵的转置，s(t)表示近场信号源向量，s(t)＝[s1(t),s2(t),…,sk(t)]t，s(t)的维数为k×1，sk(t)表示第k个近场信号源发出的信号，nx(t)表示t时刻附加在x阵列的所有阵元接收到的数据上的独立同分布加性噪声矢量，nx(t)的维数为n×1，ay表示y阵列的阵列流型矩阵，ay＝[ay(θ1,φ1,r1),ay(θ2,φ2,r2),…,ay(θk,φk,rk)]，ay的维数为m×k，ay(θk,φk,rk)＝[ay1(θk,φk,rk),ay2(θk,φk,rk),…,aym(θk,φk,rk)]t，ay(θk,φk,rk)的维数为m×1，ny(t)表示t时刻附加在y阵列的所有阵元接收到的数据上的独立同分布加性噪声矢量，ny(t)的维数为m×1；

7、步骤3：基于步骤2中的x(t)和y(t)，计算y阵列的第m个阵元接收到的数据、x阵列的第n个阵元接收到的数据、x阵列的所有阵元接收到的数据与y阵列的所有阵元接收到的数据的四阶延时互累积量，记为cmn(τ)，

8、

9、其中，cmn(τ)的维数为n×m，τ表示时延，cum{·}为累积量运算，(·)*为共轭运算，(·)h为共轭转置运算，ym(t)表示y阵列的第m个阵元接收到的数据，ym(t+τ)表示y阵列的第m个阵元经时延τ后接收到的数据，xn(t)表示x阵列的第n个阵元接收到的数据，x(t+τ)表示x阵列的所有阵元经时延τ后接收到的数据，sk(t+τ)表示对第k个近场信号源发出的信号经时延τ处理后的信号，，diag{·}表示将一个向量转换为以该向量为主对角线元素的对角矩阵，aym(θ1,φ1,r1)、aym(θ2,φ2,r2)、aym(θk,φk,rk)γ(τ)＝diag{γ1(τ),γ2(τ),…,γk(τ)}，γk(τ)表示第k个近场信号源发出的信号与经时延τ处理后的信号的累积量，

10、步骤4：对步骤3得到的cmn(τ)进行矢量化操作，得到cmn(τ)的矢量化结果，记为cmn(τ)，

11、其中，cmn(τ)的维数为(n×m)×1，vec(·)为矢量化运算，⊙为khatri-rao积运算，b表示导向矩阵，b的维数为(n×m)×k，b的第k列为bk，为kronecker积运算，

12、步骤5：基于步骤4得到的cmn(τ)，m从1至m依次取整数值，n从1至n依次取整数值，排列得到矩阵c(τ)，c(τ)＝[c11(τ),…,c1n(τ),…,c21(τ),…,cmn(τ)]；然后根据c(τ)中的元素的排列，将c(τ)表示为c(τ)＝bγ(τ)bh；

13、其中，c(τ)的维数为(n×m)×(n×m)；

14、步骤6：对步骤5得到的c(τ)进行矢量化操作，将c(τ)的矢量化结果作为虚拟阵列接收数据，并记为w(τ)，w(τ)＝vec(c(τ))＝(b*⊙b)γ(τ)＝dγ(τ)；

15、其中，w(τ)的维数为(n×m)2×1，d表示虚拟导向矩阵，d＝b*⊙b，d的维数为(n×m)2×k，d的第k列为dk，γ(τ)＝[γ1(τ),γ2(τ),…,γk(τ)]t；

16、步骤7：基于步骤6得到的w(τ)，取τ＝ts,2ts,…,npts，排列得到虚拟接收数据矩阵w，w＝[w(ts),w(2ts),…,w(npts)]＝dγs；

17、其中，ts表示虚拟采样周期，np表示伪快拍的数量，np＞1，w的维数为(n×m)2×np，γs＝[γ(ts),γ(2ts),…,γ(npts)]；

18、步骤8：计算步骤7得到的w的协方差，记为r，然后通过对r进行特征值分解，得到所有特征值；再对所有特征值从小到大排序，取前(n×m)2-k个特征值，将这些特征值对应的特征向量构成的矩阵作为噪声子空间特征向量矩阵en；

19、其中，r的维数为(n×m)2×(n×m)2，r的秩为k，λss＝γsγsh，es表示信号子空间特征向量矩阵，σs表示信号子空间特征值对角矩阵，en表示噪声子空间特征向量矩阵，σn表示噪声子空间特征值对角矩阵；

20、步骤9：基于步骤8得到的en，并利用信号子空间与噪声子空间正交的性质，定义伪功率函数p(θ,φ,r)，然后利用三维谱峰搜索，在区间(-180°,180°]内搜索θ的取值，在区间[0°,90°]内搜索φ的取值，在区间内搜索r的取值，p(θ,φ,r)存在k个极大值，第k个极大值对应的θ的取值、φ的取值、r的取值的组合即为第k个近场信号源的三维位置参数；

21、其中，θ表示方向角变量，φ表示俯仰角变量，r表示与原点的距离变量，ax(θ,φ,r)为ax(θk,φk,rk)中的θk,φk,rk对应替换为θ,φ,r，ay(θ,φ,r)为ay(θk,φk,rk)中的θk,φk,rk对应替换为θ,φ,r。

22、与现有技术相比，本发明的优点在于：

23、采用l型阵列，计算y阵列的第m个阵元接收到的数据、x阵列的第n个阵元接收到的数据、x阵列的所有阵元接收到的数据与y阵列的所有阵元接收到的数据的四阶延时互累积量，通过两次矢量化操作得到虚拟接收数据矩阵，计算虚拟接收数据矩阵的协方差并进行特征值分解提取噪声子空间特征向量矩阵，经过伪功率函数的三维谱峰搜索得到近场信号源的三维位置参数，同时，在参数估计过程中，利用数据互累积量和空时特性进行虚拟阵元扩展，能够实现近场信号源的欠定估计。与直接使用物理阵元接收数据进行估计相比，经过本发明方法虚拟阵元扩展后，具备更高的测向精度且能够实现近场信号源的欠定估计。