基于阵列划分的近场DOA估计算法

时间：2024-08-31

陈光辉，曾孝平，焦爽

（1.重庆大学微电子与通信工程学院，重庆 400044；2.中国大唐集团科学技术研究院有限公司中南电力试验研究院，河南郑州 450000）

1 引言

信源定位是阵列信号处理领域中的一个关键技术，其在雷达、声纳、电子对抗、医疗电子等领域中有着广泛的应用［1～3］.根据信源与接收阵列之间的距离，信源定位可以分为远场和近场信源定位［4，5］.对于远场信源定位，信源的传输模型被假设为平面波，仅用到达方向（Direction Of Arrival，DOA）就可以完成信源定位的任务.在过去的几十年中，学者们提出了很多针对远场DOA估计的经典算法.例如最大似然法［6］，多重信号分类（Multiple Signal Classification，MUSIC）算法［7，8］，旋转不变技术（Estimation of Signal Parameters via Rotational Invariance Technique，ESPRIT）方法［9，10］.对于近场信源定位，信源的传输模型被假设为球面波，需要估计信源的DOA和距离才可以完成信源定位的任务.针对近场DOA估计，学者们提出了二维近场MUSIC算法［11，12］，ESPRIT高阶近场估计算法［13，14］，以及一些在MUSIC和ESPRIT算法基础上改进的算法［15～19］.但这些算法在低信噪比（Signal-to-Noise Ratio，SNR）下，因测量数据的微小变化，将出现DOA估计的精度和分辨率严重下降的问题.针对这一问题，本文提出了一种基于阵列划分的近场DOA估计算法.

该算法首先利用四阶累积量和阵列划分构造只包含DOA信息的托普利兹矩阵，并利用一维MUSIC算法对托普利兹矩阵进行特征值分解和构造空间谱.然后，由于MUSIC算法的空间谱在谱峰处是一个断点、且趋近较大值.因此，为提高DOA估计分辨率，该算法求解空间谱的一阶导数来构造一个新“空间谱”.最后，利用协方差矩阵建立二维近场MUSIC空间谱，并将第二步得到的DOA估计值代入二维近场MUSIC空间谱，然后利用本文构造的新“空间谱”逐个估计近场信源的距离.此外，本文理论分析了算法的计算复杂度，并利用仿真实验验证了该算法在低信噪比下有效的提高了DOA的估计精度.

2 近场模型

近场模型如图1所示.其中，2M+1表示接收阵列的总天线个数，N表示接收阵列划分的第一个子阵列的天线个数，信源的数量为K，噪声为高斯噪声.第m个天线的输出表示为：

图1 近场模型的接收阵列

其中，sk(t)表示第k个近场窄带信号源，nm(t)表示加性高斯噪声，τmk表示第m个接收天线与第0个接收天线之间的传输延迟.定义rm表示信源到第m个接收天线的距离.由余弦定理可知：

由图1的近场模型可知，信源传输模型是球面波，因此，τmk被表示如下：

rm被带入到τmk，可得：

其中，λ、d、θk、rk分别表示信源的波长、天线间距、DOA和距离.为简化式（4），f(d/rk)被定义为：

式（5）在变量为0的点被麦克劳林展开，f(d/rk)被重写为：

其中，f1(0)、f2(0)、f3(0)分别表示f(d/rk)在零点处的一阶、二阶、三阶导数.通常d/rk的值比较小，故可以被忽略［20］.因此，τmk被重写为：

其中，wk和φk分别表示为：

由式（8）和式（9）可知，wk仅与近场信源的DOA有关，而φk既与近场信源的DOA有关，也与近场信源的距离有关.

3 本文算法的原理

3.1 DOA估计

对于近场信源的DOA估计，本文提出的算法主要包括：阵元选择、构造托普利兹矩阵和新“空间谱”.阵元选择消除四阶累积量中与距离有关的信息，避免近场信源估计中的二维谱峰搜索，降低算法的计算量；构造托普利兹矩阵可以使矩阵经过特征值分解之后其特征向量相互正交.在这种情况下，特征向量组成的矩阵仅包含一个方向向量相关的信息；构造新“空间谱”的目的是利用MUSIC算法的空间谱在谱峰处是一个断点、且趋近较大值的性质来提高信源估计的分辨率.

3.1.1 阵元选择

其次，近场模型的接收阵列被划分为三个子阵列：[-N，N]、[N+1，M]、[-M，-N-1]，且三个子阵列互不重叠.此外，m，n，p，q的取值是为了使四阶累积量中φk的系数为0，这样四阶累积量只与近场信源的DOA有关，可以有效的避免近场信源估计中的二维谱峰搜索，降低算法的计算量.具体的阵元选择方法如下：

