CN110837791A - 一种基于过完备字典的声速剖面反演方法 - Google Patents

Abstract

本发明公开了一种基于过完备字典的声速剖面反演方法，1）依据历史数据，基于K‑SVD算法，构建声速剖面的非正交且过完备字典D；2）利用一阶泰勒展开式对正演模型的非线性响应线性化，创建测量矩阵Q；3）在待测海域布放单声源和垂直阵，进行信号发射和接收；4）基于OMP算法，反演过完备字典原子的稀疏系数向量X；5）根据反演得到的时变非正交字典原子系数，计算声速剖面。相比于基于经验正交函数的方法，本发明基于过完备字典，可以通过冗余特性更好的捕捉声速剖面内在的本质特征，极大提高了声速剖面反演的精度。

Description

一种基于过完备字典的声速剖面反演方法

技术领域

本发明涉及一种利用声波的声速剖面反演方法，尤其是一种基于过完备字典的声速剖面反演方法，属于海洋环境参数反演技术领域。

背景技术

声波是在海洋中传播信息的有效载体也是探测海洋环境信息的重要手段。如何利用声信号探测海洋环境信息被人们广泛研究。声波在海洋中的传播与海洋中的物理过程紧密相关，声传播的过程中包含了丰富的海洋温度场、流场分布信息。利用海洋声学技术提取这些参数信息可实现对海洋环境进行长时间、大范围的动态实时监测。

和其他的反演方法相比，采用声波反演海洋环境参数是一种高效、便捷的获取海洋环境的方法。该方法通常只需布放若干个声源和水听器阵列便可以反演出大部分感兴趣的海洋环境信息，且可实现长时间大范围实时监测。而反演得到的海洋环境信息可以提供空间与时间平均意义上的估计,这也是传统直接测量方法难以得到的。

海水声速剖面可以由观测的声压场信号估计得到，但采样信号相对于广袤的海洋来说是有限的。此类欠定问题通常需要正则化才能取得满足物理意义的最优解。将声速剖面建模为已知的参考声速剖面和声速剖面扰动之和，通常声速剖面的扰动在一组基上是稀疏的，目前的方法大都是将声速剖面扰动在正交基上展开，然而在利用经验正交函数对未知剖面进行反演时，尤其是对样本数据之外的剖面外插处理，必须考虑已有的经验正交函数的适用范围。经验正交函数的获取受到数据完备性和数据测量时间的制约，不同季节、不同数据量获得的经验正交函数可能会差异很大，因此在使用经验正交函数时就要考虑经验正交函数的准确性和适用性。

采用超完备的冗余字典代替传统的正交基函数，为信号自适应的稀疏扩展提供了极大的灵活性。稀疏分解可以实现数据压缩的高效性，更重要的是可以利用字典的冗余特性捕捉信号内在的本质特征。使用超完备的冗余字典作为基函数，过完备基不受正交约束，字典的选择可以尽可能地符合被逼近信号的结构，从而可以体现更好的信号压缩性能，也就是说，基于过完备基的声速剖面具有更优异的稀疏性。

压缩感知是一种寻找欠定线性系统稀疏解的技术，突破了传统的奈奎斯特采样定律，用很少的采样点就可以准确重建恢复稀疏信号。海水声速剖面的估计问题可以采用稀疏分解理论进行描述，这是由于声速剖面在正交基或是过完备字典上展开后是满足稀疏条件的，因此可以看成是一个稀疏表示问题。稀疏分解理论考察的是信号的能量分布，分解处理的同时可以实现信号的降噪处理，利用少量数据即可实现高精度估计。由于声速剖面在正交基上展开受到数据完备性和数据测量时间的制约，如果能够使用更具稀疏性的过完备字典反演声速剖面，那么将为声速剖面反演精度带来极大提高。因此，发明一种基于过完备字典的声速剖面反演方法对于海洋测绘领域来说显得尤为必要。

发明内容

本发明的目的在于，为了克服现有技术中存在的不足，提供了一种基于过完备字典的声速剖面反演方法，用于解决现有方法得到的声速剖面精度差的技术问题。

为了解决背景技术所存在的问题，本发明的基于过完备字典的声速剖面反演方法，包括以下步骤如下：

步骤1：依据历史数据，基于K-SVD算法，构建声速剖面的非正交且过完备字典D。

步骤2：利用一阶泰勒展开式对正演模型的非线性响应线性化，创建测量矩阵Q。

步骤3：在待测海域布放单声源和垂直阵，进行信号发射和接收。

