被动声纳方位估计的宽带盲源分离算法

范文涛，章新华，康春玉，蒋 飚

（1.海军大连舰艇学院 博士生队，大连 116018；2.海军大连舰艇学院 信息与通信工程系，大连 116018）

方位估计，或称波达方向估计（Direction of Arrival，DOA）是被动声纳表现目标态势的方式之一［1］，其中又以方位估计精度和方位分辨率作为衡量DOA估计性能的重要指标。但是在浅海环境下，尤其是航道附近，往往存在各种强干扰源，严重影响了被动声纳对于弱目标的方位分辨性能和估计精度［1］，因此改善现有被动声纳方位估计性能显得尤为迫切。

常规波束形成［2］是被动声纳中用于目标方位估计的常用方法，其优势在于与环境无关，性能稳健。但是受到孔径的限制，使得常规波束形成的方位分辨率始终无法突破瑞利限［2］。Capon［2］和 Schmidt［3］先后提出的基于阵列观测数据协方差矩阵的高分辨方位估计算法，即 MVDR［2］和 MUSIC［3］算法，弥补了常规波束形成在方位估计方面的缺点。然而受到低信噪比、阵形失配或者源数目的限制，使得高分辨类算法在处理实际被动声纳数据时表现不够稳健［4］。

近些年发展起来的盲信号处理［5］技术得到了研究者的广泛关注，其中尤以对盲源分离技术［6-8］的研究居多。由于盲源分离算法所具有的特殊优点［6］，一些研究者将其引入到了水声信号处理中［7］。文献［6-9］成功地从多个辐射噪声混合后的信号中分离出了所需源信号，验证了盲源分离用于分离水下目标源的有效性，同时指出盲源分离算法可以分离空间上几乎重叠的多目标源［9］。然而，这些方法同时存在对于海洋多途信道影响下的宽带卷积混合源分离［7］研究较少的局限性。日本Sawatari课题组［10］针对室内多途混响下的源分离问题，提出了多种波束形成与盲源分离融合的处理结构，但是这些方法中的目标DOA估计均是为解决宽带频点间次序模糊而服务的。范文涛［11］、康春玉［12］针对被动声纳DOA估计与波形恢复，分别提出了将高分辨DOA估计算法与基于卷积混合的频域盲源分离算法进行融合，改善了传统高分辨算法的DOA估计性能。可是对于宽带频域盲源分离算法存在的频点间次序模糊和频点内尺度模糊问题［9］，依旧没有较好的解决办法。

基于时间结构的盲源分离算法［11］使用多个时间延迟二阶相关矩阵来完成目标源的提取。但是这类算法同时存在时滞量选取难，对没有时间结构特征的声源的分离效果差等问题［11］。近期提出的非参数化独立分量分析算法［13］可以对目标源的统计特性没有任何先验知识的前提下完成目标源的提取，是一种真正意义上的盲处理算法。考虑将两者融合，既利用了目标源的时间结构特性，又可以在不需预知目标源统计特性的前提下完成源信号的分离，为了仿真比较方便，本文算法简称为时间结构非参数盲分离算法（Blind Temporal Structure Non-parametric Separation，BTSNPSEP）。本文根据宽带卷积混合［14］，在每个子带内将基于时间结构与非参数化特性的两种盲源分离算法进行融合，建立新的解混矩阵代价函数，然后根据非线性最优化得到的解混矩阵得到子带内方位能量谱，最后仿照宽带频域波束形成对所有子带方位能量谱累加，得到总的宽带方位谱。通过计算机仿真与实际海试数据，分别同经典MUSIC算法和MVDR方法进行了比较，结果表明本文方法能够用于宽带被动方位估计，尤其是强干扰背景下的弱目标检测。

1 宽带盲源分离算法

1.1 问题描述

考虑满足半波长且阵元间距为d的均匀线列阵（Uniform Iinear Array，ULA），阵元数为N，各阵元之间无互耦。远场宽带目标源数目小于阵列阵元数目，即（Q＜N）。由此可以得到第i个阵元的卷积混合为［12］：

