基于高频雷达波的煤岩层位识别与追踪方法研究

许献磊,王一丹,朱鹏桥,马 正

(1.中国矿业大学(北京) 煤炭资源与安全开采国家重点实验室，北京 100083；2.中国矿业大学(北京) 地球科学与测绘工程学院，北京 100083)

0 引 言

因此，笔者提出了一种基于高频雷达波的煤岩层位识别方法和以“窗口算子”为核心的追踪算法，根据“空气-煤-岩”层的反射波特征，提出了煤岩层位种子点识别与定位算法，然后以三级“窗口算子”为核心实现煤岩界面的实时跟踪。首先阐述煤岩层位探测工作原理，接着介绍煤岩层位的识别与定位方法，讲述煤岩层位追踪算法并进行物理模型试验，随后进行矿井煤岩探测试验，最后对全文进行总结。

1 煤岩层位探测工作原理

基于高频雷达波的煤岩层位探测原理是通过向煤层发射脉冲形式的高频宽带电磁波，电磁波在煤层传播过程中遇到与煤存在电性差异的岩石层时发生反射，反射信号到达煤层上表面时由Rx(接收天线)感应接收，并经过超宽带接收机的转换、放大、量化后形成电磁回波(图1)。根据电磁回波极性、幅度、走时等参数，采用相应的算法寻找到煤层与岩石层的分界面，从而计算出煤层的厚度，其中Tx为发射天线。

图1 超宽带电磁波煤岩识别原理示意Fig.1 Schematic of the principle of ultra-wideband electromagnetic wave coal and rock identification

地质雷达所发射电磁波的传播取决于介质的电性，介质的电性主要有电导率μ和介电常数ε，前者主要影响电磁波的穿透(探测)深度，在电导率适中的情况下，后者决定电磁波在该物体中的传播速度，因此，所谓电性界面也就是电磁波传播的速度界面。不同的地质体(物体)具有不同的电性，因此，在不同电性的地质体的分界面上，都会产生能量较强的回波。

2 煤岩层位的识别与定位

2.1 雷达波层位特征分析

要实现对煤岩层的精准识别，首先要对雷达波层位特征进行分析，准确定位层位点。首先搭建理想的煤岩结构模型，如图2所示，模型煤层厚度0.2 m，岩石层厚度为0.5 m，通过GPRMax2D进行仿真模拟(图2)，单道波形如图2b所示。

图2 数值模型及仿真结果Fig.2 Numerical model and simulation results

通过对比三个界面的波形图可以发现：A为直达波，直达波是从探地雷达的发射天线出发不经过反射和折射直接到达接收天线的电磁波，其特点是信号振幅值最大；B为“空气-煤”层位回波，空气-煤界面振幅值仅次于直达波；C为“煤-岩”层位回波，因受电磁波衰减的影响，煤-岩界面振幅值最小。

预防措施：在实际操作中，黏锤现象完全可以采取一定措施进行规避，比如在冲孔工作暂停时将锤头提到起至离开孔底一定高度就可以有效减少黏锤发生的概率，考虑到冲浆清渣的效果，此高度以20-100cm为宜。如在土质黏性较大的地质区域，可以填入一定量废钢丝绳或者石块才减少孔底土壤对锤头的粘黏力。

2.2 “煤-岩”层的定位

2.2.1 “煤-岩”层的初步定位

1)直达波位置及其样点数N0。在预处理后的雷达数据中选取任一道数据序列，并确定所选数据序列的道号；其中，所选数据序列的道号记为T0(图3)。从时间窗口0点开始第一个振幅起跳点位置即为直达波的位置，其样点数记为N0。

图3 第T0道雷达回波Fig.3 T0 radar echo

2)“空气-煤”层位的位置及其样点数N1。应用已经建立的煤岩分类模型实现空气-煤层位的准确识别，进而确定空气-煤层位的位置及其样点数N1。

N1=[(2H0(N-1))/V0T]+N0

(1)

式中：H0为雷达天线悬空高度；T为采样时窗；N为采样率；V0为电磁波在空气中的传播速度。

3)“煤-岩”层位的位置及其样点数N2。对于种子道数据序列集(即预处理后的雷达数据)，煤-岩层位的位置位于空气-煤层位的下面，从空气-煤层位的样点数N1开始，应用已经建立的煤岩分类模型实现煤-岩层位的准确识别，进而确定煤-岩层位的位置及其样点数N2。

N2=[(2H1(N-1))/V1T]+N1

(2)

式中：H1为煤层高度；Vt为电磁波在煤层中的传播速度。

根据煤-岩层位的样点数和所选数据序列的道号，确定种子层位点的坐标信息，即预处理后的雷达数据中所对应的位置点(T0，N2)为煤岩层初步种子层位点。

2.2.2 “煤-岩”层的精确定位

