CN109427042B - 一种提取局部海域沉积层的分层结构和空间分布的方法 - Google Patents

Abstract

本发明提供一种提取局部海域沉积层的分层结构和空间分布的方法，包括：步骤1、数据I＝Δr×Δd作为待处理的浅剖图像；步骤2、利用二维多尺度线条滤波器对待处理的浅剖图像中的线条信息进行多尺度线条滤波，获得最大响应图像M；步骤3、对最大响应图像M在不同方向进行分解，将分解得到的低频子图像进行阈值滤波，将分解得到的水平细节子图像进行阈值增强，利用逆二维小波变换重建回原始图像，获得水平线条对比度的图像；步骤4、对步骤3中的水平线条的区域信息的图像的沉积层的界面结构和背景进行二值化处理，获得水平线条的区域边界；步骤5、水平采样距离变为r₀+Δr，重复步骤1到步骤4；步骤6、绘制出局部海域沉积层空间分层结构图。

Description

一种提取局部海域沉积层的分层结构和空间分布的方法

技术领域

本发明涉及海洋测量、海底遥感的技术领域，特别涉及一种提取局部海域沉积层的分层结构和空间分布的方法。

背景技术

获取准确的海洋环境信息是海洋资源开发以及国防建设的基础，而局部海域的海底沉积层空间分层结构和空间分布是研究海底地质变迁和声纳设备应用的重要条件。传统的提取海底沉积层分层结构和空间分布的方法是采用柱状原状取样设备，利用机械压力或者重力定点取样，通过实验室对底质进行分析得到沉积层参数信息，进而得到沉积层的分层结构信息。但是，这种提取方法的时间成本和经济成本昂贵，而且相对于广阔的海洋，取样范围极其有限，不能获得广阔海域的沉积层空间分层结构信息。

因而，根据不同沉积物分层界面处的阻抗差异，应用声学遥感的方法可以快速、准确的获得海底地形和沉积层等效的结构信息。但是，利用浅地层剖面仪，即浅剖进行海底探测时，通常会受到海面噪声、海底混响和界面多次反射波的干扰；其中，海面噪声主要包括船舶自噪声和尾流噪声；造成获得的海底沉积层反射回波数据的信噪比较低、干扰严重。这种结构模糊的数据通常需要具有相关专业知识的研究人员对其进行人工的识别和解释，其识别效率低而且由于探测海域广阔，获得的数据量大，只能给出某一测线的沉积层信息，不能实时给出局部海域的分层结构信息，不利于海洋设备自动化的发展。

发明内容

本发明的目的在于，为了解决现有的提取方法存在上述缺陷，提供了一种提取局部海域沉积层的分层结构和空间分布的方法。首先，根据海水深度对浅剖数据进行实时自适应分块图像处理；然后，综合多尺度线条滤波、二维小波分解、自适应阈值增强和滤波以及二维小波重构的多种图像处理方法对浅剖图像进行处理。该方法可以有效去除海面舰船噪声和尾流的干扰，避免海底海面多次反射回波的影响，在海底混响干扰的数据中，同步连续获取局部海域准确的沉积层空间分层界面结构，进而可以自动获得沉积层的分层数目以及各层的位置和厚度，解决了在多系统联合海底测量中，浅剖数据需要人工处理和准确性低的问题。本发明可以及时给出局部海域准确的海底空间分层模型，为实时局部海域的地声参数反演，岸基声纳设备的应用，海底声学遥感提供及时准确的海底分层结构数据，同时局部海域海底沉积层空间分层结构也为水下考古和海底地质变迁提供了分析依据。

为了实现上述目的，本发明提出了针对浅剖数据的一种提取局部海域沉积层的分层结构和空间分布的方法，主要是根据海水深度的不同，去除靠近海面的船舶自噪声和尾流噪声以及海底海面多次反射的干扰；同时自适应的截取待处理图像；对截取的待处理图像进行多尺度线条滤波处理，通过选择合适的空间尺度范围，一般尺度范围取值为2到16，抑制随机高斯噪声，获得对沉积层分层结构的最佳响应；利用二维小波变换将最佳响应图像分解为水平、垂直、对角方向细节子图和低频子图，通过对含有沉积层分层结构信息的水平方向的细节子图像进行阈值增强，和对含有界面厚度信息的低频子图进行阈值滤波，随后将处理后的图像进行重建；则重建后的图像抑制了垂直方向的干扰，增强了较深位置的分层界面；最后在适当的阈值(灰度值100)下，将界面和背景进行二值化处理，可以连续自动提取出沉积层的分层数目以及各层的厚度和位置等结构参数信息。

