CN107966677A - 一种基于空间稀疏约束的圆阵模态域方位估计方法 - Google Patents

Abstract

本发明涉及水声工程技术领域，具体涉及一种基于空间稀疏约束的圆阵模态域方位估计方法，是为了解决现有的模态域方位估计方法需要大量多快拍采样数据构造协方差矩阵，且需要对协方差矩阵求逆或特征值分解，算法执行效率较低的缺点而提出的，包括：获得均匀圆阵声压信号；构造相位模态域变换矩阵；将均匀圆阵接收到的声压信号通过相位模态域变换矩阵变换成相位模态域信号；将360°全空间方位离散化，并构造空域稀疏变换基；对均匀圆阵接收信号进行时域采样，构造相位模态域多快拍采样数据，采用l₁范数对多快拍稀疏信号进行优化求解；计算不同方位角度上的信号能量，绘制空间谱图，由最大谱峰位置得到目标方位。本发明适用于水声工程中的信号估计。

Description

一种基于空间稀疏约束的圆阵模态域方位估计方法

技术领域

本发明涉及水声工程技术领域，具体涉及一种基于空间稀疏约束的圆阵模态域方位估计方法。

背景技术

在阵列信号处理领域，空间目标定位技术备受广大研究学者关注，广泛应用于雷达、声纳、医学成像等国民经济及国防建设领域。在过去的几十年里，各式各样不同功能的阵型结构层出不穷，其中最简单实用的当属均匀线列阵(Uniform Linear Array：ULA)。以均匀线列阵为基础，许多经典的方位估计算法(例如，最大熵谱法(Maximum Entropy：ME)、最大似然法(Maximum Likelihood：ML)以及特征子空间法等)也孕育而生。然而，均匀线列阵仅能覆盖扫描水平180°范围内的方位角，并且当扫描角度临近端射角度时，扫描波束明显变宽。

相比之下，均匀圆阵(Uniform Circular Array：UCA)具有明显的优势：(1)能够对空间目标在水平范围内实现360°的全方位覆盖；(2)波束指向可任意调节且不改变波束形状。这些优点激发了广大专家学者在圆阵方位估计技术领域的浓厚兴趣，促进了该领域技术的长足发展。现阶段，该领域技术主要分为两大类：(1)传统阵列信号处理技术，该类技术直接处理阵元接收信号；(2)模态域阵列信号处理技术，该类技术基于相位模态的概念，先将阵元域信号转换成模态域信号再进行相应处理。然而，由于圆阵阵列流型不具有Vandermonde结构，许多线列阵算法无法直接移植到圆阵，这就限制了传统阵列信号处理技术在圆阵上的直接应用。相比之下，模态域处理技术有效地解决了该问题，近年来得到迅猛发展。

Meyer和Elko在2002年国际会议IEEE International Conference onAcoustics,Speech,and Signal Processing上发表的“A highly scalable sphericalmicrophone array based on an orthonormal decomposition of the soundfield”将麦克风阵列安装在刚性球体表面，提出基于球面谐波声场分解的模态域方位估计方法。Teutsch和Kellermann在2006年第120卷第5期Journal of the Acoustical Society ofAmerica上发表的“Acoustic source detection and localization based on wavefielddecomposition using circular microphone arrays”将麦克风阵列安装在刚性圆柱体表面，利用声场柱面波分解理论，提出一种多目标宽带声源方位估计方法。杨德森等人在2014年第39卷第1期声学学报上发表的“声矢量圆阵相位模态域目标方位估计”基于模态域变换思想，利用矢量水听器中声压振速联合处理的抗噪能力，结合MUSIC算法实现水下声源目标的高分辨方位估计。但上述圆阵模态域方位估计方法均需要大量多快拍采样数据构造协方差矩阵，且需要对协方差矩阵求逆或特征值分解，算法执行效率较低。此外，上述方法仅利用了阵列接收信号的时-空统计信息，没有利用空间信源与生俱来的空间稀疏性信息。

发明内容