上述确定的N0、N1和N2相当于是煤岩层位探测的单道波形中主瓣极值点所在采样点序号(图4)。而实际上主瓣极值点附近的点都可以作为我们要确定的层位点，为实现对煤岩层位的精确定位(即电磁反射波的位置)搭建天线不同悬空高度和煤层厚度的误差分析试验。模型设置参数如下：二维网格单元尺寸为0.005 m×0.005 m，时窗为15 ns，子波中心频率为1.2 GHz。在数值模型中，煤层相对介电常数ε1为6，岩石相对介电常数ε2为9。将天线的悬空高度d0设置为20 cm并保持不变，煤层厚度d1变化范围为30～50 cm。道间距设置为0.02 m，利用GPRMax2D进行探测仿真(图5)。

图4 煤岩层位探测单道波形示意Fig.4 Schematic of single-channel waveform of coal and rock layer detection

图5 物理仿真模型和其雷达剖面Fig.5 Physical simulation model and its radar profile

设在单道波形上直达波、“空气-煤”“煤-岩”层位反射波的可取样点序号的范围分别为[i1,i2]、[j1,j2]和[k1,k2],i′、j′和k′为相对于其主瓣极值点号的样点距离，称为相对序号。悬空高度计算值的相对误差和煤层厚度的相对误差分别为Er0和Er1。

在d1=30、40、50 cm处的单道波形(图6)会各自对应一个相对误差二维图谱(图7)。在每个二维图谱中设定其阈值Erm，得到其对应的Er1

图6 不同煤层厚度处的单道波形图Fig.6 Single-channel waveform diagrams at different coal seam thicknesses

图7 二维图谱Fig.7 2D map

图8 相对序号i′和j″、j′和k′之间的关系Fig.8 Relationship between relative sequence numbers i and j″, j′and k′

根据上述关系可得：预处理后的雷达数据中如果选定N0为直达波位置的样点数，那么“空气-煤”种子层位点的位置为(T0，N1+2)，“煤-岩”种子层位点的位置为(T0，N2+2)。

2.3 煤层厚度计算

因复杂地质条件，一般情况下煤矿巷道顶板并非一条直线，在煤岩层位动态探测过程中，雷达天线悬空高度H0随着煤层的起伏会有一定的变化。雷达天线悬空高度H0通过以下方法计算：根据2.1节和2.2节中提到的矿井煤岩层位识别方法，可以确定雷达数据中任意一道单道波中直达波的位置及其对应样点数N0，以及“空气-煤”层位的位置及其样点数N1+2，根据“空气-煤”层位与直达波之间的时间差和电磁波在空气中的传播速度V0即可获取雷达天线在该道处的悬空高度H0，其计算公式为：

H0=(N1+2-N0)V0T/2(N-1)

(3)

“煤-岩”层位的位置及其样点数N2+2，根据“空气-煤”层位与“煤-岩”之间的时间差和电磁波在煤中的传播速度V1即可获取煤层厚度H1。其计算公式为：

H1=(N2-N1)V1T/2(N-1)

(4)

3 煤岩层位追踪算法

3.1 三级“窗口算子”

现有的层位追踪方法可实现针对简单的煤岩层位或者理想情况下煤岩层位的追踪，而对于干扰源尚未建立判断和纠错机制。当某个位置出现追踪错误时，导致后续错误的层位追踪结果，这种状况称之为“串层”。提出设立三级窗口算子，从而减少串层现象的发生。

“窗口算子”由道相关算法改进而来，是把参考道种子点作为中心，调整上下界线范围生成三级窗口，设x(m)是参考道的数据序列，而y(n)是被追踪道的数据序列，数据序列的长度是m和n1；m1、m2和m3分别为参考道3个数据序列窗口的长度，且n1>m3>m2>m1。则有：

(5)

式中，ρxy为相关系数；N为采样点数。

根据公式(5)进行相关性分析(图9)，分析步骤如下：

图9 “三级窗口”相关分析示意Fig.9 “Three-level window” correlation analysis

1)首先判断一级窗口：在被追踪道上选定一个序列，从它的上界线首个点起始，计算参考道与被追踪道相关性，至本窗口下界线末位点结束，直至计算到最大的相关系数值才停止，接下来对二级窗口做相似的操作，直到所计算的相关系数值最大为止(计为ρ2)，如果ρ1和ρ2所对应的点号相同，那么将该点位置信息作为被追踪道的层位坐标信息。

2)如果ρ1和ρ2所对应的点号不同，即要判断3级窗口，直到所计算的相关系数值最大为止(计为ρ3)，如果ρ3和ρ1、ρ2所对应的任一点号相同，那么就将该相同点位置信息作为被追踪道的层位坐标信息。

