CN111562620A - 一种探地雷达层位自动追踪方法及系统 - Google Patents

Abstract

本发明公开了一种探地雷达层位自动追踪方法，包括：确定前一道种子点在相应道的第一位置和在预设的相关系数计算范围内的第二位置；确定当前道的映射点并在预设滑动范围限制下，将以映射点为中心且以计算范围为半径的计算区域分别向上、下滑动，确定每个计算区域与前一道种子点在计算范围内对应的区域间的相关系数，得到相应集合；计算所有相关系数的方差，判断最大相关系数是否超过预设的系数阈值且方差是否超过相关系数方差阈值，若满足，则进入下一步；基于第二位置确定最优相关区域内同位置数据点，基于此，利用预设的搜索范围，确定当前道种子点，以追踪下一道层位。本发明优化了现有层位自动追踪方法，提高了层位线的连续性和稳定性。

Description

一种探地雷达层位自动追踪方法及系统

技术领域

本发明涉及层位追踪技术领域，尤其是涉及一种探地雷达层位自动追踪方法及系统。

背景技术

探地雷达通过发射天线发射高频电磁波，遇到不同介质的分界面时发生反射或折射，反射波由接收天线接收，通过分析接收的波形信息获取地下介质信息。由于其具备快速、无损、准确的优势，在铁路基础设施检测等领域得到广泛的应用。探地雷达图像通常分为3种形式：A-SCAN、B-SCAN、C-SCAN。其中，A-SCAN指的是由各个采样点组成的一维的单道波形图，可观察一道数据的波形；B-SCAN是由若干个A-SCAN组成的二维的剖面图，可以通过不同的色标文件来呈现，以灰度图最为常见，也是探地雷达图像的最主要的形式；C-SCAN是由若干个B-SCAN进行插值组成的三维的立体图。

层位追踪是探地雷达B-SCAN图像解释的重要内容，如道床厚度检测、基床平整度检测、隧道衬砌厚度检测等。目前，层位追踪主要依赖人工完成，通过肉眼视觉观察灰度剖面图的同相轴连续性，拾取相应的层位线，但这种方式精确率较低、且效率低下，不能满足海量检测数据的处理需求。

在现有层位自动追踪算法中，通常包括如下两种方法：其一，基于波形特征的自动追踪，只利用搜索时窗内特征点的相似结构形态，但在相邻道之间不对相关系数进行计算，搜索每道波形结构相似的特征点，这种方法虽然计算量较小，但可靠性较低；其二，基于相关系数特征的自动追踪，是以前一道的种子点为中心，对当前数据与前一道相关时窗范围内的数据做相关性分析，若找到满足条件的特征点，便将当前特征点作为当前道的种子点，以继续下一道的层位追踪，这种方法虽然计算量大，但可靠性较高。

由于地下介质结构的复杂性，以及结构或含水属性的变化，导致探地雷达图像中的层位特征往往较为模糊，且连续性较差。在现有技术中，无论采用上述哪种层位自动追踪方法，都未达到较高的准确度，容易出现因前一道特征点选取误差而导致误差延续的现象发生，需要后续进行大量的人工修正工作。

发明内容

为了解决上述技术问题，本发明提供了一种探地雷达层位自动追踪方法，包括：从探地雷达二维剖面图中确定前一道种子点在相应道的第一位置和在预设的相关系数计算范围内的第二位置；基于所述第一位置确定当前道的映射点，并在预设的滑动范围的限制下，将以所述映射点为中心且以所述计算范围为半径所对应的计算区域分别向上、下滑动，确定每个计算区域与所述前一道种子点在所述计算范围内对应的区域间的相关系数，得到相应的系数集合；计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值并且所述方差是否超过预设的相关系数方差阈值，若满足，则将所述最大相关系数对应的计算区域作为最优相关区域，进入下一步；基于所述第二位置确定所述最优相关区域内的同位置数据点，基于此，利用预设的搜索范围，确定当前道的种子点，以追踪下一道层位。

优选地，所述方法，还包括：按照预设的特征点检验条件，对所述当前道的种子点进行有效性验证，其中，判断所述当前道的种子点的层位深度数据与所述前一道种子点的层位深度数据的差值是否小于预设的深度差值阈值，若满足，则进入下一步；判断所述当前道的种子点的层位深度数据是否处于预设的深度范围阈值之内，若处于，则进入下一步；确定以所述当前道的种子点为中心点且以所述搜索范围为半径所对应的区域内的最大振幅极值点和最小振幅极值点，并计算两个极值点的第一极值深度差和第一极值振幅差；判断所述第一极值深度差是否大于预设的极值深度差阈值，并且判断所述第一极值振幅差是否处于极值幅度差范围阈值之内，若满足，则确定当前种子点有效。

优选地，所述特征点的类型选自谷值、峰值和零交叉点中的一种。

优选地，所述方法还包括：对完成层位追踪的所有道的层位线进行层位置信度评价，其中，按照预设间距，将所有道分成若干组；计算每组的层位特征指标；对所述每组的层位特征指标进行加权打分处理，得到关于层位特征的置信度。

优选地，在计算每组的层位特征指标步骤中，包括：计算每组中的不连续道所占比例；和/或按照预设的第一数量，以同组内每道种子点为中心，计算每道种子点向上所述第一数量个数据点的探地雷达数据、以及向下所述第一数量个数据点的探地雷达数据的第一方差，并计算同组内所有道的所述第一方差的均值；和/或按照预设的第二数量，以同组内每道种子点为中心，计算每道种子点向上所述第二数量个数据点的探地雷达数据、以及向下所述第二数量个数据点的探地雷达数据的能量平均，并计算同组内所有道的所述能量平均的均值。

优选地，所述方法还包括：将属于不同置信度评价区间的组所对应的层位线，按照预设的展示类型，以不同形式显示在所述二维剖面图上，其中，所述展示类型选自线宽、颜色和线型差异中的一种。