所述提取局部海域沉积层的分层结构和空间分布的方法具体包括：

步骤1)给定准确处理的海底沉积层深度h₁和海底以上待处理深度h₂；在水平采样距离r₀处开始，给定步距Δr根据海水深度的不同，自适应的截取采样深度Δd，数据I＝Δr×Δd作为待处理的浅剖图像；

步骤2)选定需要提取边界特征的尺度范围，在所述尺度范围内，在不同尺度s下，将二维多尺度线条滤波器对步骤1)中所述待处理的浅剖图像中的线条信息进行多尺度线条滤波，提取所述待处理的浅剖图像中的能量对比度的结构特征，小尺度的随机高斯噪声将会被滤除掉，获得二维多尺度线条滤波器对所述待处理的浅剖图像的分层界面结构的最大响应图像M；

步骤3)采用双正交小波滤波器组的二维小波变换，对步骤2中的最大响应图像M在不同方向进行分解，将分解得到的含有界面厚度信息的低频子图像进行阈值滤波，将分解得到的含有沉积层分层结构信息的水平细节子图像进行阈值增强，然后将处理完的低频子图像和水平细节子图像利用逆二维小波变换重建回原始图像M_c，最终获得包含完整的水平线条的区域信息、滤除垂直干扰和增强水平线条对比度的图像；重建后的原始图像避免了出现将单条边界错误分解为两条的情况，获得更准确的分层边界的区域范围；而且也能有效增强对比度较低的水平分层边界结构。

步骤4)对步骤3)中的最后获得的完整的水平线条的区域信息的图像的沉积层的界面结构和背景进行二值化处理，对二值化图像逐个采样点通过列向量错位相减，获得水平线条的区域边界；最终获得海底沉积层的分层数、各层的深度和各层的厚度信息。

步骤5)水平采样距离变为r₀+Δr，重复步骤1)到步骤4)；获得不同采样点的海底沉积层的层数和各层的厚度信息；

步骤6)结合逐个采样点处的海底沉积层的层数、各层的深度信息、各层的厚度信息、每个采样点处的位置信息，即GPS信息，绘制出局部海域沉积层空间分层结构图。

在步骤1)中，在水平采样距离r₀处开始，给定步距Δr根据海水深度的不同，自适应的截取采样深度Δd，数据I＝Δr×Δd作为待处理的浅剖图像的具体步骤包括：

步骤1-1)给定准确处理的海底沉积层深度h₁和海底以上待处理深度h₂；由回波信号确定水平采样距离r₀处最大反射回波强度的深度作为海底位置h₃；

步骤1-2)如果h₁大于h₃，那么选择海底以上h₂深度和海底以下h₁深度范围作为待处理截取采样深度Δd；如果h₁小于h₃，那么选择海底以上h₂深度和海底以下h₃深度作为待处理截取深度Δd；

步骤1-3)输出待处理的分块浅剖图像I＝Δr×Δd，并将其作为待处理的浅剖图像。

在步骤2)中，获得二维多尺度线条滤波器对所述的待处理浅剖图像的分层界面结构的最大响应图像M的具体过程如下：

步骤2-1)将步骤1)中的所述待处理的浅剖图像表示为一个二维的输入矩阵L，采用高斯核函数，通过公式(1)，获得输入矩阵在尺度s下的海森矩阵H，即Hessian矩阵H；

其中，

其中，“*”表示卷积运算；

步骤2-2)通过公式(2)和(3)，获得Hessian矩阵H的本征值λ₁和λ₂；

再将本征值λ₁和λ₂带入公式(4)和(5)中，获得二维多尺度线条滤波器的尺度参数R和S；

其中，λ₁，λ₂是海森矩阵(Hessian)的本征值，而且λ₁>λ₂；A是用来描述结构信息的最大的椭圆截面积，l是椭圆长轴距离。R是海森矩阵H的本征值之比，描述截取后的浅剖图像局部特征的椭圆的长轴和短轴之比，S是海森矩阵H的本征值平方和的平方根，描述所述截取后的浅剖图像结构信息和所述截取后的浅剖图像背景噪声的比值。

步骤2-3)通过公式(6)和计算Hessian矩阵，根据在不同尺度s下的本征值来匹配图像结构的响应M(s),

其中，α是用来调节二维多尺度线条滤波器对R的敏感度，β是用来调节二维多尺度线条滤波器对S的敏感度；

步骤2-4)再通过公式(7)，获得二维多尺度线条滤波器对所述待处理的浅剖图像的分层界面结构的最大响应图像M

