基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法与流程

本发明属于信号处理技术领域，尤其涉及对雷达信号、声学信号及电磁信号的波达方向估计，具体是一种基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法，可用于无源定位和目标探测。

背景技术：

波达方向(direction-of-arrival,doa)估计是阵列信号处理领域的一个重要分支，它是指利用阵列天线接收空域信号，并通过现代信号处理技术和各类优化方法实现对接收信号统计量的有效处理，从而实现信号的doa估计，在雷达、声呐、语音、无线通信等领域有着重要的应用价值。

doa估计方法的自由度是指其能够估计的入射信号源的个数。现有的doa估计方法通常采用均匀线性阵列进行信号的接收与建模，但是基于均匀线性阵列方法的自由度受限于实际天线阵元个数。具体而言，对于一个包含l个天线阵元的均匀线性阵列，其自由度为l-1。因此，当某个空域范围内入射信号源的个数大于或等于阵列中天线阵元的个数时，现有采用均匀线性阵列的方法将无法进行有效的doa估计。

互质阵列能够在天线阵元个数一定的前提下增加doa估计的自由度，因而受到了学术界的广泛关注。作为互质采样技术在空间域上的一个典型表现形式，互质阵列提供了一个系统化的稀疏阵列架构方案，并能够突破传统均匀线性阵列自由度受限的瓶颈，实现doa估计方法自由度性能的提升。现有的基于互质阵列的doa估计方法主要通过利用质数的性质将互质阵列推导到虚拟域，并形成等价虚拟均匀线性阵列接收信号以实现doa估计。由于虚拟阵列中包含的虚拟阵元数大于实际的天线阵元数，自由度因此得到了有效的提升。但是由于从互质阵列推导而来的虚拟阵列属于非均匀阵列，因此很多现有基于均匀线性阵列的信号处理方法无法直接应用于虚拟阵列等价接收信号以实现doa估计。当前采用互质阵列的doa估计方法常用的一个解决方案是，仅利用虚拟阵列中连续的阵元部分形成一个虚拟均匀线阵以进行doa估计，但是这造成了部分原始信息的丢失和相关估计性能的降低。

同时，目前众多doa估计方法在优化问题的设计过程中，需要预先设置信号假定波达方向的空间网格点。随着对波达方向估计结果精度要求的提高，这些doa估计方法需要预先设置的空间网格点将变得越来越密集，这导致了计算复杂度的急剧增加。不仅如此，在实际情况中，难免会有一些信号的波达方向无法完全落在预先设置的网格点上，从而造成了固有的模型失配误差。

技术实现要素：

本发明的目的在于针对上述现有技术存在的不足，提出一种基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法，充分利用了非均匀虚拟阵列所提供的全部信息，并保证了无网格化的波达方向估计，从而提高了doa估计的自由度与分辨率，并在一定程度上降低了doa估计的计算复杂度。

本发明的目的是通过以下技术方案来实现的：一种基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法，包含以下步骤：

(1)接收端使用m+n-1个天线，并按照互质阵列结构进行架构；其中m与n为互质整数；

(2)假设有k个来自θ1,θ2,…,θk方向的远场窄带非相干信号源，则(m+n-1)×1维互质阵列接收信号x(t)可建模为：

其中，sk(t)为信号波形，n(t)为与各信号源相互独立的噪声分量，a(θk)为θk方向的导引矢量，表示为：

其中，pid,i＝1,2,…,m+n-1表示互质阵列中第i个物理天线阵元的实际位置，且p1＝0；d为入射窄带信号波长λ的一半，即d＝λ/2，[·]^t表示转置操作。共采集t个采样快拍，得到互质阵列接收信号的采样协方差矩阵

这里，(·)^h表示共轭转置操作；

(3)计算互质阵列接收信号所对应的等价虚拟信号：矢量化互质阵列接收信号的采样协方差矩阵获得虚拟阵列等价接收信号v：

