CN109270483A - 三维阵虚拟扩展相干源二维波达方向估计方法 - Google Patents

Abstract

三维阵虚拟扩展相干源二维波达方向估计方法，利用均匀三维阵列接收K个同频率的远场相干窄带信号，x、y和z轴方向的子阵接收数据相关矩阵特征分解得到最大特征值对应的特征矢量，利用圆信号的特性将阵列进行虚拟扩展得到x、y和z轴方向的经过Toeplitz修正后的协方差矩阵矩阵，从而实现了解相干；分别利用ESPRIT算法得到x、y和z轴方向的方向余弦矩阵，根据同一信号x、y和z轴方向的方向余弦的等量关系进行配对运算，利用配对后的x、y和z轴方向的方向余弦得到二维到达角估计结果，本发明方法克服了空间平滑算法的孔径损失，在相同快拍数下本发明算法的性能明显优于空间平滑算法，本发明方法提出的配对算法简单有效且节省了近似一半的阵元。

Description

三维阵虚拟扩展相干源二维波达方向估计方法

技术领域

本发明属于信号处理技术领域，尤其涉及一种三维虚拟扩展阵解相干二维波达方向估计方法。

背景技术

波达方向(Direction of arrival，DOA)估计是阵列信号处理的重要方向之一，广泛应用于通信、雷达、声呐和地震探测中。在实际环境中，由于多径传播与人为干扰的存在，导致入射信号中存在相干信号，相干信号的存在使得数据相关矩阵的秩出现亏损，从而导致基于数据相关矩阵特征分解的多重分类算法(Multiple Signal Classification，MUSIC算法)、空间旋转不变算法 (Estimating Signal Parameter via Rotational InvarianceTechniques ESPRIT)为代表的经典超分辨DOA估计算法不能有效估计相干信源的波达方向。针对相干信源的DOA估计问题，研究者也提出了许多解相干算法，以空间平滑类算法和基于矢量重构算法为代表的降维类解相干算法，虽然可以解相干，但是其解相干的性能是通过降低自由度即矩阵的维数获得的，所以阵列的孔径没有得到有效利用。以Toeplitz矩阵重构算法、修正MUSIC算法为代表的非降维类算法，虽然阵列孔径得到了有效利用，但在对接收数据协方差矩阵进行修正时存在较大误差，导致DOA估计误差增大。随着技术发展与实际应用的要求，二维DOA估计方法也开始成为研究的热点。另外，目前的解相干算法主要研究了线阵情况下的一维DOA估计，在实际通信中，一维DOA估计只适合于宏小区情况，当用户处于微小区或者微微小区内，如果仅仅知道用户的方位角，还是不能精确的进行定向波束形成，显然还要知道用户的仰角，利用二维DOA 能够更有效地开发信道的冗余度，可以精确地定位。因此，准确地估计各个用户在三维空间内的二维DOA也是一个需要深入研究的方向。本发明在Toeplitz 矩阵重构算法的基础上，结合虚拟扩展的思想，研究了基于阵列虚拟扩展的相干源二维DOA估计算法，将均匀分布在x轴、y轴和z轴正半周方向阵的三维阵虚拟扩展为均匀分布在x轴、y轴和z轴的对称阵，在三个坐标轴方向上分别构造Toeplitz矩阵，并分别进行特征分解得到x轴、y轴和z轴方向的方向余弦估计，通过对x轴、y轴和z轴方向的方向余弦进行配对运算得到到达角的估计。

发明内容

本发明的目的是提供一种三维阵虚拟扩展相干源二维波达方向估计方法。

为了实现上述目的，本发明采取如下的技术解决方案：

三维阵虚拟扩展相干源二维波达方向估计方法，均匀分布在坐标原点及x 轴、y轴和z轴正半轴方向的三维阵接收K个相干源窄带圆对称信号，所述阵列由坐标原点处的阵元以及M个等间隔布置于x轴正半轴，M个等间隔布置于y轴正半轴和M个等间隔布置于z轴正半轴的3M+1个阵元构成的均匀三维阵，阵元间隔分别为d_x，d_y和d_z，d_x≤λ/2，d_y≤λ/2和d_z≤λ/2，λ为入射信号的波长，坐标原点的阵元三轴共用，所以x轴子阵、y轴子阵和z轴子阵都是由M+1 个阵元构成的；

参数估计方法的步骤如下：

步骤一、利用3M+1个阵元构成的均匀三维阵作为接收阵列，获取阵列的P 次快拍接收数据，通过数据相关矩阵特征分解获取最大特征值对应的特征矢量利用特征矢量构造Toeplitz矩阵；