其中，s_min是最小的二维多尺度线条滤波器的尺度。s_max是最大的二维多尺度线条滤波器的尺度；

步骤2-5)选择尺度s范围，并重复步骤2-1)至2-5)，将不同尺度s下的二维多尺度线条滤波器响应的最大值赋值给M，即为提取得到的对沉积层分层界面的最佳匹配特征的响应，输出二维多尺度线条滤波器对沉积层分层界面结构的最大响应图像M。

在步骤3)中，获得包含完整的水平线条的区域信息、滤除垂直干扰和增强水平线条对比度的图像的具体步骤包括：

步骤3-1)对步骤2)中输出的最大响应图像M进行二维的小波分解处理，即在每一行进行一维小波变换，获得最大响应图像M的高频分量和低频分量，然后在每一列进行插值；然后，将每一行的一维高频分量图像和一维低频分量图像，在每一列进行一维小波变换，再对每一行进行插值；

步骤3-2)利用双正交小波滤波器组的二维小波变换，可以将步骤2)中输出的最大响应图像M分解为平滑的低频子图像和互相独立的方向细节子图，即水平细节子图、垂直细节子图、对角细节子图；其中，所述平滑的低频子图像中包含了水平线条的宽度信息，即界面的厚度信息；互相独立的方向细节子图像，即水平细节子图、垂直细节子图、对角细节子图中含有图像的边界信息。而水平线条的结构在水平方向的细节子图像M_h中保留了下来；所述水平细节子图像包含有沉积层水平分层结构信息；

步骤3-3)对所述水平细节子图像进行阈值增强处理，阈值T₁＝c₁×P₁；其中，c₁为阈值系数；P₁为水平细节子图像的最大灰度值；c₁大于0或小于1；逐点将水平细节子图中灰度值大于T₁的点设置为P₁；

步骤3-4)对所述含有界面厚度信息的低频子图像进行阈值滤波处理，阈值T₂＝c₂×P₂，其中，c₂为阈值系数；P₂为低频子图像的最大灰度值，c₂大于0小于1；逐点将低频子图像中灰度值小于T₂的点设置为0；

步骤3-5)将通过步骤3-3)的阈值增强处理后的水平细节子图像和步骤3-4)的阈值滤波处理后的低频子图像进行二维小波重建，得到滤除垂直干扰，增强界面对比度的二维小波重建的水平分层界面图像M_c；

步骤3-6)将得到的二维小波重建的水平分层界面图像M_c在垂直方向上作归一化处理；

步骤3-7)归一化处理后的二维小波重建的水平分层界面图像，即M_cn，将其重新输入到二维多尺度线条滤波器中，并重复步骤2)，获得光滑的、连续的水平线条结构的区域信息的图像。

在步骤4)中，获得水平线条的区域边界、海底沉积层的分层数、各层的深度和各层的厚度信息的具体步骤包括：

步骤4-1)在步骤3)输出的水平线条结构的区域信息的图像中，水平分层结构和浅剖背景已经有了明显的灰度对比，提取的结构和背景响应相差很大，那么将步骤3)输出的水平线条的结构的区域信息的图像表示为统计直方图的形式，选定灰度阈值T_h，将步骤3)输出的水平线条的结构的区域信息的图像表示为二值图像；

步骤4-2)在每一个水平距离r₀上，依据公式(8)，可以通过列向量错位相减，得到水平线条区域的边界列向量：

b＝|a(1:N-1)-a(2:N)| (8)

其中，a是某一距离处深度方向向量，即列向量，N是深度方向上的最大数据量，即垂直方向像素的最大个数；边界列向量b中的非零值所在采样点的位置代表水平线条边界的位置；

步骤4-3)通过公式(9)、(10)、(11)，输出沉积层在每一个距离上的分层数、各层的深度以及各层的厚度；

沉积层在每一个水平距离r₀上的分层数目n为：

沉积层中第i层分界面的深度d(i)为：

沉积层中第i层分界面的的厚度t(i)为：

t(i)＝d(i+1)-d(1) (11)

式中，d(1)为海底深度。

在步骤6)中，将获得水平线条的区域边界、海底沉积层的分层数、各层的深度和各层的厚度信息进行区域绘图的具体步骤包括：

步骤6-1)结合逐个采样点处的海底沉积层的层数、各层的深度信息、各层的厚度信息、每个采样点处的位置信息，即GPS信息，绘制出初步的局部海域沉积层空间分层结构图。

步骤6-2)将每条测线的海底沉积层的层数信息进行离散化处理，然后将局部数据根据临近相等的原则在线性空间中插值，绘制出最终的局部海域沉积层空间分层结构图。