其中：htq（l）表示从目标源q到阵元i之间的冲激响应；P表示冲激响应的长度，ni（t）为叠加在每个阵元上的噪声。

由于大部分盲源分离算法是基于瞬时混合模型的［7］，因此通过短时傅里叶变换 （Short-time Fourier transform，STFT）将时域卷积混合xi（t）变换为频域信号xi（k，fj），可以表示为：

其中：fj表示第j个频点，j=1，…，J；k=1，…，K表示频域快拍数。aiq（fj）为从目标源q到阵元i之间的频率响应。sq（k，fj），ni（k，fj）分别为sq（t）和ni（t）的短时傅里叶变换。由此可以得到第j个频点的阵列频域观测数据快拍矩阵［12］：

其中：X（k，fj）=［x1（k，fj），…，xN（k，fj）］T为阵列观测数据矩阵，S（k，fj）=［s1（k，fj），…，sQ（k，fj）］T为源信号矩阵，N（k，fj）=［n1（k，fj），…，nN（k，fj）］T为加性噪声矩阵。A（fj，Θ）=［aθ1（fj），…，aθQ（fj）］T表示阵列流形矩阵，其中aθ（fj）为：

1.2 基于时间结构的非参数盲源分离算法

根据文献［13］和式（3）可知，传统瞬时盲源分离算法通过求取解混矩阵来完成对目标源的分离和混合系统的辨识，即：

其中：Y（k，fj）=［y1（k，fj），…，yN（k，fj）］T为恢复出的目标源信号。W（fj）表示解混矩阵。

图1为本文所提算法框架。如图1所示，利用每个子带内估计的解混矩阵W（fj）和目标源信号Y（k，fj），求取目标方位能量谱，最后将所有子带内方位能量谱累加，即可得到总的宽带方位能量谱。

图1 算法框架Fig.1 The framework of the algorithm

文献［11］介绍了利用观测数据的时间结构特性进行盲源分离的思想，其代价函数：

其中：τb，b=1，…，B表示频点fj处的第b个时间延迟量，又称时滞量。

根据文献［13］，利用得到基于信息论准则的盲源分离代价函数：

其中：H（·）表示信息熵，H（X）为关于解混矩阵W的常量。

将式（6）与式（7）所表示的代价函数融合，可以得到基于时间延迟结构的盲源分离最优化问题：

其中：I（yi（k，fj）；yi（k+τb，fj））为：

由式（7）和（8）易得式（8）所对应的代价函数：

其中：h为内核带宽［13］，φ（·）为高斯核［13］：

Yim为内核质心［13］：

其中：X（m）（fj）表示每个子带内观测数据矩阵的第m个频域快拍。Wi（fj）表示解混矩阵W（fj）的第i行向量。

根据式（5），则yi（k，fj）=Wi（fj）X（fj），将式（13）代入式（11），可得：

目标函数（15）是解混矩阵W的非线性函数，约束（16）将最优化问题的可能解空间限定到一个有限集内。目标函数（15）可以通过牛顿法对代价函数（10）的非线性最优化解混矩阵Wopt，最优化求解过程同文献［12］类似。由此估计出的源信号Y（k，fj）=WoptX（k，fj）。

对于盲源分离恢复的结果，每一路输出的方位都能估计出来：

由式（17）在每个频点fj处估计出的N个θi（fj）值，将每个θi（fj）值在θ=-90°，…，90°方位空间中进行搜索匹配，找到其对应的方位位置，该方位位置的能量值为该路恢复出的源信号的能量。因此可以在每个子带内构造出类似于波束形成的方位能量谱：

对于式（18）需要说明的是，每个子带内盲源分离算法估计的N个θi（fj）值有可能会出现重复的情况，本文采取累加的办法。频域卷积盲源分离算法普遍存在频点间次序和频点尺度模糊的问题。然而式（18）并不受这两个问题的影响，即使是相邻频点间盲源分离输出次序不一致，对于最后的方位能量谱的构建不构成影响。仿照频域宽带波束形成的非相干累加的思想，最后总的方位能量谱为：