x轴子阵信号为x(t)＝A_x(α)S(t)+N_x(t)，其中A_x(α)是x轴子阵导向矢量矩阵A_x(α)＝[a_x(α₁)，a_x(α₂)，…，a_x(α_K)]，是第k 个信号对应的x轴子阵导向矢量，α_k＝sinθ_kcosφ_k是第k个信号对应的x轴方向的方向余弦，λ是入射信号的波长，y轴子阵信号为y(t)＝A_y(β)S(t)+N_y(t)其中A_y(β)是y轴子阵导向矢量矩阵，A_y(β)＝[a_y(β₁)，a_y(β₂)，…，a_y(β_K)]，是第k个信号对应的y轴子阵导向矢量，β_k＝sinθ_ksinφ_k是第k个信号对应的y轴方向的方向余弦，z轴子阵信号为 z(t)＝A_z(ω)S(t)+N_z(t)其中A_z(ω)是z轴子阵导向矢量矩阵，A_z(ω)＝[a_z(ω₁)，a_z(ω₂)，…，a_z(ω_K)]，是第k个信号对应的z轴子阵导向矢量，ω_k＝cosθ_k是第k个信号对应的z轴方向的方向余弦，x轴子阵、y轴子阵和z轴子阵P次快拍数据分别为Z_x、Z_y和Z_z，利用Z_x获得x轴子阵对应的数据相关矩阵利用Z_y获取y轴子阵对应的数据相关矩阵利用Z_z获取z轴子阵对应的数据相关矩阵对数据相关矩阵R_x、R_y和R_z分别进行特征分解，得到最大特征值对应的特征向量e_sx，e_sy和e_sz；利用特征矢量e_sx，e_sy和e_sz构造Toeplitz矩阵

R_Tx、R_Ty和R_Tz是经过Toeplitz修正后的协方差矩阵，为M×M维的HermitionToeplitz矩阵，因此实现了解相干处理；

步骤二、对R_Tx、R_Ty和R_Tz分别进行特征分解，利用ESPRIT算法得到信号 x轴、y轴方向和z轴方向的方向余弦估计；

对R_Tx进行特征分解得到信号子空间E_sx，E_sx的前M行元素组成矩阵E_x1， E_sx的后M行元素组成矩阵E_x2，令对Ψ_x进行特征分解，由特征值构成其中是矩阵E_x1的伪逆矩阵，是K×K的对角矩阵；对R_Ty进行特征分解得到信号子空间E_sy，E_sy的前M行元素组成矩阵E_y1，E_sy的后M行元素组成矩阵E_y2，令对Ψ_y进行特征分解，特征值构成其中是矩阵E_y1的伪逆矩阵，

是K×K的对角矩阵；对R_Tz进行特征分解得到信号子空间E_sz，E_sz的前M行元素组成矩阵E_z1， E_sz的后M行元素组成矩阵E_z2，令对Ψ_z进行特征分解，由特征值构成其中是矩阵E_z1的伪逆矩阵，Ф_z是K×K的对角矩阵，，

从而得到x轴、y轴和z轴方向的方向余弦方向的估计值：

因为x轴、y轴和z轴方向的三个方向是由三个独立的特征分解完成的，所以信号的排列顺序不同，不一定是同一个信号的对应的三个方向余弦，不能直接用来计算信号的到达角，必须进行配对运算，才能使三个方向余弦的排列顺序一致；

步骤三、根据同一个信号x轴、y轴和z轴方向的三个方向余弦之间的等量关系，对三个方向余弦进行配对运算；

对一个信号(θ_k，φ_k)来说，x轴和y轴方向余弦的和

此时由z轴方向的方向余弦ω_k可得到 sin²(arccos(ω_k))＝sin²θ_k，本发明就是利用同一个信号三个方向余弦之间满足的这个等量关系进行参数配对运算，步骤二中得到的x 轴的方向余弦与K个y轴的方向余弦进行平方求和得

T_s＝|T₁-T_z|，矩阵T_s中元素的最小值对应于配对成功的三个方向余弦分别写作对x轴方向余弦都做同样的配对处理，得到对应的的配对成功的x轴方向余弦y轴方向余弦和z轴三个方向余弦

步骤四、利用配对成功后的方向余弦和得到信号二维到达角的估计值

利用配对成功后的方向余弦从而得到到达角的估计：

前述步骤中的k＝1，...，K，rk∈[1，2，…，K]，pk∈[1，2，…，K]，qk∈[1，2，…，K]，j是虚数单位；

