大规模天线阵列信号的波达方向单快拍快速估计方法及装置与流程

本发明涉及移动通信技术领域，具体涉及一种大规模天线阵列信号的波达方向单快拍快速估计方法及装置。

背景技术：

阵列信号处理是指对在空间按一定方式排列的传感器阵列所接收的信号进行处理，增强有用信号，抑制无用信号和噪声，并从接收信号中提取有用的信号参数和特征信息(频率、角度、时延和极化等等)。相比传统的单定向传感器，传感器阵列可获得灵活的波束控制、较高的信号增益、极强的干扰抑制能力和更高的空间分辨力等。作为现代信号处理的一个重要分支，阵列信号处理已在通信、雷达、声呐、地震勘探和射电天文等领域取得了广泛应用和迅速发展。信号的波达方向(Angle of Arrival,AOA)估计是阵列信号处理的一个基本问题，也是雷达、声呐、通信等领域的重要任务之一。

长期以来，学者针对阵列信号的AOA估计展开了广泛而深入的研究，其中尤以以多重信号分类算法(Multiple Signal Classification,MUSIC)和借助旋转不变性的信号参数估计算法(Estimation of Signal Parameters via Rotational Invariance Techniques,ESPRIT)为代表的子空间算法最为引人瞩目。子空间算法建立在这样一个基本观察之上：若传感器个数大于信源个数，则阵列接收数据的信号分量一定位于一个低秩的子空间；在一定条件下，这个子空间将唯一确定信号的波达方向，并可以使用数值稳定的奇异值分解精确地估计波达方向。

目前大规模阵列下的AOA估计研究已经引起了学术界的广泛重视。顾名思义，大规模阵列是指阵元数极大的阵列，该方向的研究热潮起始于移动通信中大规模多输入多输出(Massive Multiple-input-multiple-Output,Massive MIMO)系统概念的提出。Massive MIMO系统在基站侧配置数十根甚至数百根的天线组成天线阵列,形成多个并行的数据传输通道,进而实现在相同的时频资源内同时为多个用户提供服务。相比较现有的小规模MIMO技术，大规模MIMO技术在涵盖其所有优点的同时，得益于空间自由度和阵列孔径的增加，大规模MIMO技术实现了对空间维度资源的深度挖掘,大幅度提高了频谱效率。其利用最简单的线性预编码和线性检测技术就可获得最优信号处理技术所得到的性能,消除噪声和同小区内的干扰。同时可提高天线阵列估计信号波达方向的精度,波束可集中在很窄的范围内,从而明显地降低干扰；采用一定的信号处理技术后,基站侧的每根天线的发射功率可正比于天线数平方根的倒数,即天线数量越多，总的发射功率就越小,从而节约能耗,实现绿色通信。此外，Massive MIMO技术在硬件实现上更为方便。

传统的信道估计方法通过估计信道矩阵来实现，由于Massive MIMO系统信道矩阵维数巨大，传统的信道估计算法难以实现实时处理。针对此现状，通过估计传播路径的角度信息来实现信道估计的方式被学术界视作一种可行的Massive MIMO系统下的信道估计方案。现有的结果表明，这种方式在大大减少了待估计量的同时，在误差界方面能提供最佳的估计性能。正因为如此，大规模阵列下的AOA估计问题受到了学术界的广泛重视。

传统的子空间算法需要通过特征值分解的形式获得对子空间的估计，而其中特征值分解的复杂度与阵列阵元数的三次方成正比，在大规模阵列的条件下，特征值分解将耗费巨大的计算复杂度。同时，子空间算法需要采集大量的快拍数据才能获得对子空间的准确估计，因此该类算法在少快拍、低信噪比情形下的鲁棒性欠缺。综合考虑以上两点因素，即使如ESPRIT这类具有闭式解(因而具有低复杂度)的子空间算法，其在大规模阵列下的实际应用能力仍不理想。此外，传统AOA估计算法并未考虑大规模阵列相比较传统小规模阵列可获得的新性质。

基于以上背景，适用于大规模阵列下的、低复杂度的、高精度AOA估计算法具有很高的应用价值。

技术实现要素：