令m∈[-M，-N-1]，n∈[N+1，M]，p∈[-N，0]，q∈[0，N]，且m=-n，p=-q，可以得到一个具有(MN)(N+1)个元素的向量c1.其中，c1的第(M-m)(N+1)+N+1-p个元素表示为：

令m∈[N+1，M]，n∈[-M，-N-1]，p∈[0，N]，q∈[-N，0]，且m=-n，p=-q，可以得到一个具有(MN)(N+1)个元素的向量c2.其中，c2的第(m-N-1)(N+1)+p+1个元素表示为：

令p∈[-N，N]，q∈[N，-N]，m=n=0，p=-q，可以得到一个具有2N+1个元素的向量c3，其中c3的第p个元素表示为：

令p=-N-1，q=N+1，m=-M，n=-N，可以得到一个标量c4，其被表示为：

令p=N+1，q=-N-1，m=M，n=N，可以得到一个标量c5，其被表示为：

由式（12）和（13）可知，c1与c2的倒序共轭对称.由式（14）可知，c3的2N+1个元素共轭对称.由式（15）和（16）可知，c4和c5也共轭.c1，c2，c3，c4，c5的共轭对称关系是构造托普利兹矩阵的基础.

3.1.2 托普利兹矩阵构造

构造c1，c2，c3，c4，c5之后，需要构造一个托普利兹矩阵.因为托普利兹矩阵经过特征值分解之后其特征向量相互正交.在这种情况下，特征向量组成的矩阵仅包含一个与方向向量相关的信息.具体的托普利兹矩阵构造过程见下文.

首先，一个(2(M-N+1)(N+1)+1)维的向量c被定义为：

其中，T表示矩阵转置.

其次，定义一个矩阵CH，把向量c中第((M-N+1)(N+1)+2-m)到(2(M-N+1)(N+1)+2-m)的元素组成CH的第m列.其中，CH的具体表达式如下所示：

由c1，c2，c3，c4，c5的共轭对称关系知，托普利兹矩阵CH的维度是((M-N+1)(N+1)+1)×((M-N+1)(N+1)+1).为简化式（18），定义Q=(M-N)(N+1).此外，式（18）可以被表示为矩阵形式，如下所示：

其中，H表示矩阵共轭转置.A(θ)和CS公式如下所示：

基于式（19）～（21），A(θ)可以看作是K个远场信源的虚拟方向矩阵，且其对应(M-N+1)(N+1)+1天线（虚拟天线）.另外，因为sk(t)统计独立，且CS和A(θ)均为满秩矩阵，因此CH也是满秩矩阵.所以，可以采用一维远场MUSIC算法估计近场信源的DOA.其中，基于一维远场MUSIC算法，a(θ)与Un之间的关系可以描述为［8］：

其中，Un表示由CH矩阵特征值分解之后得到的噪声空间.因此，MUSIC算法的空间谱函数定义为：

3.1.3 构造新“空间谱”

MUSIC算法空间谱的物理意义表示为：P(θ)不是真实的谱，可以理解为扫描向量a(θ)到噪声子空间Un的欧几里德距离的倒数.由MUSIC算法空间谱的物理意义和aH(θ)Un=0可知，P(θ)函数在信源的DOA处存在尖峰（也可以认为是断点）.利用此性质，对P(θ)求一阶导数，锐化谱峰，可以提高DOA估计的分辨率.因此，P(θ)的一阶导数d(θ)定义为：

为方便计算，选择合适的步长h可以用一阶差分函数y(θ)代替一阶导数d(θ)，其中一阶差分函数y(θ)被定义为：

其中，式（25）即是本文算法定义的新“空间谱”的表达式.

3.2 距离估计

为了估计近场信源的距离，本文算法又构造一个协方差矩阵R=E{y(t)y*(t)}，对其特征值分解可以构造一个二维空间谱，具体的表达式如下：

其中，U0表示协方差矩阵R=E{y(t)y*(t)}经过特征值分解之后的噪声子空间.在3.1节中获得的DOA估计值{θk，k=1，…，K}被代入式（26）中，近场信源的距离可以被一一估计.因此，式（26）可以被重写为：

在估计近场信源的距离时，同样利用本文构造的新“空间谱”估计近场信源的距离，公式如下：

3.3 本文算法的步骤

本文算法的步骤如表1所示：

表1 本文算法的步骤

4 仿真分析

4.1 DOA估计分辨率理论分析和仿真

