CN104237843A - 一种分布式信源二维中心波达角的估计方法 - Google Patents

Abstract

本发明属于无线移动通信技术领域。本发明通过设置一个由两个平行的均匀线阵组成的立体天线阵列，建立各分布式信源的数据矩阵，并采用两阵列接收数据的互相关矩阵，有效地消除了加性噪声对中心波达角估计性能的影响，提高了分布式信源的二维波达角在低信噪比时的估计精度；另外，利用特征值和特征向量一一对应的关系，对中心方位角和中心俯仰角分别进行估计，无需其它配对程序即可完成自动配对，而且配对成功率更高。本发明方法具有估计方法简单、可靠，参数配对成功率和精确度高等特点。

Description

一种分布式信源二维中心波达角的估计方法

技术领域

本发明属于无线移动通信技术领域，特别是涉及一种采用阵列天线估计多径分布式信源二维中心波达角的方法。

背景技术

在现代无线移动通信系统中，信号的波达方向(direction of arrival,简称DOA)是一个非常重要的参数。DOA参数的准确估计是实现通信系统其它功能(智能天线，无线定位等)的关键和基础。因而，DOA参数的估计精度对无线通信系统性能有着至关重要的影响。

在城区等密集通信环境中，移动端发射的信号会在空间上形成散射，尤其是当散射发生在一个小区域内，周围高大建筑物对信号也会形成绕射，这就形成了无线移动通信的多径传播现象，导致信号的能量在空间中产生一定程度的分布扩展，形成空间分布式多径信号。这时目标信号已不能被简单地等效为一个点源，而应用能够反应信号空间分布特性的数学模型描述它。这种数学模型被称为分布式信源。分布式信源参数估计方法具有重要的实际意义和广阔的应用前景，吸引着无数的国内外学者投身于此项技术的研究。

目前，已有几种低计算复杂度的分布式信源的二维中心DOA估计技术，如MP方法(Modified Propagator method)(Zhi Zheng,Guangjun Li,Yunlong Teng.2D DOA estimator for multiplecoherently distributed sources using modified propagator.Circuits,Systems&Signal Processing,2012,31(1):255-270.)以及TLS-ESPRIT方法(TLS-ESPRIT method)(Zhi Zheng,Guangjun Li,Yunlong Teng.Simplified estimation of 2D DOA for coherently distributed sources.Wireless Personal Communications,2012,62:907-922.)。其中，MP方法是利用三个平行的均匀线阵组成的立体天线阵列，建立数据矢量矩阵及传播算子矩阵，然后利用旋转矩阵特征值的配对程序，确定分布式信源的中心方位角和中心俯仰角，该技术虽然计算复杂度低，但是所需阵元数多，低信噪比时参数估计性能低。TLS-ESPRIT方法是利用双平行的均匀线阵的接收数据构建数据矢量矩阵，然后建立旋转矩阵，并通过总体最小二乘法得到包含分布式信源的二维中心DOA信息的特征值，当存在多个信源时，需要采用额外的配对程序完成中心俯仰角和中心方位角的估计。该方法所需阵元数较少，但是其它估计性能不如MP方法。尽管这两种现有方法各有优缺点，但都需要额外的参数配对程序，而且配对成功率较低。

发明内容

本发明的目的是提供一种分布式信源的二维中心波达角的估计方法，以达到降低成本、提高分布式信源二维DOA的估计精度和自动配对的成功率以及应用于智能天线系统中可提高智能天线系统在多径环境中的性能及实际应用的价值等目的。

本发明方法的流程如图1所示，具体采用如下技术方案：

步骤1：设置立体天线阵列，如图2所示；该立体天线阵列由P阵列和Q阵列构成的双平行阵列结构组成，P阵列设于z轴，Q阵列设于x-z平面且平行于z轴；所述P阵列和Q阵列均为M元均匀线阵，且同一阵列中的相邻阵元间隔以及阵列P与阵列Q之间的距离相同并均为d；