本发明提出了基于虚拟扩展的三维均匀阵列的相干信号二维到达角估计方法，该方法利用正半轴的阵元数据结果推出负半轴的数据结果，通过这种虚拟扩展提高阵元的利用率；利用同一个信号x、y、z轴方向的方向余弦关系进行配对，配对算法简单有效；

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中需要使用的附图做简单介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明三维阵列的结构示意图；

图2为本发明方法的流程图；

图3为本发明方法到达角估计散点图；

图4本空间平滑方法到达角估计散点图

图5为单次快拍情况下本发明方法和空间平滑算法到达角估计均方根误差随信噪比的变化曲线图；

图6为本发明方法的单次快拍和空间平滑算法500次快拍到达角估计均方根误差随信噪比的变化曲线图；

具体实施方式

为了让本发明的上述和其它目的、特征及优点能更明显，下文特举本发明实施例，并配合所附图示，做详细说明如下。

参照图1，本发明的二维到达角估计方法的步骤如下：

3M+1三维均匀阵列接收K个同频率远场相干信号，K为入射声源信号的数量，K≤M，

步骤一、利用3M+1个阵元构成的均匀三维阵作为接收阵列，获取阵列的P 次快拍接收数据，通过数据相关矩阵特征分解获取最大特征值对应的特征矢量利用特征矢量构造Toeplitz矩阵；

x轴子阵信号为x(t)＝A_x(α)S(t)+N_x(t)，其中A_x(α)是x轴子阵导向矢量矩阵A_x(α)＝[a_x(α₁)，a_x(α₂)，…，a_x(α_K)]，是第k 个信号对应的x轴子阵导向矢量，α_k＝sinθ_kcosφ_k是第k个信号对应的x轴方向的方向余弦，λ是入射信号的波长，y轴子阵信号为y(t)＝A_y(β)S(t)+N_y(t)其中A_y(β)是y轴子阵导向矢量矩阵，A_y(β)＝[a_y(β₁)，a_y(β₂)，…，a_y(β_K)]，是第k个信号对应的y轴子阵导向矢量，β_k＝sinθ_ksinφ_k是第k个信号对应的y轴方向的方向余弦，z轴子阵信号为 z(t)＝A_z(ω)S(t)+N_z(t)其中A_z(ω)是z轴子阵导向矢量矩阵， A_z(ω)＝[a_z(ω₁)，a_z(ω₂)，…，a_z(ω_K)]，是第k个信号对应的z轴子阵导向矢量，ω_k＝cosθ_k是第k个信号对应的z轴方向的方向余弦，x轴子阵、y轴子阵和z轴子阵P次快拍数据分别为Z_x、Z_y和Z_z，利用Z_x获得x轴子阵对应的数据相关矩阵利用Z_y获取y轴子阵对应的数据相关矩阵利用Z_z获取z轴子阵对应的数据相关矩阵对数据相关矩阵R_x、R_y和R_z分别进行特征分解，得到最大特征值对应的特征向量e_sx，e_sy和e_sz；利用特征矢量e_sx，e_sy和e_sz构造Toeplitz矩阵

R_Tx、R_Ty和R_Tz是经过Toeplitz修正后的协方差矩阵，为M×M维的HermitionToeplitz矩阵，因此实现了解相干处理；

步骤二、对R_Tx、R_Ty和R_Tz分别进行特征分解，利用ESPRIT算法得到信号 x轴、y轴方向和z轴方向的方向余弦估计；

对R_Tx进行特征分解得到信号子空间E_sx，E_sx的前M行元素组成矩阵E_x1， E_sx的后M行元素组成矩阵E_x2，令对Ψ_x进行特征分解，由特征值构成其中是矩阵E_x1的伪逆矩阵，是K×K的对角矩阵；对R_Ty进行特征分解得到信号子空间E_sy，E_sy的前M行元素组成矩阵E_y1，E_sy的后M行元素组成矩阵E_y2，令对Ψ_y进行特征分解，特征值构成其中是矩阵E_y1的伪逆矩阵，

是K×K的对角矩阵；对R_Tz进行特征分解得到信号子空间E_sz，E_sz的前M行元素组成矩阵E_z1， E_sz的后M行元素组成矩阵E_z2，令对Ψ_z进行特征分解，由特征值构成其中是矩阵E_z1的伪逆矩阵，Ф_z是K×K的对角矩阵，，

从而得到x轴、y轴和z轴方向的方向余弦方向的估计值：