本发明的优点在于：本发明的方法可以自动提取局部海域的沉积层分层结构，包括海底地形、海底沉积层分层数以及各层沉积层的厚度，同步提取，自动截取待处理的浅剖图像，提取速度快，可以通过改变滤波器参数，环境适应能力强，可以有效应用到大量浅剖数据的快速处理和分析，海洋海底遥感，局部海域多系统综合海底环境测量及反演等领域。另外，本发明的方法有效分段截取和处理浅剖数据，在一定尺度范围内对沉积层分层结构特征进行提取和分析，解决了大量浅剖数据的自动处理问题；本发明针对以往利用二维小波变换提取水平分层结构的缺点，对重建后的子图像进行处理后利用逆二维小波进行重建，重建后的图像避免了出现将单条边界错误分解为两条的情况，而且也能有效增强对比度较低的水平分层边界结构。本发明的方法可以应用与海底沉积层分层结构的快速获取，也可给出局部海域沉积层分层结构信息，为海底考古，声纳设备应用等提供直观、准确的海底结构信息。

附图说明

图1是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法的流程图；

图2是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法在步骤3中的二维小波分解的处理流程图；

图3(1)是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法的步骤1的输出自适应截取的待处理的浅剖图像；

图3(2)是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法的步骤2的输出多尺度线条滤波器的最佳响应图像；

图4(1)是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法的步骤3的输出利用逆二维小波变换将处理后的低频子图像和水平方向的细节子图像进行二维小波重建后的浅剖图像，含有单层沉积层；

图4(2)是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法的步骤4的输出的有单层沉积层的浅剖图像；

图5是本发明的一种提取局部海域沉积层的分层结构和空间分布的方法的步骤6的输出的局部海域沉积层空间分层结构数据散点图像；

图6(1)是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的沉积层水平分层数区域分布图；

图6(2)是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的第一层沉积层(不存在第一层沉积层记厚度为0)厚度的分布图；

图6(3)是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的第二层沉积层(不存在第二层沉积层记厚度为0)厚度的分布图。

具体实施方式

以下结合附图对本发明作进一步的详细说明。

本发明的方法给出区域浅剖数据处理的技术流程，和具体的获得局部海域沉积层分层信息的实施步骤。另外，本发明的方法增加了对最后区域浅剖数据的综合处理，给出了两种局部海域沉积层分层结构信息的图像表示方法。此外，本发明的方法在利用二维小波变换将多尺度滤波处理后的浅剖图像分解后，增加了在低频细节子图像中的阈值滤波和水平细节子图像中的阈值增强，之后利用逆二维小波变换将阈值滤波处理后的低频细节子图像和阈值增强处理后的水平细节子图像重建回浅剖图像。重建后的图像避免了出现将单条边界错误分解为两条的情况，获得更准确的分层边界的区域范围；而且也能有效增强对比度较低的水平分层边界结构。

如图1所示，本发明提供了一种提取局部海域沉积层的分层结构和空间分布的方法，其具体步骤包括：

步骤1)设置浅剖输出文件路径和大小；给定准确处理的海底沉积层深度h₁＝30m和海底以上待处理深度h₂＝10m；在水平采样距离r₀处开始，给定步距Δr＝600(Ping)根据海水深度的不同，自适应的截取采样深度Δd，数据I＝Δr×Δd作为待处理的浅剖图像，如图3(1)所示；读取并保存水平采样点对应的位置数据即GPS数据信息；

步骤2)选定需要提取边界特征的尺度范围，在所述尺度范围内，尺度范围取值为2到16，在不同尺度s下，将二维多尺度线条滤波器对步骤1)中所述待处理的浅剖图像中的线条信息进行多尺度线条滤波，提取所述待处理的浅剖图像中的能量对比度的结构特征，小尺度s的随机高斯噪声将会被滤除掉，获得二维多尺度线条滤波器对所述待处理的浅剖图像的分层界面结构的最大响应图像M，如图3(2)所示；

步骤3)采用双正交小波滤波器组的二维小波变换，对步骤2)中的最大响应图像M在不同方向进行分解，分解步骤如图2所示，将分解得到的含有界面厚度信息的低频子图像进行阈值滤波；将分解得到的含有水平分层信息的水平细节子图像进行阈值增强，然后将处理过的低频子图像和水平细节子图像利用逆二维小波变换重建回原始图像，最终获得包含完整的水平线条的区域信息、滤除垂直干扰和增强水平线条对比度的图像，如图4(1)所示；