设有K(K＝1～M-1)个信号源发射的信号在多径环境中传播，经由立体天线阵列中的阵列P和阵列Q接收到的数据矢量可以分别表示为x(t)和z(t)，具体如下所述：

其中，θ_k是中心方位角，φ_k是中心俯仰角；相对于中心方位角的角偏差和相对于中心俯仰角的角偏差的取值范围都为0°～10°，且相对于中心方位角的角偏差和相对于中心俯仰角的角偏差均服从对称分布(如高斯分布，均匀分布，柯西分布和指数分布等)，需要说明的是：二者可服从不同的对称分布；将原点处的阵元作为空间相位的参考点，则有M×1维P阵列导向矢量β_kl(t)为服从高斯分布的随机复增益，L为多径信号的总数，s_k(t)为第k个随机信号，n_X(t)和n_Z(t)为是加性高斯白噪声，且噪声与信号不相关，T为快拍数；

步骤2：建立参数估计信号模型；令 $\stackrel{\‾}{b}({\mathrm{\θ}}_k,{\mathrm{\φ}}_k)\≈e^{-j2\mathrm{\π}d\cos {\mathrm{\θ}}_k/\mathrm{\λ}}\stackrel{\‾}{a}\left({\mathrm{\φ}}_k\right),$ 其中， ${\left[\stackrel{\‾}{a}\left({\mathrm{\φ}}_k\right)\right]}_m\≈e^{-j2\mathrm{\π}d\cos {\mathrm{\φ}}_k/\mathrm{\λ}}{\left[w_k\right]}_m$ 和可以得到参数估计信号模型的矩阵形式为

$X=\stackrel{\‾}{A}\left(\mathrm{\φ}\right)S+N_X$

$Z=\stackrel{\‾}{B}(\mathrm{\θ},\mathrm{\φ})S+N_Z=\stackrel{\‾}{A}\left(\mathrm{\φ}\right)\mathrm{\ΦS}+N_Z$

其中，X＝[x(1),x(2),...,x(T)],Z＝[z(1),z(2),...,z(T)],信号矩阵S＝[s(1),s(2),...,s(T)]，t时刻的信号矢量s(t)＝[s₁(t),s₂(t),...,s_K(t)]^T，Φ＝diag(η₁,η₂,...,η_K)，其中Φ、η_k均为中间变量，N_X和N_Z均为高斯矩阵，

步骤3：由立体天线阵列的P阵列和Q阵列各接收T次数据矢量，可得到M×M维互相关矩阵R_xz＝E[XZ^H]，其中E[·]是期望运算；

步骤4：对互相关矩阵R_xz进行奇异值分解，分别得到K个非零特征值对应的M×K维的左、右奇异向量矩阵和以及M-K个零奇异值对应的M×(M-K)维左、右奇异向量矩阵和将向量矩阵和合并成一个2M×K维矩阵W₁，将向量矩阵和合并成一个2M×(M-K)维矩阵W₂，且满足根据广义导向矢量与噪声子空间是正交的关系，可知 $\stackrel{\‾}{A}{\left(\mathrm{\φ}\right)}^H{\tilde{U}}_2=0$ 和 $\stackrel{\‾}{B}{(\mathrm{\θ},\mathrm{\φ})}^H{\tilde{V}}_2=0,$ 将和合并成一个2M×K维矩阵C，可以得到 ${\mathrm{CW}}_2^H=0;$ 进一步地，可得span(C)＝span(W₁)，其中，span(·)表示生成的列空间；所述的2M×K维矩阵W₁、2M×(M-K)维矩阵W₂、2M×K维矩阵C的具体形式分别如下：

$W_1=\left[\begin{array}{c}{\tilde{U}}_1\\ {\tilde{V}}_1\end{array}\right],W_2=\left[\begin{array}{c}{\tilde{U}}_2\\ {\tilde{V}}_2\end{array}\right],C=\left[\begin{array}{c}\stackrel{\‾}{A}\left(\mathrm{\φ}\right)\\ \stackrel{\‾}{B}(\mathrm{\θ},\mathrm{\φ})\end{array}\right];$