优选地，还包括：从探地雷达二维剖面图中确定前一道的种子点，其中，对所述二维剖面图进行去干扰预处理，所述预处理包括背景去噪、滤波和增益处理；从经过预处理的所述二维剖面图中，选取前一道的初始特征点；以所述初始特征点为中心在所述搜索范围内，选取相应的种子点。

优选地，在计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值、并且所述方差是否超过预设的相关系数方差阈值步骤中，若不满足任一标准，则停止追踪，并标记所述当前道与前一道为不连续道。

另一方面，本发明还提供了一种探地雷达层位自动追踪系统，包括：预处理模块，其配置为从探地雷达二维剖面图中确定前一道种子点在相应道的第一位置和在预设的相关系数计算范围内的第二位置；相关集合生成模块，其配置为基于所述第一位置确定当前道的映射点，并在预设的滑动范围的限制下，将以所述映射点为中心且以所述计算范围为半径所对应的计算区域分别向上、下滑动，确定每个计算区域与所述前一道种子点在所述计算范围内对应的区域间的相关系数，得到相应的系数集合；最优相关处搜索模块，其配置为计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值并且所述方差是否超过预设的相关系数方差阈值，若满足，则将所述最大相关系数对应的计算区域作为最优相关区域，执行种子点追踪模块；所述种子点追踪模块，其配置为基于所述第二位置确定所述最优相关区域内的同位置数据点，基于此，利用预设的搜索范围，确定当前道的种子点，以追踪下一道层位。

优选地，所述系统，还包括：种子点检验模块，其配置为按照预设的特征点检验条件，对所述当前道的种子点进行有效性验证，其中，所述种子点检验模块，具备如下单元：种子点深度差判断单元，其配置为判断所述当前道的种子点的层位深度数据与所述前一道种子点的层位深度数据的差值是否小于预设的深度差值阈值，若满足，则执行种子点深度范围判断单元；所述种子点深度范围判断单元，其配置为判断所述当前道的种子点的层位深度数据是否处于预设的深度范围阈值之内，若处于，则执行极值搜索单元；所述极值搜索单元，其配置为确定以所述当前道的种子点为中心点且以所述搜索范围为半径所对应的区域内的最大振幅极值点和最小振幅极值点，并计算两个极值点的第一极值深度差和第一极值振幅差；极值条件诊断单元，其配置为判断所述第一极值深度差是否大于预设的极值深度差阈值，并且判断所述第一极值振幅差是否处于极值幅度差范围阈值之内，若满足，则确定当前种子点有效。

与现有技术相比，上述方案中的一个或多个实施例可以具有如下优点或有益效果：

本发明提供一种探地雷达层位自动追踪方法及系统。该方法及系统通过加强对种子点待搜索区域的搜索条件、以及待搜索区域内种子点搜索条件的约束性，一方面滤除掉虽满足搜索条件，但当前道的种子点待搜索区域与前一道相关性程度不强的非突出连续性区域的干扰；还消除了属于极值点，但不具有强层位代表性的特征点作为种子点所造成的误差，从而对现有基于相关的层位自动追踪方法进行了优化，提高了现有方法得到的层位线的连续性和稳定性。另外，本发明还提供了一种对层位线内不同曲线段之间连续性程度的置信区间的评价方法，使得本发明得到的层位线展示了含有不同曲线段之间对应的清晰度和显著性程度等级的信息，更加直观的展现出层位的连续性和稳定性程度，为基于层位线进行后续图像解释提供了丰富的依据，进一步提高了后续地质解释的效率。

虽然在下文中将结合一些示例性实施及使用方法来描述本发明，但本领域技术人员应当理解，为并不旨在将本发明限制于这些实施例。反之，旨在覆盖包含在所附的权利要求书所定义的本发明的精神与范围内的所有替代品、修正及等效物。

本发明的其他优点、目标，和特征在某种程度上将在随后的说明书中进行阐述，并且在某种程度上，基于对下文的考察研究对本领域技术人员而言将是显而易见的，或者可以从本发明的实践中得到教导。本发明的目标和其他优点可以通过下面的说明书，权利要求书，以及附图中所特别指出的结构来实现和获得。

附图说明

附图用来提供对本发明的进一步理解，并且构成说明书的一部分，与本发明的实施例共同用于解释本发明，并不构成对本发明的限制。在附图中：

图1为人工拾取探地雷达层位线的示意图。

图2为本申请实施例的探地雷达层位自动追踪方法的步骤图。

图3为本申请实施例的探地雷达层位自动追踪方法中最大相关系数诊断过程的流程图。

图4为本申请实施例的探地雷达层位自动追踪方法中相关系数集合生成过程的原理示意图。

图5为本申请实施例的探地雷达层位自动追踪方法中种子点检验过程的流程图。

图6为现有层位自动追踪算法中特征点检验方法的缺陷示意图。

图7为应用本申请实施例的探地雷达层位自动追踪方法得到的层位线与现有层位自动追踪算法得到的层位线的对比效果示意图。

图8为本申请实施例的探地雷达层位自动追踪系统的模块框图。

具体实施方式

以下将结合附图及实施例来详细说明本发明的实施方式，借此对本发明如何应用技术手段来解决技术问题，并达成技术效果的实现过程能充分理解并据以实施。需要说明的是，只要不构成冲突，本发明中的各个实施例以及各实施例中的各个特征可以相互结合，所形成的技术方案均在本发明的保护范围之内。

另外，在附图的流程图示出的步骤可以在诸如一组计算机可执行指令的计算机系统中执行，并且，虽然在流程图中示出了逻辑顺序，但是在某些情况下，可以以不同于此处的顺序执行所示出或描述的步骤。

探地雷达通过发射天线发射高频电磁波，遇到不同介质的分界面时发生反射或折射，反射波由接收天线接收，通过分析接收的波形信息获取地下介质信息。由于其具备快速、无损、准确的优势，在铁路基础设施检测等领域得到广泛的应用。探地雷达图像通常分为3种形式：A-SCAN、B-SCAN、C-SCAN。其中，A-SCAN指的是由各个采样点组成的一维的单道波形图，可观察一道数据的波形；B-SCAN是由若干个A-SCAN组成的二维的剖面图，可以通过不同的色标文件来呈现，以灰度图最为常见，也是探地雷达图像的最主要的形式；C-SCAN是由若干个B-SCAN进行插值组成的三维的立体图。