本发明的目的是为了解决现有的模态域方位估计方法需要大量多快拍采样数据构造协方差矩阵，且需要对协方差矩阵求逆或特征值分解，算法执行效率较低，进而使得预测时间过长，准确率低的缺点，而提出一种基于空间稀疏约束的圆阵模态域方位估计方法，包括：

步骤一、在空间待测区域布放均匀圆阵，获得均匀圆阵声压信号P(t)；

步骤二、构造相位模态域变换矩阵T；

步骤三、将均匀圆阵接收到的声压信号P(t)通过相位模态域变换矩阵T变换成相位模态域信号X(t)；

步骤四、将360°全空间方位离散化，并构造空域稀疏变换基

步骤五、对均匀圆阵接收信号进行时域采样，构造相位模态域多快拍采样数据X，根据相位模态域多快拍采样数据X和空域稀疏变换基采用l₁范数对多快拍稀疏信号进行优化求解；

步骤六、根据步骤五中得到的多快拍稀疏信号计算不同方位角度上的信号能量，以方位角度为横轴，信号能量为纵轴，绘制空间谱图，由最大谱峰位置得到目标方位。

本发明的有益效果为：1)利用空间信源与生俱来的空间稀疏性，将传统模态域方位估计方法与空间稀疏约束技术有机结合起来，明显地提高了算法的分辨能力；2)对相干源信号具有适用性，不需要解相干处理，可直接对相干源信号进行方位估计；3)不需要构造协方差矩阵或对协方差矩阵求逆，明显降低了算法的运算复杂度；4)可在小快拍数条件下对目标进行准确的方位估计；5)本发明的单信源方位估计结果相较于其它方法具有更尖锐的谱峰和最低的旁瓣级(-25dB)，即在单信源入射的条件下，本发明内容所提出的方法。具有更好的方位估计性能。6)在双相干信源方位估计分析中，本发明具有较尖锐的谱峰和较低的旁瓣级。即在双相干源入射的条件下，本发明内容所提出的方法EB-CS具有更好的方位估计性能。

附图说明

图1为本发明的基于空间稀疏约束的圆阵模态域方位估计方法的流程图；

图2为具体实施方式二中均匀圆阵接收信号模型的示意图；

图3为实施例一的单信源方位估计结果；其中图3(a)为直角坐标系空间谱图，图3(b)为极坐标系空间谱图；

图4为实施例二的双相干信源方位估计结果；其中图4(a)为直角坐标系空间谱图，图4(b)为极坐标系空间谱图。

具体实施方式

具体实施方式一：本实施方式的基于空间稀疏约束的圆阵模态域方位估计方法，如图1所示，包括：

步骤一、在空间待测区域布放均匀圆阵，获得均匀圆阵声压信号P(t)；

步骤二、构造相位模态域变换矩阵T；

步骤三、将均匀圆阵接收到的声压信号P(t)通过相位模态域变换矩阵T变换成相位模态域信号X(t)；

步骤四、将360°全空间方位离散化，并构造空域稀疏变换基

步骤五、对均匀圆阵接收信号进行时域采样，构造相位模态域多快拍采样数据X，根据相位模态域多快拍采样数据X和空域稀疏变换基采用l₁范数对多快拍稀疏信号进行优化求解；

步骤六、根据步骤五中得到的多快拍稀疏信号计算不同方位角度上的信号能量，以方位角度为横轴，信号能量为纵轴，绘制空间谱图，由最大谱峰位置得到目标方位。

具体实施方式二：本实施方式与具体实施方式一不同的是：步骤一具体为：

将半径为a的N元均匀圆阵置于水平面x-y平面内，设空间存在Q个窄带信号s_q(t)(q＝1,…,Q)且与均匀圆阵位于同一平面，入射角度分别为θ_q(q＝1,…,Q)，如图2所示。

则第n个阵元接收到的声压信号可以表示为：

其中，为波数，f为信号频率，c为声速；M为圆阵可激发的最大相位模态数，取值为比ka大的最小整数；为虚数单位；J_m(·)为m阶贝塞尔函数；e为数学常数。

进一步，均匀圆阵接收到的声压信号P(t)可以表示为：

P(t)＝FBA(θ)S(t)+n(t) (2)

