适用于浅海多途环境的水下噪声源高稳健性聚焦定位方法

摘要

本发明提供的是适用于浅海多途环境的水下噪声源高稳健性聚焦定位方法。利用相干多途信道特点，在射线理论基础上建立了符合水声传播特点的阵列信号模型，生成与实际声传播特性相匹配的空间聚焦导向矢量，从而有效克服了多途效应的影响，提高了定位精度；利用最差性能最优的稳健聚焦处理器设计方法，对空间聚焦导向矢量实施约束，并通过二阶锥规划方法求解最优权矢量，最终得到优化后的空间谱形式，从而克服了多种海洋环境声学参数失配误差对MVDR高分辨算法空间分辨率及背景起伏抑制能力的影响，可以获得更大的动态范围、更尖锐的聚焦峰尺度以及更强的背景噪声级抑制能力，可以满足在浅海多途条件下，水下噪声源定位对于高分辨算法稳健性的迫切需求。

说明书

技术领域

本发明涉及一种水声定位方法，具体地说是一种适用于浅海多途环境的水下噪声源高稳健性聚焦定位方法。

背景技术

在浅海环境下，声信道是一个包括海面、海底和海水介质的复杂环境，存在明显的多途效应。对于浅海近程噪声源定位问题而言，多途效应的存在会导致通常使用的常规及高分辨类聚焦定位方法存在模型失配，产生定位偏差，同时在空间谱图像上出现较高的旁瓣级起伏，对聚焦定位性能产生严重影响。

针对浅海环境下的信号处理问题，基于模型的信号处理思想（Model Based Processing：MBP）受到了广大学者的青睐，该思想是设计一种融合海洋传播模型，并且可以用来完成各种各样的信号处理功能的处理器。Candy等人将MBP方法广泛的应用到水下声纳被动定位，海洋环境参数反演，水下目标检测、估计、识别等方面，取得了一定的成果。目前，在浅海噪声源聚焦定位方面，对浅海多途效应的影响以及可行性方法的研究普遍集中于时间反转镜方法（惠娟，胡丹，惠俊英等．聚焦波束形成声图测量原理研究．声学学报,2007，32(4)：356-361；罗方方,生雪莉,梅继丹,郭咏.基于MVDR高分辨算法的时反定位技术研究.哈尔滨工程大学学报，2010,31(7)：945-950），但仍未出现基于多途模型匹配思想的噪声源聚焦定位方法研究结果。

发明内容

本发明的目的在于提供了一种可有效提高MVDR高分辨聚焦定位方法在浅海多途条件下的定位精度以及存在失配误差及海洋环境声学参数失配误差下的稳健性的适用于浅海多途环境的水下噪声源高稳健性聚焦定位方法。

本发明所采用的技术方案包括以下步骤：

（a）在近程声传播条件下采用射线理论中的虚源法对阵列信号进行建模，建立浅海多途条件下的水平声压线阵接收信号模型，得到基阵接收信号矩阵X及采样数据协方差矩阵

（b）在与声源等深的水平面S上进行逐点扫描，在扫描点坐标及分析频率f上生成空间聚焦导向矢量A'(f)；

（c）为消除由于不同扫描点到达基阵参考阵元的距离不同带来的对匹配聚焦的影响，对空间聚焦导向矢量进行二范数下的归一化处理，得到归一化后的空间聚焦导向矢量

（d）利用最差性能最优方法对空间聚焦导向矢量实施约束优化；

（e）利用二阶锥规划对该优化问题进行求解，得到最优化权矢量w_opt；

（f）将最优权矢量w_opt代入目标函数中，得到最优权矢量下的阵列输出功率P；

（g）根据运算量要求，设置合适扫描步长，重复（b）至（f）的步骤，进行完整平面搜索并比较输出功率谱图，由谱峰位置确定噪声源所在位置。

本发明将多途信道模型与最差性能最优稳健约束优化方法相结合，提出了具有高稳健性的浅海噪声源聚焦定位方法。

本发明的有益效果是：该方法一方面综合了多途模型匹配思想，利用相干多途信道特点，在射线理论基础上建立了符合水声传播特点的阵列信号模型，生成与实际声传播特性相匹配的空间聚焦导向矢量，从而有效克服了多途效应的影响，提高了定位精度。另一方面，利用最差性能最优的稳健聚焦处理器设计方法，对空间聚焦导向矢量实施约束，并通过二阶锥规划方法求解最优权矢量，最终得到优化后的空间谱形式，从而克服了多种海洋环境声学参数失配误差对MVDR高分辨算法空间分辨率及背景起伏抑制能力的影响，可以获得更大的动态范围、更尖锐的聚焦峰尺度以及更强的背景噪声级抑制能力，可以满足在浅海多途条件下，水下噪声源定位对于高分辨算法稳健性的迫切需求。