层位追踪是探地雷达B-SCAN图像解释的重要内容，如道床厚度检测、基床平整度检测、隧道衬砌厚度检测等。目前，层位追踪主要依赖人工完成，通过肉眼视觉观察灰度剖面图的同相轴连续性，拾取相应的层位线，但这种方式精确率较低、且效率低下，不能满足海量检测数据的处理需求。

探地雷达资料的地质解释是探地雷达测量的一项重要目的，对雷达数据剖面进行准确的地质解释，需要比较准确地提取反射层，即“层位追踪”。探地雷达资料中具有能够展现出二维剖面结构的GPR图像，GPR图像为二维矩阵的灰度图，每一列表示1道数据，其中，1道数据中包括自上而下排列的多个数据点，每个数据点对应一组探地雷达数据，每组探地雷达数据包括：当前数据点对应的层位深度数据、当前数据点的振幅数据(反射强度)等。层位追踪的目标是每1道找到一个特征点(种子点)，从而由多条数据道内的特征点形成为层位线，如图1所示。

在现有层位自动追踪算法中，通常包括如下两种方法：其一，基于波形特征的自动追踪，只利用搜索时窗内特征点的相似结构形态，但在相邻道之间不对相关系数进行计算，搜索每道波形结构相似的特征点，这种方法虽然计算量较小，但可靠性较低；其二，基于相关系数特征的自动追踪，是以前一道的种子点为中心，对当前数据与前一道相关时窗范围内的数据做相关性分析，若找到满足条件的特征点，便将当前特征点作为当前道的种子点，以继续下一道的层位追踪，这种方法虽然计算量大，但可靠性较高。

由于地下介质结构的复杂性，以及结构或含水属性的变化，导致探地雷达图像中的层位特征往往较为模糊，且连续性较差。在现有技术中，无论采用上述哪种层位自动追踪方法，都未达到较高的准确度，容易出现因前一道特征点选取误差而导致误差延续的现象发生，需要后续进行大量的人工修正工作。又由于地下地质界面在局部区域具有一定的连续性、稳定性，由此，反映到探地雷达二维剖面图上，则表现为相邻道上的反射层的同相轴具有在振幅上具有稳定性和渐变性，并在波形上具有连续性和相似性的特征。在实际的层位追踪方法中所形成的层位线本应在更接近表面的位置，但由于层位追踪误差，使得特征点逐渐向下跑偏，进而需要手动重新设置种子点，才能继续向右追踪，由此，参考图7(a)中的层位线并未体现出层位结构的连续性和相似性特征。

层位自动追踪是，在探地雷达二维剖面图的不同道(测线)上放置“种子点”，用来约束控制层位自动追踪的计算，然后根据种子点的特征准则，计算相邻道之间的相似性，找出相邻道之间的相似特征种子点，若在指定的条件下搜索到了满足条件的特征种子点就取出来，以这个特征点为新的特征种子点，再搜索计算下一道。在层位自动追踪的时候，一般需指定自动追踪时自动追踪的特征、振幅范围的大小以及自动追踪时的控制时窗范围(相关时窗范围)，在层位自动追踪过程中，若没能够搜索到满足指定条件的特征种子点，自动追踪就停止在这一道上，或者跳过这一道，继续追踪下一道，直到结束。

因此，本发明为了解决现有层位自动追踪方法的技术问题，提出了一种对现有方法进行优化的探地雷达层位自动追踪方法。该方法在得到关于当前道与前一道(当前道的左侧道)的相关系数集合后对最大相关系数的选取过程(当前道的最佳相关处)进行了条件优化限制，不但需要将系数集合内的方差与预设的相关系数方差阈值进行对比，以对当前道与前一道的整体相关性进行检验，还需要将最大相关系数与预设系数阈值进行对比，以寻找出与前一道具有较强相关性的当前道内的最佳相关处，从而使得优选出的当前最大相关系数所属的待搜索区域与前一道在层位结构连续性上的相关性程度较为明显、突出。这样，本发明能够在于前一道具有突出的相关性数据的待搜索区域内来追踪当前道的种子点，提高了当前道种子点定位的准确性，更为以当前道种子点为基础的针对下一道进行层位追踪的后续操作打下了良好的基础，从而提高了这一后续操作的准确性。另一方面，本发明利用这种改进的当前道最佳相关区域的选取方法，还能够在当前道与前一道整体相关性较差时，自动停止或跳过，以避免产生误差延续的问题。

进一步，本发明对在优选出的最大相关系数所属的搜索区域内寻找特征点的条件进行了优化，通过增加极值差阈值的方式，提高了特征点搜索条件的约束性，改善了选取出的特征点无法准确代表当前道层位的问题。

另外，本发明基于对现有层位自动追踪方法改进后得到的层位线，对层位线设定评价指标(层位特征指标)、并基于层位特征指标归纳了置信度评价区间，使得本发明实施例不仅能够得到准确性和可靠性较高的层位线，还能够更加直观的知晓不同段层位线的数据指标所属等级，大大提高地质解释的直观性、准确性及效率，为基于层位线进行后续图像解释提供了丰富的依据。

图2为本申请实施例的探地雷达层位自动追踪方法的步骤图。下面参考图2，对本发明中的探地雷达层位自动追踪方法的步骤进行说明。首先，步骤S210，从探地雷达二维剖面图中确定前一道种子点，并确定出该种子点在相应道的第一位置、以及该种子点在预设的相关系数计算范围内的第二位置。需要说明的是，此处的前一道指的是当前待追踪范围内的第一道。在步骤S210中，需要先获取到关于待追踪区域的原始二维剖面图，在从当前二维剖面图中选取出待追踪道的范围，确定待追踪范围内的最左边道(第一道)，即为前一(数据)道，将第二(数据)道确定为当前段(数据)道；而后，根据每道所包含的数据点的数量，确定当前追踪算法的初始化参数，基于此，确定前一道内的种子点的位置。