针对现有技术的不足，本发明提供一种大规模天线阵列信号的波达方向单快拍快速估计方法及装置，通过DFT这一简单有效的方法，仅需要单次快拍的数据就能实现对于波达方向的准确估计。

本发明提出一种大规模天线阵列信号的波达方向单快拍快速估计方法，包括：

步骤一：对阵列接收信号进行离散傅里叶变换，搜索变换后的信号功率谱，由谱峰位置得到信号波达方向的初始估计；

步骤二：将所述阵列接收信号进行相位旋转，对所述旋转后的信号进行离散傅里叶变换，并获取相位最佳旋转因子；

步骤三：根据所述波达方向的初始估计和所述相位最佳旋转因子，计算出所述波达方向的最终估计值。

所述相位最佳旋转因子使得所述离散傅里叶变换后功率谱最集中。

所述天线阵列为均匀线阵或均匀面阵，当为均匀面阵时，所述离散傅里叶变换为二维离散傅里叶变换。

所述相位旋转具体包括：

步骤一：定义长度为(C+1)×1的矢量Λ，矢量Λ第c个元素为其中，矢量Λ中的元素为旋转因子，所述旋转因子在内均匀分布，M为全向天线阵元的个数；

步骤二：对旋转因子[Λ]_c构建对角矩阵得到相位旋转矩阵，将所述阵列接收信号右乘所述相位旋转矩阵得到相位旋转后的信号。

所述获取相位最佳旋转因子具体包括：

步骤一：对旋转后的信号离散傅里叶变换功率谱内每一点，设置一定长度的滑动窗，设置一定长度的滑动窗，并计算该点功率谱与所述滑动窗内所有点功率谱和的比值，将所述比值定义为该点的稀疏指数；

步骤二：对所述阵列接收信号进行C+1次旋转，对每一次旋转后的离散傅里叶变换功率谱，计算每一点的稀疏指数并记录，得到一(C+1)×M稀疏指数矩阵Q；

步骤三：根据所述稀疏指数矩阵Q中初始估计值谱峰位置附近的最大值，获取对应的相位最佳旋转因子。

本发明还提供了一种大规模天线阵列信号的波达方向单快拍快速估计装置，包括：

初始估计模块，用于对阵列接收信号进行离散傅里叶变换，搜索变换后的信号功率谱，由谱峰位置得到信号波达方向的初始估计；

相位获取模块，用于将所述阵列接收信号进行相位旋转，对所述旋转后的信号进行离散傅里叶变换，并获取相位最佳旋转因子，；

最终估计模块，用于根据所述波达方向的初始估计和所述相位最佳旋转因子，计算出所述波达方向的最终估计值。

所述相位获取模块通过如下方式进行所述相位旋转：

定义长度为(C+1)×1的矢量Λ，矢量Λ第c个元素为其中，矢量Λ中的元素为旋转因子，所述旋转因子在内均匀分布，M为全向天线阵元的个数；

对旋转因子[Λ]_c构建对角矩阵得到相位旋转矩阵，将所述阵列接收信号右乘所述相位旋转矩阵得到相位旋转后的信号。

所述相位获取模块通过如下方式获取所述相位最佳旋转因子：

对旋转后的信号离散傅里叶变换功率谱内每一点，设置一定长度的滑动窗，设置一定长度的滑动窗，并计算该点功率谱与所述滑动窗内所有点功率谱和的比值，将所述比值定义为该点的稀疏指数；

对所述阵列接收信号进行C+1次旋转，对每一次旋转后的离散傅里叶变换功率谱，计算每一点的稀疏指数并记录，得到一(C+1)×M稀疏指数矩阵Q；

根据所述稀疏指数矩阵Q中初始估计值谱峰位置附近的最大值，获取对应的相位最佳旋转因子。

本发明提供一种大规模天线阵列信号的波达方向单快拍快速估计方法及装置，充分考虑了大阵列孔径的特性，利用均匀大规模天线阵列下离散傅里叶变换DFT功率谱的精确性，通过DFT这一简单有效的方法，仅需要单次快拍的数据就能实现对于波达方向的准确估计，其性能逼近克拉美罗界这一性能上限。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明的阵列模型图；