式中：

P(t)＝[P₁(t),P₂(t),…,P_N(t)]^T为均匀圆阵接收到的声压信号，其中(·)^T表示转置运算；

F＝[F_-M,F_-M+1,…,F_M]^T为空间傅里叶变换矩阵，其中

B＝diag[b_-M,b_-M+1,…,b_M]为传输函数对角矩阵，其中b_m＝j^mJ_m(ka)(m＝-M,…,M)表示传输函数，diag[·]表示对角矩阵；

A(θ)＝[a(θ₁),a(θ₂),…,a(θ_Q)]为模态域阵列流型矩阵；

为第q个信号的模态域导向矢量；

S(t)＝[s₁(t),s₂(t),…,s_Q(t)]^T为源信号矩阵；

n(t)＝[n₁(t),n₂(t),…,n_N(t)]^T为噪声矩阵，其中n_n(t)(n＝1,…,N)表示第n个阵元通道噪声。

其它步骤及参数与具体实施方式一相同。

具体实施方式三：本实施方式与具体实施方式一或二不同的是：步骤二具体为：

将传输函数对角矩阵B求逆，然后乘以空间傅里叶变换矩阵F的复共轭转置，最后除以阵元个数N，则构造得相位模态域变换矩阵T，即

其中，(·)^-1表示矩阵求逆运算，(·)^H表示矩阵复共轭转置。

其它步骤及参数与具体实施方式一或二相同。

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是：步骤三具体为：

将相位模态域变换矩阵T乘以声压信号P(t)，得到相位模态域信号X(t)，即

X(t)＝TP(t) (4)

将式(2)代入式(4)，进一步得到：

其中，N(t)＝Tn(t)代表相位模态域噪声场。

其它步骤及参数与具体实施方式一至三之一相同。

具体实施方式五：本实施方式与具体实施方式一至四之一不同的是：步骤四具体为：

将360°全空间方位离散化为个方位，设这些方位依次编号为：则离散化的个方位可以表示为

在此基础之上，构造空域稀疏变换基

其中，为第个离散方位角所对应的模态域导向矢量。

另一方面，Q×1维源信号矩阵S(t)＝[s₁(t),s₂(t),…,s_Q(t)]^T可以推广成维稀疏表示形式其中，

此时，相位模态域信号X(t)可以通过信号的稀疏形式和空域稀疏变换基进行表示

其它步骤及参数与具体实施方式一至四之一相同。

具体实施方式六：本实施方式与具体实施方式一至五之一不同的是：步骤五具体为：

对均匀圆阵接收信号进行时域连续采样，设采样时刻(每一个采样时刻称为一个快拍)依次为：t＝t₁,…,t_T，则表达式(8)可拓展为多快拍表达形式

表达式(9)可进一步写成矩阵形式

其中，X＝[X(t₁),X(t₂),…,X(t_T)]表示相位模态域多快拍采样数据，

表示多快拍稀疏信号，N＝[N(t₁),N(t₂),…,N(t_T)]表示多快拍噪声。

此时，为维矩阵，其具体形式为记矩阵第行向量为并计算该行向量的l₂范数

将矩阵所有行向量的l₂范数(即)重新组成一个新的维行向量

在已获得相位模态域多快拍采样数据X和空域稀疏变换基的基础之上，根据表达式(10)，对向量采用l₁范数进行稀疏约束，并考虑噪声的影响，可得到求解的凸优化问题

其中，min表示求最小，||·||₁表示向量的l₁范数，s.t.表示约束条件，||·||₂表示向量的l₂范数，σ是对当前噪声功率的一个估值。

表达式(13)属于凸优化问题，可利用凸优化数学工具箱CVX进行有效求解。即，可利用凸优化数学工具箱CVX求得多快拍稀疏信号

其它步骤及参数与具体实施方式一至五之一相同。

具体实施方式七：本实施方式与具体实施方式一至六之一不同的是：步骤六具体为：

第(5)步骤中所求得的多快拍稀疏信号为

记矩阵第行向量为将向量中所有元素取平方运算并累加求和，即计算第个方位角度上的信号能量