2 数据检验

2.1 仿真数据

宽带仿真实验考虑均匀线列阵位于两宽带目标源的远场，两目标源均用海上实录舰船辐射噪声代替。两目标源的原始时域波形与频谱如图2所示。宽带分析频段为800～1 600 Hz，阵元数为32，阵元间距满足分析频段最高频率对应的半波长。采样频率为25 kHz。环境噪声添加方式使用 MATLAB中 Awgn函数［12］对宽带阵列数据整体添加加性高斯白噪声。

图2 两目标源原始时域波形与频谱Fig.2 The original time domain waveform and frequency spectrum of the two targets

2.1.1 不同信噪比下方位分辨率比较

两个目标等功率，分析数据长度为1 s，变换不同环境噪声信噪比（-20 dB～10 dB），采用两种方位间隔，测试BTSNPSEP与MVDR的方位分辨性能，其中每个信噪比处的估计运行100次蒙特卡洛仿真［12］。图3为8°目标与11°目标（3°间隔）在不同信噪比条件下的方位分辨性能。图3（a）为BTSNPSEP方法结果，图3（b）为MVDR方法结果。比较图2中两个子图可以看出，BTSNPSEP方法在低信噪比时的分辨率要好于MVDR方法。

图4为8°目标与13°目标（5°间隔）在不同信噪比条件下的方位分辨性能。图4（a）为BTSNPSEP方法结果，图4（b）为MVDR方法结果。比较图3中两个子图可见当目标间距拉大以后，MVDR的性能变得较好，BTSNPSEP能够在低信噪比时能够区分两个目标，但对于13°目标的估计误差了1°。

2.1.2 不同目标源数目下方位分辨率比较

考虑到阵元数较多时的BTSNPSEP方法运算速度问题，通过加入白化预处理步骤［7］，使得BTSNPSEP方法也涉及到源数目设置的问题，由此同MUSIC方法比较不同源数目设置时的方位分辨率。两个目标依旧等功率，同时环境噪声信噪比固定为-5 dB，其中每个源数目设置处的估计运行100次蒙特卡洛仿真［12］。图5为8°目标与11°目标（3°间隔）在不同源数目设置下的方位分辨性能。图5（a）为BTSNPSEP方法结果，图5（b）为MVDR方法结果。图6为8°目标与13°目标（5°间隔）在不同源数目设置下的方位分辨性能。图6（a）为BTSNPSEP方法结果，图6（b）为MVDR方法结果。从图5（a）可以看出，当方位间隔3°时，源数目设置的变化对于BTSNPSEP方法影响还是非常大的，对11°目标的估计不如MUSIC算法。从图6可以看出当方位间距扩大（5°间隔），BTSNPSEP方法能够在各种源数目设置下区分两个目标，尤其是当源数目为1个时（欠定情况）要好于MUSIC。

2.1.3 信干比为-6 dB下弱目标检测能力

设定8°声源为干扰源，11°声源为所需目标源，信干比［12］为 -6 dB，环境噪声信噪比固定为 -5 dB，BTSNPSEP采用预白化处理源数目设置为3，MUSIC方法同之。共进行20次仿真，每次为100次蒙特卡洛仿真结果。图7为分别为BTSNPSEP、MVDR、MUSIC三种方法估计结果。比较图7三个图可以看出BTSNPSEP对于弱目标的检测效果要好于MVDR与MUSIC。

图7 信干比为-6 dB下弱目标检测能力比较Fig.7 The comparison of weak target detection capability under SIR=-6dB.

2.2 海上实测数据

试验海区位于东经121°33＇～121°42＇，北 38°46＇～38°52＇的主航道附近，海区平均水深46 m，声传播速度约1 500 m/s。海区地形复杂，水流急，过往船只多，离岸近、干扰源多，属典型的复杂海区［11］。接收阵列由28个水听器等间隔组成，阵元间距0.225 m，接收阵列深度约为20 m。水听器接收信号经信号调理机后到Sony sir1000i录音机磁带记录，整个记录数据期间，接收船一直辅机发电，有较大噪声。试验时目标船在大约5.3 km的距离沿正横经过接收船［11］。同时，目标船运动时在视觉范围内还发现有渔船目标运动，而且渔船先于目标船经过接收船正横。