步骤4：基于OMP算法，反演过完备字典原子的稀疏系数向量X。

步骤5：根据反演得到的时变非正交字典原子系数，计算声速剖面。

本发明基于过完备字典的声速剖面反演方法，包括如下具体步骤：

步骤1：在待测海域选取M个声速剖面，每个声速剖面在全海深标准层上进行差值，生成K个垂直层，声速剖面训练信号的集合表示为

其中训练信号为剔除平均声速剖面之后的均值为0的矩阵。令

分别代表字典、训练信号、训练信号的稀疏表示稀疏向量，

为C的解向量的集合。从线性组合角度看，KSVD训练算法的目标方程可表示为：

其中，T₀为稀疏表示系数中非零分量数目的上限，即系数向量中的最大差异度。式(1)的求解是一个迭代过程。使用OMP算法可以得到字典D上，C的稀疏表示的系数矩阵X，然后根据系数矩阵X更新字典。字典的更新是逐列进行的，将要更新字典的第k列d_k，系数矩阵相应的第k行，令其为则抽离了

后，样本矩阵与字典逼近的误差为：

得到当前误差矩阵E_k后，调整d_k和

使其乘积与E_k的误差尽可能的小。如果此时就用奇异值分解(SVD)更新的d_k和

SVD能找到距离E_k最近的秩为1的矩阵，这能有效地减小式(2)代表的误差，但是，如此得到的

将不满足稀疏性，所以只要更新

的非0项即可。

步骤2：y(X)为声压信号，是声速为c(X)时通过简正波模型计算得到的，由于声压与声速剖面为非线性关系，在扰动相对于参考声速梯度c₀很小时，利用一阶泰勒展开式将y(X)线性化。

其中，y(0)是声速剖面为平均声速c₀时的声压，这里令

其中Q＝[q₁…q_P]∈R^N×P，通过计算字典原子系数微小变化时的有限差分获得。比如q_i是通过改变第i个原子的系数，从而进行有限差分而获取的：

步骤3：在待测海域布放单声源和垂直阵，进行信号发射和接收；

由N个水听器组成的垂直阵列测量得到声压p_obs。

步骤4：基于OMP算法，反演过完备字典原子的稀疏系数向量X；

通过压缩感知方法，可以利用实测的声压估计声速剖面，具体算法如下：

其中ε是误差容限，T是非零系数的个数，上式可以通过正交匹配追踪算法OMP计算得到稀疏稀疏矩阵X。

正交匹配追踪算法的基本流程如下：

4.1：用x表示信号，初始化残差e₀＝x。

4.2：选择与e₀内积绝对值最大的原子，表示为

4.3：将选择的原子作为列组成矩阵

定义

列空间的正交投影算子为

通过从e₀减去其在所张成空间上的正交投影得到残差e₁e₁＝e₀-Pe₀＝(I-P)e₀。

4.4：对残差迭代执行4.2、4.3步；e_m+1＝e_m-Pe_m＝(I-P)e_m其中I为单位阵。需要注意的是在迭代过程中为所有被选择过的原子组成的矩阵，因此每次都是不同的，所以由它生成的正交投影算子矩阵P每次都是不同的。

4.5：直到达到某个指定的停止准则后停止算法。

步骤5：根据步骤4估计得到的系数矩阵X，利用下面的公式(6)就可以得到声速剖面：

c(X)≈c₀+DX (6)

与现有技术相比，本发明的有益效果是：相比较与利用经验正交函数(EOFs)生成声速剖面的正交基，本方法利用K-SVD算法生成声速剖面的过完备字典，提出基于过完备字典的声速剖面压缩感知方法。由于K-SVD算法不需要强制使用正交条件，对于训练数据更加灵活，从而可以使用更少数的原子组合达到更高的重构精度。过完备字典可以通过冗余特性更好的捕捉声速剖面内在的本质特征，从而能够极大提高声速剖面反演的精度。

附图说明

图1是本发明的工作流程图。

图2是声速剖面所在的位置。

图3是训练集声速剖面的取值。

图4是训练集声速剖面减去均值以后的扰动。

图5是前10阶经验正交原子。

图6是前10阶非正交过完备字典原子。

图7海洋声速剖面反演环境。

图8测试集声速剖面的取值。

图9基于过完备字典(LDs)估计得到的声速剖面与真实值的比较。

图10LDs和EOFs对31条测试集内声速剖面估计误差进行比较。

图11基于过完备字典(LDs)对测试集声速剖面估计的标准偏差随深度的分布。

附图说明：图2中星号为训练集，圆圈为测试集，数字表示的是测试集中样本的编号。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

