最小方差无失真响应波束形成解卷积后处理算法

宋其岩 马晓川 李 璇 詹 飞

（1.中国科学院声学研究所中科院水下航行器信息技术重点实验室，北京 100190；2.中国科学院大学，北京 100049）

1 引言

来波方向（Direction of Arrival，DOA）估计一直是阵列信号处理的重要内容。声波频率越高，介质吸收损失越大，为了使声纳的作用距离增大，需要尽可能降低探测信号的频率，根据半波长布阵原则，频率越低阵列的阵元间距往往越大，阵元数量越少，波束主瓣越大，旁瓣造成的能量泄露越严重，阵列的分辨能力越差［1］。因此，如何在尺寸固定的声纳平台，尽可能提高方位估计的分辨力和估计精度，一直是阵列信号处理的难点问题。另外，随着主动传感器阵列作用范围的增大，搜索范围内很可能会出现多个目标，例如海上航母护航编队、拖曳诱饵、水下航行器集群，如何对多目标进行准确、快速的方位估计亟待解决。在空气声领域，基于麦克风阵列的室内说话人声源定位、语音增强同样需要应用到多目标方位估计技术［2-6］。

经过几十年的发展，众多学者提出了众多方位估计算法［7-10］。1948 年，Bartlett 提出传统波束形成技术（Conventional Beamforming，CBF）［11］，CBF算法将各个阵元的接收数据进行延时求和（Delay and Sum，DAS），对于不同的波达方向，阵元的加权值不同，然后求波束输出功率，功率谱的峰值对应的角度就是真实的波达方向。CBF 算法计算速度快，在工程中得到了普遍的应用，但是由于CBF算法的波束主瓣较宽，分辨能力较差，旁瓣较高引起能量泄露，功率强的信号容易掩蔽功率低的信号，最终导致模糊效应［12-13］。1969 年Capon 提出了高分辨的最小方差无失真响应（Minimum Variance Distortionless Response，MVDR）波束形成器［14］。1986 年Schmidt 提出了高分辨多重信号分类算法（Multiple Signal Classification，MUSIC）［15］，MUSIC算法需要对数据采样协方差矩阵执行特征分解，大特征值对应的特征向量是信号主特征向量，通过小特征值对应的特征向量构造噪声子空间，利用信号子空间与噪声子空间的正交性实现波达方向高分辨，但是MUSIC 估计的空间功率谱属于伪谱，无法准确估计出声源的功率。1986 年Roy 等人利用规则阵列的旋转不变性提出了旋转不变参数估计技术（Estimating Signal Parameter via Rotational Invariance Techniques，ESPRIT）［16］，计算复杂度比MUSIC 算法低，但是分辨能力低于MUSIC。继解卷积后处理技术在图像恢复领域取得了成功的应用后，有学者将解卷积技术扩展应用到了参数估计领域，2018 年，T.C.Yang 将常规波束形成表示为卷积形式，并将观察方向为0°的波束图作为点扩散函数（Point Spreading Function，PSF），将Richardson-Lucy（RL）［17］算法应用于CBF 算法进行解卷积后处理，提高了分辨能力，2019 年谢磊等人对MUSIC 的方位谱进行后处理，获得了背景级更低的方位谱［18］。Ma 等人提出了声源图像的快速解卷积方法（Deconvolution Approach for the Mapping of Acoustic Source，DAMAS）［19］，2017 年，Lylloff 将快速迭代收缩阈值算法（Fast Iterative Shrinkage-Thresholding Algorithm，FISTA）解卷积技术应用于平面麦克风阵列声源定位中［20］，并指出如果阵列具有平移不变性则可以采用快速傅里叶变换（Fast Fourier Transform，FFT）执行解卷积计算，提高了解卷过程中的计算效率，这类后处理算法的显著优势是能够继承被解卷积算法的优势，在付出少量的计算代价下，能够提高原始算法的分辨性能。

MVDR 算法分辨能力高，估计出的方位谱为功率真实谱，然而低信噪比（Signal-to-Noise Ratio，SNR）下，方位谱背景级较高，估计精度较差。本文对MVDR 算法的方位谱重新表述为卷积模型，并运用两种解卷积算法对MVDR 方位谱进行解卷积，算法在付出少量的计算代价下，能够提高MVDR 算法的方位分辨能力。仿真实验显示提出的算法具有高分辨能力。

2 信号模型

假设有K个远场窄带信号以角度θk，k=1，2，...，K入射到均匀直线阵列，阵元间距d=λ2，阵元m的输出可表示为