步骤5：利用步骤4所得的结果，可以确定一个K×K维的非奇异矩阵Η使得其中，(·)⁺表示广义逆，(·)^-1表示逆；令旋转矩阵并对旋转矩阵Ψ进行特征分解，得到旋转矩阵Ψ的特征值及其特征向量矩阵

步骤6：由旋转矩阵Ψ的特征值的相位信息参数确定各分布式信源对应的中心方位角；

步骤7：利用旋转矩阵Ψ的特征向量矩阵对导向矢量矩阵进行估计得到其估计值 $\hat{A}\left(\mathrm{\φ}\right)={\tilde{U}}_1\hat{U},$ 其中， $\hat{A}\left(\mathrm{\φ}\right)=[\hat{a}\left({\mathrm{\φ}}_1\right),\hat{a}\left({\mathrm{\φ}}_2\right),...,\hat{a}\left({\mathrm{\φ}}_K\right)],$ 表示第k个信源导向矢量的估计值。

步骤8：利用导向矢量矩阵估计值的前M-1行和后M-1行构造近似旋转不变关系，确定各分布式信源对应的中心俯仰角。

本发明技术方案带来的有益效果：

本发明的阵列配置的阵元数较少，降低了现有技术所需的成本；另外本发明利用两阵列接收数据的互相关矩阵，有效地消除了加性噪声对中心DOA估计性能的影响，提高了分布式信源二维DOA在低信噪比时的估计精度；最后，本发明利用特征值和特征向量一一对应的关系，对中心方位角和中心俯仰角进行估计，无需其它配对程序即可完成自动配对，而且配对成功率更高。

附图说明

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

图2为本发明的天线阵列设置示意图；

图3为现有的MP方法的中心DOA估计散布图；

图4为本发明方法的中心DOA估计散布图；

图5为现有的MP方法、TLS-ESPRIT方法和本发明提供的方法的中心DOA均方根误差随信噪比变化的曲线图。

具体实施方式

下面结合两个实施例对本发明作进一步说明。

实施例1：本发明方法对多个分布式信源的中心DOA的估计配对性能。

本实施例设置六个分布式信源，其中心方位角和中心俯仰角分别为：(θ₁,φ₁)＝(40°,50°)，(θ₂,φ₂)＝(50°,70°)，(θ₃,φ₃)＝(70°,90°)，(θ₄,φ₄)＝(80°,40°)，(θ₅,φ₅)＝(60°,80°)，(θ₆,φ₆)＝(30°,60°)，在信噪比为10dB，采样数T＝200时进行500次独立的Monte Carlo仿真。

如图1所示，本发明的二维分布式信源中心DOA估计方法步骤如下：

步骤1：如图2所示设置立体天线阵列，该立体天线阵列由P阵列和Q阵列构成双平行阵列结构，P阵列设于z轴，Q阵列设于x-z平面且平行于z轴；每个阵列均为M＝16元均匀线阵，且相邻阵元间隔及阵列P与阵列Q之间的距离相同均为半波长，根据上述初始设定有K＝6个远场窄带不相关的分布式信源从不同方向入射到立体天线阵列的两个阵列P和Q上，阵列P和阵列Q接收数据矢量可以表示为x(t)和z(t)具体如下：

其中， 相对于中心方位角的角偏差和相对于中心俯仰角的角偏差的取值范围都为0°～10°，且和均服从高斯分布；将原点处的阵元作为空间相位的参考点，有16×1维的P阵列导向矢量β_kl(t)为服从高斯分布的随机复增益，多径信号的总数L＝200，n_X(t)和n_Z(t)为是加性高斯白噪声，且噪声与信号不相关，设置快拍数T＝200；

