一种水下强干扰环境下的宽带稀疏空间谱估计方法

本发明涉及水声目标入射噪声信号特征提取与识别，具体为一种水下强干扰环境下的宽带稀疏空间谱估计方法。

背景技术：

1、水中目标辐射噪声是指水面船只或水下航行器在航行时由于机械振动或螺旋桨旋转向周围水域辐射的宽带噪声，它可以被一定距离的水听器阵列接收。因此，水声目标辐射噪声是被动声呐探测目标的主要信息源。根据香农公式，水中目标辐射噪声的带宽越宽，所携带的目标信息越多，因此其对应的抗干扰能力也越强。目前，基于宽带目标辐射噪声的到达方向(doa)估计已在声呐上得到广泛应用。

2、然而，随着现代减震降噪技术的迅速发展，水中目标辐射噪声的源级在不断降低(甚至“隐身”)，此外，在复杂的水下环境中，弱目标的辐射噪声会被强干扰所掩盖，这些都会使得水中目标的探测变得困难。常用的高分辨率doa估计算法，例如多信号分类(music)算法、最小方差无失真响应(mvdr)算法和稀疏谱拟合(spsf)算法，都会变得失效。因此，对强干扰的有效抑制对于弱目标的doa估计有着重要意义。

3、矩阵滤波器是一种有效的强干扰抑制方法，它在抑制阻带强干扰的同时尽可能低失真地通过通带弱目标信号。传统的矩阵滤波器(cmf)在阻带提供了一个固定的抑制，当干扰的强度远超cmf提供的抑制时，残留的干扰将会破坏矩阵滤波后的协方差矩阵，这将会给后续的doa估计精度带来严重的影响。此外，cmf使用内点法去求解它的二阶锥规划(socp)问题会导致长的运算时间，这种长的运算时间在处理宽带信号时体现尤其明显。

4、申请号：201810677070.9，专利名称：一种基于稀疏分解的电磁目标空间谱估计方法，本发明专利主要在于解决水下强干扰环境下宽带声信号稀疏空间谱的估计问题，主要包含强干扰抑制和宽带稀疏空间谱估计两个部分。其中，强干扰抑制是通过本专利中的tnnr-mf空间矩阵滤波器来实现的；接着，基于矩阵滤波后的协方差矩阵，利用修改的稀疏谱拟合算法进行空间谱估计。

5、而对比专利的主要应用领域为电磁目标空间谱估计，而无强干扰抑制处理步骤；且其空间谱估计方法与本专利的空间谱估计方法(本专利的步骤(10))有着明显的差异，详见后续的不同点对比。

6、其他区别如下：

7、1)步骤(2)中的标校矩阵c(θ)用于对导向矩阵进行校准，而本专利步骤(10)中的矩阵滤波器g(fl)用于空间滤波。目的不同。

8、2)步骤(3)、(4)中将稀疏变换基用于信号模型构建，而本专利的步骤(10)中未使用稀疏变换基。

9、3)步骤(6)构建了一个二次型规划(socp)问题，本专利的步骤(10)利用矩阵滤波后的协方差矩阵构建了一个半定规划(sdp)问题。构建的优化问题有着本质不同。

10、本发明提出了一种水下强干扰环境下的宽带稀疏空间谱估计方法，基于截断核范数规则化的快速自适应矩阵滤波器能够在自适应地抑制强干扰的同时抑制矩阵滤波器引入的误差，保证了滤波后的协方差矩阵的准确性，从而保证了后续空间谱估计的精度。

技术实现思路

1、为了针对水下强干扰环境下弱目标的doa估计问题，本发明提出了一种水下强干扰环境下的宽带稀疏空间谱估计方法，该方法基于截断的核范数规则化构建了自适应矩阵滤波器优化问题，其能够自适应地在强干扰方向形成抑制零点并抑制矩阵滤波器引入的误差，构建的矩阵滤波器优化问题通过一个快速的两层apgl算法进行求解。基于apgl算法求解后同时获得滤波后的准确协方差矩阵和矩阵滤波器操作符，接着一个修改的spsf算法被给出用来进行空间谱估计。