其中，为(m+n-1)²×k维虚拟阵列导引矩阵，包含k个入射信号源的功率，为噪声功率，iv＝vec(im+n-1)。这里，vec(·)表示矢量化操作，即把矩阵中的各列依次堆叠以形成一个新的矢量，(·)^*表示共轭操作，表示克罗内克积，im+n-1表示(m+n-1)×(m+n-1)维单位矩阵。矢量v对应的虚拟阵列中各虚拟阵元的位置为

去除集合中重复的元素所对应位置上的重复虚拟阵元，得到一个非均匀的虚拟阵列其对应的等价虚拟信号vc可通过选取矢量v中相对应位置上的元素获得；

(4)构造内插虚拟阵列及其接收信号并建模：首先对于非均匀的虚拟阵列在保留其原有虚拟阵元位置不变的前提下，向其中非连续的位置插入若干虚拟阵元，从而将非均匀虚拟阵列转化为间距为d、阵列孔径与互质阵列相同、且虚拟阵元数目增加的均匀虚拟阵列该内插均匀虚拟阵列共包含个虚拟阵元，其中|·|表示集合的势，其对应的等价虚拟信号vi可通过往矢量vc中插入0获得，插入0的位置与中插入的虚拟阵元的位置相对应；

(5)构造内插虚拟阵列多采样快拍信号及其采样协方差矩阵：将切割为li个长度为li的连续子阵列，其中

相应地，内插虚拟阵列的多采样快拍信号可通过截取矢量vi中对应的元素获得，即：由vi中第li+1-l到第2li-l个元素组成。接着，vi的采样协方差矩阵rv可以由如下方式得到：

其中，<vi>i表示位置为id的虚拟阵元所对应的等价接收信号；

(6)构造投影矩阵并定义投影运算：投影矩阵p的维度与rv相同，如果矩阵rv中某个元素为0，则投影矩阵p中相同位置的元素值也为0；反之投影矩阵p中相应位置的元素值为1。定义为投影运算，其中括号内变量为与p维度相同的矩阵，投影运算通过变量矩阵中每一个元素与投影矩阵p中相应位置上的元素一一相乘实现，得到一个与矩阵p维度相同的矩阵；

(7)设计基于内插虚拟阵列信号原子范数最小化的优化问题并求解：假设内插虚拟阵列信号的理想协方差矩阵为表示以矢量z为第一列的厄米特对称toeplitz矩阵。定义矢量z的原子范数：

其中α为非零实数，inf表示下确界。对应于每一个不同的θk，称r(θk)为一个用于描述矢量z的原子。利用(5)得到的协方差矩阵rv作为参考值，寻找一个由最少的原子表示，且由其构建的toeplitz矩阵与rv差异小于某一阈值的矢量z。可构建如下以矢量z为变量的优化问题：

其中，∈为阈值常数，用于约束协方差矩阵的重建误差；保证了重建的协方差矩阵满足半正定的条件；‖·‖f表示frobenius范数。求解上述凸优化问题可得到最优化值进而可以得到内插虚拟阵列协方差矩阵

(8)根据重建的内插虚拟阵列协方差矩阵进行波达方向估计。

进一步地，步骤(1)所述的互质阵列结构可具体描述为：首先选取一对互质整数m、n；然后，构造一对稀疏均匀线性子阵列，其中第一个子阵列包含m个间距为nd的天线阵元，其位置为0,nd,…,(m-1)nd，第二个子阵列包含n个间距为md的天线阵元，其位置为0,md,…,(n-1)md；接着，将两个子阵列按照首个阵元重叠的方式进行子阵列组合，获得实际包含m+n-1个天线阵元的非均匀互质阵列架构。

进一步地，步骤(5)所构建的vi的采样协方差矩阵rv也可以由下述方法等价得到：

进一步地，步骤(7)中的凸优化问题可转化为如下以矢量z为变量的优化问题：

其中μ为正则化参数，用于在最小化过程中权衡矩阵重建误差和z的原子范数。

进一步地，步骤(8)中的波达方向估计，可采用以下方法：多重信号分类方法、旋转不变子空间方法、求根多重信号分类方法、协方差矩阵稀疏重建方法等。

