罗雪

(北京理工大学 信息与电子学院，北京 100081)

多信源测向在雷达、声纳、移动通信等领域受到了极大关注. 以多重信号分类(multiple signal classification, MUSIC)[1]和基于旋转不变技术的信号参数估计(estimating signal parameter via rotational invariance techniques,ESPRIT)[2]为代表的超分辨算法，能够利用子空间理论识别入射角度相距很近的信源，因此在过去几十年间被广泛应用[3]. 尽管特征结构方法有很多潜在优势，但由于它们对阵列流形信息的精度有极高要求，将其直接应用在实际系统中是很困难的. 换句话说，超分辨算法的性能严重依赖于阵列流形的精确性. 实际上，阵列流形不可避免地会受到互耦、位置变动和阵列增益或相位不确定性的影响，这将造成阵列接收信号明显的幅度与相位失真. 直接用基于特征结构的方法将会造成测向性能严重恶化[4-7]. 阵列校准是一个能有效缓解波达方向(direction of arrival, DOA)估计偏差的方法.

事实上，在过去几十年间，人们已经研究了多种校准方法. 包括插值法[8-9]，已知位置校准信源的最大似然法[10]，最小均方法[11]，自校准方法[12-14]以及最大后验估计器[15-16]. 在近些年，多个基于互耦矩阵(mutual coupling matrix, MCM)特征的算法被发展. YE[17]通过在阵列两侧设置传感器获得了DOA的初始估计. 因此，运用原始阵列数据能够直接执行空间谱的一维搜索. 此外，还可以将估计的互耦参数用于补偿，搜索结果能够通过迭代进一步细化. 从中心子阵收到的信号被直接运用到传统的MUSIC，因而MCM被假设有一个复对称的Toeplitz矩阵[18]. 为了提升DOA估计的性能，LIAO等[19]中提出利用MCM的特殊结构参数化了方向向量，从而用整个阵列来实现DOA估计，并通过互耦补偿提升了结果的精确性. 此外，这些校准方法也被应用于远场信源[20]和非圆信号的DOA问题中[21].

本文假设MCM是有有限个非零元素的带对称Toeplitz矩阵. 通过运用文献[19]中参数化的方向向量，能够用整个阵列获得初步的DOA估计. 之后，估计角度依赖系数(angularly dependent coefficients , ADCs)，利用子空间的正交性和有全部ADCs信息的互耦系数进行进一步的互耦补偿.

1 信号模型与问题分析

本文考虑了由M个传感器组成的均匀线阵(uniform linear array, ULA)收到了K个窄带远场信源s1(t),s2(t),…,sK(t). 信源分别从与阵列法线相隔θ1,θ2,…,θK角度的方向入射. 假设相邻传感器间的距离为d，载波波长为λ，由于传感器间的相互作用，互耦影响不能够被忽略.M×1维的阵列输出向量可以表示为

x(t)=CAs(t)+n(t)

(1)

式中：x(t)=[x1(t)x2(t)…xM(t)]T；s(t)=[s1(t)s2(t)…sK(t)]T；n(t)=[n1(t)n2(t)…nM(t)]T分别表示了接收信号向量，源信号向量和噪声向量;符号[·]T表示转置;A=[a(θ1)a(θ2)…a(θK)]是阵列流形矩阵，其中的变量a(θk)=[1β(θk)β2(θk)…βM-1(θk)]T是方向向量，β(θk)=exp{-j2πdsinθk/λ};C∈M×M是MCM，通常被认为独立于信号DOA[22]，代表了阵列中任意两个传感器间的相互作用. 假设加性传感器噪声ni(t),i=1,2,…,M为独立同分布的高斯白噪声，有相同的方差σ2.si(t)和ni(t)为零均值广义平稳随机过程.

关于耦合模型的几个研究[12,23]表明一对传感器间的耦合几乎是相同的. 因此，ULA的MCM是一个带状对称的Toeplitz矩阵. 此外，两个传感器间的互耦与它们的距离成反比，如果两个传感器相距几个波长，系数将变为0. 用cij=cji=c|i-j|表示ULA的第i个阵元与第j个阵元间的互耦系数，并假设MCM中共有P个不同的非零元素，其中自耦c0被归一化为1. 那么，M×M维的矩阵C能够被表示为