附图说明

图1浅海水平线阵接收信号模型示意图。

图2(a)-图2(d)空间聚焦导向矢量失配下的空间谱对比效果图，其中：图2(a)多途常规；图2(b)匹配常规；图2(c)匹配MVDR；图2(d)匹配稳健MVDR。

图3(a)-图3(b)空间聚焦导向矢量失配下的空间谱切片对比效果图（N=21，d=0.75m），其中：图3(a)y向；图3(b)x向。

图4(a)-图4(b)空间聚焦导向矢量失配下的空间谱切片对比效果图（N=11，d=1.5m），其中：图4(a)y向；图4(b)x向。

图5(a)-图5(b)空间聚焦导向矢量失配下的空间谱切片对比效果图（N=7，d=2.5m），其中：图5(a)y向；图5(b)x向。

图6(a)-图6(b)空间聚焦导向矢量失配下的空间谱切片对比效果图（N=5，d=3.75m），其中：图6(a)y向；图6(b)x向。

图7(a)-图7(b)水深存在5%误差下的空间谱切片（ε=0.25），其中：图7(a)y向；图7(b)x向。

图8(a)-图8(b)水中声速取1480m/s下的空间谱切片（ε=0.03），其中：图8(a)y向；图8(b)x向。

图9(a)-图9(b)海水和海底介质的密度之比m_ρ存在30%的误差下的空间谱切片（ε=0.03），其中：图9(a)y向；图9(b)x向。

图10(a)-图10(b)海水与海底中的声速之比n_c存在30%的误差下的空间谱切片（ε=0.2），其中：图10(a)y向；图10(b)x向。

具体实施方式

下面结合附图举例对本发明进一步说明。

（a）在近程声传播条件下采用射线理论中的虚源法对阵列信号进行建模，建立浅海多途条件下的水平声压线阵接收信号模型，得到基阵接收信号矩阵X及采样数据协方差矩阵

以浅海水平声压线阵测试系统为例（如图1所示）。设海水深度为H，一个N元均匀线阵入水深度为z_a，阵元间距为d，设1号阵元为参考阵元。空间中共存在M个点声源，入水深度均为z_s，第m个单频点声源的空间位置坐标为(x_m,y_m,z_s)。

略去时间因子，基阵第n号参考阵元的接收到的由第m个声源发射的声压信号可表示为：

$>p^{(m,n)}=\underset{l=0}{\overset{\infty }{\Sigma }}[{\left({\mathrm{VV}}_{l1}^{(m,n)}\right)}^l\frac{e^{\mathrm{jk}{R}_{l1}^{(m,n)}}}{R_{l1}^{(m,n)}}+{\left({\mathrm{VV}}_{l2}^{(m,n)}\right)}^l{V}_{l2}^{(m,n)}\frac{e^{\mathrm{jk}{R}_{l2}^{(m,n)}}}{R_{l2}^{(m,n)}}+{\left({\mathrm{VV}}_{l3}^{(m,n)}\right)}^{l}V\frac{e^{\mathrm{jk}{R}_{l3}^{(m,n)}}}{R_{l3}^{(m,n)}}+{\left({\mathrm{VV}}_{l4}^{(m,n)}\right)}^{l+1}\frac{e^{\mathrm{jk}{R}_{l4}^{(m,n)}}}{R_{l4}^{(m,n)}}]---\left(1\right)$>

$>\left(\begin{array}{c}{Z}_{l1}=2\mathrm{Hl}+z_s-z_a\\ Z_{l2}=2H(l+1)-z_s-z_a\\ Z_{l3}=-2\mathrm{Hl}-z_s-z_a\\ Z_{l4}=-2H(l+1)+z_s-z_a\end{array}\right)---\left(2\right)$>

$>R_{\mathrm{li}}^{\left(\text{m,n}\right)}=\sqrt{{(x_m-x_a\left(n\right))}^2+{(y_m-y_a\left(n\right))}^2+{\left(Z_{\mathrm{li}}\right)}^2}---\left(3\right)$>

其中，k=2πf/c为波数。表示第l阶虚源第i条声线所对应的虚源与观察点之间的距离，Z_li为其垂直距离，i＝1~4即每增加一阶虚源所增加的4条声线。当l＝0时，Z₀₁代表直达声，Z₀₂代表声线经过一次海底反射，Z₀₃代表声线经过一次海面反射，Z₀₄代表声线先后经过一次海底和一次海面反射；V表示海面反射系数，表示第l阶虚源第i条声线的海底反射系数，满足瑞利反射，其表达式如下：