步骤4)对步骤3)中的最后获得的完整的水平线条的区域信息的图像的界面结构和背景进行二值化处理，对二值化图像逐个采样点通过列向量错位相减，获得水平线条的区域边界，如图4(2)所示；最终获得海底沉积层的层数、各层的深度和各层的厚度信息；

步骤5)水平采样距离变为r₀+Δr，重复步骤1)到步骤4)；

步骤6)结合逐个采样点处的海底沉积层的层数、各层的深度信息、各层的厚度信息、每个采样点处的位置信息，即GPS信息，绘制出局部海域沉积层空间分层结构图。

在步骤1)中，在水平采样距离r₀处开始，给定步距Δr根据海水深度的不同，自适应的截取采样深度Δd，数据I＝Δr×Δd作为待处理的浅剖图像的具体步骤包括：

步骤1-1)给定准确处理的海底沉积层深度h₁＝30m和海底以上待处理深度h₂＝10m；由回波信号确定水平采样距离r₀处最大反射回波强度的深度作为海底位置h₃；

步骤1-2)如果h₁大于h₃，那么选择海底以上h₂深度和海底以下h₁深度范围作为待处理采样深度Δd；如果h₁小于h₃，那么选择海底以上h₂深度和海底以下h₃深度作为待处理采样深度Δd；

步骤1-3)输出待处理的分块浅剖图像I＝Δr×Δd，并将其作为待处理浅剖图像；

在步骤2)中，获得二维多尺度线条滤波器对所述待处理的浅剖图像的分层界面结构的最大响应图像M的具体过程如下：

步骤2-1)将步骤1中的所述待处理的浅剖图像表示为一个二维的输入矩阵L，采用高斯核函数，通过公式(1)，获得二维的输入矩阵L在尺度s下的海森矩阵H，即Hessian矩阵H；

其中，

其中，“*”表示卷积运算；

步骤2-2)通过公式(2)和(3)，获得Hessian矩阵H的本征值λ₁和λ₂；

再将本征值λ₁和λ₂带入公式(4)和(5)中，获得二维多尺度线条滤波器的尺度参数R和S；

其中，λ₁，λ₂是海森矩阵(Hessian)的本征值，而且λ₁>λ₂；A是用来描述结构信息的最大的椭圆截面积，l是椭圆长轴距离。R是海森矩阵H的本征值之比，描述截取后的浅剖图像局部特征的椭圆的长轴和短轴之比，S是海森矩阵H的本征值平方和的平方根，描述所述截取后的浅剖图像结构信息和所述截取后的浅剖图像背景噪声的比值；

步骤2-3)通过公式(6)和计算Hessian矩阵，根据在不同尺度s下的本征值来匹配图像结构的响应M(s)，

其中，α是用来调节二维多尺度线条滤波器对R的敏感度，β是用来调节二维多尺度线条滤波器对S的敏感度；

步骤2-4)再通过公式(7)，获得二维多尺度线条滤波器对沉积层的分层界面结构的最大响应图像M：

其中，s_min是最小的二维多尺度线条滤波器的尺度；s_max是最大的二维多尺度线条滤波器的尺度；

步骤2-5)选择尺度s范围，并重复步骤2-1)至2-4)，将不同尺度s下的二维多尺度线条滤波器响应的最大值赋值给M，即为提取得到的对沉积层分层界面的最佳匹配特征的响应，输出二维多尺度线条滤波器对沉积层分层界面结构的最大响应图像M。

在步骤3)中，获得包含完整的水平线条的区域信息、滤除垂直干扰和增强水平线条对比度的图像的具体步骤包括：

步骤3-1)对步骤2)中输出的最大响应图像M进行二维的小波分解处理，即在每一行进行一维小波变换，获得所述最大响应图像M的高频分量和低频分量，然后在每一列进行插值；然后，将每一行的一维高频分量图像和一维低频分量图像，在每一列进行一维小波变换，再对每一行进行插值；

步骤3-2)利用双正交小波滤波器组的二维小波变换，可以将步骤2)中输出的最大响应图像M分解为平滑的低频子图像和互相独立的方向细节子图，即水平细节子图、垂直细节子图、对角细节子图；其中，所述平滑的低频子图像中包含了水平线条的宽度信息，即界面的厚度信息；互相独立的方向细节子图像，即水平细节子图、垂直细节子图、对角细节子图中含有图像的边界信息。而水平线条的结构在水平方向的细节子图像M_h中保留了下来；所述水平细节子图像包含有沉积层水平分层结构信息；