其中，在本发明一实施例中，初始化参数包括：(种子点)搜索范围、相关系数计算范围和滑动半径范围。参考图6，在本发明实施例中，在对种子点进行搜索时需要在指定搜索区域内进行确定，其中，指定搜索区域是以特定数据点为中心并以搜索范围r1为半径所限定的区域；在选取当前道内与前一道具有最强相关性的最优相关区域时，需要先确定出多个待评价相关区域(计算区域)，每个待评价相关区域的范围即为以特定数据点为中心且以相关系数计算范围r3为半径所限定的区域，而后从多个待评价相关区域筛选出最优相关区域，从而在最优相关区域内选取出指定搜索区域，使得在对满足最优相关条件的道内的最优位置范围内选取出准确代表当前到层位的种子点；在筛选最优相关区域时所对应的滑动范围边界由滑动半径范围r2所限定。例如：若每道的数据点数量为512，可设置r1＝5，r2＝5，r3＝30。

进一步，由于第一道的层位追踪过程并没有前一道的种子点作为基础进行追踪，因此，在这一过程中需要解释人员的辅助，以对当前追踪范围的第二道的层位进行追踪。需要说明的是，在按照下述步骤S220～步骤S240确定出第二道甚至是后续道的种子点后，都可以基于前一(相邻)连续道的种子点，再次利用步骤S220～步骤S240对下一道的种子点进行选取，而无需解释人员的参与，便能够完成追踪范围内后续所有道的层位追踪。

进一步，在确定前一道(第一道)的种子点的位置时，首先，需要对上述原始二维剖面图进行去干扰预处理。在本发明实施例中，去干扰预处理的目的是压制规则的和随机的干扰信号，提高地质雷达图像的分辨率，突出电磁波振幅、波形及其随时间、里程变化等有用的异常信息。其中，去干扰预处理包括但不限于：背景去噪、水平滤波、垂直带通滤波和增益。而后，对需要从经过预处理的二维剖面图中，采用人工看图法，选取前一道(第一道)的初始特征点A0。接着，根据已设定好的种子点搜索范围r1，以初始特征点A0为中心，在种子点搜索范围为半径的区域内选取出更准确的种子点A1。

进一步，在本发明一实施例中，特征点的类型选自：谷值、峰值和零交叉点中的一种。需要说明的是，本发明对特征点的类型不作具体限定，本领域技术人员可根据实际需求进行选择。当确定出上述初始特征点A0后，以初始特征点A0为中心、并且以种子点搜索范围为半径，选取出当前搜索范围内的特征点，此处的选取方法可以按照现有的波形特征追踪方法或其他合理方法得到，本发明对此不作限定。

在完成上述整个预处理过程后，结束步骤S210，进入到步骤S220，开始对第二道的层位进行追踪。步骤S220基于前一道种子点在前一道的(第一)位置，确定当前道的映射点，并在预设的滑动半径范围边界的限制下，将以当前道映射点为中心、且在上述相关系数计算范围所对应的计算区域，分别向上、向下进行滑动，确定每个计算区域与以前一道种子点在相关系数计算范围内对应的区域(在本发明实施例中，以前一道种子点在相关系数计算范围对应的区域，优选为，以前一道种子点为中心点、且以相关系数计算范围为半径所对应的区域)之间的相关系数，得到相应的相关系数集合。

更具体地说，在上述步骤S220中，首先，需要以前一道种子点A1为中心、并且以相关系数计算范围r3为半径，将当前层位搜索范围内的所有数据点组成为第一向量L1。而后，基于前一道种子点的深度位置，在当前道上的相同深度位置处确定出前一道种子点的映射点。接着，限定出以当前映射点为中心、且以相关系数计算范围r3为半径对应的计算区域(第二向量L2)，将当前计算区域内的所有数据点组成为第二向量L2(例如，若r3为30表示当前相关系数计算范围包括30个数据点，那么，以当前映射点为中心所对应的第二向量L2，包括：中心点、中心点向上依次排列的30个数据点、以及中心点向下依次排列的30个数据点)，从而根据第一向量L1内的所有数据点的探地雷达数据和第二向量L2内的所有数据点的探地雷达数据，计算出两个向量之间的相关系数(R值)。

接下来，以当前映射点为原点，在上述滑动半径范围作为滑动边界范围的限定下，将当前道内的以当前映射点为中心对应的第二向量L2分别基于原点整体向上滑动N次、基于原点整体向下滑动N次，一次滑动一个数据点(例如，若r2＝N＝5表示当前滑动半径范围包括5个数据点，在整体向上滑动时，依次以当前原点向上相邻排列的5个数据点作为新的中心点所分别构成的5个第二向量L2；在整体向下滑动时，依次以当前原点向下相邻排列的5个数据点作为新的中心点所分别构成的5个第二向量L2；结合以当前映射点为中心所对应的第二向量L2，从而得到11个第二向量L2)，同时，计算每次移动一个数据点所构成的第二向量L2与第一向量L1之间的相关系数，从而得到2N+1个相关系数(R值)，继而将这些相关系数数据形成为相关系数集合。其中，集合内的每个相关系数对应有相应的(第二道)计算区域(即第二向量L2)。

在得到相关系数集合后，进入到步骤S230中，来搜索与前一道内的第一向量具有明显相关性的种子点待搜索区域(最优相关区域)。现有技术中，在对计算区域进行选取时，常常仅利用预设的相关系数阈值来与集合中的最大相关系数进行对比，在集合中的最大相关系数超过相关系数阈值，便将当前最大相关系数所对应的区域作为与前一道具有连续性特征的种子点待搜索区域，但这种方法得到的种子点待搜索区域并不是与前一道具有突出的连续性特征的待搜索区域，往往难以获取到准确代表当前道层位的种子点。由于实际情况中，除层面极为规整的区域外，大部分数据道的连续性或者相关性并不好，故相关系数阈值通常不会设的很大，以免自动追踪的结果过于断续或者频繁拾取新的种子点，而相关系数阈值设的小的弊端就是，会在某一道错误的位置一直延续下去。