$>V_{\mathrm{li}}^{(m,n)}=\frac{m_{\rho }\cos {\theta }_{\mathrm{li}}^{(m,n)}-\sqrt{n_c^2-{\sin }^2{\theta }_{\mathrm{li}}^{(m,n)}}}{m_{\rho }\cos {\theta }_{\mathrm{li}}^{(m,n)}+\sqrt{n_c^2-{\sin }^2{\theta }_{\mathrm{li}}^{(m,n)}}}---\left(4\right)$>

式中，m_ρ为海水和海底介质的密度之比，n_c为海水与海底中的声速之比，为第l阶虚源第i条声线对应的俯仰角。可看出，随着海底入射角度的变化而变化，一般来讲，它是一个具有实部和虚部的复数，这是因为在反射时产生相移。由于测量时间较短，可以认为介质的密度、声速和边界的反射系数不随时间变化。

将基阵接收信号写为如下矩阵形式：

$>X=\left(\begin{array}{c}{x}_1\\ x_2\\ ·\\ ·\\ ·\\ x_N\end{array}\right)=\underset{m=1}{\overset{M}{\Sigma }}\left(\begin{array}{c}{P}^{(m,1)}\\ P^{(m,2)}\\ ·\\ ·\\ ·\\ P^{(m,N)}\end{array}\right)s^{\left(m\right)}=\underset{m=1}{\overset{M}{\Sigma }}{A}^{\left(m\right)}{s}^{\left(m\right)}---\left(5\right)$>

A^(m)为第m个声源对应的浅海多途条件下的空间聚焦导向矢量。

由(6)式可知，空间聚焦导向矢量的构造与传统方位估计中导向矢量的生成存在本质上的不同，空间聚焦导向矢量是多根本征声线贡献的叠加，具有明显的相干结构，不同本征声线在空间聚焦导向矢量中的贡献将会对定位结果产生严重影响。

对于序列长度L，采样数据协方差矩阵可表示为：

$>\hat{R}={\mathrm{XX}}^H/L---\left(7\right)$>

（b）在与声源等深的水平面S上进行逐点扫描，在扫描点坐标及分析频率f上生成空间聚焦导向矢量A'(f)。

$>R_{\mathrm{li}}^{\prime \left(n\right)}=\sqrt{{(\hat{x}-x_a\left(n\right))}^2+{(\hat{y}-y_a\left(n\right))}^2+{\left(Z_{\mathrm{li}}\right)}^2}---\left(9\right)$>

$>V_{\mathrm{li}}^{\prime \left(n\right)}=\frac{m_{\rho }\cos {\theta }_{\mathrm{li}}^{\prime \left(n\right)}-\sqrt{n_c^2-{\sin }^2{\theta }_{\mathrm{li}}^{\prime \left(n\right)}}}{m_{\rho }\cos {\theta }_{\mathrm{li}}^{\prime \left(n\right)}+\sqrt{n_c^2-{\sin }^2{\theta }_{\mathrm{li}}^{\prime \left(n\right)}}}---\left(10\right)$>

其中，为扫描点至第n号阵元的距离，表示扫描点对应的第l阶虚源第i条声线的海底反射系数。

（c）为消除由于不同扫描点到达基阵参考阵元的距离不同带来的对匹配聚焦的影响，对空间聚焦导向矢量进行二范数下的归一化处理，得到归一化后的空间聚焦导向矢量

$>{\hat{A}}^{\prime }\left(f\right)=\frac{A^{\prime }\left(f\right)}{{\left|\right|A^{\prime }\left(f\right)\left|\right|}_2}---\left(11\right)$>

空间聚焦导向矢量的生成是至关重要的一步，不同聚焦算法的差异仅体现在使用了不同的处理器形式。本专利中将采用常规处理器及MVDR处理器的算法分别称为匹配常规与匹配MVDR算法。

（d）利用最差性能最优方法对空间聚焦导向矢量实施约束优化。

空间聚焦导向矢量误差ΔA的范数可以由常数ε＞0进行约束：

||ΔA||≤ε                               (12)

则实际的空间聚焦导向矢量属于下面的集合：

$>A\left(ϵ\right)=\left\{{\hat{A}}_{\mathrm{true}}^{\prime }\left(f\right)\right|{\hat{A}}_{\mathrm{true}}^{\prime }\left(f\right)={\hat{A}}^{\prime }\left(f\right)+\mathrm{\Delta A},\left|\right|\mathrm{\Delta A}\left|\right|\le ϵ\}---\left(13\right)$>