进一步地，步骤8中，通过多重信号分类方法进行波达方向估计，具体为：画出虚拟域空间谱pmusic(θ)：

其中d(θ)是li×1维内插虚拟阵列导引矢量，对应于位置为由0到(li-1)d的一段虚拟均匀阵列；en是li×(li-k)维矩阵，表示内插虚拟阵列协方差矩阵的噪声子空间；θ是假定的信号波达方向；通过谱峰搜索寻找空间谱pmusic(θ)上的峰值，并将这些峰值所对应的响应值从大到小排列，取前k个峰值所对应的角度方向，即为波达方向估计结果。

本发明与现有技术相比具有以下优点：

(1)本发明在互质阵列等价虚拟域上引入了阵列内插的思想，充分利用了虚拟阵列提供的全部信息。通过在非均匀虚拟阵列中内插虚拟阵元的方式构建出均匀线性虚拟阵列，在保留了由原始非均匀虚拟阵列接收到的全部信息的同时，使得构建的虚拟域信号模型满足奈奎斯特采样定律；

(2)本发明基于内插虚拟阵列信号原子范数最小化的思想设计优化问题，在优化问题设计的过程中无需预先定义空间网格点，保证了无网格化的波达方向估计。由于优化问题是从最小化虚拟域等价信号的原子范数角度进行协方差矩阵的重建，因此利用该优化问题重建内插虚拟阵列协方差矩阵的过程中不需要预先设置信号波达方向的空间网格点，即该波达方向估计方法是无网格化的，保证了波达方向估计的分辨率以及计算效率；

(3)本发明的优化问题在重建内插虚拟阵列协方差矩阵时保证了优化求解结果为厄米特对称的toeplitz矩阵，使得最优解与理论协方差矩阵之间的误差更小。由于均匀线性阵列的理论协方差矩阵满足toeplitz结构，因此利用它的toeplitz特性作为先验约束条件进行协方差矩阵的重建，可以使得重建结果与真实值差异更小，从而提高doa估计方法的性能。

附图说明

图1是本发明的方法总体流程框图。

图2是本发明中组成互质阵列的一对稀疏均匀子阵列结构示意图。

图3是本发明中互质阵列的结构示意图。

图4是本发明中内插虚拟阵列的结构示意图。

图5是本发明中内插虚拟阵列分割方法的示意图。

图6是用于体现本发明所提方法自由度性能的空间功率谱示意图。

图7是用于体现本发明所提方法分辨率性能的归一化空间谱示意图。

具体实施方式

以下参照附图，对本发明的技术方案和效果作进一步的详细说明。

对于doa估计方法在实际系统中的应用，互质阵列由于其可以通过等价虚拟阵列信号的计算和统计信号处理，突破物理阵元数量对自由度的限制而备受关注。但是受限于虚拟阵列的非均匀性，目前很多方法都会选择利用其中连续的虚拟阵元部分进行doa估计，从而造成了信息损失。同时，很多方法在进行doa估计之前会预先设置假定波达信号方向的空间网格点，这造成了固有失配误差以及计算复杂度与估计精度之间的矛盾。为了充分利用非均匀虚拟阵列中所包含的所有信息，并避免由于预定义空间网格点所造成的估计分辨率下降，本发明提供了一种基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法，参照图1，本发明的实现步骤如下：

步骤一：在接收端使用m+n-1个天线阵元架构互质阵列；首先，选取一组互质整数m、n；然后，参照图2，构造一对稀疏均匀线性子阵列，其中第一个子阵列包含m个间距为nd的天线阵元，其位置为0,nd,…,(m-1)nd；第二个子阵列包含n个间距为md的天线阵元，其位置为0,md,…,(n-1)md；单位间距d取为入射窄带信号波长λ的一半，即d＝λ/2；接着，将两个子阵列的首个天线阵元视为参考阵元，参照图3，将两个子阵的参考阵元重叠以实现子阵列组合，获得实际包含m+n-1个天线阵元的非均匀互质阵列架构。