图4为本申请实施例的探地雷达层位自动追踪方法中最大相关系数选取原则的示意图。图4中的两条曲线展示了两种不同类型的相关系数集合内系数数据的分布情况。第一曲线的弧度较为平滑，虽有数据点超过阈值，但超过阈值的数据点的数值并不高，实际反映出两道数据的相关性并不好，也反映出设置的阈值过低的问题，在遇到这种现象时，这一道应该标记为不连续道或直接跳过，若从第一曲线中所选取的最优相关处来确定种子点，就会导致出现后续道将种子点误差延续的问题。第二曲线的弧度变化较大，具有明显的拐点，超过相关系数阈值的数据点较多，为了选取到与前一道相关性更强的待搜索区域，而滤除掉像第一曲线中那些超出相关系数阈值较少的多个非最佳数据点的影响。

为了解决上述问题，本发明实施例在选取最优相关区域(即种子点所属待搜索区域)的过程中，还加入了方差阈值的对比操作。步骤S230计算上述相关系数集合内所有相关系数的方差，而后，判断当前相关系数集合内的最大相关系数是否超过预设的相关系数阈值、并且判断上述方差是否超过预设的相关系数方差阈值，若同时满足，则将当前最大相关系数对应的计算区域作为最优相关区域。

图3为本申请实施例的探地雷达层位自动追踪方法中最大相关系数诊断过程的流程图。如图3所示，首先，步骤S301需要计算步骤S220得到的相关系数集合内所有相关系数数据的(相关系数)方差，而后进入到步骤S302中，以该方差值作为一种待评价因素，对最大相关系数进行评估。步骤S302确定出当前相关系数集合的最大相关系数数据，步骤S303先将当前最大相关系数与预设的相关系数阈值进行对比，当超过时，进入到步骤S304中；当最大相关系数小于或等于上述相关系数阈值时，进入到步骤S306中。步骤S304继续将步骤S301得到的当前相关系数集合内的所有相关系数数据的方差与预设的相关系数方差阈值(相关系数方差阈值表示能够体现出集合内所有超过相关系数阈值的数据点数量及超过程度具有显著差距特征的最小方差阈值)进行对比，在当前集合方差超过相关系数方差阈值时，进入到步骤S305中；在当前集合方差小于或等于相关系数方差阈值时，进入到步骤S306中。步骤S305将当前最大相关系数所属的计算区域确定为特征点搜索目标区域，即最优相关区域，步骤S230结束，以等待进入到步骤S240中进行种子点搜索。

步骤S306当前最大相关关系无法表征与前一道具有突出的连续性特征，立即停止追踪，并标记当前道与前一道为不连续道关系。

进一步，在完成对当前最大相关系数所属的计算区域的评估并得到特征点搜索目标区域后，进入到步骤S240中。步骤S240基于前一道种子点在相关系数计算范围内的第二位置，确定最优相关区域内的同位置数据点，而后，基于当前同位置数据点，利用预设的搜索范围，确定当前道的种子点，以将当前种子点作为追踪下一道层位的基础。

具体地，在步骤S240中，首先，需要确定前一道种子点在相关系数计算范围对应的区域，即为，以前一道种子点为中心点、并以相关系数计算范围r3为半径所对应的第一区域(上述第一向量L1)，并记录下前一道种子点在第一向量L1内所处的数据点编号位置(也就是说，确定出前一道种子点在第一向量L1中属于第几个编号位置的数据点)。而后，确定出以前一道种子点在当前道的映射点为中心点、且以相关系数计算范围r3为半径所对应的第二区域(上述以前一道种子点在当前道的映射点为中心对应的第二向量L2)，并确定出第二区域中具有上述相同数据点编号位置的数据点，记为初级种子点(例如，前一道种子点在第一向量L1中为第10个数据点，那么当前道上的相同位置处的数据点，指的是第二向量L2中的第10个数据点)。接着，根据已设定好的种子点搜索范围r1，以当前初级种子点为中心、且以种子点搜索范围为半径所对应的种子点最优搜索区域(最优相关区域)内选取出更准确的种子点，作为当前道的种子点(例如：在当前特征点类型为波谷时，需要从种子点最优搜索区域内选取出振动幅度最大的波谷值对应的数据点)。

进一步，由于现有技术中，在对待搜索区域内的种子点进行选取后时，常常还需要对选取出的当前道种子点进行有效性验证，从而检验当前种子点是否能够准确代表当前道层位。具体地，常常将待搜索区域内种子点的探地雷达数据与前一道种子点的探地雷达数据的差值直接利用深度差值阈值比较、和/或将待搜索区域内种子点的探地雷达数据与深度范围阈值(包括最大深度和最小深度)来进行比较的方式进行筛选检验，但这种方法得到的种子点并不是真正的层位点，该点虽在所属的最优搜索区域内，但并没有落入发生在振幅反相且相位差(需要说明的是，此处的相位差指的是最大极值点与最小极值点的深度之差)较大的最优极值点深度范围内的点。在实际层位追踪过程中，真正要找到的分界面(层位线)两侧的介质往往反射系数较大，进而导致雷达波形相应位置处发生振幅反相且相位差较大。因此，往往选择出的种子点需要是落入上述最优极值点深度范围内的点，才能准确且突出的代表当前道的层位点。这样，利用这一特点对当前种子点进行验证，能够进一步提升层位追踪的精度。

因此，为了解决这一问题，本发明所述的探地雷达层位自动追踪方法还包括：步骤S250。步骤S250按照预设的特征点检验条件，对当前道的种子点进行有效性验证。本发明在步骤S250中提高了特征点验证条件的约束性，验证出最优的种子点。

图6为现有层位自动追踪算法中特征点检验方法的缺陷示意图。如图6所示，该曲线为某一搜索区域内所有数据点的探地雷达数据幅值分布情况，横坐标为(深度)位置，纵坐标为探地雷达数据的幅值。其中，d1表示应用现有特征点检验方法没有滤除掉的种子点，该种子点并没有落入发生在振幅反相且反向振幅之间的相位差较大的最优极值点深度范围内，并不是最能够显著且准确代表当前道层位的特征点；d2展示了一种落入发生在振幅反相且反向振幅之间的相位差较大的最优极值点深度范围内的种子点，该类种子点是能够相较于d1更具有层位代表性。

