导向矢量失配情况下基于稀疏表示的波达方向估计方法与流程

本发明属于阵列信号处理领域，具体涉及到提出了一种在非理想高斯噪声环境中，传感器导向矢量适配情况下，利用稀疏表示方法进行信号源波达方向估计的算法。

背景技术：

信号源波达方向(DOA)估计一直是目标定位和追踪领域的重要研究方向之一，在雷达、声纳、灾害学、通信等方面都有着极大的应用。近些年来随着人们在航天、海洋、灾害预测等领域的研究日益深入，目标DOA估计受到更加广泛的关注，并被实际应用在空间目标追踪、水声通信、管线保护和地震预测等领域。

对DOA估计的研究经历了很长的发展历程，产生了一系列经典的的估计算法，比如Capon，MUSIC和ESPRIT等，但这些算法在信噪比较低、快拍数较少的情况下信号的分辨率和测向精度下降严重。近些年来，基于信号在空间分布的稀疏性特点，研究人员将稀疏重构的思想应用在DOA估计方面，并提出了一些基于稀疏表示的DOA估计算法。早期稀疏表示方法是通过构建仅与信号源DOA相关的累积量域数据，利用加权l₁范数最小化方法进行DOA估计，并在此基础上对DOA估计值进行稀疏化用于估计距离参数，该方法能够有效地避免信号源物理距离相近对DOA估计产生的影响。最近研究人员提出利用聚焦方法思想提出了一种基于稀疏表示的宽带信号源DOA估计算法(FSP)，用单一频点的基矩阵代替频率和角度联合构建的基矩阵，减少了基矩阵的列数，解决了传统稀疏表示方法中由于基矩阵维数过大而使得存储量大以及计算复杂等难题。已有研究人员在实际环境中利用稀疏表示方法验证DOA估计的性能，根据嵌套阵和互素阵两种阵列的阵型结构推导出稀疏阵列信号处理模型，采用空域等角度稀疏和等正弦空间稀疏两种模型进行DOA估计，这种方法能够用较少的阵元达到较大的自由度，从而提高目标分辨的精度，但并未考虑实际环境中不稳定因素对接收信号产生的影响。

在一些实际环境中，由于噪声分量重尾现象比较严重，因此需要选择一个合适的噪声模型用于噪声信号处理，通常采用的是冲击噪声或者SαS噪声进行模拟，而本发明采用合成圆对称广义高斯分布来模拟噪声环境，它的衰减速率更加缓慢，能够更加真实地描述实际环境噪声分布。同时在实际应用中，由于受到自身以及外界环境因素的影响，传感器的导向矢量可能会产生波动，生成一个非理想增益，使之与真实值之间存在误差，这对DOA估计的精确度会产生很大的影响。

技术实现要素：

为克服现有技术的不足，本发明旨在实现在传感器导向矢量失配的情况下利用稀疏表示方法对信号源波达方向进行估计。本发明采用的技术方案是，导向矢量失配情况下基于稀疏表示的波达方向估计方法，首先构建传感器阵列的接收信号模型；然后采用合成圆对称广义高斯分布对重尾严重的实际环境噪声信号进行模拟，利用分数低阶矩方法对接收信号模型进行处理，并对处理后的信号模型中存在的未知增益值即导向矢量失配而生成的不确定参数增益值进行优化；最后，采用基于稀疏表示方法进行信号源波达方向DOA估计，得到最优的信号源波达方向DOA估计值。

首先构建传感器阵列的接收信号模型，假设处于远场环境下的Q个窄带随机信号源发射的信号波到达由P个传感器组成的均匀线性阵列ULA，其中Q个信号源之间相互独立，并且噪声与信号互不相关；第q个信号源的到达角记为θ_q，在t时刻，第p个传感器处的接收信号为式中，a(θ_q)为第q个信号源s_q(t)在其波达方向θ_q上第p个传感器处的导向矢量，表达式为d为传感器之间的距离，λ为信号波长，g_p为实际环境中第p个传感器处的导向矢量增益值，s_q(t)为空间中随机分布的窄带信号，n_p(t)为空间噪声信号；

根据每个传感器处的接收信号求得传感器阵列的接收信号模型：

x(t)＝GAs(t)+n(t)

其中

x(t)＝[x₁(t),x₂(t),...,x_P(t)]^T