步骤3-3)对水平细节子图像进行阈值增强处理，阈值T₁＝c₁×P₁，其中，c₁为阈值系数；P₁为水平细节子图像的最大灰度值；c₁大于0或小于1；逐点将水平细节子图像中的灰度值大于T₁的点设置为P₁；

步骤3-4)对低频子图像进行阈值滤波处理，阈值T₂＝c₂×P₂，其中，c₂为阈值系数；P₂为低频子图像的最大灰度值，c₂大于0或小于1；由于低频子图像中的低灰度值中含有其它方向的干扰，因而需要阈值滤波，因此，逐点将低频子图像中的灰度值小于T₂的点设置为0；

步骤3-5)将阈值增强处理后的水平细节子图像和阈值滤波处理后的低频子图像进行二维小波重建，得到滤除垂直干扰，增强界面对比度的二维小波重建的水平分层界面图像M_c；

步骤3-6)将得到的二维小波重建的水平分层界面图像M_c在垂直方向上作归一化处理；

步骤3-7)把归一化处理后的水平分层界面图像M_cn重新输入到二维多尺度线条滤波器中，并重复步骤2)，获得光滑的、连续的水平线条结构的区域信息的图像。

在步骤4)中，获得水平线条的区域边界、海底沉积层的分层数、各层的深度和各层的厚度信息的具体过程如下：

步骤4-1)在步骤3)输出的水平线条的结构的区域信息的图像中，水平分层结构和浅剖背景有明显的灰度对比，提取的结构和背景有不同的响应，那么将步骤3)输出的水平线条的结构的区域信息的图像表示为统计直方图的形式，选定灰度阈值T_h，将步骤3)输出的水平线条的结构的区域信息的图像表示为二值图像；

步骤4-2)在每一个水平距离r₀上，依据公式(9)，可以通过列向量错位相减，得到水平线条区域的边界：

b＝|a(1:N-1)-a(2:N)| (9)

其中，a是某一距离处深度方向向量，即列向量，N是深度方向上的最大数据量，即垂直方向像素的最大个数；b中的非零值所在采样点的位置代表水平线条边界的位置；

步骤4-3)通过公式(10)、(11)、(12)，输出沉积层在每一个距离上的分层数以及各层的厚度；

沉积层在每一个水平距离r₀上的分层数目n为：

沉积层中第i层分界面的深度d(i)为：

沉积层中第i层分界面的的厚度t(i)为：

t(i)＝d(i+1)-d(1) (12)

其中，d(1)为海底深度。

在步骤6)中，将获得水平线条的区域边界、海底沉积层的层数和各层的厚度信息进行区域绘图的步骤包括：

步骤6-1)结合逐个采样点处的海底沉积层的层数、各层的深度信息、各层的厚度信息、每个采样点处的位置信息，即GPS信息，如图5所示，绘制出初步的局部海域沉积层空间分层结构图。

步骤6-2)将每条测线的海底沉积层的层数信息进行离散化处理，然后将局部数据根据临近相等的原则在线性空间中插值，如图6(1)所示，绘制出最终的局部海域沉积层空间分层结构图；

步骤6-3)将每条测线的第一层和第二层海底沉积层厚度信息进行均值平滑处理，每500个采样点取平均值，然后将局部数据根据临近相等的原则在线性空间中插值，如图6(2)和图6(3)所示，绘制出局部海域的第一层沉积层厚度分布图和局部海域的第二层沉积层厚度分布图。

图3(1)是浅剖输出的图像，有单层沉积层，沉积层厚度约为8m，图中从上到下第一条水平线条是海底表面深度约为54m，第二条水平线条是沉积层和海底底质的分界面由于阻抗差异产生的能量回波。

图3(2)是对图3(1)进行多尺度线条滤波的最佳响应图像。图3(2)中可以看出，在非线条结构区域，随机的噪声被抑制；在水平线条区域，水平线条的基本结构被提取出来。

图4(1)是对图3(2)进行二维小波分解和重建的结果。图4(1)中可以看出，位置较深的水平线条的能量对比度得到了增强。

图4(2)是对图4(1)进行再次的多尺度线条滤波的最佳响应图像。图4(2)中可以看出，水平线条区域已经被完全提取出来，而且能量对比度不强的深层沉积层分层界面也得到了加强。因而，仅需要简单的数学运算就可以将沉积层的层数及对应的厚度提取出来，达到了本发明所要实现的效果。

图5是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的沉积层水平分层数分布散点图。从图5中可以看出该区域大部分海底没有明显的沉积层分层结构。

图6(1)是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的沉积层水平分层数区域分布图。从图6(1)中可以清晰的看出该局部海域中沉积层分层结构的分布情况。