图5为本申请实施例的探地雷达层位自动追踪方法中种子点检验过程的流程图。如图5所示，首先，步骤S501确定以当前道的种子点为中心且以上述搜索范围为半径所限定的检验区域内各数据点的探地雷达数据，步骤S502计算当前种子点的探地雷达数据中的层位深度数据与前一道种子点的探地雷达数据中的层位深度数据的(相邻道种子点)深度差值，步骤S503将步骤S502得到的(相邻道种子点)深度差值与预设的深度差值阈值进行对比，并判断当前(相邻道种子点)深度差值是否小于深度差值阈值，若小于，则进入到步骤S504中。另外，在步骤S503中，若当前(相邻道种子点)深度差值大于或等于深度差值阈值，则进入到步骤S510中。

步骤S504继续将当前道的种子点的探地雷达数据中的层位深度数据与预设的深度范围阈值(即最大深度阈值和最小深度阈值)进行比较，并判断当前道的种子点的层位深度数据是否处于上述深度范围阈值之内(即是否位于最大深度阈值和最小深度阈值之间)，若满足这一条件，则进入到步骤S505中。另外，在步骤S504中，若当前道的种子点的层位深度数据超出上述深度范围阈值之外(即当前道的种子点的层位深度数据大于或等于最大深度阈值，或者小于或等于最小深度阈值)，则进入到步骤S510中。

步骤S505确定上述检验区域内的最大振幅极值点和最小振幅极值点，步骤S506计算当前最大振幅极值点和当前最小振幅极值点的深度差、及振幅差，记为第一极值深度差和第一极值振幅差，从而进入到步骤S507中。

步骤S507将上述步骤S506得到的第一极值深度差与预设的极值深度差阈值进行对比，判断当前第一极值深度差是否大于预设的极值深度差阈值，若大于，则进入到步骤S508中。另外，在步骤S507中，若当前第一极值深度差小于或等于上述极值深度差阈值，则进入到步骤S510中。

进一步，步骤S508继续将上述S506得到的第一极值振幅差与预设的极值幅度差范围阈值(包括极值幅度差最大阈值和极值幅度差最小阈值)进行对比，判断当前第一极值振幅差是否处于极值幅度差范围阈值之内，若满足这一条件，则进入到步骤S509中。步骤S509标记当前种子点为当前到的有效种子点。另外，在步骤S508中，若当前第一极值振幅差超出上述极值幅度差范围阈值之外(即当前第一极值振幅差大于或等于极值幅度差最大阈值，或者小于或等于极值幅度差最小阈值)，则进入到步骤S510中。

步骤S510当前待搜索区域内没有符合特征点搜索条件的种子点，立即停止追踪，并标记当前道与前一道为不连续道关系。

需要说明的是，本发明对上述深度差值阈值、深度范围阈值、极值深度差阈值和极值幅度差范围阈值的数值不作具体限定，本领域技术人员可根据层位介质性质(例如介质电性参数)等因素进行设定并调整。

由此按照上述步骤S220～步骤S240所述的方法追踪到了当前道的准确的种子点，消除了属于极值点但不具有较强层位代表性的特征点作为种子点所造成的误差，避免了将这种误差传递下去的现象，从而以当前道的更为准确的种子点作为基础，对当前道的下一相邻道的种子点(层位)进行追踪。另外，在通过加强种子点检验条件的约束性，一方面提高了形成更为准确且明显的层位线，更精准的发现没有层位点或者层位特征不明显的数据道，从而跳过这一道或在这一道停止，另一方面，更精准的定位最优相关处和特征点。由此，完成针对整个追踪区域的层位追踪过程，提高了最终得到的层位线的准确性、连续性和稳定性。

进一步，在本发明中，在完成针对整个追踪区域的层位追踪过程后，还能够对完成层位追踪的所有道的层位线进行层位置信度评价。具体地，首先，需要按照预设间距，将所有道分成若干组。例如：每10道一组。而后，计算每组的层位特征指标。接着，对每组的层位特征指标进行加权打分处理，从而得到关于当前层位线在层位特征指标上的(不同等级的)置信度区间。最后，将代表不同等级程度的置信度区间进行归一化处理。优选地，将每个置信度区间都归一化至(0，1)区间范围内。

需要说明的是，层位特征指标为描述层位线内不同道之间清晰度及明显性程度等级的特征类型，该特征类型可根据层位解释需求进行自定义。在本发明实施例中，层位特征指标可定义为如下类型中的一种或几种的组合，也可以定义成其他合理类型。在本发明一实施例中，层位特征指标优选为不连续道所占比例；和/或，在本发明一实施例中，层位特征指标优选为第一方差的均值；和/或，在本发明一实施例中，层位特征指标优选为能量平均的均值。其中，上述第一方差的均值按照如下方式计算：按照预设的第一数量，以同组内每道种子点为中心，根据每道种子点向上第一数量个数据点的探地雷达数据、以及向下第一数量个数据点的探地雷达数据，计算每道的探地雷达数据的方差，记为第一方差，而后计算同组内所有道的第一方差的平均值。上述能量平均的均值按照如下方式计算：按照预设的第二数量，以同组内每道种子点为中心，根据每道种子点向上第二数量个数据点的探地雷达数据(中的能量值)、以及向下第二数量个数据点的探地雷达数据(中的能量值)，计算每道的能量值的平均值，记为能量平均，而后计算同组内所有道的能量平均的平均值。进一步，本发明对上述预设数量的大小不作具体限定，本领域技术人员可根据实际需求进行设定。