2、为实现上述目的，本发明采取的技术方案是：

3、一种水下强干扰环境下的宽带稀疏空间谱估计方法，包括如下步骤：

4、(1)读入阵列水声目标辐射噪声数据，将其分段，并将每段信号分别变换到频域；

5、(2)计算给定频点fl处的平均样本协方差矩阵；

6、(3)将方位角空间[-90,90]°离散化成q个候选方向，利用样本协方差矩阵构建水下强干扰环境下基于截断核范数规则化的矩阵滤波器tnnr-mf优化问题，并初始化其参数；

7、步骤(3)具体包括如下步骤：

8、(3.1)将方位角空间[-90,90]°离散化成q个候选方向θ＝[θ1,θ2,...,θq]t，(·)t表示转置操作符；

9、(3.2)根据步骤(2)得到的平均样本协方差矩阵在频点fl处构建水下强干扰环境下tnnr-mf的优化问题：

10、

11、其中，(·)h表示共轭转置操作符，tr(·)、||·||*和||·||f分别表示矩阵的迹、核范数和frobenius范数，ε、δ和σ为容限参数，τ和β为规则化参数，r为截断核范数参数，和分别为通带和阻带区间，a(fl,θp)和a(fl,θs)分别为θp和θs对应的导向矩阵，y(fl)和g(fl)分别为频率fl处的滤波后的协方差矩阵和矩阵滤波器操作符，和分别为如下所示的一阶导数算子和选择矩阵：

12、

13、(3.3)将步骤(3.2)中的优化问题约束移动到目标函数，得到如下优化问题：

14、

15、(3.4)设定tnnr-mf优化问题规则化参数τ、β、λ、γ和η，设定截断核范数参数r，设定内部迭代容限εin和外部迭代容限εout；

16、(4)初始化外部迭代参数，并对进行奇异值分解，计算获得对应的前r个左、右奇异向量矩阵，接着开始外部迭代，在第i次外部迭代中，将对应的前r个左、右奇异向量矩阵的共轭转置分别赋值给ai、bi；

17、(5)初始化内部迭代变量；

18、(6)开始内部迭代，利用加速近端梯度法apgl，在给定的点γ(fl)处构建关于优化问题的一个近似，接着利用构建的优化问题近似来更新y(fl)和矩阵滤波器操作符g(fl)；

19、步骤(6)具体包括如下步骤：

20、(6.1)设定z(fl):＝(y(fl),g(fl))，令：

21、g(z(fl))＝||y(fl)||*

22、以及：

23、

24、利用g(z(fl))和hi(z(fl))，将步骤(3.3)中的优化问题转换为如下优化问题：

25、

26、(6.2)利用apgl算法求解，在给定的点γ(fl)＝(h1(fl),h2(fl))处构建关于函数f(z(fl))的近似，即：

27、

28、上式转换为下列形式：

29、

30、(6.3)通过迭代的更新z(fl)、γ(fl)和t可完成步骤(6.1)中优化问题的求解，在第k次内部迭代中，使用q(z(fl),γ(fl))的唯一极小值来更新：

31、

32、其中：

33、

34、

35、式中，h1(fl)和h2(fl)的初始值分别设定为包含以下组成部分：

36、

37、

38、上述式中：

39、

40、t1＝a(fl,θs)ah(fl,θs)

41、t2＝a(fl,θp)ah(fl,θp)

42、

43、(6.4)基于步骤(6.3)中的优化问题，和通过下面的方式进行更新：

44、

45、

46、上述优化问题中，通过奇异值阈值法svt进行求解，如下所示：

47、

48、式中，σm,m＝1,...,m、u和v分别是的奇异值、左和右奇异值向量矩阵；

49、(7)若优化问题函数下降速度小于设定的阈值或内部迭代次数大于设定的最大值，则退出内部迭代，否则返回步骤(6)继续迭代；

50、(8)将步骤(7)的迭代结果赋值给第i次外部迭代结果；