步骤二：采用互质阵列接收信号并建模。假设有k个来自θ1,θ2,…,θk方向的远场窄带非相干信号源，采用步骤一架构的非均匀互质阵列接收入射信号，得到(m+n-1)×1维互质阵列接收信号x(t)，可建模为：

其中，sk(t)为信号波形，n(t)为与各信号源相互独立的噪声分量，a(θk)为θk方向的互质阵列导引矢量，表示为

其中，pid,i＝1,2,…,m+n-1表示互质阵列中第i个物理天线阵元的实际位置，且p1＝0，[·]^t表示转置操作。共采集t个采样快拍，得到互质阵列接收信号的采样协方差矩阵

其中，(·)^h表示共轭转置操作。

步骤三：计算互质阵列接收信号所对应的等价虚拟信号。矢量化互质阵列接收信号的采样协方差矩阵获得虚拟阵列等价接收信号v：

其中，为(m+n-1)²×k维虚拟阵列导引矩阵，包含k个入射信号源的功率，为噪声功率，iv＝vec(im+n-1)。这里，vec(·)表示矢量化操作，即把矩阵中的各列依次堆叠以形成一个新的矢量，(·)^*表示共轭操作，表示克罗内克积，im+n-1表示(m+n-1)×(m+n-1)维单位矩阵。矢量v对应的虚拟阵列中各虚拟阵元的位置为其中

去除集合中重复的元素所对应位置上的重复虚拟阵元，得到一个非均匀的虚拟阵列其对应的等价虚拟信号vc可通过选取矢量v中相对应位置上的元素获得。

步骤四：构造内插虚拟阵列及其接收信号建模。参照图4，对于非均匀的虚拟阵列在保留其原有虚拟阵元位置不变的前提下，向其中存在孔洞的位置插入若干虚拟阵元(如图4中的空心圆所示)，从而将非均匀虚拟阵列转化为间距为d、阵列孔径与互质阵列相同、且虚拟阵元数目增加的均匀虚拟阵列内插虚拟阵列共包含个虚拟阵元，其中|·|表示集合的势。内插虚拟阵列对应的等价虚拟信号vi可通过往矢量vc中孔洞的相应位置填入0获得。

步骤五：构造内插虚拟阵列多采样快拍信号及其采样协方差矩阵。参照图5，将内插虚拟阵列切割为li个长度为li的连续子阵列，其中

由于中的虚拟阵元以零位对称分布，始终为奇数，故li为整数。相应地，内插虚拟阵列的多采样快拍信号可通过截取矢量vi中对应的元素获得，即：其中vi,l,l＝1,2,…,li由vi中第li+1-l到第2li-l个元素组成。接着，vi的采样协方差矩阵rv可以由如下方式得到：

其中，〈vi〉i表示位置为id的虚拟阵元所对应的等价接收信号。由于内插虚拟阵列中虚拟阵元关于零位对称分布，因此其上的等价虚拟接收信号对应于零位呈共轭关系，所以上述采样协方差矩阵也可以通过如下方式等价得到：

步骤六：构造投影矩阵并定义投影运算。由于步骤五所得的协方差矩阵rv中包含有在步骤四中插入的0，因此其相应位置对角线上的元素全部为0。根据这样的结构定义一个与rv维度相同的投影矩阵p，如果rv中某一位置上的元素是0，则投影矩阵p中相同位置的元素值也为0；反之则投影矩阵p中相应位置的元素值为1。定义为投影运算，其中括号内变量为与p维度相同的矩阵，投影运算通过变量矩阵的每一个元素与投影矩阵p中相应位置上的元素一一相乘实现，得到一个与矩阵p维度相同的矩阵。

步骤七：设计基于内插虚拟阵列信号原子范数最小化的优化问题并求解。根据均匀阵列非相干接收信号理想协方差矩阵的toeplitz特性，假设内插虚拟阵列信号的理想协方差矩阵为表示以矢量z为第一列的厄米特对称toeplitz矩阵。定义矢量z的原子范数：