G＝diag[g₁,g₂,...,g_P]

A＝[a(θ₁),a(θ₂),...,a(θ_Q)]

s(t)＝[s₁(t),s₂(t),...,s_q(t)]^T

n(t)＝[n₁(t),n₂(t,...,n_p(t))]^T

式中，x(t)为传感器阵列接收信号模型，[·]^T表示阵列转置，G是各个传感器处的增益值组成的对角阵，A为所有导向矢量组成的阵列流形，a(θ_q)为各个传感器在信号源q处的导向矢量，s(t)为信号源矩阵，s_q(t)为第q个信号源的发射信号，n(t)为所有传感器处的噪声信号，n_p(t)为第p个传感器处的噪声信号。

采用合成圆对称广义高斯分布来模拟实际噪声环境，噪声信号关于随机变量w的概率密度函数表示为：

式中，Γ(·)表示伽马函数，α为稳定分布系数，σ为函数方差，并且当α＝2时，噪声呈现高斯分布。

采用基于分数低阶矩矩阵方法FLOM对其进行处理，具体步骤如下：均匀线性传感器阵列的接收信号x(t)分数低阶矩矩阵是一个P×P的矩阵，它的第(i,k)个分量为

式中，E(·)表示期望值，x^*(t)表示其共轭函数，x(t)中各个信号之间互不相关，ε为小于α的指数参数；

将信号时域模型代入到传感器阵列接收信号x(t)的分数低阶矩矩阵表达式中

将R_ik分解为两部分C_ik和D_ik，其中

式中，为克罗内克符号，i，k为分数低阶矩矩阵分量的下标，θ_q,θ_r分别为第q个和第r个信号源的到达角，s_q(t),s_r(t)分别为第q和第r个信号源的发射信号，且在上式中

因此，传感器阵列的分数低阶矩矩阵表示为

R＝GAΓA^H+γI

式中，I为P×P阶的单位矩阵，A为阵列流形矩阵，H表示共轭转置。

由于上式中含有传感器自身运动以及外界环境干扰产生的增益值G，利用加权最小二乘方法对其进行最优化估计，求出其最优解：首先将得到的传感器阵列接收信号x(t)按照时间间隔T进行K次采样得到采样信号x(n)，则得到其协方差矩阵的表达式

根据获取到的采样数据，利用加权最小二乘法得到增益函数，再通过求取增益函数的最小值来获得无定向增益的最优解

式中，W是权重系数，通这里先利用W＝Ι对其进行初始化，之后用对每次的采样数据进行修改，则上式中的增益函数改写为

由于传感器阵列是由P个传感器组成的，则第p个传感器处的增益函数记为：

式中，w_n＝ζ(W_n)且z＝ζ_p(GAΓA^H)/g_p，其中操作符ζ_p(·)表示从一个矩阵中选择其第p列作为列向量，并去掉向量中第p个元素。

这样最优增益值求解问题就转化为求解P个线性最小二乘问题，然后对增益函数求最小值来求得最终结果，在这里引入公式z_w＝z⊙w_n，则求出的最优增益值为

最后，基于稀疏表示方法进行DOA估计，具体步骤是：首先对增益值优化后的信号模型进行矢量化处理；在得到矢量化的分数低阶矩矩阵之后，对DOA角度进行稀疏化处理，构造冗余的角度集合，并采用奇异值分解SVD的方法得到信号子空间，将角度求解转化为二阶锥规划问题。

对增益值优化后的分数低阶矩矩阵进行矢量化

式中，vec()表示矢量化运算，和⊙分别表示矩阵的Kronecker积和Khatri-Rao积。

则K个采样帧内所有的分数低阶矩矩阵矢量化后得

式中，y_k为每个采样帧内的矢量化分数低阶矩矩阵，定义由于响应矩阵矩阵Y的维数变为P²，可测得的信号源数目Q与阵元数P需满足Q≤2(P-1)。

在得到矢量化的分数低阶矩矩阵之后，采用稀疏表示方法对DOA角度进行冗余化处理，稀疏表示DOA估计方法的主要思想是用角度集合覆盖所有可能的目标信号入射方向，且满足N＞＞Q，N是划分的空间角度的个数，基于此可构成一个过完备的冗余字典它包含源的所有可能的位置信息，采取某些方法从这些位置信息中提取到所需的角度信息，采用奇异值分解SVD的方法对K帧内的矢量化分数低阶矩矩阵进行处理