由此，本发明实施例在按照上述方法对当前追踪范围内的层位线进行了层位置信度评价，并得到了适应于当前层位线的层位置信度评价结果，还需要将当前二维剖面图中的层位线，按照上述层位置信度评价结果，确定出当前层位线内不同组所属的置信度评价区间。进一步，按照预设的展示类型，将属于不同区间的层位线以不同的形式显示在二维剖面图上，从而直观的向解释人员展示出不同组之间的层位线的明显及清晰度的量化程度。在层位追踪过程中，层位线明显程度或清晰程度越高的组段，表示对应道的层位反射能量越强。其中，展示类型选自线宽、颜色和线型差异等中的一种。也就是说，在本发明实施例中，最终所展现的层位线，可以以不同线宽的形式表示出属于不同置信区间的层位线段，可以以不同线的颜色的形式表示出属于不同置信区间的层位线段，也可以以不同线的类型(点划线、虚线、实线等)的形式表示出属于不同置信区间的层位线段，还可以以其他合理形式来区别出层位线中属于不同置信区间的层位线段。

这样，本发明不仅提高了追踪到的层位线的准确性和连续性，还提供了一种对层位线进行置信区间评价的确定方法，将层位线内展现出的不同道之间的层位的明显程度或清晰程度进行了量化表示，减轻了后续对层位线进行评价时的计算负担，为对基于层位线进行后续地质解释提供了更多的依据，进一步提高了后续地质解释的效率。

图7为应用本申请实施例的探地雷达层位自动追踪方法得到的层位线与现有层位自动追踪算法得到的层位线的对比效果示意图。图7(a)中的方框内展示了应用现有层位自动追踪算法得到的层位线，图7(b)中的方框内展示了应用本发明实施例的探地雷达层位自动追踪方法得到的层位线，在图7(b)中，用实心点表示追踪到的层位点，虚心点表示与前一相邻道为不连续道的位置，由实心点和虚心点将层位线中的断续处进行了区分。从图7中可以看出，采用本发明所述的探地雷达层位自动追踪方法得到的层位线，准确地表现出了层位在局部区域内具有连续性和稳定性的特征。

另一方面，基于上述探地雷达层位自动追踪方法，本发明还提出了一种探地雷达层位自动追踪系统。图8为本申请实施例的探地雷达层位自动追踪系统的模块框图。如图8所示，本发明实施例所述的探地雷达层位自动追踪系统，包括：预处理模块81、相关集合生成模块82、最优相关处搜索模块83和种子点追踪模块84。其中，预处理模块81按照上述步骤S210所述的方法实施，配置为从探地雷达二维剖面图中确定前一道种子点在相应道的第一位置和在预设的相关系数计算范围内的第二位置。相关集合生成模块82按照上述步骤S220所述的方法实施，配置为基于上述第一位置确定当前道的映射点，并在预设的滑动范围的限制下，将以映射点为中心且以相关系数计算范围为半径所对应的计算区域分别向上、向下滑动，确定每个计算区域与前一道种子点在相关系数计算范围内对应的区域之间的相关系数，得到相应的系数集合。最优相关处搜索模块83按照上述步骤S230所述的方法实施，配置为计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值、并且当前方差是否超过预设的相关系数方差阈值，若同时满足，则将最大相关系数对应的计算区域作为最优相关区域，而后执行种子点追踪模块84。种子点追踪模块84按照上述步骤S240所述的方法实施，配置为基于上述第二位置，确定最优相关区域内的同位置数据点，而后，基于该数据点，利用预设的搜索范围，确定当前道的种子点，以追踪下一道层位。

另外，本发明实施例所述的探地雷达层位自动追踪系统，还包括：种子点检验模块85。种子点检验模块85按照上述步骤S250所述的方法实施，配置为按照预设的特征点检验条件，对当前道的种子点进行有效性验证。

其中，种子点深度差判断单元851，配置为判断当前道的种子点的层位深度数据与前一道种子点的层位深度数据的差值是否小于预设的深度差值阈值，若满足，则执行种子点深度范围判断单元852。种子点深度范围判断单元852，配置为判断当前道的种子点的层位深度数据是否处于预设的深度范围阈值之内，若处于，则执行极值搜索单元853。极值搜索单元853，配置为确定以当前道的种子点为中心点且以搜索范围为半径所对应的区域内的最大振幅极值点和最小振幅极值点，并计算最大振幅极值点和最小振幅极值点的第一极值深度差和第一极值振幅差。极值条件诊断单元954，配置为判断第一极值深度差是否大于预设的极值深度差阈值，并且判断第一极值振幅差是否处于极值幅度差范围阈值之内，若满足，则确定当前种子点有效。

本发明提供一种探地雷达层位自动追踪方法及系统。该方法及系统通过加强对种子点待搜索区域的搜索条件、以及待搜索区域内种子点搜索条件的约束性，一方面滤除掉虽满足搜索条件，但当前道的种子点待搜索区域与前一道相关性程度不强的非突出连续性区域的干扰；还消除了属于极值点，但不具有强层位代表性的特征点作为种子点所造成的误差，从而对现有基于相关的层位自动追踪方法进行了优化，提高了现有方法得到的层位线的连续性和稳定性。另外，本发明还提供了一种对层位线内不同曲线段之间连续性程度的置信区间的评价方法，使得本发明得到的层位线展示了含有不同曲线段之间对应的清晰度和显著性程度等级的信息，更加直观的展现出层位的连续性和稳定性程度，为基于层位线进行后续图像解释提供了丰富的依据，进一步提高了后续地质解释的效率。例如：在铁路路基检测中，道床与基床界面的连续性是基床密实度、平整度的重要反映，道床内部界面的连续性和稳定性是道床脏污、含水程度的重要反映。

以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉该技术的人员在本发明所揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

应该理解的是，本发明所公开的实施例不限于这里所公开的特定结构、处理步骤或材料，而应当延伸到相关领域的普通技术人员所理解的这些特征的等同替代。还应当理解的是，在此使用的术语仅用于描述特定实施例的目的，而并不意味着限制。

说明书中提到的“一个实施例”或“实施例”意指结合实施例描述的特定特征、结构或特性包括在本发明的至少一个实施例中。因此，说明书通篇各个地方出现的短语“一个实施例”或“实施例”并不一定均指同一个实施例。