其中α是一个非零的实数，inf表示下确界。对应于每一个不同的θk，我们称r(θk)为一个用于描述矢量z的原子。我们取α＝1，利用步骤五得到的协方差矩阵rv作为参考值，寻找一个由最少原子表示，且由其构建的toeplitz矩阵与rv差异小于某一阈值的矢量z，可构建如下以矢量z为变量的优化问题：

其中，∈为阈值常数，用于约束协方差矩阵的重建误差；保证了重建的协方差矩阵满足半正定的条件；‖·‖f表示frobenius范数。

求解上述凸优化问题可得到最优化值上述凸优化问题转化为以下以矢量z为变量的优化问题：

其中μ为正则化参数，用于在最小化过程中权衡矩阵重建误差和z的原子范数。求解上述优化问题可得到最优化值相应地，重建的toeplitz矩阵为内插虚拟阵列协方差矩阵。

步骤八：根据重建的内插虚拟阵列协方差矩阵进行波达方向估计。通过引入经典的方法，如多重信号分类方法、旋转不变子空间方法、求根多重信号分类方法、协方差矩阵稀疏重建方法等对重建的内插虚拟阵列协方差矩阵进行操作，可以求得波达方向估计结果。以多重信号分类方法为例，画出虚拟域空间谱pmusic(θ)：

其中d(θ)是li×1维内插虚拟阵列导引矢量，对应于位置为由0到(li-1)d的一段虚拟均匀阵列；en是li×(li-k)维矩阵，表示内插虚拟阵列协方差矩阵的噪声子空间；θ是假定的信号波达方向；通过谱峰搜索寻找空间谱pmusic上的峰值，并将这些峰值所对应的响应值从大到小排列，取前k个峰值所对应的角度方向，即为波达方向估计结果。

本发明一方面引入虚拟阵列内插的思想，在推导的原始虚拟阵列基础上内插入虚拟阵元，从而将原始的非均匀虚拟阵列转化为虚拟均匀阵列，同时保留了原始非均匀虚拟阵列上的所有信息，避免了因原始虚拟阵列的非均匀性所导致的统计信号处理模型失配及传统方法截取虚拟均匀子阵所导致的信息损失问题；另一方面，引入了基于虚拟阵列信号原子范数最小化的思想来设计优化问题，以重建内插虚拟阵列的协方差矩阵，在虚拟域上进行信号处理以实现doa估计。

下面结合仿真实例对本发明的效果做进一步的描述。

仿真实例1：采用互质阵列接收入射信号，其参数选取为m＝3，n＝5，即架构的互质阵列共包含m+n-1＝7个物理阵元。假定入射窄带信号个数为9，且入射方向均匀分布于-50°至50°这一空间角度域范围内；信噪比设置为30db，采样快拍数t＝500；正则化参数μ设置为0.25。

本发明所提出的基于内插虚拟阵列信号原子范数最小化的互质阵列波达方向估计方法空间功率谱如图6所示，其中垂直虚线代表入射信号源的实际方向。可以看出，本发明所提方法能够有效分辨这9个入射信号源。而对于传统采用均匀线性阵列的方法，利用7个物理天线阵元最多只能分辨6个入射信号，以上结果体现了本发明所提方法实现了自由度的增加。

仿真实例2：采用互质阵列接收入射信号，其参数同样选取为m＝3，n＝5，即架构的互质阵列共包含m+n-1＝7个物理天线阵元；假定入射窄带信号个数为2，且入射方向为-0.5°至0.5°，其余参数设置与仿真实例1保持一致。由图7所示的归一化空间谱可以看出，本发明所提方法可以有效地分辨出这两个近距离信号源的波达方向，说明了本方法良好的分辨率性能。

综上所述，本发明所提方法充分利用了非均匀虚拟阵列上的全部信息，能够在信号源个数大于等于物理天线个数的情况下实现入射信号的有效估计，增加了doa估计的自由度和分辨率。此外，与传统采用均匀线性阵列的方法相比，本发明所提方法在实际应用中所需的物理天线阵元及射频模块也能够相应减少，体现了经济性和高效性。