式中，和分别为左奇异矩阵和右奇异矩阵，为对角阵；

得到信号子空间

式中，I为Q×Q阶的单位矩阵，

因此，将信号源DOA求解问题转化为

式中，||·||₁和||·||₂分别表示l₁范数和l₂范数，规则化参数β用来平衡信号稀疏性和噪声的影响；

根据阵列输出数据的有限快拍数只能得到的估计结果并且它们近似相等把上式的最优化问题转化为

式中，参数λ用来平衡l₁范数和l₂范数，目标函数的前一项反映失配程度，后一项反映稀疏性要求，上式是一个二阶锥规划问题，利用内点法进行求解，进而根据非零元素的位置求得信号源的DOA估计值。

本发明的特点及有益效果是：

1.本发明考虑实际环境中噪声重尾严重时利用合成圆对称广义高斯分布进行模拟，并采用分数低阶矩矩阵对信号模型进行处理，克服了此噪声环境下，噪声二阶以上矩不存在时传统算法不能进行准确估计的缺点；

2.本发明对实际环境中传感器导向矢量非理想的情况进行处理，利用加权最小二乘方法对导向矢量增益值进行优化估计；

3.本发明克服了传统DOA估计算法中当传感器数目少于源数目时不能进行精确DOA估计的问题；

4.本发明利用改进的稀疏表示算法进行信号源DOA估计，提高了估计精度，并采用奇异值分解降低了运算量。

附图说明：

附图1是本发明算法的方向流程图；

附图2是本发明算法及其它经典算法的DOA估计结果；

附图3是信号DOA估计的均方误差与信噪比的关系曲线；

附图4是信号DOA估计的均方误差与快拍数的关系曲线；

附图5是信号DOA估计的检测概率与信噪比的关系曲线。

具体实施方式：

针对现有技术的问题，本发明提出了一种新的基于稀疏表示的DOA估计算法。该算法选用合成圆对称广义高斯分布模拟实际环境噪声，但由于其二阶以上矩都不存在，传统的基于协方差矩阵分解的算法不再适用，因此本发明算法采用基于分数低阶矩阵的方法对噪声信号进行处理。然后对导向矢量失配而生成的增益值，利用加权最小二乘法进行优化估计。最后采用稀疏表示方法构建过完备的冗余字典，并利用奇异值分解降低维数，将其转化为二阶锥规划问题进行求解，利用内点法得到最优的DOA估计值。

本发明主要解决的问题是在传感器导向矢量失配的情况下利用稀疏表示方法对信号源波达方向进行估计。本发明的实现过程如下：

步骤一：构建信号模型。

假设处于远场环境下的Q个窄带随机信号源发射的信号波到达由P个传感器组成的均匀线性阵列(ULA)，其中这Q个信号源之间相互独立，并且噪声与信号互不相关。第q个信号源的到达角记为θ_q，在t时刻，第p个传感器处的接收信号为式中，a(θ_q)为信号s_q(t)在其波达方向θ_q上第p个传感器处的导向矢量，表达式为d为传感器之间的距离，λ为信号波长。g_p为实际环境中第p个传感器处的导向矢量增益值，s_q(t)为空间中随机分布的窄带信号，n_p(t)为空间噪声信号。

根据每个传感器处的接收信号可以求得传感器阵列的接收信号模型

x(t)＝GAs(t)+n(t)

其中

x(t)＝[x₁(t),x₂(t),...,x_P(t)]^T

G＝diag[g₁,g₂,...,g_P]

A＝[a(θ₁),a(θ₂),...,a(θ_Q)]

s(t)＝[s₁(t),s₂(t),...,s_Q(t)]^T

n(t)＝[n₁(t),n₂(t,...,n_P(t))]^T

式中，x(t)为传感器阵列接收信号模型，G是各个传感器处的增益值组成的对角阵，A为所有导向矢量组成的阵列流形，a(θ_q)为各个传感器在信号源q处的导向矢量，s(t)为信号源矩阵，s_q(t)为第q个信号源的发射信号，n(t)为所有传感器处的噪声信号，n_p(t)为第p个传感器处的噪声信号。