步骤2：令 $\stackrel{\‾}{b}({\mathrm{\θ}}_k,{\mathrm{\φ}}_k)\≈e^{-j\mathrm{\π}\cos {\mathrm{\θ}}_k}\stackrel{\‾}{a}\left({\mathrm{\φ}}_k\right),$ 其中， ${\left[\stackrel{\‾}{a}\left({\mathrm{\φ}}_k\right)\right]}_m\≈e^{-j\mathrm{\π}m\cos {\mathrm{\φ}}_k}{\left[w_k\right]}_m$ 和 $w_k=\underset{l=1}{\overset{L}{\mathrm{\Σ}}}{\mathrm{\β}}_{k1}\left(t\right)e^{\mathrm{j\πm}{\widetilde{\mathrm{\φ}}}_{k1}\left(t\right)\sin {\mathrm{\φ}}_k},$ 可以得到参数估计信号模型的矩阵形式为

$X=\stackrel{\‾}{A}\left(\mathrm{\φ}\right)S+N_X$

$Z=\stackrel{\‾}{B}(\mathrm{\θ},\mathrm{\φ})S+N_Z=\stackrel{\‾}{A}\left(\mathrm{\φ}\right)\mathrm{\ΦS}+N_Z$

其中，X和Z都是16×200维矩阵，X＝[x(1),x(2),...,x(200)]，Z＝[z(1),z(2),...,z(200)]，S＝[s(1),s(2),...,s(200)]，s(t)＝[s₁(t),s₂(t),...,s₆(t)]^T，Φ＝diag(η₁,η₂,...,η_K)， $\stackrel{\‾}{A}\left(\mathrm{\φ}\right)=[\stackrel{\‾}{a}\left({\mathrm{\φ}}_1\right),\stackrel{\‾}{a}\left({\mathrm{\φ}}_2\right),...,\stackrel{\‾}{a}\left({\mathrm{\φ}}_6\right)],$ $\stackrel{\‾}{B}(\mathrm{\θ},\mathrm{\φ})=[\stackrel{\‾}{b}({\mathrm{\θ}}_1,{\mathrm{\φ}}_1),\stackrel{\‾}{b}({\mathrm{\θ}}_2,{\mathrm{\φ}}_2),...,\stackrel{\‾}{b}({\mathrm{\θ}}_6,{\mathrm{\φ}}_6)];$

步骤3：P阵列和Q阵列接收数据的互相关矩阵R_xz的最大似然估计由下式确定：

${\hat{R}}_{\mathrm{xz}}=(1/200){\mathrm{\Σ}}_{t=1}^{200}x\left(t\right)z{\left(t\right)}^H$

步骤4：对互相关矩阵进行奇异值分解，得到获得样本数据的左右奇异向量矩阵和和分别是对应非零奇异值的16×6维的左右奇异向量矩阵；

步骤5：令旋转矩阵对令旋转矩阵进行特征值分解，其中，(·)⁺表示广义逆，(·)^-1表示逆，旋转矩阵的特征值向量和特征向量矩阵分别为和其中，特征值向量和特征向量矩阵是一一对应的；

步骤6：利用旋转矩阵的特征值的相位信息获得中心方位角的估计值可以表示为：

${\widehat{\mathrm{\θ}}}_k=\arccos \left(\frac{-\mathrm{angle}\left({\widehat{\mathrm{\γ}}}_k\right)}{\mathrm{\π}}\right),k=\mathrm{1,2},\·\·\·,6$

其中，angle(·)表示取相位运算。

步骤7：对P阵列的导向矢量矩阵进行估计，得到中心俯仰角的估计估，具体表示为：

$\hat{A}\left(\mathrm{\φ}\right)={\tilde{U}}_1\hat{U}$

其中， $\hat{A}\left(\mathrm{\φ}\right)=[\hat{a}\left({\mathrm{\φ}}_1\right),\hat{a}\left({\mathrm{\φ}}_2\right),...,\hat{a}\left({\mathrm{\φ}}_K\right)],$ 表示第k个信源导向矢量的估计值。

步骤8：利用的前15个元素和后15个元素构造近似旋转不变关系，得到中心俯仰角的估计值：