C=Toeplitz{1,c1,…,cP-1,0M-P}

(2)

其中，Toeplitz{·}表示了一个对称的Toeplitz矩阵，由P×1维向量c=[1c1…cP-1]T构造. 定义方向θ的等效方向向量为

am(θ)=Ca(θ)

(3)

从公式(1)中可见，接收信号的协方差矩阵为

Rx=E[x(t)xH(t)]=CARsAHCH+σ2I=

(4)

式中：(·)H表示了共轭转置操作；Rs=E[s(t)sH(t)]为s(t)的协方差矩阵；Am=[am(θ1)am(θ2)…am(θk)]是等效的阵列流形矩阵. 如文献[2]中所证明的，运用Rx的特征值分解(eigenvalue decomposition, EVD)，式(4)能够被重新写为

Rx=[Es,En]Λ[Es,En]H

(5)

其中Λ=diag{Λs}为一个矩阵，M个特征值降序排列在对角线上，其余位置为0.Es∈M×K和En∈M×(M-K)分别由K个最大特征值与M-K个较小特征值对应的特征向量组成. 根据著名的子空间的算法MUSIC[1]，Am=CA的列张成的信号子空间与En张成的噪声子空间是正交的. 因此下式成立

(6)

2 互耦补偿与波达方向估计

首先通过参数化MCM重新制定了等效的方向向量. 之后，在存在未知互耦的情况下，通过应用整个阵列采样数据解决了DOA估计问题，获得更好的性能[19,24]. 最后，为了进一步细化，本文提出了一种新的互耦补偿算法，该算法利用了MCM的特殊结构，没有信息损失.

2.1 利用全部阵列采样的DOA估计

定义rk为M×1维等效方向向量am(θ)的第k个元素. 结合公式(2)(3)，am(θ)能够被表示为

am(θ)=[r1r2…rM]T

(7)

其中

am(θ)=g0Γ(θ)a(θ)

(8)

其中，Γ(θ)=diag{μ1,…,μP-1,1,…,α1,…,αP-1}，中间有M-2(P-1)个1，且

μk=rkβ1-k(θ)/g0,k=1,2,…,P-1

(9)

αk=rk+(M-P+1)β1-k(θ)/g0,k=1,2,…,P-1

(10)

是ADCs，它们为ck和β(θ)的函数. 式(8)～(10)表明，与角度无关的MCM被转换成角度有关的表达式g0Γ(θ). 因为Γ(θ)是有2(P-1)个未知变量的对角矩阵，a(θ)为一个向量，交换他们的元素并将Γ(θ)中1的数量减少为1个[19]，可得

am(θ)=g0T(θ)v(θ)

(11)

v(θ)=[μ1…μP-11α1…αP-1]T为一个(2P-1)×1维向量，并且

(12)

是一个M×(2P-1)维矩阵，由一个(M-2P+2)×1维向量和两个(P-1)×(P-1)维对角矩阵组成.

考虑式(6)中的子空间特性和式(11)中am(θ)的定义，能够在g0≠0的条件下得到如下的等式：

vH(θ)Q(θ)v(θ)=0

(13)

(14)

要注意的是，Q(θ)并不包含c的任何信息. 因此即使在互耦未知的情况下也能使用这个谱函数. 与Svantesson提出的算法[17]相比，这个谱函数利用了全部的阵列采样，不需要提取中间阵列进行互耦消除. 因此谱估计没有任何信息损失. 此方法仅适用于如下条件被满足的情况：

K≤M-2P+1

(15)

在这个情况下，(M-K)×(2P-1)维矩阵EHT(θ)是满列秩的.

2.2 互耦补偿

(16)

(17)

Q′(∶,2:2P-1)v′(2∶P)=-qP

(18)

因此，通过求解上述方程能够得到μk和αk的估计：

v′(2∶P)=-Q′(∶,2∶2P-1)†qP