其中，阵元数为M，阵元号m=0，1，...，M-1，快拍数为N，n=1，2，...，N。

将式（1）写成矢量形式

3 方位估计算法

3.1 常规波束形成

经典波束形成法（Conventional Beamforming，CBF）又称延迟相加法（Delay and Sum，DAS），是DOA 领域最基础的方法，通过权系数对各个传感器记录的数据执行加权求和，从而达到调整阵列的波束指向的目的，令波束主瓣指向所期望的角度方向。

各阵元的权矢量，

阵列输出为，

则整个阵列输出的平均功率为

对波达方向在[-90°，90°]范围内进行角度扫描，计算阵列输出功率，最大功率点所对应的角度θs即为声源的DOA。

3.2 Richardson-Lucy解卷积常规波束形成

由于传统方位估计CBF 的分辨率受到阵列阵元数量的限制，即存在瑞利限。因此对CBF 的波束输出进行后处理是提高CBF 分辨能力的有效技术途径，文献［17］采用Richardson-Lucy 算法对CBF 方位谱进行解卷积后处理（CBF Richardson-Lucy，CBF-RL），提高了CBF算法的性能。

重新考虑常规波束形成方位谱公式（6），

将式（7）重新表示为卷积的形式，

PSFCBF(sinθ)是CBF 算法的PSF，PSF 的物理意义是空间中的点光源经过成像系统冲激响应函数卷积后的扩散后的图像，如图1 所示，点光源经过采样函数卷积后形成涟漪状的扩散圆环图像。在图像恢复领域，Richardson-Lucy（RL）是一种最流行的解卷积技术，RL 迭代解卷积过程如式（10）所示［21］，

其中，S=S(sinθ)，sinθ∈[-1，1]，Si代表第i次估计 结果，PCBF=P(sinθ)，sinθ∈[-1，1]，PSFCBF是CBF 算法的波束图，波束形状与sinθ=0 处的单目标的功率谱图相同。式（10）的物理意义是PSFCBF⊗Si是第i次迭代对PCBF的拟合，即预测的模糊图像，Perror=是PCBF与其预测图像之间的误差的度量，是以点扩散函数PSFCBF为加权值对估计误差Perror的加权求和，当前估计Si小于最优解时，下一次迭代Si+1会变大，当前估计Si大于最优解时，下一次迭代Si+1会变小，以趋近最优解。

3.3 最小方差无失真响应波束形成

传统波束形成聚焦于目标信号的相干累加，没有考虑其他方位信号（视为干扰）和噪声的信息，当不存在其他方位信号且环境噪声为白噪声时，则DAS 波束形成具有最大的输出信噪比。当存在其他方位信号或环境噪声非白噪声时，如能有效利用干扰和噪声的统计信息则会实现比DAS更高的输出信干噪比。Capon 法或MVDR 的角度分辨率好于CBF，该算法将阵列自由度分成两部分，其中一部分用于在期望角度处形成波束，另一部分自由度用于排除干扰的影响，即在非期望信号方向设计较深的凹陷。Capon 法的准则是让阵列的输出功率达到最小同时又保证在期望角度方向上的增益维持不变，因此可最大程度地排除干扰的影响［14］。即

运用拉格朗日算子可将上式转化为无约束最优化模型，其解为：

从上式可以看出，权值W依赖于数据协方差矩阵R的逆，因此MVDR属于数据依赖的自适应算法。

3.4 解卷积MVDR

对现有算法进行解卷积后处理是提升算法性能的有效技术，现把式（12）表示为

将式（14）代入式（13），

其中，σ2是噪声功率δ(sinθ-sinθk)是多目标角度的真实空间分布，PSFMVDR=|WH(sinθ)a(0)|2是MVDR的PSF。

根据PSF 的物理意义可知，一种可行的确定PSFMVDR的方法是：假设感兴趣区间的中心位置有一个点目标，目标的功率为MVDR 方位谱谱峰的平均值，即根据式PMVDR(sinθ)=计 算PSFMVDR，其 中R=P0a(0)a(0)H+I，P0是MVDR 估计出的声源的功率的平均值，〚·〛1表示归一化操作。图3直观显示了CBF算法与MVDR算法的PSF，显然，CBF 的方位谱主瓣更宽，旁瓣更高，分辨能力差，MVDR 算法的PSF 更加尖锐，分辨能力更高。

3.4.1 Richardson-Lucy解卷积MVDR

在图像处理领域RL 算法是一种已被理论证实和工程实际证明有效的后处理解卷积算法，RL解卷积MVDR迭代公式如下所示，