最后，以方位角度为横轴，信号能量为纵轴，绘制空间谱图，由最大谱峰位置即可获得目标入射方位角度。

其它步骤及参数与具体实施方式一至六之一相同。

上面对发明内容各部分的具体实施方式进行了说明。基于空间稀疏约束的圆阵模态域方位估计方法，利用空间信源与生俱来的空间稀疏性，将传统模态域方位估计方法与空间稀疏约束技术有机结合起来，可获得准确的高分辨方位估计结果。下面通过仿真实例对本发明具体实施做进一步描述。

实施例一：单信源方位估计分析

实例参数设置如下：设空间存在一均匀圆阵，阵元个数为12且圆阵半径为a＝0.75m。空间信源频率为1.5kHz，入射方位角为125°。接收阵元采样快拍数为24，信噪比为5dB。

图3给出单信源条件下，模态域常规波束形成处理器(简记为EBF)、模态域最小方差无畸变处理器(简记为EB-MVDR)，模态域多重信号分类处理器(简记为EB-MUSIC)以及本发明基于空间稀疏约束的圆阵模态域方位估计方法(EB-CS)的方位估计结果图。

分析图3可知：由于阵列物理孔径瑞利限的限制，EBF空间谱图的主瓣最宽且旁瓣最高。EB-MVDR和EB-MUSIC虽然采用了高分辨处理算法，突破了瑞利限限制，但旁瓣级仍然较高(小于-20dB)。相比之下，由于利用了空间信源与生俱来的空间稀疏性，本发明EB-CS方法较其它方法具有最尖锐的谱峰和最低的旁瓣级(-25dB)。由此可见，在单信源入射的条件下，本发明内容所提出的方法EB-CS具有更好的方位估计性能。

实施例二：双相干信源方位估计分析

为了更好的表明本发明方法EB-CS的优良性能，接下来对空间双相干信源进行方位估计分析。

实例参数设置如下：设空间存在一均匀圆阵，阵元个数为12且圆阵半径为a＝0.75m。空间存在两相干信源，频率均为1.5kHz，入射方位角分别为90°和120°。接收阵元采样快拍数为24，信噪比为20dB。图4给出双相干信源条件下，EBF、EB-MVDR、EB-MUSIC以及本发明方法EB-CS的方位估计结果图。

分析图3可知：在此参数条件下，由于阵列物理孔径瑞利限的限制，EBF空间谱图仅显示出一个主峰，根本无法分辨出空间两个目标。与此同时，EB-MVDR和EB-MUSIC的空间谱图虽然能够辨认出两个谱峰，但两个谱峰中间的凹槽非常浅，表明这两种方法的估计性能明显下降。相比之下，本发明EB-CS方法能够直接处理相干源信号，且具有较尖锐的谱峰和较低的旁瓣级。由此可见，在双相干源入射的条件下，本发明内容所提出的方法EB-CS具有更好的方位估计性能。

本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，本领域技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。

Claims (7)

1.一种基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，包括：

步骤一、在空间待测区域布放均匀圆阵，获得均匀圆阵声压信号P(t)；

步骤二、构造相位模态域变换矩阵T；

步骤三、将均匀圆阵接收到的声压信号P(t)通过相位模态域变换矩阵T变换成相位模态域信号X(t)；

步骤四、将360°全空间方位离散化，并构造空域稀疏变换基

步骤五、对均匀圆阵接收信号进行时域采样，构造相位模态域多快拍采样数据X，根据相位模态域多快拍采样数据X和空域稀疏变换基采用l₁范数对多快拍稀疏信号进行优化求解；

步骤六、根据步骤五中得到的多快拍稀疏信号计算不同方位角度上的信号能量，以方位角度为横轴，信号能量为纵轴，绘制空间谱图，由最大谱峰位置得到目标方位。

2.根据权利要求1所述的基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，步骤一具体包括：

将半径为a的N元均匀圆阵置于水平面内，设空间存在Q个窄带信号s_q(t)，q＝1,…,Q，且与均匀圆阵位于同一平面，入射角度分别为θ_q；