(19)

(·)†代表了伪逆操作.

vg=

(20)

那么式(21)成立：

(21)

其中：

(22)

(23)

结合式(20)、(23)，互耦参数向量c′能够通过求解式(21)得到

(24)

(25)

通过最小二乘法求解公式(25)，可以得到c′的精确估计为

(26)

上面的方法提供了一个互耦补偿的方式，即当向量c′被确定后，矩阵C可以通过将其元素定位在相应的对角线上形成. 因此，DOA估计能够被获得通过搜索谱函数的峰值：

(27)

DOA估计性能可以被进一步的提升通过重复上述过程. 文中提出的DOA估计与互耦补偿算法总结如下.

(1) 在t=t1,t2,…,tL时刻阵列收到信号x(t)的L个快拍，得到数据矩阵：

X=[x(t1)x(t2)…x(tL)]

(28)

(2) 运用上述的数据矩阵生成协方差矩阵：

Rx=XXH/L

(29)

(3) 对Rx进行特征值分解，得到噪声子空间En.

(7) 用估计的C和式(27)提升DOA估计性能.

(8) 重复第(5)～(7)步可以得到更精确的方向估计.

3 波达方向估计性能对比

考虑两个独立信源从远场入射到ULA，有入射角度θ1=-10°和θ2=20°. 阵列中的传感器间隔为d=λ/2. 有效互耦参数的数量为P=3. 设置c=[10.433 01-0.25i0.141 42-0.141 42i]，它被用在了YE[17]中的第2个仿真，保证了在任意方向θ有g0≠0.

3.1 无互耦补偿下的估计性能比较

第一个实验评估了没有互耦补偿的情况下第(4)步中的DOA估计性能. 假设阵列中阵元数为M=7，快拍数为500. 均方根误差(RMSE)被用于比较不同算法的DOA估计精度，它通常能够被计算为

(30)

RMSE作为信噪比的函数如图1(a)所示. 在低信噪比时，由于用到了整个阵列，本文提出的算法优于文献[7]中的方法. 随着信噪比增加，所有方法DOA估计的RMSE都逐步降低. 当信噪比高于10 dB时，用于矩阵C未知情况下的两种方法有几乎相同的估计准确性. 图1(b)显示了阵列尺寸对于RMSE的影响，信源数M的选择应该满足式(15). 能够看到在M≤10时本文的方法稍优于文献 [17]中的方法.

图1 DOA估计的RMSE作为SNR的函数和DOA估计的RMSE作为阵元数M的函数Fig.1 RMSE of DOA versus SNR and RMSE of DOA versus the array number M

运用初始DOA估计可以得到互耦向量c. 在第二个仿真中，考虑与第一个仿真相同的情况，定义c的RMSE为

(31)

图2 互耦矩阵(MCM)估计的RMSE作为SNR和快拍数的函数Fig.2 RMSE of mutual coupling matrix(MCM) estimates and RMSE of MCM estimates versus snapshots

3.2 有互耦补偿下的估计性能比较

考虑两个信源分别从角度θ1和θ2=θ1+Δθ入射到阵元数为M=7的ULA上，RSN=0 dB，快拍数L=500. 定义：

|θ1-θ1|+|θ2-θ2|<|θ2-θ1|

图3 DOA的可识别率作为角度间隔的函数以及DOA估计的RMSE分别作为SNR、快拍数、阵列尺寸的函数Fig.3 Probability of DOA identification versus the angle interval RMSE of DOA versus SNR,snapshots and the array number

4 结 论

本文解决了存在未知互耦时的DOA估计问题. 该方法基于子空间理论，在没有校准源和辅助传感器的情况下首先用整个阵列估计初始DOA. 假设在空间中不存在盲角，该方法通过从角度相关系数间接估计耦合系数来进行互耦补偿. 最后，在不丢弃任何ADCs的情况下，通过求解最小二乘问题来确定互耦系数. 仿真结果表明，提出的算法在低SNR和小尺寸阵列的情况下可获得更好的性能. 同时，提出方法的鲁棒性也得到了验证.