其中，PMVDR=PMVDR(sinθ)，sinθ∈[-1，1]。

文献［21］应用矢量外推算法提出了一种快速RL算法，本文仿真计算均是基于快速RL算法，详细执行过程请查看文献［21］。

RL解卷积MVDR（MVDR-RL）的执行流程框图如下所示：

3.4.2 FISTA解卷积MVDR

Oliver 等人提出可以采用快速迭代收缩阈值算法（Fast Iterative Shrinkage-Thresholding Algorithm，FISTA）完成解卷积过程［20］，FISTA 解卷积MVDR（MVDR-FISTA）算法执行流程图如下所示。注意RL算法与FISTA 算法的使用方式完全一样，算法的输入和输出完全一样，均是首先计算方位谱，然后采用解卷积算法利用PSF对方位谱解卷。

4 仿真实验

4.1 单次仿真实验

本节通过仿真分析提出的两种解卷积MVDR算法的性能。假定阵列为半波长布阵的均匀直线阵列，阵元数为M=10，两个不等强度的非相干声源波达方向为sinθ1=-0.08，sinθ2=0.1，SNR 分别为5 dB 和3 dB。快拍数为100（对于窄带信号一个快拍数据是一个时域采样点，采样率为100 kHz时，100 个快拍数据对应1 ms 时长），为防止卷绕误差，将sinθ∈[-1，1]扩展为[-2，2]。其中SNR 的定义为

其中，代表信号能量，代表噪声能量。

CBF 算法是阵列信号处理领域的常用算法，其波束主瓣的3 dB 宽度约为110o，由于CBF 算法主瓣宽、旁瓣高，当多个目标的角度间隔小于等于波束主瓣时，CBF 算法无法分辨［1］。通过图4 可以看出CBF 算法无法准确估计这两个声源信号的DOA。CBF-RL 算法通过解卷积运算，降低了CBF算法宽主瓣，高旁瓣的不利影响，能够成功分辨出两个信号，但是估计的声源位置有偏移。MVDR 算法能够分辨出两个声源，但是由于功率谱图的谱峰和谱谷的差值为6 dB，没有将两个谱峰完全分离。CBF-RL，MVDR-RL，MVDR-FISTA 三种算法能够准确分辨两个声源，由于采用解卷积后处理，两个声源的谱峰完全分离。其中MVDR-RL 和MVDRFISTA算法的方位谱主峰宽度一致，为CBF-RL算法主峰宽度的16。

将上述实验进行1000次仿真，统计不同算法的迭代次数和计算时间如表1所示。

表1 不同算法计算效率（Intel（R）CPU：i5-1035G1@3.4 GHz）Tab.1 Computational efficiency of different algorithms

表1 表明，CBF-RL 算法的计算时间最短，但是对比图4，CBF-RL 算法的并没有准确估计目标方位，因此没有寻找到全局最优值。MVDR-RL 的迭代次数少于MVDR-FISTA 算法，两种算法的迭代策略不同，RL 算法采用矢量加速的算法进行迭代，更容易寻找到最优值，FISTA 算法采用梯度下降思想进行迭代计算，收敛速度较慢。但是由于RL 算法需要计算加速因子，因此平均每次迭代的计算时间比FISTA算法的计算时间更长。

4.2 时间方位历程谱图仿真实验

假设有两个不等强度的非相干目标，其中目标1 的角度sinθ=-0.05，SNR=5 dB，目标2 的角度sinθ从-0.5 到0.5 变化，SNR=3 dB，每次方位估计取100 个快拍数据，以下是不同算法的时间方位谱估计结果。

从图5可以看出，与其他算法相比较，CBF 算法的功率谱图的主瓣宽，当运动目标角度在sinθ=±0.2范围内时，两个目标不可分，而且旁瓣较高较多。图6 表明CBF-RL 算法能够提高CBF 算法的分辨能力，这是由于解卷积后处理可以降低CBF 算法波束主瓣宽，旁瓣高的不利影响，CBF-RL 主瓣更窄，能够去除大部分旁瓣，功率谱的背景能量更低。当运动目标角度在sinθ=±0.15 范围内，CBF-RL 算法难以分辨两个目标。图7表明MVDR算法的估计结果好于CBF 算法，这是由于MVDR 的优势是抑制干扰，属于高分辨算法。图8 和图9 是分别采用RL 算法和FISTA 算法对MVDR 的功率谱进行解卷积后处理，两者差别不大，都能够提高MVDR 的性能，其分辨能力高于CBF 算法、CBF-RL 算法和MVDR 算法。解卷积MVDR 的背景级更低，几乎没有旁瓣。另外，在方位交汇附近，由于两个目标太近，此时所有算法均无法准确分辨两个声源，因此会形成模糊，并最终两个谱峰合并为一个谱峰。