因为x轴、y轴和z轴方向的三个方向是由三个独立的特征分解完成的，所以信号的排列顺序不同，不一定是同一个信号的对应的三个方向余弦，不能直接用来计算信号的到达角，必须进行配对运算，才能使三个方向余弦的排列顺序一致；

步骤三、根据同一个信号x轴、y轴和z轴方向的三个方向余弦之间的等量关系，对三个方向余弦进行配对运算；

对一个信号(θ_k，φ_k)来说，x轴和y轴方向余弦的和

此时由z轴方向的方向余弦ω_k可得到 sin²(arccos(ω_k))＝sin²θ_k，本发明就是利用同一个信号三个方向余弦之间满足的这个等量关系进行参数配对运算，步骤二中得到的x 轴的方向余弦与K个y轴的方向余弦进行平方求和得

T_s＝|T₁-T_z|，矩阵T_s中元素的最小值对应于配对成功的三个方向余弦分别写作对x轴方向余弦都做同样的配对处理，得到对应的的配对成功的x轴方向余弦y轴方向余弦和z轴三个方向余弦

步骤四、利用配对成功后的方向余弦和得到信号二维到达角的估计值

利用配对成功后的方向余弦从而得到到达角的估计：

前述步骤中的k＝1，...，K，rk∈[1，2，…，K]，pk∈[1，2，…，K]，qk∈[1，2，…，K]，j是虚数单位；

本发明引入利用虚拟扩展Toeplitz矩阵重构方法构造x、y、z方向的方向三个矩阵，通过矩阵重构解相干，利用三个方向余弦之间的关系进行配对运算，该配对算法简单有效，利用一半的阵元即可以得到二维到达角的估计结果，本方法不存在阵列孔径的损失，估计性能优于空间平滑算法；

本发明的效果可以通过以下的仿真结果进一步说明：坐标原点放1阵元；

仿真实验条件如下：

图1所示为本发明实施例的阵列结构示意图，本发明的阵列为三维均匀十字阵列，坐标原点放置1个阵元；在各轴正半轴放置M＝4个阵元，阵元间距均为半波长；在每一轴负半轴通过虚拟构造4个阵元，那么每轴上总共9个阵元，两个相干信源的入射方向为[10°，20°]，[80°，75°]，信噪比为20dB，本发明方法的阵元数是13，空间平滑算法的阵元数是25；图中阵元用“●”表示，虚拟阵元用“○”表示；

图3为单次快拍虚拟重构算法执行100次得到的散点图，图4为单次快拍空间平滑算法执行100次得到的散点图，对比图3和图4可以发现，本发明方法的估计结果集中在真实到达角附近，而空间平滑算法的估计结果较为散乱；图5为单次快拍空间平滑算法和单次快拍虚拟重构算法均方根误差对比图，图6 为500次快拍空间平滑算法和单次快拍虚拟重构算法均方根误差对比图。从图5 可以看出，在单次快拍数据的情况下，虚拟重构算法的性能不仅要优于空间平滑算法，而且所需的传感器阵元数远远小于空间平滑算法；从图6可以看出，只有将空间平滑所需快拍数提升至500，才能达到和单快拍数据虚拟重构算法相近的性能；所以本发明算法可以估计二维相干源下的到达角，是一个十分有效并且易于工程实现的相干源到达角估计算法。

以上所述，仅是本发明的较佳实施例而已，并非对本发明做任何形式上的限制，虽然本发明已以较佳实施例揭露如上，然而并非用以限定本发明，任何熟悉本专业的技术人员，在不脱离本发明技术方案范围内，当可利用上述揭示的技术内容做出些许更动或修饰为等同变化的等效实施例，但凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属于本发明技术方案的范围内。

Claims (1)

1.三维阵虚拟扩展相干源二维波达方向估计方法，其特征在于：

本发明方法所用的接收阵列由坐标原点处的阵元以及M个等间隔布置于x轴正半轴，M个等间隔布置于y轴正半轴和M个等间隔布置于z轴正半轴的3M+1个阵元构成的均匀三维阵，阵元间隔分别为d_x，d_y和d_x，d_x≤λ/2，d_y≤λ/2和d_z≤λ/2，λ为入射信号的波长；

多参数联合估计方法的步骤如下：阵列接收K个同频率的远场相干窄带入射信号，

步骤一、利用3M+1个阵元构成的均匀三维阵作为接收阵列，获取阵列的P次快拍接收数据，通过数据相关矩阵特征分解获取最大特征值对应的特征矢量利用特征矢量构造Toeplitz矩阵；

