基于改进TLS-ESPRIT算法的相干定位信号波达方向空间平滑

胡茂厅，杨福兴，邓中亮，刘 鋆

(北京邮电大学智能通信、导航与微纳系统实验室，北京 100876)

0 引言

在波达方向(Direction of Arrival，DOA)估计算法中，常见的算法有很多种，例如多重信号分类(Multiple Signal Classification,MUSIC)算法、最大似然(Maximum Likelihood，ML)算法、最大熵(Maximum Entropy，ME)算法、旋转不变信号参数估计技术(Estimating Signal Parameter via Rotational Invariance Techniques，ESPRTIT)算法等，算法的优良性决定了DOA估计精度，进一步决定了角度定位的误差。

其中ESPRIT算法[1]最早由Roy、Paulra和Kailath等提出。由于在参数估计等方面的优势，ESPRIT算法应用广泛。该算法主要研究由阵列的移不变特性而引起的信号子空间的旋转不变性，利用旋转因子估计DOA。与DOA估计算法中常见的MUSIC算法相比，ESPRIT算法不需要构造空间谱，减少了计算量，鲁棒性强，分辨率高；但在实际应用中只能用阵列输出采样值的相关函数来估计信号输出的协方差，不能用在根据阵列的移不变性分解得到的子阵长成的相似的信号子空间计算中[2]。后来学者在ESPRIT算法的基础上，又提出了基于最小二乘(Least-Squares，LS)准则，引入一个与子阵张成子信号空间正交的矩阵以解决上述问题。后来在LS-ESPRIT算法的基础上提出了基于总体最小二乘(Total Least-Squares，TLS)准则的TLS-ESPRIT算法,提高了估计准确性[3]。

考虑由于传播环境的复杂性，入射到实际阵列中的强相关或相干源信号是普遍存在的。而且对于强相关或者相干信号，这些基于特征分解类算法都无法有效进行DOA估计，因此研究有效的解相干算法是当前阵列信号处理领域的一个重要的研究方向[4-5]。

常见的解相干处理方法有最大似然方法、Toeplitz矩阵重构类算法及空间平滑类预处理方法等。最大似然方法[6]是通过概率密度模型来进行相干源信号处理的，具有较好的解相干效果，但是它最终需要非线性多维搜索来实现DOA估计，计算量巨大。Toeplitz矩阵重构类算法[7]的最大优点是没有阵列孔径的损失，但是它在矩阵重构过程中没有充分考虑到信号的先验信息，所以在信号源功率不相同的场合，其估计精度会相对变差。空间平滑类方法最早是由Evans[8]等提出，随后经Shan等[9]改进，成为一种常用且有效的解相干预处理方法。这种解相干方法是以牺牲阵列的有效孔径为代价，为了减少阵列孔径的损失，Pillai等[10]在前人的基础上提出了一种双向空间平滑算法，即前后向平滑算法。旺晋宽等[11]在ESPRIT算法的基础上利用空间平滑进行DOA估计(M-ESPRIT算法)，解决了相关信号源无法进行DOA估计问题，但其估计DOA的准确性受外界参数影响较大，准确性较低，且平均误差较大。

针对上述问题，本文将双向空间平滑方法与TLS-ESPRIT算法相结合，提出了一种基于改进TLS-ESPRIT算法的空间平滑方法(SS-TLSESPRIT)。该方法首先对输出相干信号源进行空间平滑，然后把计算得到的子阵协方差应用在TLS-ESPRIT算法，使其最大限度利用信号子空间的信息。

1 阵列信号数学模型

设有一个天线阵列，它由M个具有任意方向性的阵元按等间距排列，构成K个具有相同中心频率ω0、波长为λ、阵元间距为d=λ/2的空间窄带平面波(M>K)，分别以角θ1,θ2,θ3…，θk入射到该阵列，如图1所示，则第M个阵元的输出为

nM(t)σ2

(1)