图2为本发明一实施例提供的大规模天线阵列信号的波达方向单快拍快速估计方法的流程示意图；

图3为本发明一实施例提供的理想无泄漏情形下阵列流型的DFT和DTFT功率谱示意图；

图4为本发明一实施例提供的存在功率泄漏情形下阵列流型的DFT和DTFT功率谱示意图；

图5为本发明一实施例提供的大规模天线阵列信号的波达方向单快拍快速估计方法的算法流程示意图；

图6为本发明一实施例提供的未经相位旋转直接DFT变换得到的功率谱示意图；

图7为本发明一实施例提供的相位旋转后信源一获得最佳旋转时的DFT功率谱示意图；

图8为本发明一实施例提供的相位旋转后信源二获得最佳旋转时的DFT功率谱示意图；

图9为本发明实施例中的提供的大规模天线阵列信号的波达方向单快拍快速估计装置的结构示意图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

图1示出了本发明的阵列模型图。如图1所示，以M个全向天线阵元构成的均匀线阵(Uniform Linear Array,ULA)为例,其阵元间隔为入射信号波长的一半。然而所述方法可以被直接移植到均匀面阵(Uniform Rectangular Array,URA)。M的典型取值大于100。假设阵元彼此独立，不存在互耦现象。空间内有K个远场、窄带、非相干信源，其发射的信号s＝[s₁，…，s_K]^T以平面波的形式照射到阵列，入射角度分别为θ＝[θ₁，…，θ_K]，其中θ_k∈[-π/2，π/2]，k＝1，…，K。取最左侧的阵元为参考点，此时阵列接收信号y可以表示为

y＝As+n

其中，A为阵列流型矩阵，其具体表达式为

n为零均值的加性高斯白噪声。

图2示出了本发明一实施例提供的大规模天线阵列信号的波达方向单快拍快速估计方法的流程示意图，如图2所示，本实施例的大规模天线阵列信号的波达方向单快拍快速估计方法包括如下步骤。

101、对阵列接收信号进行离散傅里叶变换，搜索变换后的信号功率谱，由谱峰位置得到信号波达方向的初始估计。

本实施例中，初始估计的得到是基于以下一基本发现，即大规模阵列的接收信号经DFT变换后，其功率谱具有稀疏特性，谱峰位置与信号的入射角度有关。下面将就这一发现做具体说明。

定义归一化的DFT矩阵F，其第(p，q)个元素为对第k个信源的方向矢量a(θ_k)进行DFT变换可得其第q个元素的表达式为：

取的第q个元素，其幅值满足以下关系式：

在理想大规模阵列的前提下，天线数趋向于无穷大，显然，总存在某个q(称其为)满足此时对应于某信源的入射角θ_k，取所有元素的幅值，其将仅在处具有非零值，即反映在DFT功率谱上，DFT功率谱将在点形成明显的尖峰。在理想情形下，呈现高度稀疏特性，所有的功率都集中在第点，不存在功率泄漏现象。如图3所示，我们通过设置一特殊情形展现了这一高度系数特性：天线数M＝128，空间仅存在一个信源θ＝30°。图5绘制了该情形下阵列流型的DFT和离散事件傅里叶变换(Discrete Time Fourier Transform,DTFT)谱，DTFT谱的尖峰值即是入射信源AOA对应的理论位置，其仅在32点处具有非零幅值。

在实际情形中，我们需对接收信号进行DFT变换，定义其为x＝Fy，其第q个元素为

结合大规模下天线数极大，下列关系式仍然成立

因此，DFT功率谱包含AOA空域信息的结论仍然成立。通过寻找DFT功率谱中K个较大的峰值的位置我们可获得AOA的初始估计

然而，在实际情形下，尽管我们考虑的是大规模阵列情形，但阵元数仍只能是有限值，而不是无限大。因此，等式很可能因不是整数而不成立。如图4所示，此时|x|仍然是高度稀疏，大部分功率集中在但功率会由点泄露到其余点，其中表示取最接近的整数。易得功率泄露的程度与阵元数M成反比而与的值成正比。

102、将所述阵列接收信号进行相位旋转，对所述旋转后的信号进行离散傅里叶变换，并获取相位最佳旋转因子。