则第n个阵元接收到的声压信号可以表示为：

其中，为波数，f为信号频率，c为声速；M为圆阵可激发的最大相位模态数，取值为比ka大的最小整数；为虚数单位；J_m(·)为m阶贝塞尔函数；e为数学常数；

均匀圆阵接收到的声压信号P(t)可以表示为：

P(t)＝FBA(θ)S(t)+n(t) (2)

式中：

P(t)＝[P₁(t),P₂(t),…,P_N(t)]^T为均匀圆阵接收到的声压信号；

F＝[F_-M,F_-M+1,…,F_M]^T为空间傅里叶变换矩阵，其中 m＝-M,…,M；

B＝diag[b_-M,b_-M+1,…,b_M]为传输函数对角矩阵，其中b_m＝j^mJ_m(ka)表示传输函数，diag[·]表示对角矩阵；

A(θ)＝[a(θ₁),a(θ₂),…,a(θ_Q)]为模态域阵列流型矩阵；

为第q个信号的模态域导向矢量；

S(t)＝[s₁(t),s₂(t),…,s_Q(t)]^T为源信号矩阵；

n(t)＝[n₁(t),n₂(t),…,n_N(t)]^T为噪声矩阵，其中n_n(t)表示第n个阵元通道噪声，n＝1,…,N。

3.根据权利要求1或2所述的基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，步骤二具体为：

将传输函数对角矩阵B求逆，再乘以空间傅里叶变换矩阵F的复共轭转置，最后除以阵元个数N，则构造得相位模态域变换矩阵T，即

其中，(·)^-1表示矩阵求逆运算，(·)^H表示矩阵复共轭转置。

4.根据权利要求3所述的基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，步骤三具体为：

将相位模态域变换矩阵T乘以声压信号P(t)，得到相位模态域信号X(t)，即

X(t)＝TP(t) (4)

将式(2)代入式(4)，进一步得到：

其中，N(t)＝Tn(t)代表相位模态域噪声场。

5.根据权利要求4所述的基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，步骤四具体为：

将360°全空间方位离散化为个方位，设这些方位依次编号为：则离散化的个方位可以表示为

在此基础之上，构造空域稀疏变换基

其中，为第个离散方位角所对应的模态域导向矢量；

将Q×1维源信号矩阵S(t)＝[s₁(t),s₂(t),…,s_Q(t)]^T变换为维稀疏表示形式

其中，

此时，相位模态域信号X(t)可以通过信号的稀疏形式和空域稀疏变换基进行表示

。

6.根据权利要求5所述的基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，步骤五具体为：

对均匀圆阵接收信号进行时域连续采样，设采样时刻依次为：t＝t₁,…,t_T，则表达式(8)可拓展为多快拍表达形式

将公式(9)表示成矩阵形式

其中，X＝[X(t₁),X(t₂),…,X(t_T)]表示相位模态域多快拍采样数据，

表示多快拍稀疏信号，N＝[N(t₁),N(t₂),…,N(t_T)]表示多快拍噪声；

此时，为维矩阵，其具体形式为记矩阵第行向量为并计算该行向量的l₂范数

将矩阵所有行向量的l₂范数重新组成一个新的维行向量

在已获得相位模态域多快拍采样数据X和空域稀疏变换基的基础之上，根据表达式(10)，对向量采用l₁范数进行稀疏约束，并考虑噪声的影响，可得到求解的凸优化问题

其中，min表示求最小，||·||₁表示向量的l₁范数，s.t.表示约束条件，||·||₂表示向量的l₂范数，σ是对当前噪声功率的一个估值；

利用凸优化数学工具箱CVX进行有效求解得到多快拍稀疏信号

7.根据权利要求6所述的基于空间稀疏约束的圆阵模态域方位估计方法，其特征在于，步骤六具体为：

步骤五中所求得的多快拍稀疏信号为

记矩阵第行向量为将向量中所有元素取平方运算并累加求和，即计算第个方位角度上的信号能量

最后，以方位角度为横轴，信号能量为纵轴，绘制空间谱图，最大谱峰位置即为目标入射方位角度。