其中，sk(t)为入射到阵列的第k个源信号；nm(t)为第M个阵元的加性噪声，方差为σ2；θk为第k个源信号的DOA；λ为信号源波长。定义θk的 方向向量为

图1 信号接收模型Fig.1 Signal reception model

(2)

则阵列输出信号可写成如下矩阵形式

X(t)=AS(t)+G(t)

(3)

式中

X(t)=[x1(t),x2(t),…,xM(t)]T
A=[a(θ1),a(θ2),…,a(θk)]
S(t)=[s1(t),s2(t),…,sk(t)]T
G(t)=[g1(t),g2(t),…，gM(t)]T

根据式(3)可得阵列输出信号协方差矩阵为

R=E[X(t)XH(t)]=ARsAH+σ2IM

(4)

式中：Rs=E[S(t)SH(t)] 为信号协方差矩阵；IM为M阶单位阵。

2 本文方法

2.1 空间平滑技术

空间平滑技术(Spatial Smoothing Techni-ques, SST)是处理强相关或相干信号的有效方法，其基本思想是将具有相同子阵列流型的均匀线阵分成若干个相互重叠的子阵列，然后再将各子阵列的协方差矩阵相加平均后以取代原来意义上的协方差矩阵，从而使得矩阵的秩恢复到等于信号源的个数[9-10]。

如图2所示，有一个由M个阵元构成的线阵，利用前向滑动方式划分为L个子阵列，其中每个子阵列有N个阵元，即N=M-L+1。取最左边的子阵阵列为参考阵列，则第l个前向子阵列的输出为

=ANDl-1s(t)+nl(t),1≤l≤L

(5)

图2 空间平滑示意Fig.2 Spatial smoothing representation

其中，AN为N×K维的方向矩阵，其列为N维的导向矢量aN(θi)(i=1,2,…,K)。

(6)

式中：IN表示N阶单位阵。

前向空间平滑技术是通过求取各前向子阵列的协方差矩阵的均值来实现秩恢复的，即前向空间平滑协方差矩阵为

(7)

同理，可得后向空间平滑协方差矩阵为

(8)

前向空间平滑协方差矩阵Rf和后向空间平滑协方差矩阵Rb实质上互为共轭倒序阵，利用共轭倒序不变性，可得前后向空间平滑协方差矩阵为

Rfb=(Rf+Rb)/2

(9)

2.2 改进的TLS-ESPRIT(SS-TLSESPRIT)算法

2.2.1 ESPRIT算法模型

根据上文建立的线性信号模型，把阵列分为2个平移量为Δ的子阵Zx和Zy。2个子阵阵元的输出信号可分别表示为(在这里设子阵阵元数为M个)

(10)

nyi(t),i=1,2,…,M

(11)

其中，sk(t)为子阵Zx和Zy接收到的第k个信号;θk为第k个信号到达方向;ai(θk)为第i个阵元对第k个信源的响应;c为电波的传导速度，nxi(t)和nyi(t)为子阵Zx和Zy的第i个阵元上的加性噪声。2个子阵列的每个阵元在t时刻的输出信号矢量表达式为

x(t)=As(t)+nx(t)

(12)

y(t)=AΦs(t)+ny(t)

(13)

其中

x(t)=[x1(t),x2(t),…,xM(t)]T
y(t)=[y1(t),y2(t),…,yM(t)]T
s(t)=[s1(t),s2(t),…,sK(t)]T
nx(t)=[nx1(t),nx2(t),…,nxM(t)]T
ny(t)=[ny1(t),ny2(t),…,nyM(t)]T

A=[a(θ1),a(θ2),…,a(θk)]称为方向矩阵。Φ为K×K的对角矩阵，称为旋转算子，由于子阵的移不变特性形成了2个子阵信号的旋转不变性，即可以把子阵Zx和Zy联系起来，Zy可以由Zx乘以一个旋转算子得到。其对角元素为K个信号在任意一个阵元偶之间的相位延迟，表示为