其中，即为由（10）式得到的理论上无任何失配误差下的空间聚焦导向矢量。

对所有属于集合A(ε)的导向矢量进行约束，即阵列响应的绝对值不小于1：

$>\left|w^H{\hat{A}}_{\mathrm{true}}^{\prime }\left(f\right)\right|\ge 1,{\hat{A}}_{\mathrm{true}}^{\prime }\left(f\right)\in A\left(ϵ\right)---\left(14\right)$>

则利用最差性能最优方法得到的约束优化问题可表示为：

$>\left(\begin{array}{c}\underset{w}{\min }{w}^H\hat{R}w\\ s.t.\left|w^H{\hat{A}}_{\mathrm{true}}^{\prime }\left(f\right)\right|\ge 1\mathrm{forall}{\hat{A}}_{\mathrm{true}}^{\prime }\left(f\right)\in A\left(ϵ\right)\end{array}\right)---\left(15\right)$>

（e）利用二阶锥规划对该优化问题进行求解，得到最优化权矢量w_opt。

式（15）可进一步等价为具有单一非线性约束形式的二次最小化问题：

$>\left(\begin{array}{c}\underset{w}{\min }{w}^H\hat{R}w\\ s.t.w^H{\hat{A}}^{\prime }\left(f\right)\ge ϵ\left|\right|w\left|\right|+1,\mathrm{Im}\left\{w^H{\hat{A}}^{\prime }\left(f\right)\right\}=0\end{array}\right)---\left(16\right)$>

对采样数据协方差矩阵进行Cholesky分解并将约束优化问题中的目标函数转化为：

$>\underset{w}{\min }{w}^H\mathrm{Rw}=\underset{w}{\min }{\left(\mathrm{Uw}\right)}^H\left(\mathrm{Uw}\right)=\underset{w}{\min }{\left|\right|\mathrm{Uw}\left|\right|}^2---\left(17\right)$>

引进一个非负标量c，并构造一个新的约束||Uw||≤c，则约束优化问题可进一步写为：

$>\left(\begin{array}{c}\underset{c,w}{\min }c\\ s.t.ϵ\left|\right|w\left|\right|\le w^H{\hat{A}}^{\prime }\left(f\right)-1,\mathrm{Im}\left\{w^H{\hat{A}}^{\prime }\left(f\right)\right\}=0,\left|\right|\mathrm{Uw}\left|\right|\le c\end{array}\right)---\left(18\right)$>

上式的实值形式可表示为：

其中，$>\stackrel{\cup }{w}={[\mathrm{Re}{\left\{w\right\}}^T,\mathrm{Im}{\left\{w\right\}}^T]}^T,$>$>\stackrel{\cup }{{\hat{A}}^{\prime }\left(f\right)}={[\mathrm{Re}{\left\{{\hat{A}}^{\prime }\left(f\right)\right\}}^T,\mathrm{Im}{\left\{{\hat{A}}^{\prime }\left(f\right)\right\}}^T]}^T,$>

$>\stackrel{\cup }{U}=\left(\begin{array}{cc}\mathrm{Re}\left\{U\right\}& -\mathrm{Im}\left\{U\right\}\\ \mathrm{Im}\left\{U\right\}& \mathrm{Re}\left\{U\right\}\end{array}\right)$>

上式可进一步表达为：

$>\left(\begin{array}{c}\min {\underset{y}{d}}^{T}y\\ s.t.f+F^{T}y\in \mathrm{SO}{C_1}^{2N+1}\times {{\mathrm{SOC}}_2}^{2N+1}\times \left\{0\right\}\end{array}\right)---\left(20\right)$>

其中，d＝[1,0^T]∈R^(2N+1)×1，f＝[0^T,-1,0^T]^T∈R^(4N+1)×1，

利用Sedumi软件对约束优化问题进行求解，得到的最优化权矢量w_opt：

$>w_{\mathrm{opt}}={[\stackrel{\cup }{w_1},...,\stackrel{\cup }{w_N}]}^T+j{[\stackrel{\cup }{w_{N+1}},...,\stackrel{\cup }{w_{2N}}]}^T---\left(21\right)$>

（f）将最优权矢量w_opt代入目标函数中，得到最优权矢量下的阵列输出功率P。

$>P=w_{\mathrm{opt}}^H{\hat{R}w}_{\mathrm{opt}}---\left(22\right)$>

将本专利中的新方法称为匹配稳健MVDR算法。仿真实例中将对匹配常规、匹配MVDR算法以及匹配稳健MVDR算法的性能进行对比分析。

（g）根据运算量要求，设置合适扫描步长，重复（b）至（f）的步骤，进行完整平面搜索并比较输出功率谱图，由谱峰位置确定噪声源所在位置。