虽然本发明所揭露的实施方式如上，但所述的内容只是为了便于理解本发明而采用的实施方式，并非用以限定本发明。任何本发明所属技术领域内的技术人员，在不脱离本发明所揭露的精神和范围的前提下，可以在实施的形式上及细节上作任何的修改与变化，但本发明的专利保护范围，仍须以所附的权利要求书所界定的范围为准。

Claims (10)

1.一种探地雷达层位自动追踪方法，包括：

从探地雷达二维剖面图中确定前一道种子点在相应道的第一位置和在预设的相关系数计算范围内的第二位置；

基于所述第一位置确定当前道的映射点，并在预设的滑动范围的限制下，将以所述映射点为中心且以所述计算范围为半径所对应的计算区域分别向上、下滑动，确定每个计算区域与所述前一道种子点在所述计算范围内对应的区域间的相关系数，得到相应的系数集合；

计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值并且所述方差是否超过预设的相关系数方差阈值，若满足，则将所述最大相关系数对应的计算区域作为最优相关区域，进入下一步；

基于所述第二位置确定所述最优相关区域内的同位置数据点，基于此，利用预设的搜索范围，确定当前道的种子点，以追踪下一道层位。

2.根据权利要求1所述的方法，其特征在于，所述方法，还包括：按照预设的特征点检验条件，对所述当前道的种子点进行有效性验证，其中，

判断所述当前道的种子点的层位深度数据与所述前一道种子点的层位深度数据的差值是否小于预设的深度差值阈值，若满足，则进入下一步；

判断所述当前道的种子点的层位深度数据是否处于预设的深度范围阈值之内，若处于，则进入下一步；

确定以所述当前道的种子点为中心点且以所述搜索范围为半径所对应的区域内的最大振幅极值点和最小振幅极值点，并计算两个极值点的第一极值深度差和第一极值振幅差；

判断所述第一极值深度差是否大于预设的极值深度差阈值，并且判断所述第一极值振幅差是否处于极值幅度差范围阈值之内，若满足，则确定当前种子点有效。

3.根据权利要求1或2所述的方法，其特征在于，所述特征点的类型选自谷值、峰值和零交叉点中的一种。

4.根据权利要求1～3中任一项所述的方法，其特征在于，所述方法还包括：对完成层位追踪的所有道的层位线进行层位置信度评价，其中，

按照预设间距，将所有道分成若干组；

计算每组的层位特征指标；

对所述每组的层位特征指标进行加权打分处理，得到关于层位特征的置信度。

5.根据权利要求4所述的方法，其特征在于，在计算每组的层位特征指标步骤中，包括：

计算每组中的不连续道所占比例；和/或

按照预设的第一数量，以同组内每道种子点为中心，计算每道种子点向上所述第一数量个数据点的探地雷达数据、以及向下所述第一数量个数据点的探地雷达数据的第一方差，并计算同组内所有道的所述第一方差的均值；和/或

按照预设的第二数量，以同组内每道种子点为中心，计算每道种子点向上所述第二数量个数据点的探地雷达数据、以及向下所述第二数量个数据点的探地雷达数据的能量平均，并计算同组内所有道的所述能量平均的均值。

6.根据权利要求4或5所述的方法，其特征在于，所述方法还包括：将属于不同置信度评价区间的组所对应的层位线，按照预设的展示类型，以不同形式显示在所述二维剖面图上，其中，所述展示类型选自线宽、颜色和线型差异中的一种。

7.根据权利要求1～6中任一项所述的方法，其特征在于，所述方法，还包括：从探地雷达二维剖面图中确定前一道的种子点，其中，

对所述二维剖面图进行去干扰预处理，所述预处理包括背景去噪、滤波和增益处理；

从经过预处理的所述二维剖面图中，选取前一道的初始特征点；

以所述初始特征点为中心在所述搜索范围内，选取相应的种子点。

8.根据权利要求1～7中任一项所述的方法，其特征在于，在计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值、并且所述方差是否超过预设的相关系数方差阈值步骤中，若不满足任一标准，则停止追踪，并标记所述当前道与前一道为不连续道。

9.一种探地雷达层位自动追踪系统，包括：

预处理模块，其配置为从探地雷达二维剖面图中确定前一道种子点在相应道的第一位置和在预设的相关系数计算范围内的第二位置；

相关集合生成模块，其配置为基于所述第一位置确定当前道的映射点，并在预设的滑动范围的限制下，将以所述映射点为中心且以所述计算范围为半径所对应的计算区域分别向上、下滑动，确定每个计算区域与所述前一道种子点在所述计算范围内对应的区域间的相关系数，得到相应的系数集合；

最优相关处搜索模块，其配置为计算所有相关系数的方差，判断集合内的最大相关系数是否超过预设的系数阈值并且所述方差是否超过预设的相关系数方差阈值，若满足，则将所述最大相关系数对应的计算区域作为最优相关区域，执行种子点追踪模块；

所述种子点追踪模块，其配置为基于所述第二位置确定所述最优相关区域内的同位置数据点，基于此，利用预设的搜索范围，确定当前道的种子点，以追踪下一道层位。

10.根据权利要求9所述的系统，其特征在于，所述系统，还包括：种子点检验模块，其配置为按照预设的特征点检验条件，对所述当前道的种子点进行有效性验证，其中，所述种子点检验模块，具备如下单元：

种子点深度差判断单元，其配置为判断所述当前道的种子点的层位深度数据与所述前一道种子点的层位深度数据的差值是否小于预设的深度差值阈值，若满足，则执行种子点深度范围判断单元；

所述种子点深度范围判断单元，其配置为判断所述当前道的种子点的层位深度数据是否处于预设的深度范围阈值之内，若处于，则执行极值搜索单元；

所述极值搜索单元，其配置为确定以所述当前道的种子点为中心点且以所述搜索范围为半径所对应的区域内的最大振幅极值点和最小振幅极值点，并计算两个极值点的第一极值深度差和第一极值振幅差；

极值条件诊断单元，其配置为判断所述第一极值深度差是否大于预设的极值深度差阈值，并且判断所述第一极值振幅差是否处于极值幅度差范围阈值之内，若满足，则确定当前种子点有效。

Images