Φ=diag{exp(jμ1),…,exp(jμk)}

(14)

其中

μk=ω0Δsinθk/c(k=1,2,…,K)

把2个子阵加以合并得

(15)

取t=t1,t2,…,tN的N个时刻的快拍组成2M×N数据矩阵，则式(15)可表示为

(16)

其中

Z=[z(t1),z(t2),…,z(tN)]
S=[s(t1),s(t2),…,s(tN)]
Nz=[nz(t1),nz(t2),…,nz(tN)]

2.2.2 信号输出相关矩阵双向平滑

阵列输出信号矩阵Z的自相关矩阵为

(17)

自相关矩阵Rzz的特征分解为

(18)

其中，Es为信号子空间，EN为噪声子空间。由于本文采用的是相干信号，受到相干信号的影响，Rzz奇异，即非满秩矩阵。Rzz特征值中较大的特征值数K1即信号子空间的维数是小于信源数K个的。这样无法进行TLS-ESPRIT算法计算。

解决方法：采用2.1节中式(9)双向空间平滑技术，对阵列进行平滑，这里选用平滑的子阵阵元数N一定要大于K值(N>K)，即可得出平滑后的子阵自相关矩阵Rzzz，Rzzz中较大特征值数就达到K个值，即λ1≥…≥λK>λK+1=…λ2M=σ2，Rzzz实现满秩。这样就可以完成接下来的计算工作。

2.2.3 平滑后改进TLS-ESPRIT算法解算

设平滑后阵列信号输出矩阵的协方差矩阵Rzzz中K较大的特征值对应的特征向量为Ess，其中，存在唯一非奇异的K×K维满秩矩阵T，使得式(19)成立

(19)

根据信号阵列的移不变特性，可以将Ess分解为2个部分，即EX∈CN×K和EY∈CN×K，分别对应子阵Zx和Zy(这里的子阵是ESPRIT基础算法模型中的2个子阵)，即

(20)

从式(20)中可以得出

EY=EXT-1ΦT=EXΨ

(21)

其中，Ψ=T-1ΦTΨ。至此可知，EX和EY张成相似的子空间，且矩阵Φ的对角线元素为Ψ的特征值。

(22)

(23)

所以，根据TLS准则，算法通过解以下最小值问题来获得Ψ的最小二乘解。

定义矩阵F

(24)

使得V值最小

(25)

且满足

FHF=I

(26)

定义了矩阵F，则Ψ的估计由式(27)给出

(27)

2.3 算法步骤

算法SS-TLSESPRIT的步骤归纳如下：

1)建立相干定位信号模型s(t)，加入高斯白噪声，方差设为σ2，均值为0。

2)根据TLS-ESPRIT算法特性，划分2个子阵，确定子阵的阵元数要大于信源数，这是为了确保第3步的实现。

4)计算

将E分解

3 仿真验证

仿真1 对于信源为相干信号TLS-ESPRIT与平滑SS-TLSESPRIT算法的DOA估计对比

仿真中，采用阵元间距为半波长λ/2、阵元数M=12的均匀线阵，平滑子阵阵元为8，信号源为相干信号且信源数已知。信源数分别取2、3、4个，DOA分别是(0° 30°)、(-30° 0° 30°)和(-30° 0° 30° 50°)。仿真快拍数为500，信噪比取10dB。对于DOA分别用经典TLS-ESPRIT算法和SS-TLSESPRIT算法进行估计，如表1所示，为了真实看到DOA估计，取一次仿真实验结果。

表1 不同信源和不同DOA估计算法对比Tab.1 Comparison of different sources and different DOA estimation algorithms

由表1可以看出，采用没有平滑的TLS-ESPRIT算法不能对相干信号进行DOA估计；而加入平滑改进的SS-TLSESPRIT算法可以实现任意信号包含相干信号的DOA估计，具有有效性和高精度性，DOA估计误差达到要求。