图6(2)是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的第一层沉积层(不存在第一层沉积层记厚度为0)厚度的分布图。从图6(2)中可以得到该局部海域中第一层沉积层的厚度范围以及分布情况。

图6(3)是综合本发明的一种提取局部海域沉积层的分层结构和空间分布的方法得到的局部海域的第二层沉积层(不存在第二层沉积层记厚度为0)厚度的分布图。从图6(2)中可以得到该局部海域中第一层沉积层的厚度范围以及分布情况。

最后所应说明的是，以上实施例仅用以说明本发明的技术方案而非限制。尽管参照实施例对本发明进行了详细说明，本领域的普通技术人员应当理解，对本发明的技术方案进行修改或者等同替换，都不脱离本发明技术方案的精神和范围，其均应涵盖在本发明的权利要求范围当中。

Claims (5)

1.一种提取局部海域沉积层的分层结构和空间分布的方法，包括如下步骤：

步骤1)在水平采样距离r₀处开始，给定步距Δr根据海水深度的不同，自适应的截取采样深度Δd，将数据I＝Δr×Δd作为待处理的浅剖图像；

步骤2)利用二维多尺度线条滤波器对步骤1)中所述待处理的浅剖图像中的线条信息进行多尺度线条滤波，提取所述待处理的浅剖图像中的能量对比度的结构特征，获得所述待处理的浅剖图像的分层界面结构的最大响应图像M；

在步骤2)的具体过程如下：

步骤2-1)将步骤1)中的所述待处理的浅剖图像表示为一个二维的输入矩阵L，采用高斯核函数，通过公式(1)，获得输入矩阵在尺度s下的海森矩阵H；

其中，

其中，“*”表示卷积运算；

步骤2-2)通过公式(2)和(3)，获得海森矩阵H的本征值λ₁和λ₂；

再将本征值λ₁和λ₂带入公式(4)和(5)中，获得二维多尺度线条滤波器的尺度参数R和S；

其中，λ₁，λ₂是海森矩阵H的本征值，而且λ₁＞λ₂；A是用来描述结构信息的最大的椭圆截面积，l是椭圆长轴距离；R是海森矩阵H的本征值之比，描述截取后的浅剖图像局部特征的椭圆的长轴和短轴之比，S是海森矩阵H的本征值平方和的平方根，描述所述截取后的浅剖图像结构信息和所述截取后的浅剖图像背景噪声的比值；

步骤2-3)通过公式(6)和计算海森矩阵H，根据在不同尺度s下的本征值来匹配图像结构的响应M(s),

其中，α是用来调节二维多尺度线条滤波器对R的敏感度，β是用来调节二维多尺度线条滤波器对S的敏感度；

步骤2-4)再通过公式(7)，获得二维多尺度线条滤波器对所述待处理的浅剖图像的分层界面结构的最大响应图像M：

其中，s_min是最小的二维多尺度线条滤波器的尺度；s_max是最大的二维多尺度线条滤波器的尺度；

步骤2-5)选择尺度s范围，并重复步骤2-1)至2-5)，将不同尺度s下的二维多尺度线条滤波器响应的最大值赋值给M，提取得到的对沉积层分层界面的最佳匹配特征的响应，输出二维多尺度线条滤波器对沉积层分层界面结构的最大响应图像M；

步骤3)采用双正交小波滤波器组的二维小波变换，对步骤2)中的最大响应图像M在不同方向进行分解，将分解得到的低频子图像进行阈值滤波，将分解得到的水平细节子图像进行阈值增强，然后将处理完的低频子图像和水平细节子图像利用逆二维小波变换重建回原始图像，获得包含完整的水平线条的区域信息、滤除垂直干扰和增强水平线条对比度的图像；

步骤4)对步骤3)中的最后获得的完整的水平线条的区域信息的图像的沉积层的界面结构和背景进行二值化处理，对二值化图像逐个采样点通过列向量错位相减，获得水平线条的区域边界；最终获得海底沉积层的分层数、各层的深度和各层的厚度信息；

步骤5)水平采样距离变为r₀+Δr，重复步骤1)到步骤4)；获得不同采样点的海底沉积层的层数和各层的厚度信息；

步骤6)结合逐个采样点处的海底沉积层的层数、各层的深度信息、各层的厚度信息、每个采样点处的位置信息，绘制出局部海域沉积层空间分层结构图。

2.根据权利要求1所述的提取局部海域沉积层的分层结构和空间分布的方法，其特征在于，所述的步骤1)的具体步骤包括：

