一种部分稀疏l阵及其二维doa估计方法
【专利摘要】本发明公开了一种部分稀疏L阵及其二维DOA估计方法。本发明包括两个子阵构成L阵，第一子阵的阵元间距为半倍波长，第二子阵的阵元间距为n倍波长；并在第二子阵上设置一个距参考阵元的距离为半倍波长的辅助阵元。在进行DOA估计处理时，将子阵分别置于x轴和z轴，利用互相关矩阵不受噪声影响的特性，求接收数据的互相关并提取信号子空间；利用ULA的平移不变性和信号子空间，求取z轴阵列流型矩阵的旋转矩阵，并通过对旋转矩阵进行特征值分解，得到可能的俯仰角估计值，再基于辅助阵元对估计结果去模糊处理；再基于其进行信源波形的估计及x轴的阵列流型矩阵，求取相应的方位角。本发明用于雷达，声呐等，其实现低成本，低运算量，测向精度高。
【专利说明】
-种部分稀疏L阵及其二维DOA估计方法
技术领域
[0001] 本发明属于无线移动通信技术领域，特别是设及一种利用线阵构造的部分稀疏L 阵列及其二维波达方向(DOA)估计方法。
【背景技术】
[0002] 空分多址技术是继时分复用、频分复用和码分复用之后又一个可提高网络容量的 关键技术，空分多址的关键在于精确地测量信号源方位和形成高定向波束，因此，WDOA估 计和数字波束形成为核屯、的阵列信号处理技术越来越得到广大研究者的研究和关注。
[0003] 现有的DOA估计算法大多是基于传统满阵，即阵列中相邻阵元间距不大于入射信 号半波长的天线阵列。但满阵由于阵元间距的限制，要增加阵列孔径，提升测向精度和分辨 率就必须增加阵元数，造成系统过于复杂和系统成本增加。鉴于传统满阵存在的上述问题， 人们提出了稀疏阵，即阵元间距大于半波长的均匀或非均匀阵列系统。和常规满阵相比，阵 元数相同时它拥有更大的阵列孔径，减小了阵元间的互禪，改善了测向精度、测向分辨率和 最大可处理信源数等性能;此外，阵列孔径相同时稀疏阵所需阵元数更少，运就意味着更小 规模的接收系统和信号处理系统等，极大地降低了设备成本。
[0004] 目前基于稀疏阵列的DOA估计主要分为两大类:一类是基于稀疏线阵的一维DOA估 计；另一类是基于由几个稀疏线阵组成的简化面阵的二维DOA估计。其中，基于简化面阵的 二维DOA估计中，L阵由于有更大的有效孔径、更小的运算量、更易实现、更强的方法适用性 等优点得到了广泛的关注和应用。但现有的L阵大多是由常规ULA(均匀线阵）构成，利用稀 疏线阵构成的L阵还比较少，已有的基于稀疏线阵的L阵的表达式不够简洁，且在计算二维 DOA时，要么需要配对算法，要么就是需要谱峰捜索等，不能充分利用稀疏阵列的优势。
【发明内容】
[0005] 本发明的发明目的在于:针对上述存在的问题，提供一种构造简单，基于线阵的部 分稀疏L阵及其能获得较好性能且计算简单的DOA估计方法。
[0006] 本发明的一种部分稀疏的L阵列，包括阵元数相同的线性第一子阵、第二子阵构成 L字形的阵列，第一子阵和第二子阵的交点定义为参考阵元，其特征在于，所述第一子阵的 阵元间距等于半倍波长，第二子阵的阵元间距等于n倍波长，其中整数n的取值范围为1~ 10;并在第二子阵上设置一个辅助阵元，所述辅助阵元距离参考阵元的距离为半倍波长。其 中L阵列的第一子阵、第二子阵可分别放置于X轴和Z轴，X轴对应方位角，Z轴对应俯仰角;或 者将L阵列的第一子阵、第二子阵分别放置于Z轴和X轴，X轴对应方位角，Z轴对应俯仰角。
[0007] 同时，本发明还公开了一种用于所述的L阵列的二维波大方向的估计方法，包括下 列步骤：
[000引步骤1:将L阵列的第一子阵、第二子阵分别放置于X轴和Z轴；
[0009]步骤2:L阵列接收K个不相关信源的入射信号，得到X轴、Z轴各阵元的接收数据，其 中K小于第一子阵的陈元数M;
[0010] 步骤3:计算各入射信号的估计俯仰角；
[0011] 步骤301:计算接收数据(t)和(t)在N(N的取值为大于等于10的整数)次采样 下的互相关矩阵RzV，并从Rz'X'中任意提取K列构造信号子空间Uz,其中接收数据(t)为X 轴上除参考阵元外的所有陈元的接收数据，接收数据(t)为Z轴上除辅助阵元外的所有阵 元的接收数据；
[001^ 将信号子空间Uz划分为上下两个(M-I)体维的信号子空间Uzl和Uz2,基于Uzl和Uz2求 取旋转矩阵Q Z，令矩阵
I矩阵Tz为K X K的非奇异矩阵；
[001引对矩阵Fz进行特征值分解，得到相应的K个特征值4，其中k= 1，2，...，K;
[0014] 步骤302:计算每个入射信号的两个备选估计俯仰角：
[0015] 根据
计算K个第一备选估计俯仰角g，其中dz = nA，A表示信 号波长，符号angle( ?)表不取相位角；
[0016] 根据公式
计算估计空间相位差卽，若0,.+2;r>2;r，则第二备选估 计俯仰角
；否化
[0017] 步骤4:对每个入射信号的备选估计俯仰角去模糊角处理，得到估计俯仰角％;
[001引步骤401:计算接收数据z(t)在N次采样下的自相关矩阵Rzz,其中接收数据z(t)为Z 轴上所有阵元的接收数据；
[0019] 步骤402:对Rzz进行特征值分解：
，其中Us和Un分别为 (M+1) XK和(M+1) X (M+1-K)维的信号子空间和噪声子空间；Xs为从Rzz的M+1个特征值中， 取前K个最大特征值组成对角矩阵Xs,剩余的M+1-K个特征值组成对角矩阵Sn,符号（?）H 表示矩阵共辆；
[0020] 步骤403:根据公;
计算每个入射信号的备选估计俯仰角4 在Z轴的导向矢量设），其中远包括第一备选估计俯仰角殘和第二备选估计俯仰角耸，0, 为对应4的估计空间相位差，其中e为自然常数，符号（.）T表示矩阵转置；
[0021] 将导向矢量SJ运)与噪声子空间Un正交的备选估计俯仰角爲作为入射信号k的估 计俯仰角 4，k=l，2，...，K;
[0022] 步骤5:基于K个估计俯仰角4，得到Z轴的阵列流型矩阵4 ;
[0023] 根据接收数据z(t)、阵列流型矩阵4计算估计信号源
.其中符 号（?）+表示M-P广义逆；
[0024] 步骤6:估计各入射信号的方位角是：
[0025] 计算接收数据^(t)和X(t)在N次采样下的互相关矩阵Rxz",其中接收数据^(t)为 Z轴上除参考阵元的所有阵元的接收数据;接收数据X(t)为X轴上所有阵元的接收数据；
[0026] 根据公式
开算X轴的阵列流型矩阵的估计值1,，其中估计值
，氏为估计信号源如)在N次采样下的自相关矩阵；
[0027] 根据公式
//了]得到每个入射信号的方位角4 :，其中k =1，2,…，K A，A表示矩阵At的第k列的前M-I行元素的子向量，&一:表示矩阵At的第k列的 后M-I行元素的子向量。
[0028] 综上所述，由于采用了上述技术方案，本发明的有益效果是:可W实现低成本，低 运算量，高测向精度的二维DOA估计，该方法可W应用于雷达，声响及无线通信等领域，解决 现有常规整理所不能解决的问题。
【附图说明】
[0029] 图1是本发明提出的由稀疏ULA和常规ULA构造的部分稀疏L阵的阵列结构；
[0030] 图2是本发明所提部分稀疏L阵的Z轴阵元配置；
[0031 ]图3是利用本发明所提算法估计二维DOA的角度散布图；
[0032] 图4是本发明所提算法和现有的算法的性能随信噪比变化的对比图；
[0033] 图5是本发明所提算法和现有的算法的性能随采样快拍数变化的对比图。
【具体实施方式】
[0034] 为使本发明的目的、技术方案和优点更加清楚，下面结合实施方式和附图，对本发 明作进一步地详细描述。
[0035] 本发明为了利用稀疏阵有更大的阵列孔径和均匀线阵的平移不变性，提出一种基 于稀疏ULA和常规ULA的L阵结构，然后利用互相关矩阵不受噪声影响的特性，ULA的平移不 变性，解模糊算法和LS(最小二乘法)技术求解二维D0A。本发明可用于雷达，声响，无线通信 等领域，为了更好地理解本发明，W下结合附图对本申请作进一步描述。
[0036] 步骤1:构建稀疏L阵：
[0037] 本实施方式中，设第一、二线阵的阵元个数均为M，将第一、二线阵的交点阵元定义 为参考阵元，其中，第一子阵的陈元间距di = V2,即常规的ULA;第二子阵的陈元间距cb = n A，其中A表示信号波长，本实施例中取n = l，即稀疏ULA。参见图1，将第一子阵放置于X轴、第 二子阵放置于Z轴(疏化A和常规ULA在坐标轴上的位置可W互换），参考阵元所在位置定义 为原点。且在第二子阵设置辅助阵元，其与原点处的参考阵元间距为di。其中Z轴的系数ULA 用来估计俯仰角，X轴的常规ULA用来估计方位角。
[0038] 步骤2:构建系统模型，从而得到本发明的DOA估计的理论支持：
[0039] 设K个化<M)个远场窄带不相关信号5(〇 = [31(0，...，31((0^从1(个方向（01， (60 4 = 1, ...，K入射到本发明的发明L阵上，参见图2,其X轴和Z轴的阵列流型矩阵Az和Ax 分别为:Ax=[ax( 4 1)，. . .，ax( 4K)]、Az = [az(目1)，. . .，az(目K)]，其中，
分别为第k个信号在X轴和Z轴的导向矢量，其中dx = di，dz = d2。考虑系统噪声nx(t) = [nx,o (t),. . . ,nx,M-i(t)]T和nz(t) = [nz,o(t),. . . ,nz,M(t)]T分别为X轴和Z轴的MX I维和(M+1) X I 维的独立同分布的加性高斯白噪声。从而，可将X轴和Z轴的阵元在t时刻的接收数据矢量分 别表示为：
[0040]
[0041]
[0042] 其中，L阵在原点处的参考阵元是X轴和Z轴所共有的，因此有xo(t) = zo(t)，nx,o(t) = nz,o(t) O
[0043] 步骤3:求俯仰角：
[0044] 在得到本发明的L阵得到接收数据X(t)和Z(t)之后，为了充分利用ULA的平移不变 性和互相关矩阵不受加性噪声影响的特性，选取图2所示的X轴第1到M-I个阵元的接收数据 (t) = [Xl(t), . . .，XM-l(t)]哺Z轴除辅助阵元外其余阵元的接收数据z/ (t) = [Z0(t)，Z2 (t)，. . .，ZM(t)]T，求N(取值范围为20~1000)次采样下，(t)和(t)的互相关矩阵RzV， 即：
[0045]
[0046] 由上式可知:互相关矩阵RzV的每一列都是Az的行的线性组合，符号"A(a:b，：r表 示对应矩阵A的第a至第b行。本发明中选取RzV的前K列来构造信号子空间Uz,则存在一个K X K的非奇异矩阵Tz,满足:Uz =AzTzd
[0047] 在得到信号子空间Uz后，将Uz划分为上下两个(M-I) XK维的信号子空间Uzi和Uz2， 即：Uzi = UzU :M-1，：）=AziTz=Az(1:M-1，：）Tz、Uz2 = Uz(2:M，：）=Az2Tz=Az(2:M，：）Tz。
[004引根据Uzi和Uz2，利用Z轴稀疏ULA的平移不变性，有 ....
巧中，旋转矩阵Q Z表示包含K个入射信号 俯仰角的所有信息的对角阵。令
，则可得Uz2 = Uz1Fz。从而：
，由于矩阵 Fz和旋转矩阵Q Z相似，因此，Fz的特征值是Q Z的对角线元素的排列。通过对Fz进行特征值分 解，得到相应的K个特征值Al，A2，…，Ak，再计算每个入射信号的两个备选估计俯仰角：
[0049] 根据
计算K个第一备选估计俯仰角g，其中dz = nA，A表示信 号波长，符号angle( ?)表示取相位角；
[0化0]根据公式
计算估计空间相位差贫I，若
，则第二备选估 计俯仰角
[0051]步骤4:利用辅助阵元解模糊，即对每个入射信号的备选估计俯仰角去模糊角处 理，得到估计俯仰角爲：
[0052]由于信号子空间和噪声子空间正交，因此，本发明中，先利用Z轴所有的接收数据，求其 自相关矩^
I然后通过特征值分解得到：
其中，Us和Un分别为(M+1) XK和(M+1) X (M+1-K)维的信号子空间和噪声子空间；X S为从Rzz 的M+1个特征值中，取前K个最大特征值组成对角矩阵Xs,剩余的M+1-K个特征值组成对角 矩阵Sn;
[0化3]根据公式
计算每个入射信号的备选估计俯仰角轴的 导向矢量5:(马），其中&包括马'和耸，卽为对应每的估计空间相位差，其中e为自然常数； [0化4]将导向矢量与噪声子空间Un正交的备选估计俯仰角4作为入射信号k的估计俯仰 角4，k=l，2,…，K，即将每个信号的备选估计俯仰角中的不正交删除。
[0055] 步骤5:利用已得的俯仰角与LS技术估计信源波形。
[0056] 本发明中通过求解信源波形自的，建立俯仰角和方位角之间一一对应的关系，可W实现 俯仰角和方位角的自动配对，不需要额外的配对过程。在求解信源波形时，首先基于K个估计俯仰 角4，得到Z轴的阵列流型矩阵4，即
， 然后将求解鼠肖的问题转化为最小二乘问题
，从而可得：
[0化7]步骤6:估计方位角。
[005引先求解 X(tWP^(t) = [Zl(t)，Z2(t)，...，ZM(t)]T 的互相关矩阵 Rxz",即：
[0059 ]
，然后利用最小二乘法求A X R S的估计 值，即：
[0060]
。由于方位角的所有信息都包含在Ax中，因此，求方位角 的问题就转变为求解Ax的估计值的问题，因而本发明在得到后，先求解其自相关矩阵
'然后利用和I；，求X轴阵列流型矩阵Ax的估计值，即： ^ 一列应用平移不变性，可W求得每一列的旋转因子，从而可W 得到相应的方位角4 1，即：！
[0061] 另外，本发明的第一子阵和第二子阵所在的位置可W互换，即将L阵列的第一子 阵、第二子阵分别放置于Z轴和X轴，X轴对应方位角，Z轴对应俯仰角。则在DOA估计的处理过 程中，采用上述相同的原理，先求取第二子阵所在的X轴的信号子空间山，然后再估计每个 入射信号的两个备选估计方位角，进行去模糊处理，基于每个信号唯一对应的方位角得到 其在X轴的阵列流型矩阵，再基于信号波形的估计求取Z轴的阵列流型矩阵，进而基于已估 计的方位角得到估计俯仰角，实现DOA估计。
[0062] 图3是取K = 2，M = 5，N = 200，A = 0.8m，SNR = 5dB，（目l，(6l) = (55°，60°)，（目2，(62) = (65°，70°)时的仿真结果图，由该角度散布图可知，本发明的俯仰角和方位角比较集中的 分布于真实值附近。
[0063] 图4、5是本发明和现有方式（其中ULA-CCM是文献"D0A Estimation using cross? correlation matrix" 中所 采用的 方案， JSVD 是文献 "Joint SVd Of two crosscorrelation matrices to achieve automatic pairing in 2-D angle estimation problems"中所采用的方案）的性能随信噪比变化的对比图、W及随采样快拍 数变化的对比图。图4分别是第一个俯仰角和方位角的估计误差随信噪比的变化曲线，其 中，K = 2，M = 5，N=200 ，A = O.8m，SNR=(O~4O)dB，（0i，<l)i) = (5O°，6O°)，（02，<l)2) = (7O°， 50° )，由图4可知，采用本发明提出的稀疏L阵和相应的二维DOA估计算法之后，俯仰角的性 能可W比对比算法提升大约15地，而方位角可W提升大约10地。图5分别是第一个俯仰角和 方位角的估计误差随采样快拍数N的变化曲线，其中，K = 2 ,M= 5，N= 200，A = 0.8m，SNR = 5dB，（0i，(61) = (50°，60°)，（目2, (62) = (70°，50°)，由图5可知，采用本发明提出的部分稀疏 L阵和相应的二维DOA估计方法之后，俯仰角和方位角的性能均可W比对比算法提升至少 12地。
[0064] 因此，本发明所提的新的L阵列及其相应的二维DOA估计方法能够很好的提高二维 DOA估计的测向精度，和已有的方式相比，本发明W较少的成本和计算量实现了比已有方式 更好的性能。
【主权项】
1. 一种部分稀疏L阵列，包括阵元数相同的线性第一子阵、第二子阵构成L字形的阵列， 第一子阵和第二子阵的交点定义为参考阵元，其特征在于，所述第一子阵的阵元间距等于 半倍波长，第二子阵的阵元间距等于η倍波长，其中整数η的取值范围为1~10;并在第二子 阵上设置一个辅助阵元，所述辅助阵元距离参考阵元的距离为半倍波长。2. 如权利要求1所述的L阵列，将L阵列的第一子阵、第二子阵分别放置于X轴和ζ轴，X轴 对应方位角，ζ轴对应俯仰角。3. 如权利要求1所述的L阵列，将L阵列的第一子阵、第二子阵分别放置于ζ轴和X轴，X轴 对应方位角，ζ轴对应俯仰角。4. 一种用于权利要求1所述的L阵列的二维波大方向的估计方法，其特征在于，包括下 列步骤： 步骤1:将L阵列的第一子阵、第二子阵分别放置于X轴和ζ轴； 步骤2:L阵列接收K个不相关信源的入射信号，得到X轴、ζ轴各阵元的接收数据，其中K 小于第一子阵的陈元数M; 步骤3:计算各入射信号的估计俯仰角； 步骤301:计算接收数据ζ'(t)和Y (t)在N次采样下的互相关矩阵RzV，并从RzV中任意 提取K列构造信号子空间Uz，其中接收数据YU)为X轴上除参考阵元外的所有陈元的接收 数据，接收数据y(t)为ζ轴上除辅助阵元外的所有阵元的接收数据； 将信号子空间Uz划分为上下两个(M-I) *K维的信号子空间UzI和UZ2,基于UzI和UZ2求取旋 转矩阵Ω z，令矩阵€ = Γ. U .T；，矩阵T^K X K的非奇异矩阵； 对矩阵?2进行特征值分解，得到相应的K个特征值疋，其中k=l，2,…，Κ; 步骤302:计算每个入射信号的两个备选估计俯仰角： 根据计算K个第一备选估计俯仰角(?，:其中dz = ηλ，λ表示信号波 长，符号angle( ·)表不取相位角； 根据公式合,_ = 〇?s7e(X)计算估计空间相位差，若+._2:^?2τγ，贝1J第二备选估计俯仰步骤4:对每个入射信号的备选估计俯仰角去模糊角处理，得到估计俯仰角身： 步骤401:计算接收数据ζ (t)在N次采样下的自相关矩阵Rzz，其中接收数据ζ (t)为ζ轴上 所有阵元的接收数据； 步骤402:对Rzz进行特征值分解：，其中Us和仏分别为(M+1) XK和(M+l) X (M+1-K)维的信号子空间和噪声子空间；Xs为从Rzz的M+1个特征值中，取前K 个最大特征值组成对角矩阵Σ s，剩余的M+1-K个特征值组成对角矩阵Σ n; 步骤403:根据公式\(堯）=[1,2,..^_^计算每个入射信号的备选估计俯仰角易在2轴 的导向矢量屺(?)，其中4包括第一备选估计俯仰角豸和第二备选估计俯仰角老，<为对 应氣的估计空间相位差，其中e为自然常数； 将导向矢量?成）与噪声子空间U1UH交的备选估计俯仰角:?作为入射信号k的估计俯仰 角 4:,k=l，2,…，K; 步骤5:基于K个估计俯仰角堯，得到z轴的阵列流型矩阵& ; 根据接收数据z(t)、阵列流型矩阵iz计算估计信号源鉍>，即§(/) = λ_Ζ(〇,其中符号 (· )+表示M-P广义逆； 步骤6:估计各入射信号的方位角表； 计算接收数据Z〃（t)和X(t)在N次采样下的互相关矩阵Rxz 〃，其中接收数据Z〃（t)为ζ轴 上除参考阵元的所有阵元的接收数据;接收数据X(t)为X轴上所有阵元的接收数据； 根据公式=為式疫;1计算X轴的阵列流型矩阵的估计值Ajc，其中估计值? §⑴在N次采样下的自相关矩阵； 4得到每个入射信号的方位角i，其中k=l， 2,…，,表示矩阵A的第k列的前M-I行元素的子向量χ-2表示矩阵人t的第k列的后Μ-? 行元素的子向量。5.如权利要求4所述的方法，其特征在于，采样数N的取值为大于等于10的整数。
【文档编号】G01S3/14GK106019213SQ201610300904
【公开日】2016年10月12日
【申请日】2016年5月9日
【发明人】郑植, 杨雨轩, 刘柯宏, 闫波, 杨海芬, 林水生, 葛琰
【申请人】电子科技大学