本节首先分析本文算法定义的新“空间谱”如何确定近场信源的DOA，然后分析本文算法的DOA估计分辨率.假设θ0是近场信源的DOA，由式（22），（23）和MUSIC算法的空间谱的物理意义可知，P(θ0)是一个较大的正值，P(θ0-h)，P(θ0+h)，P(θ0+2h)远小于P(θ0).因此，可以计算出y(θ0-h)是一个较大的正值，y(θ0+h)是一个较大的负值.所以，根据式（25）中y(θ)函数是否先有一个较大的正值，后有一个较大的负值来确定近场信源的DOA.此外，本文算法通过对空间谱函数求一阶导数，一定会先产生一个极大值，然后产生一个极小值.这反映到谱峰上是先产生一个正的谱峰，后又产生一个负的谱峰.由于本文算法在求空间谱一阶导数的过程中，采用从-90°到90°的递增顺序，所以在谱峰搜索时会舍弃搜索到的负谱峰，这样即不会对DOA估计结果造成影响，也不会造成DOA估计分辨率下降.当θ1和θ2是两个近场信源的DOA时，MUSIC算法的DOA估计分辨率定义为［7］：

其中，

基于式（29），本文算法的DOA估计分辨率为：

由式（31）可知，本文算法DOA估计分辨率近似等于MUSIC算法的一阶导数.由MUSIC算法的空间谱函数的物理意义可知，空间谱函数在DOA处是断点，且其在DOA处的值远大于相邻角度的值.因此，本文算法DOA估计分辨率φ(θ)高于MUSIC的估计分辨率φmusic(θ).此外，本文又通过实验验证了上述分析，实验结果如图2和图3所示.

由图2（a）可知，在DOA=（10°，13°），SNR=0 dB条件下，传统的MUSIC方法未估计出两个DOA，仅仅在角度等于10°处估计出一个谱峰，且不明显.通过对5°～20°范围的空间谱放大处理，显示在图2（a）的子图中，发现在角度等于13°处的空间谱的谱峰也不明显.但由图2（b）可知，在DOA=（10°，13°），SNR=0 dB条件下，本文算法可以在角度等于10°和13°处估计出两个DOA，且空间谱的谱峰明显.因此，本文算法在SNR=0 dB条件下，可以有效的提高DOA的估计分辨率.

由图3（a）可知，在DOA=（10°，13°），SNR=20 dB条件下，传统的MUSIC方法估计出两个DOA，且空间谱的谱峰较明显.由图3（b）可知，在DOA=（10°，13°），SNR=20 dB条件下，本文算法也可估计出两个DOA，且空间谱的谱峰更明显.但由于MUSIC算法在高SNR条件下，DOA估计分辨率已经较高，所以本文算法在SNR=20 dB条件下，DOA估计分辨率和MUSIC算法是相似的.最终，由图2和图3可知，本文算法提高了低信噪比下的DOA估计分辨率.

图2 DOA=(10°,13°),信噪比SNR=0 dB的对比图

图3 DOA=(10°,13°),信噪比SNR=20 dB的对比图

4.2 DOA和距离估计分辨率仿真

仿真参数如表2所示，噪声为高斯白噪声，分别利用3个不同的近场信源对比了本文算法和文献［19］中的算法.其中，3个近场信源分别由K1、K2、K3表示，K1是（［0°，3λ］，［1°，4λ］），K2是（［10°，3λ］，［13°，4λ］），K3是（［20°，3λ］，［25°，4λ］），DOA角度差Δθ分别为1°、3°、5°.DOA和距离估计结果如图4～6所示.

②Stufflebeam，D．L．，“The CIPP Model for Evaluation”，International Handbook of Educational Evaluation，Springer Netherlands，2003．

表2 仿真参数表

在图4中，当Δθ为1°和3°时，文献［19］的算法只估计出一个DOA.当Δθ为5°时，文献［19］算法估计两个DOA.在图5中，当Δθ为1°时，本文算法估计出一个DOA，当Δθ为3°和5°时，本文算法估计两个DOA，且空间谱的谱峰较明显.原因是当信源的DOA相近时，信源之间的空间相关性较强.但本文算法利用空间谱函数在DOA处是断点并具有较大值的物理特性和空间谱函数的一阶导数构造一个新“空间谱”来减少DOA处的空间相关性，因此提高了DOA估计分辨率.

图4 文献[19]中DOA估计结果图

图5 本文算法的DOA估计结果图

图6（a）～（c）是文献［19］的距离仿真结果，图6（d）～（f）是本文算法的距离仿真结果图.图6（a）、图6（b）、图6（d）中仅有一个空间谱，而图6（c）、图6（e）、图6（f）中有两个空间谱.原因是文献［19］对K1和K2近场信源只估计出一个DOA，对K3近场信源估计出两个DOA，而本文算法对K1近场信源也只估计出一个DOA，而对K2和K3近场信源估计出两个DOA.

图6 文献[19]和本文算法的距离估计结果图

4.3 DOA和距离估计精度仿真