为证明基于过完备字典的声速剖面反演方法对提高声速剖面反演的精度的有效性，下面选取阿拉伯东北部海域(范围在14-19°N，65-70°E)为实验海域2011年-2017年在9月份的Argo剖面序列计算声速剖面数据，其中2011年-2016年间的M＝224个声速剖面为训练集，2017年的31个声速剖面作为测试集。声速剖面的深度为0m-1000m，图4可以看出声速扰动较大，最大时可以超过10m/s。

步骤1：将实验海域224个声速剖面中的每个声速剖面在全海深标准层上进行差值，生成K个垂直层，声速剖面训练信号的集合表示为

其中训练信号为剔除平均声速剖面之后的均值为0的矩阵；令

分别代表字典、训练信号、训练信号的稀疏表示稀疏向量，

为C的解向量的集合；从线性组合角度看，KSVD训练算法的目标方程可表示为：

公式(1)中T₀为稀疏表示系数中非零分量数目的上限，即稀疏度。

步骤2：利用一阶泰勒展开式对正演模型的非线性响应线性化，利用式(4)计算测量矩阵Q；

y(X)为声压信号，是声速为c(X)时通过简正波模型计算得到的，由于声压与声速剖面为非线性关系，在扰动相对于参考声速梯度c₀很小时，利用一阶泰勒展开式将y(X)线性化。

公式(3)中y(0)是声速剖面为平均声速c₀时的声压，这里令

其中Q＝[q₁…q_P]∈R^N×P，通过计算字典原子系数微小变化时的有限差分获得；比如q_i是通过改变第i个原子的系数，从而进行有限差分而获取的：

其中声源位于海深40m，发射信号频率为300HZ，垂直阵与声源距离为10km，深度为从0m到1000m均匀分布的20个水听器。海洋声速剖面反演环境如图7所示，其中海底声速、密度和衰减系数分别为c_b＝1800m/s，ρ_b＝2.0g/cm³和α_b＝0.1dB/m。

步骤3：在待测海域布放单声源和垂直阵，进行信号发射和接收；由于没有实测数据，这里通过简正波模型生成仿真声场y(X)，求得的声压加上高斯白噪声n得到观测声压信号：

p＝y(X)+n (7)

其中计算仿真声场时的环境如图7所示，声速剖面为测试集中的声速剖面如图8所示，信噪比为30dB。

步骤4：基于OMP算法，反演过完备字典原子的稀疏系数向量X；利用式(5)，可以通过正交匹配追踪算法OMP计算得到稀疏稀疏矩阵X。通过压缩感知方法，可以利用实测的声压估计声速剖面，具体算法如下：

公式(5)中ε是误差容限，T是非零系数的个数；

步骤5：根据反演得到的时变非正交字典原子系数，计算声速剖面。由步骤4得到的系数矩阵X，利用公式(6)就可以得到声速剖面。

c(X)≈c₀+DX (6)

将基于过完备字典(LDs)估计得到的声速剖面与真实值进行比较，从图9图中可以看出，声速剖面的估计值与真实值吻合程度较好，因此说明了本方法的切实可行。

为了将本方法与传统的基于经验正交函数EOFs这两种方法得到的结果进行比较。

首先，将K-SVD与EOF方法分别生成的前10阶原子进行对比，这里过完备字典原子个数N＝K＝224，即两种方法的基函数个数相同。从图5可以看出，EOF方法生成的基函数只有前几阶与声速起伏相关性较大，而从图6可以看出K-SVD生成的基函数与声速起伏都较为相似。

其次，将两种方法对31条测试集内声速剖面估计误差进行比较，平均重构误差由式(8)计算得到：

通过图10可见，对所有的31个声速剖面测试样本，基于过完备字典的海洋声速剖面反演方法在一定程度上降低了声速剖面的反演误差，能显著提高声速剖面反演的精度。其中，第18，19和20号声速剖面的估计误差很小，这主要是由于这三条声速剖面位于训练集中间，从而数据的统计特性比较好；而对第21-31号声速剖面来讲，估计误差稍大，这主要是由于这些声速剖面位于训练集的边缘。

图11给出了基于本方法对测试集声速剖面估计的标准偏差随深度的分布，图中从内到外分别是1STD，2STD和3STD，其中99.7％置信区间(3STD)内，大部分深度处的误差位于1m/s之内，而在水深50m-100m之间的误差较大，这主要是由于此深度处的声速扰动最大引起的。