步骤1-1)给定准确处理的海底沉积层深度h₁和海底以上待处理深度h₂；由回波信号确定水平采样距离r₀处最大反射回波强度的深度作为海底位置h₃；

步骤1-2)如果h₁大于h₃，那么选择海底以上h₂深度和海底以下h₁深度范围作为待处理截取采样深度Δd；如果h₁小于h₃，那么选择海底以上h₂深度和海底以下h₃深度作为待处理截取深度Δd；

步骤1-3)输出待处理的分块浅剖图像I＝Δr×Δd，并将其作为待处理的浅剖图像。

3.根据权利要求1所述的一种提取局部海域沉积层的分层结构和空间分布的方法，其特征在于，在步骤3)的具体步骤包括：

步骤3-1)对步骤2)中输出的最大响应图像M进行二维的小波分解处理，即在每一行进行一维小波变换，获得最大响应图像M的高频分量和低频分量，然后在每一列进行插值；然后，将每一行的一维高频分量图像和一维低频分量图像，在每一列进行一维小波变换，再对每一行进行插值；

步骤3-2)利用双正交小波滤波器组的二维小波变换，可以将步骤2)中输出的最大响应图像M分解为平滑的低频子图像和互相独立的方向细节子图，即水平细节子图、垂直细节子图、对角细节子图；其中，所述平滑的低频子图像中包含了水平线条的界面的厚度信息；互相独立的方向细节子图像，包括：水平细节子图、垂直细节子图、对角细节子图中含有图像的边界信息；而水平线条的结构在水平方向的细节子图像M_h中保留了下来；所述水平细节子图像包含有沉积层水平分层结构信息；

步骤3-3)对所述水平细节子图像进行阈值增强处理，阈值T₁＝c₁×P₁；其中，c₁为阈值系数；P₁为水平细节子图像的最大灰度值；c₁大于0或小于1；逐点将水平细节子图中灰度值大于T₁的点设置为P₁；

步骤3-4)对所述含有界面厚度信息的低频子图像进行阈值滤波处理，阈值T₂＝c₂×P₂，其中，c₂为阈值系数；P₂为低频子图像的最大灰度值，c₂大于0小于1；逐点将低频子图像中灰度值小于T₂的点设置为0；

步骤3-5)将通过步骤3-3)的阈值增强处理后的水平细节子图像和步骤3-4)的阈值滤波处理后的低频子图像进行二维小波重建，得到滤除垂直干扰，增强界面对比度的二维小波重建的水平分层界面图像M_c；

步骤3-6)将得到的二维小波重建的水平分层界面图像M_c在垂直方向上作归一化处理；

步骤3-7)将归一化处理后的二维小波重建的水平分层界面图像M_cn，将其重新输入到二维多尺度线条滤波器中，并重复步骤2)，获得光滑的、连续的水平线条结构的区域信息的图像。

4.根据权利要求1所述的一种提取局部海域沉积层的分层结构和空间分布的方法，其特征在于，在步骤4)的具体步骤包括：

步骤4-1)将步骤3)输出的水平线条的结构的区域信息的图像表示为统计直方图的形式，选定灰度阈值T_h，将步骤3)输出的水平线条的结构的区域信息的图像表示为二值图像；

步骤4-2)在每一个水平距离r₀上，依据公式(8)，通过列向量错位相减，得到水平线条区域的边界列向量：

b＝|a(1：N-1)-a(2：N)| (8)

其中，a是某一距离处深度方向向量，记为列向量，N是深度方向上的最大数据量，记为垂直方向像素的最大个数；边界列向量b中的非零值所在采样点的位置代表水平线条边界的位置；

步骤4-3)通过公式(9)、(10)、(11)，输出沉积层在每一个距离上的分层数、各层的深度以及各层的厚度；

沉积层在每一个水平距离r₀上的分层数目n为：

沉积层中第i层分界面的深度d(i)为：

沉积层中第i层分界面的的厚度t(i)为：

t(i)＝d(i+1)-d(1) (11)

式中，d(1)为海底深度。

5.根据权利要求1所述的一种提取局部海域沉积层的分层结构和空间分布的方法，其特征在于，在步骤6)的具体步骤包括：

步骤6-1)结合逐个采样点处的海底沉积层的层数、各层的深度信息、各层的厚度信息、每个采样点处的位置信息，绘制出初步的局部海域沉积层空间分层结构图；

步骤6-2)将每条测线的海底沉积层的层数信息进行离散化处理，然后将局部数据根据临近相等的原则在线性空间中插值，绘制出最终的局部海域沉积层空间分层结构图。

Images