统计意义的均方根误差被定义为DOA的估计精度.其中，近场信源的均方根误差被定义为：

其中，T是蒙特卡洛实验的次数（T=500），K是近场信源的个数，βi，j是第j次蒙特卡洛实验中αi，j的估计值.本文的算法和文献［13，16，19］中的方法进行了比较.同时为了更好的验证本文算法的性能，文献［21］中CRB方法也被比较，结果如图7～10所示.

在图7和图8中，随着信噪比的增加，均方根误差逐渐减少，近场信源的DOA和距离估计精度被逐渐提高.对于DOA估计，文献［16，19］和本文算法的均方根误差低于文献［13］中的ESPRIT-Like算法，但都高于文献［21］中CRB方法.当信噪比低于10 dB时，本文算法的均方根误差低于文献［16，19］和ESPRIT-Like算法.当信噪比高于10 dB时，本文算法和文献［16，19］的均方根误差是相似的.因为在低信噪比时，接收阵列的统计矩阵的托普利兹结构被破坏.但本文算法首先利用阵列划分和四阶累积量构造一个新的托普利兹矩阵，同时又利用空间谱函数的一阶导数构造一个新“空间谱”来减小DOA处的空间相关性.因此，低信噪比下，本文算法提高了近场信源的DOA的估计精度.

图7 DOA的均方根误差与信噪比关系图

对于距离估计，当信噪比低于10 dB时，本文算法的均方根误差低于文献［13，19］，相似于文献［16］，当信噪比高于10 dB时，文献［16］的均方根误差低于本文算法和文献［13，19］.因为对于距离估计，文献［19］和本文算法的接收阵列的统计矩阵由协方差矩阵构造的，而文献［16］的接收阵列的统计矩阵由四阶累积量构造的.众所周知，四阶累积量矩阵可以有效的抵抗高斯噪声，而协方差矩阵不具有此特性.所以，文献［16］在高信噪比下的均方根误差低于本文算法和文献［13，19］.此外，本文算法同样采用新“空间谱”估计近场信源的距离，这可以提高低信噪比下距离的估计精度.因此，在低信噪比下，本文算法的均方根误差相似于文献［16］，低于文献［13，19］.

在图9和图10中，随着采样快拍数的增加，均方根误差逐渐减少，近场信源的DOA和距离估计精度被逐渐提高.对于DOA估计，本文算法的均方根误差低于文献［13，16，19］，但均高于文献［21］CRB方法.对于距离估计，本文算法的均方根误差低于文献［13，19］，相似于文献［16］，但均高于文献［21］中CRB方法.原因是图9和10的仿真结果是在信噪比SNR=0 dB条件下得到的，由图7可知，低信噪比下本文算法的DOA估计的均方根误差低于文献［13，16，19］，由图8可知，低信噪比下本文算法的距离估计的均方根误差相似于文献［16］，低于文献［13，19］.因此，随着采样快拍数的增加，在低信噪比下，本文算法的DOA估计精度高于文献［13，16，19］，本文算法的距离估计精度相似于文献［16］，高于文献［13，19］.

图9 DOA的均方根误差与采样快拍数关系图

图10 距离的均方根误差与采样快拍数关系图

4.4 计算复杂度分析

在计算复杂度分析中，考虑了几个主要步骤：统计矩阵的构建，特征值分解和谱峰搜索.为了更好的比较本文算法的计算复杂度，对比了本文算法和文献［11，13，15，16，19］中算法，具体的计算复杂度结果如表3所示.

在表3中，M、K、L分别表示接收阵列的天线个数，信源个数和采样快拍数.ng表示在DOA空间搜索的次数，ng=180/Δ，其中Δ是搜索步长，Δ=0.01，nl表示在距离空间搜索的次数其中D是阵列中天线的孔径，λ为波长.

表3 计算复杂度

本文算法的计算复杂度高于文献［13］的ESPRITLike方法，相近于文献［16］中的方法，低于文献［11］中的2D-MUSIC算法和文献［15，19］中的算法.原因是文献［13］中的ESPRIT-Like算法不需要谱峰搜索，所以计算量是低于其他算法的，但从图7和图8可知，ESPRITLike算法的DOA和距离估计精度是远低于其它算法的；文献［16］中算法只构造了一个四阶累积量，但增加了K次施密特正交化，所以计算量和本文相似，但从图7可知，文献［16］的DOA估计精度在低信噪比下低于本文算法；文献［11］需要二维谱峰搜索，文献［15，19］构造的四阶累计量的维度高于本文算法，而本文算法不需要二维谱峰搜索，只需要K+1次一维谱峰搜索（其中，一次谱峰搜索是对DOA估计，K次谱峰搜索是对距离估计）.