仿真2 2个相干信源SS-TLSESPRIT算法的DOA估计相关参数效果仿真

1)为了验证算法的稳定性，在不同信噪比条件下对实验进行1000次Monte Carlo仿真，并采用均方根误差(Root Mean Square Error, RMSE)来检验算法的有效性，RMSE值如式(28)所示。仿真数据：低信噪比和高信噪比都不易于工程实现，所以信噪比取-10dB～30dB之间，2个相干信源DOA取(0° 30°)，快拍数采用中等大小，取值500，如图3所示。

(28)

图3 DOA估计均方根误差Fig.3 Root mean square error of DOA estimation

2)再次改变快拍数，观察随着快拍数的改变，DOA估计RMSE随着信噪比的变化曲线，以此来验证算法的优越性。快拍数分别取250、500和750，如图4所示，可以看出，在相同信噪比情况下，随着快拍数增加，DOA估计RMSE降低，说明了在实际应用中，DOA估计精度与采样数有正向关系，增加一定量的快拍数可以降低DOA估计误差，提高定位精度。

图4 快拍数与DOA估计关系图Fig.4 Quick shot number and DOA estimation diagram

3)改变仿真信源DOA间隔，观察算法的DOA分辨能力。在不同DOA间距下，分别进行了1000次Monte Carlo仿真，观察仿真效果，如表2所示。

表2 DOA间隔与RMSE之间的关系Tab.2 Relationship between DOA interval and RMSE

由表2可以看出，随着DOA间隔增大，DOA估计精度提高，达到实验要求；但是在DOA间距在8°以下时，DOA估计误差开始增大，DOA估计成功概率降低。

值得一提的是，在实际应用中，信号多为混合信号，当信号中同时含有相干信号和独立信号时，利用本文算法进行仿真，对于混合信号DOA的估计精度有很大的提高，精度提高到0.1°。

仿真3 SS-TLSESPRIT算法与M-ESPRIT算法的DOA估计对比[11]

在仿真2基础上做算法对比，同时还是以RMSE值为评估参数，1000次Monte Carlo仿真，快拍数500，DOA取(0° 30°)，如图5所示。SS-TLSESPRIT算法与M-ESPRIT算法对比，定位精度都得到了提高，特别是在低信噪比条件下本文算法提高效果明显。

图5 SS-TLSESPRIT与M-ESPRIT对比Fig.5 SS-TLSESPRIT vs. M-ESPRIT

从以上仿真结果中可以看出，算法的DOA估计准确性达到室内角度定位要求，与其他算法相比，定位精度得到了提高，能够应用到单角度定位或者融合定位中去[12]，且DOA估计的稳定性和鲁棒性都大大提高。

4 结论

对于相干定位信号DOA估计，提出了一种基于改进TLS-ESPRIT算法的双向空间平滑方法。算法分析与实验结果表明：

1)在TLS-ESPRIT算法的基础上实现了相干定位信号的空间平滑。通过对算法中子阵阵列平滑，使输出定位信号的协方差矩阵实现满秩，平滑了秩亏损，解决了相干信号条件下基础算法无法准确估计DOA的问题，加上改进TLS-ESPRIT算法，使得DOA估计的准确性大大提高。

2)SS-TLSESPRIT算法的各项仿真结果都达到理想效果。通过RMSE值对算法进行评估，再加上不同参数改变，都体现了算法的稳定性和良好性。稳定测量到精确的DOA值，对于室内角度定位提供了良好的数据基础，为室内角度精确定位或者与TOA/TDOA/RSSI融合定位奠定了基础。

3)本文提出的平滑算法在仿真效果方面，DOA分辨性能不够，当相干信号源之间的角度较小时，算法会出现DOA估计分辨力下降；特别是复杂环境下，定位信号分辨力如果不够，则更需要其他融合定位，因此还需待进一步完善。