3)如果ρ3、ρ1和ρ2所对应的点号都不同，那就把与ρ3相应的位置信息当作被追踪道的层位坐标信息。

3.2 层位追踪

利用“窗口算子”，以种子层位点为初始参考中心，对预处理后的雷达数据进行层位追踪，确认全部被追踪道和起始参考道这两部分数据序列的层位地理坐标信息，起始参考道的数据序列是把种子层位点作为中心的数据序列。具体为：

计种子层位点为(T0，N2+2)，层位追踪是把(T0，N2+2)作为起始参考中心，使用道相关算法实时地追踪层位。其中，以(T0，N2+2)为参考道数据序列的中心，设x(m)、y(n)分别为参考道数据序列和被追踪道数据序列，m和n为数据序列的长度，且n>m。详细的操作流程如下：

1)以种子层位点为参考道数据序列的中心，将设定好的步长作为数据采择窗口，确认从参考道数据序列中提选的数据。

2)把与参考道数据序列相邻的数据序列作为被追踪道数据序列，利用“窗口算子”，根据式(5)计算参考道数据序列选取的数据与被追踪道数据序列的相关系数，并确定被追踪道数据序列的层位坐标。其中，横坐标值X为道数，纵坐标值Y为煤岩层位高度，即N2+2。

3)把参考道数据序列的中心更换为被追踪道数据序列的层位坐标，将参考道数据序列更新为被追踪道数据序列，返回确定参考道数据序列选取的数据步骤，循环迭代，直到确定出所有被追踪道数据序列的层位坐标信息。

3.3 物理模型试验

为验证本方法的有效性，在实验室搭建物理模型，模型剖面为直角梯形(图10)，其长为200 cm，短边宽22.5 cm，长边宽35 cm。模型内部分为上下两层，上层为岩石层，下层为煤层。在煤岩层位起伏平缓区选择5处进行打标，层位起伏明显区选择3处进行打标，测量并记录打标处的悬空高度和煤层实际厚度(表1)。

图10 室内煤岩模型Fig.10 Indoor coal and rock model

表1 悬空高度及煤层实际厚度

利用中国矿业大学(北京)自主研发的1.2 GHz空气耦合式煤岩识别天线，放置在煤岩物理模型下方，从左向右匀速进行探测，随天线的运动动态采集煤岩层位雷达数据。根据上述层位追踪算法，得到模型煤岩层位追踪结果(图11)。应用本文提出的煤岩层位追踪方法和煤层厚度计算方法获得打标位置点煤层厚度探测结果，并进行精度分析。

由图11a可以看出：利用上述算法能准确识别并追踪到煤岩层位，图11b中对打标点处的测量煤厚和实际煤厚进行误差分析，平均误差为±0.12 cm，最大误差百分比为4.76%，最小误差百分比为0.75%，平均误差百分比为2.18%。

图11 模型追踪结果及煤层厚度分析Fig.11 Model tracking results and coal seam thickness analysis

4 现场试验及分析

为进一步验证研究方法的效果，应用中心频率为1.2 GHz的空气耦合式煤岩识别天线在上湾煤矿12306工作面进行煤岩识别探测试验，共计采集了8组数据。为验证探测结果的准确性，在探测区域范围选取8个位置进行顶板开挖，准确量取煤层厚度值。数据处理结果如图12、图13、图14所示。

图12 彩色雷达图谱及追踪层位Fig.12 Color radar map and tracking horizon

图13 灰度雷达图谱及追踪层位Fig.13 Grayscale radar map and tracking horizon

图14 实际煤厚与探测煤厚对比Fig.14 Comparison of actual coal thickness and detected coal thickness

对所获取煤岩界面雷达数据进行识别追踪并与现场定点位置实测煤层厚度数据对照分析，上述煤岩识别与追踪算法得到的探测数据与实际数据之间的平均误差为±0.71 cm，平均误差百分比为3.53%，验证了该算法在实际工程中的可行性。

5 结 论

1)提出了一种基于高频空气耦合雷达的煤岩层位识别与追踪方法，高频空气耦合雷达的中心频率为1.2 GHz。物理模型试验和现场探测试验结果表明：高频空气耦合雷达及煤岩层识别方法，可在煤岩层连续的情况下实时动态获取到1 m范围内煤岩层位信息，探测精度达到厘米级，其中物理模型探测中平均误差为±0.12 cm，平均误差百分比为2.18%，矿井工作面探测平均误差值为±0.71 cm，平均误差百分比为3.53%。

2)研究成果在煤岩层存在明显电性差异的工作面可实现煤岩层位的快速探测、识别与追踪。下一步的工作重点是针对煤岩层中出现夹矸、不连续等情况完善煤岩层位追踪算法，同时研究在煤层顶板高度发生连续变化情况下，基于采煤机坐标系下的煤岩层位坐标转换方法。