对于本领域技术人员而言，显然本发明不限于上述示范性实施例的细节，而且在不背离本发明的精神或基本特征的情况下，能够以其他的具体形式实现本发明。因此，无论从哪一点来看，均应将实施例看作是示范性的，而且是非限制性的，本发明的范围由所附权利要求而不是上述说明限定，因此旨在将落在权利要求的等同要件的含义和范围内的所有变化囊括在本发明内。不应将权利要求中的任何附图标记视为限制所涉及的权利要求。

此外，应当理解，虽然本说明书按照实施方式加以描述，但并非每个实施方式仅包含一个独立的技术方案，说明书的这种叙述方式仅仅是为清楚起见，本领域技术人员应当将说明书作为一个整体，各实施例中的技术方案也可以经适当组合，形成本领域技术人员可以理解的其他实施方式。

Claims (3)

1.一种基于过完备字典的声速剖面反演方法，其特征在于，包含以下具体步骤：

步骤1：依据历史数据，基于K-SVD算法，构建声速剖面的非正交且过完备字典D；

在待测海域选取M个声速剖面，每个声速剖面在全海深标准层上进行差值，生成K个垂直层，声速剖面训练信号的集合表示为

其中训练信号为剔除平均声速剖面之后的均值为0的矩阵；令

分别代表字典、训练信号、训练信号的稀疏表示稀疏向量，为C的解向量的集合；从线性组合角度看，KSVD训练算法的目标方程可表示为：

公式(1)中T₀为稀疏表示系数中非零分量数目的上限，即稀疏度；

步骤2：利用一阶泰勒展开式对正演模型的非线性响应线性化，创建测量矩阵Q；

y(X)为声压信号，是声速为c(X)时通过简正波模型计算得到的，由于声压与声速剖面为非线性关系，在扰动相对于参考声速梯度c₀很小时，利用一阶泰勒展开式将y(X)线性化；

公式(3)中y(0)是声速剖面为平均声速c₀时的声压，这里令

其中Q＝[q₁…q_P]∈R^N×P，通过计算字典原子系数微小变化时的有限差分获得；比如q_i是通过改变第i个原子的系数，从而进行有限差分而获取的：

步骤3：在待测海域布放单声源和垂直阵，进行信号发射和接收；

由N个水听器组成的垂直阵列测量得到声压p_obs；

步骤4：基于OMP算法，反演过完备字典原子的稀疏系数向量X；

通过压缩感知方法，可以利用实测的声压估计声速剖面，具体算法如下：

公式(5)中ε是误差容限，T是非零系数的个数；

步骤5：根据反演得到的时变非正交字典原子系数，计算声速剖面；

c(X)≈c₀+DX (6)

根据步骤4估计得到的系数矩阵X，利用公式(6)就可以得到声速剖面。

2.根据权利要求1所述的一种基于过完备字典的声速剖面反演方法，其特征在于，所述的步骤1中公式(1)的求解是一个迭代过程，使用OMP算法可以得到字典D上，C的稀疏表示的系数矩阵X，然后根据系数矩阵X更新字典；字典的更新是逐列进行的，将要更新字典的第k列d_k，系数矩阵相应的第k行，令其为

则抽离了

后，样本矩阵与字典逼近的误差为：

得到当前误差矩阵E_k后，调整d_k和

使其乘积与E_k的误差尽可能的小；如果此时就用奇异值分解(SVD)更新的d_k和

SVD能找到距离E_k最近的秩为1的矩阵，这能有效地减小式(2)代表的误差，但是，如此得到的

将不满足稀疏性，所以只要更新

的非0项即可。

3.根据权利要求1所述的一种基于过完备字典的声速剖面反演方法，其特征在于，所述的步骤4中公式(5)可以通过正交匹配追踪算法OMP计算得到稀疏稀疏矩阵X；正交匹配追踪算法的基本流程如下：

4.1：用x表示信号，初始化残差e₀＝x；

4.2：选择与e₀内积绝对值最大的原子，表示为

4.3：将选择的原子作为列组成矩阵

定义列空间的正交投影算子为

通过从e₀减去其在

所张成空间上的正交投影得到残差e₁e₁＝e₀-Pe₀＝(I-P)e₀；

4.4：对残差迭代执行4.2、4.3步；e_m+1＝e_m-Pe_m＝(I-P)e_m其中I为单位阵；需要注意的是在迭代过程中为所有被选择过的原子组成的矩阵，因此每次都是不同的，所以由它生成的正交投影算子矩阵P每次都是不同的；

4.5：直到达到某个指定的停止准则后停止算法。

Images