${\widehat{\mathrm{\φ}}}_k=\arccos \left(\frac{-\mathrm{angle}\left({\left({\left[\hat{a}\left({\mathrm{\φ}}_k\right)\right]}_{1:15}\right)}^{+}{\left[\hat{a}\left({\mathrm{\φ}}_k\right)\right]}_{2:16}\right)}{\mathrm{\π}}\right),k=\mathrm{1,2},\·\·\·,6$

现有的MP方法和本发明提供的中心DOA估计方法的散布图分别如图3和图4所示。仿真结果表明，现有的MP方法仅完成了(80°,40°)和(70°,90°)的配对；对于(30°,60°)和(60°,80°)，现有的MP方法只完成了部分准确地配对；对于(40°,50°)和(50°,70°)，现有的MP方法不能准确地完成参数的估计。由图3和图4可知，本发明提供的分布式信源二维中心波达角的估计方法不仅能准确估计六个信源的中心DOA，而且无需额外的配对程序就能够完成准确地配对，且每一对参数的配对成功率都为100％。

实施例2：考察一个分布式信源的中心DOA均方根误差随信噪比变化的性能。

本实施例所采用的分布式信源的中心方位角和中心俯仰角为(50°,70°)，总阵元数为10，在采样数据T＝200时进行500次独立的Monte Carlo仿真。

本实施例的二维分布式信源中心DOA估计方法步骤与发明内容及实施例1的过程相同，流程如图1所示。图5是现有的MP方法、TLS-ESPRIT方法和本发明提供的方法的中心DOA均方根误差随信噪比变化的对比。从图5可以看出，在信噪比大于-8dB时，本发明方法所得的中心DOA均方根误差与现有MP方法所得的中心DOA均方根误差的曲线基本重合；而在低信噪比时，本发明方法的中心DOA估计性能更好。

Claims (6)

1.一种二维分布式信源中心波达角的估计方法，其特征在于，包括以下步骤：

步骤1：设置立体天线阵列，该立体天线阵列由P阵列和Q阵列构成的双平行阵列结构组成，P阵列设于z轴，Q阵列设于x-z平面且平行于z轴；所述P阵列和Q阵列均为M元均匀线阵，且同一阵列中的相邻阵元间隔以及阵列P与阵列Q之间的距离相同并均为d；

设有K(K＝1～M-1)个信号源发射的信号在多径环境中传播，经由立体天线阵列中的阵列P和阵列Q接收到的数据矢量分别表示为x(t)和z(t)，具体如下：

其中，θ_k是中心方位角，φ_k是中心俯仰角，分别是相对于中心方位角的角偏差、相对于中心俯仰角的角偏差；将原点处的阵元作为空间相位的参考点，则有P阵列的M×1维导向矢量β_kl(t)为服从高斯分布的随机复增益，L为多径信号的总数，s_k(t)为第k个随机信号，n_X(t)和n_Z(t)为是加性高斯白噪声，且噪声与信号不相关，T为快拍数；

步骤2：利用立体天线阵列接收到的数据矢量x(t)和建立相应的参数估计信号模型；令其中，广义导向矢量和变量w_k分别满足 ${\left[\stackrel{\‾}{a}\left({\mathrm{\φ}}_k\right)\right]}_m\≈e^{-j2\mathrm{\π}d\cos {\mathrm{\φ}}_k/\mathrm{\λ}}{\left[w_k\right]}_m$ 和 $w_k=\underset{l=1}{\overset{L}{\mathrm{\Σ}}}{\mathrm{\β}}_{k1}\left(t\right)e^{j2\mathrm{\πmd}{\widetilde{\mathrm{\φ}}}_{k1}\left(t\right)\sin {\mathrm{\φ}}_k/\mathrm{\λ}},$ 由此可得参数估计信号模型的矩阵形式，表示为矩阵X和矩阵Z，具体如下：

$X=\stackrel{\‾}{A}\left(\mathrm{\φ}\right)S+N_X$

$Z=\stackrel{\‾}{B}(\mathrm{\θ},\mathrm{\φ})S+N_Z=\stackrel{\‾}{A}\left(\mathrm{\φ}\right)\mathrm{\ΦS}+N_Z$