图8为三种方法的估计得方位历程图比较，其中BTSNPSEP与MUSIC的源数目均设置为20。从图中可以看出对于-80°～-40°区域的弱干弱目标BTSNPSEP的估计结果相对清晰。第0～100 s时间，对于20°～40°区域的目标船和渔船的估计BTSNPSEP方法也比MUSIC、MVDR方法要清晰，空间增益提高比较明显。

图8 三种方法海试数据比较Fig.8 The comparison of sea trial result during three methods

3 结论

本文针对被动声纳信号的特殊性，结合两种盲源分离算法的特点，提出了一种融合时间延迟结构和非参化特性的盲源分离算法（BTSNPSEP），并将其用于被动声纳方位估计。分别利用宽带仿真和海试数据与经典MVDR和MUSIC高分辨方位估计方法进行了比较，得到以下结论：

（1）BTSNPSEP方法能够突破瑞利限的限制分辨空间接近的多目标，且在低信噪下的方位分辨性能要好于MVDR方法。

（2）BTSNPSEP方法能够在欠定条件下分辨出多目标源。

（3）BTSNPSEP方法比起MVDR和MUSIC方法，对于强干扰背景下的弱目标检测具有一定的优势。

［1］Meng T Z，Buck J R.Rate distortion bounds on passive sonar performance［J］.IEEE Trans.on Acoustics.Speech and Signal Processing，2010，58（1）：326-336.

［2］ Capon J.High resolution frequency wavenumber spectrum analysis［C］//Proc.of the IEEE，1969，57（8）：1408-1418.

［3］Stoica P，Nehorai A.Music，maximum likelihood，and cram r-rao bound［J］.IEEE Trans.on Acoustics.Speech and Signal Processing，1989，37（5）：720-741.

［4］ Bao C Y.Adaptive beamforming for sonar audio［J］.Proc.of ACOUSTICS 2005.Australia：Australian Acoustical Society Press，2005：483-485.

［5］蔡艳平，李艾华，石林锁，等.基于盲解卷积的柴油机振动信号分离研究［J］.振动与冲击，2010，29（9）：38-41.

［6］Benchekroun N，Mansour A.A general structure for the separation of underwater acoustic signals［C］//Proc.of Symposium on Communications，Control and Signal Processing.Marrakech：IEEE Press，2006：273-278.

［7］De Moura，Seixas N N，Filho J M.Independent component analysis for optimal passive sonar signal detection［C］//Seventh International Conference on Intelligent Systems Design and Applications.Rio de Janeiro：IEEE Press，2007：671-678.

［8］雷衍斌，李舜酩，门秀花，等.基于自相关降噪的混叠转子振动信号分离［J］.振动与冲击，2011，30（1）：218-222.

［9］张安清.盲分离技术及其在水声信号中的应用研究［D］.大连：海军大连舰艇学院，2006.

［10］Saruwatari H，Kurita S，Takeda K.Blind source separation combining independent component analysis and beamforming［J］.EURASIP Journal on Applied Signal Processing，2003，1135-1146.

［11］康春玉，章新华，韩 东.盲源分离与高分辨融合的DOA估计与信号恢复方法［J］.自动化学报，2010，36（3）：443-445.

［12］章新华，范文涛.波束形成与独立分量分析融合的宽带高分辨方位估计方法［J］.声学学报，2009，34（4）：175-180.

［13］ Boscolo R，Pan H，Roychowdhury V P.Independent component analysis based on nonparametric density estimation［J］.IEEE Trans.on Neural Networks，2004，15（1）：55-65.

［14］范 涛，李志农，岳秀廷.基于变分贝叶斯算法和MLP网络的后非线性混合盲源分离方法研究［J］.振动与冲击，2010，29（6）：21-24.