为了消除功率泄露对AOA估计的影响，本发明实施例提出了相位旋转(phase rotation)技术，其可以加强信号的系数程度，即使更多的功率集中到更少的DFT点上。相位旋转技术的详细步骤为：

定义以长度为(C+1)×1的矢量Λ，其第c个元素为：矢量Λ中的元素称为“旋转因子”，旋转因子在内均匀分布。

对某一旋转因子[Λ]_c，构建对角矩阵称为相位旋转矩阵。对接收信号右乘相位旋转矩阵即得到相位旋转后的信号，对旋转后的信号进行DFT变换可得表示对应于第c个旋转因子旋转后的矢量。

根据前文分析，当C趋向于无穷大，总有某相位旋转因子可使得成立，此时取的幅值仅在某一点处具有非零值由此，即使天线数有限，对任意AOA，也可通过相位旋转的技术找到对应相位旋转因子来实现增强稀疏性的目的。

图5示出了本发明一实施例提供的大规模天线阵列信号的波达方向单快拍快速估计方法的算法流程示意图。如图5所示，对阵列接收信号进行一系列相位旋转，对旋转后的信号分别进行DFT变换。针对信号功率谱内每一点，设置某一适合长度的滑动窗并计算该点的能量集中程度，选取能量最集中的情形为最佳情形，并由此确定旋转角度。

在实际中，考虑到C的取值不可能趋向无穷大以及存在噪声的影响，DFT功率谱仍然不会呈现理想稀疏特性。因此我们设置了一滑动窗来衡量考量DFT功率谱的稀疏程度，其具体过程为：

对旋转后的DFT幅度谱某一点设置一长度为2P+1的滑动窗W_c，q，该滑动窗内包含了以下2P+1个点：

，其中mod(a，b)表示变量a对b取模。P的取值在1到M之间，考虑到复杂度，其典型值在5到10之间。

定义q点处的稀疏指数表示第c次旋转得到的DFT幅度谱中第q点幅度与其对应滑动窗W_c，q内所有点幅值和的比值。显然Q_c，q越大，该点能量越集中，稀疏程度亦越高。

对接收信号进行C+1次旋转，对每一次旋转后的DFT谱，计算每一点的稀疏指数并记录，由此将得到一(C+1)×M矩阵Q。最佳旋转可由下式估计得到：

其中表示Q的第列。

为展示相位旋转技术的效果，本发明实施例进行了如下仿真。考虑考虑一128阵元的大规模ULA阵列，空间存在两个信号源，分别称为信源一和信源二，入射角分别为10.5°和70.5°，未加噪声。图6展示了未对接收信号相位旋转而直接进行DFT变换的DFT幅度谱，右侧的两小图分别是信源一和信源二的局部放大图。可以发现，在噪声和阵列孔径的限制下，功率泄露现象明显存在。图7和图8分别展示了相位旋转后，信源一和信源二分别获得最佳旋转时的DFT幅度谱。可见当信源一和信源二分别获得最佳旋转后，相应位置的DFT谱和DTFT谱重叠程度更高，功率泄露现象得到了有效抑制，稀疏性得到了有效增强。

103、根据所述波达方向的初始估计和所述相位最佳旋转因子，计算出所述波达方向的最终估计值。

具体的，所述波达方向的最终估计值通过如下公示计算得到：

可见，上述方法，充分考虑了大阵列孔径的特性，利用均匀大规模天线阵列下离散傅里叶变换DFT功率谱的精确性，通过DFT这一简单有效的方法，大大降低了复杂度；同时，通过对时域信道进行空间相位旋转，使得功率泄露现象得到了有效抑制，加强了空间域信道的稀疏性，从而实现仅需要单次快拍的数据就能实现对于波达方向的准确估计，提高了估计精度。

图9示出了本发明一实施例提供的大规模天线阵列信号的波达方向单快拍快速估计方法的装置结构示意图，如图9所示，本实施中大规模天线阵列信号的波达方向单快拍快速估计方法的装置包括：

初始估计模块91，用于对阵列接收信号进行离散傅里叶变换，搜索变换后的信号功率谱，由谱峰位置得到信号波达方向的初始估计；