步骤二：基于分数低阶矩的信号模型处理及增益值最优化估计。

实际环境中的噪声信号不总是理想的高斯白噪声，这会对目标角度估计产生一定的影响，尤其是在噪声重尾现象比较严重时这种影响尤为突出。针对此情况我们采取以下步骤进行解决：首先采用合成圆对称广义高斯分布对重尾严重的噪声信号进行模拟，并利用分数低阶矩方法处理接收信号模型，接下来对信号模型中存在的未知增益值进行优化。

一些实际环境中的噪声分量重尾现象非常地严重，在进行噪声处理时需要选择合适的分布类型对其进行模拟，通常采用的是冲击噪声或者SαS噪声分布类型，而本发明采用合成圆对称广义高斯分布来模拟噪声环境，它的衰减速率更加缓慢，能够更加真实地描述实际环境噪声分布。噪声信号关于随机变量w的概率密度函数表示为

式中，Γ(·)表示伽马函数，α为稳定分布系数，σ为函数方差，并且当α＝2时，噪声呈现高斯分布。

由于噪声信号为合成圆对称广义高斯噪声时，接收信号协方差矩阵的二阶以上矩不存在，因此传统的MUSIC、ESPRIT等基于协方差矩阵特征分解的DOA估计算法都不再适用。这里我们采用分数低阶矩矩阵(FLOM)方法对噪声分量进行处理。

均匀线性传感器阵列的接收信号x(t)分数低阶矩矩阵是一个P×P的矩阵，它的第(i,k)个分量为式中，E(·)表示期望值，x^*(t)表示其共轭函数，x(t)中各个信号之间互不相关。

由于FLOM矩阵中存在因导向矢量失配而生成的不确定参数增益值g，这对后续数据的处理造成极大的困扰，因此我们需要采取合理的方法对其包含的g值进行估计，求出其最优解。在本发明中我们利用加权最小二乘方法迭代求出最优增益值。通过求取增益函数的最小值来获得无定向增益的最优解式中，W是权重系数，通常情况下是未知的，但由于在实际环境中，环境噪声对信号处理的影响会比较大，因此这里先利用W＝Ι对其进行初始化，之后用对每次的采样数据进行修改。

步骤三：基于稀疏表示方法进行DOA估计。

在实际应用中我们经常会遇到传感器数目少于信号源数目的情况，这时利用传统的DOA估计方法不再能够进行精确测量，针对这种情况本文采取矢量化方法进行处理，增加FLOM矩阵的维数，进而增加可测得的源数目。

在得到矢量化的分数低阶矩矩阵之后，采用稀疏表示方法对DOA角度进行冗余化处理。稀疏表示DOA估计方法的主要思想是将包含导向矢量信息的流形矩阵扩展成一个过完备的冗余字典，它包含源的所有可能的位置信息，进而从冗余的位置信息中提取所需的角度信息。随着字典冗余度的增大，计算量也急剧增加，为了减少计算量，本文采用奇异值分解(SVD)的方法进行处理。

下面结合附图和具体实施方式进一步详细说明本发明。

本发明提出了一种导向矢量失配情况下的基于稀疏表示的信号源DOA估计算法，具体的实施过程可以分为三个部分。

第一部分为信号模型的构建。

假设处于远场环境下的Q个窄带随机信号源发射的信号波到达由P个传感器组成的均匀线性阵列(ULA)，其中这Q个信号源之间相互独立，并且噪声与信号互不相关。第q个信号源的到达角记为θ_q，在t时刻，第p个传感器处的接收信号为式中，a(θ_q)为第q个信号源s_q(t)在其波达方向θ_q上第p个传感器处的导向矢量，表达式为d为传感器之间的距离，λ为信号波长。g_p为实际环境中第p个传感器处的导向矢量增益值，s_q(t)为空间中随机分布的窄带信号，n_p(t)为空间噪声信号。

第二部分为基于分数低阶矩的信号模型处理及增益值最优化估计。

本发明利用合成圆对称广义高斯噪声对实际噪声环境进行模拟，相对普通的噪声分布其衰减速率更加缓慢，能够更加真实地模拟实际环境噪声分布，其概率密度函数(probability density function，PDF)表示为

式中，Γ(·)表示伽马函数，α为稳定分布系数，σ为函数方差，并且当α＝2时，噪声呈现高斯分布。