通过上述仿真实验可以发现，CBF 算法通过旁瓣产生能量泄露，因此当两个目标相互接近时，方位估计结果会相互干扰，影响分辨性能。MVDR 算法对干扰会自适应形成凹陷，即形成空域陷波器，因此即使是两个相近的目标，也不会通过能量泄露，彼此干扰，本文从MVDR 方位估计的结果着手分析，实际上隐含地利用了MVDR 抑制干扰这一优势。仿真结果显示，解卷积MVDR 优于解卷积CBF算法。

毋庸置疑，解卷积CBF 算法理论非常完美，提高了CBF 算法的分辨能力，因此在水声方位估计领域，解卷积思想得到越来越多的重视。解卷积MVDR 尽管只是一种近似表示，但是仿真结果表明解卷积MVDR 算法提高了MVDR 分辨能力，这是由于MVDR 会抑制干扰，属于高分辨算法，对其进行解卷积，能够进一步提升分辨能力，其性能好于解卷积CBF。以下通过蒙特卡洛仿真实验进行详细说明。

4.3 蒙特卡洛统计实验

为进一步定量统计不同算法的均方根误差（Root Mean Square Error，RMSE），假设两个等强度的非相干信号，波达方向为sinθ1=-0.08，sinθ2=0.1，SNR 从-10 dB～10 dB 范围内变化，对于每个SNR 取值，进行100 次蒙特卡洛仿真实验，其中RMSE的定义为，

当方位谱的峰值与波谷相差3 dB 时认为成功分辨出了两个声源。在低SNR 条件下，算法存在无法以概率1 分辨两个声源的情况，本文将无法分辨定义为：1.两个较近的目标的波峰合为一个主峰，无法通过主峰的位置确定两个目标的波达方向，如图4 中的CBF 算法的方位谱；2.波峰与两个目标的波峰之间的波谷的差值小于3 dB。当第p次蒙特卡洛仿真实验无法分辨两个目标时，此时方位谱最高峰对应的角度为θpeak，此时令目标1，2 的估计角度

图10 是不同算法的分辨概率与SNR 的关系曲线，所有算法的分辨概率随SNR 的增大而增大。MVDR 算法的分辨概率最差，主要是由于在低SNR环境下，波谷与波峰之间的差值小于3 dB，分辨概率低。CBF-RL、MVDR-RL、MVDR-FISTA 算法由于采用解卷积后处理，消除模糊效应，因此分辨率更高。MVDR-RL 算法的分辨概率在低SNR 条件下，高于MVDR-FISTA 算法。当SNR 大于-6 dB 时，MVDR 解卷积算法的分辨概率达到1，当SNR 大于-4 dB时，CBF-RL算法的分辨概率达到1。

图11 显示了不同算法在不同SNR 条件下的RMSE，由于CBF 算法无法分辨角度间隔小于瑞利限的声源，因此没有展示CBF 的RMSE 曲线。通过图11 可以看出，本文提出的两种解卷积MVDR 算法均具有比传统MVDR 更低的RMSE，MVDR-RL和MVDR-FISTA 算法的估计性能非常接近。在低SNR 环境下，MVDR-RL 的分辨概率大于MVDRFISTA 算法，因此MVDR-RL 的估计精度明显好于MVDR-FISTA 算法。当信噪比大于2 dB 时，MVDR算法和解卷积MVDR 算法的估计精度相同。图12是不同算法的估计偏差，偏差的定义为：偏差=，由于偏差有正有负，因此正偏差和负偏差分开统计，图12纵轴的正半轴是偏差为正的部分，负半轴是偏差为负的部分，图12 结果与图11 基本一致，MVDR-RL 和MVDR-FISTA 算法的效果最好，低SNR 下，MVDR 的性能最差，估计偏差最大。

5 结论

为提高MVDR 算法的方位分辨性能，本文将MVDR 算法的输出功率谱重新表示为卷积的形式，并运用两种解卷技术对MVDR 的方位谱进行后处理，提升了算法的性能。该方法将角度空间中心位置的单个声源的MVDR 方位谱作为PSF，并利用RL算法和FISTA 算法分别对MVDR 的方位谱进行解卷积后处理，获得具有更低背景级的MVDR-RL 和MVDR-FISTA 方位谱，同时提高了分辨能力和估计精度。仿真实验显示了所提算法的良好性能。