其中，X＝[x(1),x(2),...,x(T)],Z＝[z(1),z(2),...,z(T)]，信号矩阵S＝[s(1),s(2),...,s(T)]，t时刻的信号矢量s(t)＝[s₁(t),s₂(t),...,s_K(t)]^T，Φ＝diag(η₁,η₂,...,η_K)，其中Φ、η_k均为中间变量，N_X和N_Z均为高斯矩阵，

步骤3：由立体天线阵列的P阵列和Q阵列各接收T次信号所得的数据矢量，可得到M×M维互相关矩阵R_xz＝E[XZ^H]，其中E[·]是期望运算；

步骤4：对互相关矩阵R_xz进行奇异值分解，分别得到K个非零特征值对应的M×K维的左、右奇异向量矩阵和以及M-K个零奇异值对应的M×(M-K)维左、右奇异向量矩阵和将向量矩阵和合并成一个2M×K维矩阵W₁，将向量矩阵和合并成一个2M×(M-K)维矩阵W₂，且满足由于广义导向矢量与噪声子空间是正交的关系，可知 $\stackrel{\‾}{A}{\left(\mathrm{\φ}\right)}^H{\tilde{U}}_2=0$ 和 $\stackrel{\‾}{B}{(\mathrm{\θ},\mathrm{\φ})}^H{\tilde{V}}_2=0,$ 将和合并成一个2M×K维矩阵C，且满足 ${\mathrm{CW}}_2^H=0;$ 进一步地，可得span(C)＝span(W₁)，其中，span(·)表示生成的列空间；

步骤5：利用步骤4所得的结果，可以确定一个K×K维的非奇异矩阵Η使得其中，(·)⁺表示广义逆，(·)^-1表示逆；定义旋转矩阵Ψ，令旋转矩阵并对旋转矩阵Ψ进行特征分解，得到旋转矩阵Ψ的特征值及其特征向量矩阵

步骤6：由旋转矩阵Ψ的特征值的相位信息参数确定各分布式信源对应的中心方位角；

步骤7：利用旋转矩阵Ψ的特征向量矩阵对导向矢量矩阵进行估计得到其估计值矩阵其中， 表示第k个信源导向矢量的估计值；

步骤8：利用导向矢量矩阵估计值矩阵的前M-1行和后M-1行构造近似旋转不变关系，确定各分布式信源对应的中心俯仰角。

2.根据权利要求1所述的一种分布式信源二维中心波达角的估计方法，其特征在于，所述的2M×K维矩阵W₁、2M×(M-K)维矩阵W₂、2M×K维矩阵C的具体形式分别如下：

$W_1=\left[\begin{array}{c}{\tilde{U}}_1\\ {\tilde{V}}_1\end{array}\right],W_2=\left[\begin{array}{c}{\tilde{U}}_2\\ {\tilde{V}}_2\end{array}\right],C=\left[\begin{array}{c}\stackrel{\‾}{A}\left(\mathrm{\φ}\right)\\ \stackrel{\‾}{B}(\mathrm{\θ},\mathrm{\φ})\end{array}\right].$

3.根据权利要求1所述的一种分布式信源二维中心波达角的估计方法，其特征在于，所述相对于中心方位角的角偏差和相对于中心俯仰角的角偏差的取值范围都为0°～10°，且相对于中心方位角的角偏差和相对于中心俯仰角的角偏差均服从对称分布。

4.根据权利要求1或3所述的一种分布式信源二维中心波达角的估计方法，其特征在于，所述的相对于中心方位角的角偏差的分布形式是高斯分布、均匀分布、柯西分布、指数分布中的一种。

5.根据权利要求1或3所述的一种分布式信源二维中心波达角的估计方法，其特征在于，所述的相对于中心俯仰角的角偏差的分布形式是高斯分布、均匀分布、柯西分布、指数分布中的一种。

6.根据权利要求1所述的一种分布式信源二维中心波达角的估计方法，其特征在于，所述的接收数据的采样次数T的范围是50～1000次。