51、(9)若关于y(fl)和g(fl)变化值的frobenius范数小于设定的阈值或外部迭代次数大于设定的最大值，则退出外部迭代，否则返回步骤(4)继续迭代；

52、(10)基于y(fl)和g(fl)的外部迭代结果，利用修改的稀疏谱拟合算法spsf估计频点fl处的空间谱；

53、步骤(10)具体包括如下步骤：基于步骤(9)获得的在频率fl处的滤波后的协方差矩阵y(fl)和矩阵滤波器操作符g(fl)，通过求解下面半定规划sdp问题获得在频率fl处的空间谱：

54、

55、其中，vec(·)为矩阵垂直化操作符，||·||1和||·||2分别表示向量的1范数和2范数，*和分别表示共轭操作和kronecker积，μ为规则化参数，矩阵的第q列为式中l(fl,θq)为矩阵的第q列，q＝1,...,q，其中a(fl,θ)为θ对应的导向矩阵，为估计的空间功率谱，为估计的高斯白噪声功率；

56、(11)返回步骤(2)，依次完成所有频点对应的空间谱估计；

57、步骤(11)具体包括如下步骤：重复步骤(2)到步骤(10)，依次完成所有频点对应的空间谱估计；

58、(12)将所有频点对应的空间谱估计结果累加，得到最终的宽带空间谱估计；

59、步骤(12)具体包括如下步骤：将所有频点对应的空间谱估计结果累加，得到最终的宽带空间谱估计：

60、

61、作为本发明进一步改进，步骤(1)具体包括如下步骤：

62、(1.1)读入m阵元阵列水声目标辐射噪声时域数据，将其划分成b段；

63、(1.2)设定离散傅立叶变换(dft)点数l，将b段时域信号分别变换到频域，得到阵列在频点fl处关于第b段信号的输出为xl(b)，b＝1,...,b，l＝1,...,l。

64、作为本发明进一步改进，步骤(2)具体包括如下步骤：阵列在频点fl处输出的平均样本协方差矩阵为：

65、

66、作为本发明进一步改进，步骤(4)具体包括如下步骤：

67、(4.1)初始化矩阵滤波器操作符g0(fl)＝jm、滤波后的协方差矩阵和外部迭代次数i＝1，设定最大外部迭代次数iin和最大内部迭代次数iout，其中jm是一个大小为m×m的全1矩阵；

68、(4.2)对协方差矩阵进行奇异值分解，可得在第i次外部迭代中，设定ai＝(u1,…,ur)h、bi＝(v1,…,vr)h分别为u和v的前r个奇异向量矩阵的共轭转置。

69、作为本发明进一步改进，步骤(5)具体包括如下步骤：初始化内部迭代次数k＝1、参数t1＝1，设定

70、作为本发明进一步改进，步骤(7)具体包括如下步骤：

71、(7.1)计算

72、(7.2)令k＝k+1，计算若δk＜εin或内部迭代次数大于设定的最大值iin，则退出内部迭代，否则返回步骤(6)中依次更新和tk+1、γk+1(fl)从而继续迭代。

73、作为本发明进一步改进，步骤(8)具体包括如下步骤：将步骤(7)的迭代结果赋值给第i次外部迭代结果，即

74、作为本发明进一步改进，步骤(9)具体包括如下步骤：

75、(9.1)计算

76、(9.2)令i＝i+1，若δi-1＜εout或外部迭代次数大于设定的最大值iout，则退出外部迭代并输出y(fl)＝yi-1(fl)和g(fl)＝gi-1(fl)，否则返回步骤(4.2)继续迭代。

77、与现有技术相比，本发明公开的方法具有以下优点：基于截断核范数构建的矩阵滤波器在自适应抑制强干扰的同时能够抑制矩阵滤波器引入的误差，从而获得准确的滤波后的协方差矩阵；快速的两层apgl迭代求解算法同时输出准确的滤波后的协方差矩阵和矩阵滤波器操作符，保证了后续spsf算法空间谱估计的精度。整个算法运算时间短、分辨率高，适用于快速、准确的宽带空间谱估计。