上面对发明内容各部分的具体实施方式进行了说明。综合多途模型匹配思想和稳健聚焦优化设计方法，可有效提高MVDR高分辨聚焦定位方法在浅海多途条件下的定位精度的稳健性，下面对仿真实例进行分析。

实例一：空间聚焦导向矢量扰动下的空间谱

实例参数设置如下：利用空间聚焦导向矢量扰动量来度量其存在失配误差的程度。

坐标系如图1所示，点声源入水深度4m，其直角坐标系下的坐标为(20,5,-4)m，基阵入水深度为10m，阵元个数11个，阵元间距1.5m。海水深度为15m，海水中声速为1500m/s，海水和海底介质的密度之比m_ρ为1.8，海水与海底中的声速之比n_c为0.85，设海面反射系数为-1。系统采样率为20kHz，单频信号频率为1kHz，信噪比为20dB，数据快拍数为4096，失配误差Disturb=-10dB。在x向的扫描范围0m~40m，在y向的扫描范围-10m~10m，扫描步长均为0.25m。图2给出多途条件下的常规算法（简称多途常规）与三种基于多途模型匹配类算法（即匹配常规、匹配MVDR和匹配稳健MVDR）的空间谱对比效果图。其中，针对仿真中的失配程度，匹配稳健MVDR算法的约束参数选取为ε=0.3。

保持阵列孔径不变，考虑稀疏阵列的情况，图3给出阵元个数21个，阵元间距0.75m时四种算法分别在声源位置处的x向和y向空间谱切片对比效果图。图4给出阵元个数11个，阵元间距1.5m时的空间谱切片对比效果图。图5给出阵元个数7个，阵元间距2.5m时的空间谱切片对比效果图。图6给出阵元个数5个，阵元间距3.75m时的空间谱切片对比效果图。

由存在空间聚焦导向矢量扰动误差下的空间谱结果可以看出：

（1）多途下的常规算法由于无法与真实的声信道模型向匹配，存在一定的定位偏差，且x向的定位偏差明显大于y向的定位偏差；同时，由于多途效应的存在，加剧了在垂直与基阵方向上（即x向）的背景起伏。

（2）本发明中给出的三种适用于多途条件下的基于多途模型匹配思想的算法，由于本质上利用了基于模型和数据匹配的相干处理思想，充分利用了多途信道信息，生成与实际声传播特性相“匹配”的空间聚焦导向矢量，从而有效克服了多途效应的影响，提高了定位精度。

（3）对比三种基于多途模型匹配思想的算法的空间谱图可知，匹配常规算法的背景起伏剧烈，这对于目标与背景区分十分不利；匹配MVDR算法的谱峰较为尖锐，背景起伏较小，但在存在失配误差的情况下，出现了明显的性能下降，无法体现高分辨MVDR算法的优势；匹配稳健MVDR算法由于利用了约束优化思想，有效改善了高分辨MVDR算法的稳健性，空间谱表现出更大的动态范围、更为尖锐的聚焦峰尺度以及更强的背景起伏压制能力。

（4）多途下的常规算法及匹配常规算法受稀疏阵列的影响较大，容易出现栅瓣；匹配MVDR和匹配稳健MVDR算法在一定程度上压制了栅瓣的影响，较常规算法具有更优的处理效果，这是由高分辨处理器自身的特点决定的。利用高分辨处理器可以在获取较好的处理效果的同时，简化基阵系统的复杂度。

实例二：海洋环境声学参数失配情况下的空间谱

实例参数设置如下：基本参数保持不变，设真实海水深度为15m，海水中声速为1500m/s，海水和海底介质的密度之比m_ρ为1.8，海水与海底中的声速之比n_c为0.85。信噪比为20dB，数据快拍数为4096。以下给出水深存在误差、水中声速存在误差、海水和海底介质的密度之比m_ρ存在误差以及海水与海底中的声速之比n_c存在误差等多种环境声学参数失配情况下的空间谱。图7至图10依次给出四种算法分别在声源位置处的x向和y向空间谱切片对比效果图。

由以上海洋环境声学参数失配情况下的空间谱结果可知：

（1）海洋环境声学参数失配对于匹配MVDR算法的影响较大，此时的空间谱在分辨率及背景抑制能力上均出现明显的退化，这是高分辨算法普遍存在的问题。

（2）由于水深存在误差、水中声速存在误差、海水和海底介质的密度之比m_ρ存在误差以及海水与海底中的声速之比n_c存在误差等情况均可等效视为空间聚焦导向矢量存在误差，因此分析结果呈现的规律与实例一中相同。