相位获取模块92，用于将所述阵列接收信号进行相位旋转，对所述旋转后的信号进行离散傅里叶变换，并获取相位最佳旋转因子，；

最终估计模块93，用于根据所述波达方向的初始估计和所述相位最佳旋转因子，计算出所述波达方向的最终估计值。

本发明的另一个优选的实施例中，所述相位获取模块92通过如下方式进行所述相位旋转：

定义长度为(C+1)×1的矢量Λ，矢量Λ第c个元素为其中，矢量Λ中的元素为旋转因子，所述旋转因子在内均匀分布，M为全向天线阵元的个数；

对旋转因子[Λ]_c构建对角矩阵得到相位旋转矩阵，将所述阵列接收信号右乘所述相位旋转矩阵得到相位旋转后的信号。

本发明的另一个优选的实施例中，所述相位获取模块92通过如下方式获取所述相位最佳旋转因子：

对旋转后的信号离散傅里叶变换功率谱内每一点，设置一定长度的滑动窗，设置一定长度的滑动窗，并计算该点功率谱与所述滑动窗内所有点功率谱和的比值，将所述比值定义为该点的稀疏指数；

对所述阵列接收信号进行C+1次旋转，对每一次旋转后的离散傅里叶变换功率谱，计算每一点的稀疏指数并记录，得到一(C+1)×M稀疏指数矩阵Q；

根据所述稀疏指数矩阵Q中初始估计值谱峰位置附近的最大值，获取对应的相位最佳旋转因子。

可理解的是，上述装置与上述方法是一一对应的关系，上述方法中的实施细节也适用于上述装置，因此，本实施例不再对上述装置的具体实施细节进行详细描述。

综上，上述装置，充分考虑了大阵列孔径的特性，利用均匀大规模天线阵列下离散傅里叶变换DFT功率谱的精确性，通过DFT这一简单有效的方法，大大降低了复杂度；同时，通过对时域信道进行空间相位旋转，使得功率泄露现象得到了有效抑制，加强了空间域信道的稀疏性，从而实现仅需要单次快拍的数据就能实现对于波达方向的准确估计，提高了估计精度。

本领域技术人员可以理解，可以对实施例中的设备中的模块进行自适应性地改变并且把它们设置在于该实施例不同的一个或多个设备中。可以把实施例中的模块或单元或组件组合成一个模块或单元或组件，以及此外可以把它们分成多个子模块或子单元或子组件。除了这样的特征和/或过程或者单元中的至少一些是互相排斥之处，可以采用任何组合对本说明书(包括伴随的权利要求、摘要和附图)中公开的所有特征以及如此公开的任何方法或者设备的所有过程或单元进行组合。除非另外明确陈述，本说明书(包括伴随的权利要求、摘要和附图)中公开的每个特征可以由提供相同、等同或相似目的的替代特征来代替。

此外，本领域的技术人员能够理解，尽管在此所述的一些实施例包括其它实施例中所包括的某些特征而不是其它特征，但是不同实施例的特征的组合意味着处于本发明的范围之内并且形成不同的实施例。例如，在下面的权利要求书中，所要求保护的实施例的任意之一都可以以任意的组合方式来使用。

本发明的各个部件实施例可以以硬件实现，或者以在一个或者多个处理器上运行的软件模块实现，或者以它们的组合实现。应该注意的是上述实施例对本发明进行说明而不是对本发明进行限制，并且本领域技术人员在不脱离所附权利要求的范围的情况下可设计出替换实施例。在权利要求中，不应将位于括号之间的任何参考符号构造成对权利要求的限制。单词“包含”不排除存在未列在权利要求中的元件或步骤。位于元件之前的单词“一”或“一个”不排除存在多个这样的元件。本发明可以借助于包括有若干不同元件的硬件以及借助于适当编程的计算机来实现。在列举了若干装置的单元权利要求中，这些装置中的若干个可以是通过同一个硬件项来具体体现。单词第一、第二、以及第三等的使用不表示任何顺序。可将这些单词解释为名称。

最后应说明的是：本领域普通技术人员可以理解：以上各实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述各实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分或者全部技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明权利要求所限定的范围。