在本发明中我们采用基于分数低阶矩(FLOM)的方法对这种噪声进行处理。均匀线性传感器阵列的接收信号x(t)分数低阶矩矩阵是一个P×P的矩阵，它的第(i,k)个分量为

式中，E(·)表示期望值，x^*(t)表示其共轭函数，x(t)中各个信号之间互不相关，ε为小于α的指数参数。

将信号时域模型代入到传感器阵列接收信号x(t)的分数低阶矩矩阵表达式中

将R_ik分解为两部分C_ik和D_ik，其中

式中，为克罗内克符号，i，k为分数低阶矩矩阵分量的下标，θ_q,θ_r分别为第q个和第r个信号源的到达角，s_q(t),s_r(t)分别为第q和第r个信号源的发射信号，且在上式中

因此，传感器阵列的分数低阶矩矩阵可表示为

R＝GAΓA^H+γI

式中，I为P×P阶的单位矩阵，A为阵列流形矩阵，H表示共轭转置。

将得到的传感器阵列接收信号x(t)按照时间间隔T进行K次采样得到采样信号x(n)，则可以得到其协方差矩阵的表达式

根据获取到的采样数据，我们可以利用加权最小二乘法得到增益函数，再通过求取增益函数的最小值来获得无定向增益的最优解

式中，W是权重系数，通常情况下是未知的，但由于在浅海环境中，环境噪声对信号处理的影响会比较大，因此这里先利用W＝Ι对其进行初始化，之后用对每次的采样数据进行修改，则上式中的增益函数可改写为

由于传感器阵列是由P个传感器组成的，则每个传感器处的增益函数可记为

式中，w_n＝ζ(W_n)且z＝ζ_p(GAΓA^H)/g_p，其中操作符ζ_p(·)表示从一个矩阵中选择其第p列作为列向量，并去掉向量中第p个元素。

这样最优增益值求解问题就转化为求解P个线性最小二乘问题，然后对增益函数求最小值来求得最终结果，在这里我们引入公式z_w＝z⊙w_n，则求出的最优增益值为

第三部分为基于稀疏表示的DOA估计。

为了增加传感器阵列可测得的信号源数目，对增益值优化后的分数低阶矩矩阵进行矢量化

式中，vec()表示矢量化运算，和⊙分别表示矩阵的Kronecker积和Khatri-Rao积。

则K个采样帧内所有的分数低阶矩矩阵矢量化后可得

式中，y_k为每个采样帧内的矢量化分数低阶矩矩阵，定义我们可以看到由于响应矩阵矩阵Y的维数变为P²，大于通常所采用的物理维数P(P>1时)，因此能够使用少于目标源数目的传感器对其进行测量，但可测得的信号源数目Q与阵元数P需满足Q≤2(P-1)。

下面利用得到的矢量化FLOM矩阵采用基于稀疏表示的方法进行DOA估计。在得到矢量化的分数低阶矩矩阵之后，采用稀疏表示方法对DOA角度进行冗余化处理。稀疏表示DOA估计方法的主要思想是用角度集合覆盖所有可能的目标信号入射方向，且满足N＞＞Q，N是划分的空间角度的个数，基于此可构成一个过完备的冗余字典它包含源的所有可能的位置信息，进而采取某些方法从冗余的位置信息中提取出所需的角度信息。随着字典冗余度的增大，计算量也急剧增加，为了减少计算量，本文采用奇异值分解(SVD)的方法对K帧内的矢量化分数低阶矩矩阵进行处理。

式中，和分别为左奇异矩阵和右奇异矩阵，为对角阵。

可以得到信号子空间

式中，I为Q×Q阶的单位矩阵，

因此，可以将信号源DOA求解问题转化为

式中，||·||₁和||·||₂分别表示l₁范数和l₂范数，规则化参数β用来平衡信号稀疏性和噪声的影响。

根据阵列输出数据的有限快拍数我们只能得到的估计结果并且它们近似相等因此我们把上式的最优化问题转化为

式中，参数λ用来平衡l₁范数和l₂范数，目标函数的前一项反映失配程度，后一项反映稀疏性要求。

上式是一个二阶锥规划问题，我们可以利用内点法进行求解，进而根据非零元素的位置求得信号源的DOA估计值，结果如图2所示。