x轴子阵信号为x(t)＝A_x(α)S(t)+N_x(t)，其中A_x(α)是x轴子阵导向矢量矩阵A_x(α)＝[a_x(α₁)，a_x(α₂)，…，a_x(α_K)]，是第k个信号对应的x轴子阵导向矢量，α_k＝sinθ_kcosφ_k是第k个信号对应的x轴方向的方向余弦，λ是入射信号的波长，y轴子阵信号为y(t)＝A_y(β)S(t)+N_y(t)其中A_y(β)是y轴子阵导向矢量矩阵，A_y(β)＝[a_y(β₁)，a_y(β₂)，…，a_y(β_K)]，是第k个信号对应的y轴子阵导向矢量，β_k＝sinθ_ksinφ_k是第k个信号对应的y轴方向的方向余弦，z轴子阵信号为z(t)＝A_z(ω)S(t)+N_z(t)其中A_z(ω)是z轴子阵导向矢量矩阵，A_z(ω)＝[a_z(ω₁)，a_z(ω₂)，…，a_z(ω_K)]，是第k个信号对应的z轴子阵导向矢量，ω_k＝cosθ_k是第k个信号对应的z轴方向的方向余弦，x轴子阵、y轴子阵和z轴子阵P次快拍数据分别为Z_x、Z_y和Z_z，利用Z_x获得x轴子阵对应的数据相关矩阵利用Z_y获取y轴子阵对应的数据相关矩阵利用Z_z获取z轴子阵对应的数据相关矩阵对数据相关矩阵R_x、R_y和R_z分别进行特征分解，得到最大特征值对应的特征向量e_sx，e_sy和e_sz；利用特征矢量e_sx，e_sy和e_sz构造Toeplitz矩阵R_Tx、R_Ty和R_Tz是经过Toeplitz修正后的协方差矩阵，为M×M维的Hermition Toeplitz矩阵，因此实现了解相干处理；

步骤二、对R_Tx、R_Ty和R_Tz分别进行特征分解，利用ESPRIT算法得到信号x轴、y轴方向和z轴方向的方向余弦估计；

对R_Tx进行特征分解得到信号子空间E_sx，E_sx的前M行元素组成矩阵E_x1，E_sx的后M行元素组成矩阵E_x2，令对Ψ_x进行特征分解，由特征值构成其中是矩阵E_x1的伪逆矩阵，是K×K的对角矩阵；对R_Ty进行特征分解得到信号子空间E_sy，E_sy的前M行元素组成矩阵E_y1，E_sy的后M行元素组成矩阵E_y2，令对Ψ_y进行特征分解，特征值构成其中是矩阵E_y1的伪逆矩阵，

是K×K的对角矩阵；对R_Tz进行特征分解得到信号子空间E_sz，E_sz的前M行元素组成矩阵E_z1，E_sz的后M行元素组成矩阵E_z2，令对Ψ_z进行特征分解，由特征值构成其中是矩阵E_z1的伪逆矩阵，Φ_z是K×K的对角矩阵，

从而得到x轴、y轴和z轴方向的方向余弦方向的估计值：

因为x轴、y轴和z轴方向的三个方向是由三个独立的特征分解完成的，所以信号的排列顺序不同，不一定是同一个信号的对应的三个方向余弦，不能直接用来计算信号的到达角，必须进行配对运算，才能使三个方向余弦的排列顺序一致；

步骤三、根据同一个信号x轴、y轴和z轴方向的三个方向余弦之间的等量关系，对三个方向余弦进行配对运算；

对一个信号(θ_k，φ_k)来说，x轴和y轴方向余弦的和此时由z轴方向的方向余弦ω_k可得到sin²(arccos(ω_k))＝sin²θ_k，本发明就是利用同一个信号三个方向余弦之间满足的这个等量关系进行参数配对运算，步骤二中得到的x轴的方向余弦与K个y轴的方向余弦进行平方求和得 T_s＝|T₁-T_z|，矩阵T_s中元素的最小值对应于配对成功的三个方向余弦分别写作对x轴方向余弦都做同样的配对处理，得到对应的的配对成功的x轴方向余弦y轴方向余弦和z轴三个方向余弦

步骤四、利用配对成功后的方向余弦和得到信号二维到达角的估计值

利用配对成功后的方向余弦从而得到到达角的估计：

前述步骤中的k＝1，...，K，rk∈[1，2，…，K]，pk∈[1，2，…，K]，qk∈[1，2，…，K]，j为虚数单位。