煤层数字高程模型构建与动态修正方法

李娟莉，杜文勇，谢嘉成，王学文

(1.太原理工大学 机械与运载工程学院，山西 太原 030024;2.煤矿综采装备山西省重点实验室，山西 太原 030024)

0 引 言

地质信息透明化水平是煤矿智能化领域无人开采技术取得突破的重要基础，如何为煤矿开采提供高精度煤层三维地质模型是无人开采技术发展过程中亟待解决的瓶颈问题[1-3]。在国家能源局、国家矿山安全监察局印发的《煤矿智能化建设指南(2021年版)》[4]中明确将智能地质保障系统作为井工煤矿智能化系统基本建设内容之一，井工煤矿中煤层地质精确建模技术作为其中需重点研究关键技术，其重要性不言而喻。

DEM是基于GIS的数字地形分析的理论、方法与应用[5]，其能够充分表达地质体高程变化特征，具有二次开发与应用方便、模型数据储存与变换容易等优点[6]。与构建地表地质体DEM不同，限于当前地下煤层勘测技术水平，一次性获取全部所需的建模数据源存在高成本、高难度、低效率等难以克服的困难[7-9]。有学者结合煤矿生产特点提出“一次构建，逐步修正”的煤层DEM构建思路，即将煤层开采过程中暴露出的地质信息逐步引入到原有煤层地质模型中，以提高煤层模型的精度。贾庆仁等[10]提出在开采过程中获得新的地质数据后，首先确定新数据对原有模型的影响范围，并以此为基础对原有模型进行动态修正；刘万里等[11]提出利用上一回采阶段新揭露的煤层信息，基于克里格空间曲面插值法动态修正下一回采阶段煤层地质模型；李鹏等[12]提出融合地质、钻探、物探数据构建回采工作面煤层地质模型等。

虽然新思路下的建模方法能够显著提升煤层DEM建模精度，但是目前这些方法在实际生产中依然难以得到广泛应用[13-14]。笔者认为，究其原因是这些建模方法没有充分考虑地下煤层地质建模中数据来源特点，即获取采区中部的地质数据是十分困难的，而开采前往往能够依托巷道获取到可靠的、精确的、集中在采区边缘的煤层地质数据。现有建模方法在建立初始煤层DEM时过于依赖采区中部的煤层地质数据[15-17]，这无疑大幅增加了煤层DEM建模实现难度。

为此，笔者提出一种煤层数字高程模型构建与动态修正方法。在构建初始煤层DEM时，基于双轨扫掠、加权融合原理最大限度利用采区边缘的煤层地质数据；在向模型中引入采区内分散钻孔勘测所得的煤层地质数据时，通过构造残差曲面的方法实现煤层DEM宏观形态修正；在向模型中引入从综采工作面勘测得到的连续煤层地质数据时，通过分段建模，局部覆盖的方法实现对煤层DEM动态精细修正，为提升煤矿地质信息透明化水平提供了切实可行的途径。

1 总体技术架构

在煤层开采前获取充分的煤层地质数据，一次性构建精准的煤层DEM是较为困难的，因此，随着采煤工作面的推进，需要将各类新获取的煤层地质数据源源不断地引入煤层DEM中，以实现煤层DEM的动态修正[18-19]。笔者针对这一过程，提出基于采区边缘煤层地质数据的初始煤层DEM构建方法、引入采区内分散煤层地质数据修正采区煤层DEM方法及引入工作面上连续煤层地质数据动态修正采区煤层DEM方法，总体技术架构如图1所示。

图1 总体技术架构Fig.1 Overall technical architecture

1)基于采区边缘煤层地质数据构建初始煤层DEM。基于采区四周巷道中获取的煤层地质数据，通过沿不同方向分组双轨扫掠的方法得到多个煤层DEM，然后依据每个DEM在不同节点处高程值的可靠性，通过加权融合方法对所有由单方向扫掠构建的模型进行融合，最终得到总体上较为可靠的初始煤层DEM。

2)引入采区内分散煤层地质数据修正煤层DEM。引入通过钻孔勘测所得的一些采区内煤层地质数据，在原有模型的基础上计算煤层DEM中相应节点处高程误差值，然后在这些误差值的基础上通过高次拟合的方法得到残差曲面，最后利用残差曲面对DEM中节点进行整体性偏移，最终实现煤层DEM宏观形态修正。

3)引入工作面上连续煤层地质数据修正煤层DEM。引入某一工作面上一系列连续的煤层地质数据，继而得到煤层DEM中相应一系列节点处的高程值，然后在该工作面处产生新的模型分段面，将分段面两侧本属于同一区域的DEM分为不同的区域，在各自区域中重构更加精准的区域性煤层DEM，最终用这些新构建的区域性煤层DEM对原有煤层DEM进行局部覆盖，实现煤层DEM的动态精细修正。

2 构建初始煤层DEM

2.1 建模基础

如图2所示，在采区内适当位置处，选取巷道方向或垂直方向等容易确定的方向作为坐标轴方向，建立XY面与地面平行的三轴直角坐标系XYZ。将此坐标系所对应的三轴与坐标原点作为真实地质体与DEM模型的匹配基准，DEM模型中节点坐标便可基于这一基准与真实地质体中的节点相互映射。将采区煤层在XY面上的矩形投影规则网格化后，依托这个网格建立煤层DEM，其煤层顶板与底板面中的高程信息以数字高程矩阵H的形式储存：

图2 坐标系建立示意Fig.2 Schematic of coordinate system establishment

(1)

式中，m和n分别为网格中交叉点的总行数和总列数；tij为煤层顶板面上与网格中第i行j列交叉点相垂直对应节点的高程值；bij为煤层底板上与网格中第i行j列交叉点相垂直对应点的高程值，m。

随着网格节点密度的提升，网格中所有节点的水平坐标与数字高程矩阵H中能够包含更多的煤层地质信息，煤层DEM精度也会得到相应的提升。

根据巷道与煤层相对位置关系的不同分别采用直接测量法、向上或向下钻孔法，可以从巷道中获得煤层顶板面与煤层底板面上延巷道分布的一系列节点坐标[20]，如图3所示。用三次样条插值的方法对获得的边缘节点进行插值处理后，便可以得到沿巷道方向延伸的曲线，如图4所示，这些曲线能够很好地反映煤层顶底板面在采区边缘的变化情况，并且容易得到数字高程矩阵H中所有处于边缘位置元素的值。

1—巷道顶部；2—巷道底部；3—待开采煤层；4—向上钻孔勘测点；5—直接测量勘测测点；6—向下钻孔勘测点

图4 煤层DEM中的边缘节点Fig.4 Edge nodes in coal seam DEM

2.2 煤层DEM构建

以求解采区煤层底板面DEM所对应的数字高程矩阵B为例，首先基于规则网格以及已有的煤层数据构建初始底板数字高程矩阵B0，B0中所有元素值未知的位置用0填充：

(2)

式中，bij为煤层底板面上与网格中第i行j列交叉点相垂直对应节点的高程。

如图5所示，从AB出发以AD与BC为导轨进行双轨扫掠，得到数字高程矩阵BAB的过程如下：

图5 双轨扫掠示意Fig.5 Schematic of dual-track sweep

1) 计算矩阵CAB

设：

(3)

(4)

则：

(5)

2)基于CAB将B0补充完整得到BAB

设：

(6)

式中：cij为CAB中第i行j列元素。

(7)

同理，可以分别求解从BC、CD、AD出发，沿对应不同方向进行双轨扫掠后所构建的矩阵BBC、BCD、BAD。

由煤层变化的连续性易知，在上述方法构建数字高程矩阵时，距扫掠起始线段越近的节点处，所求解出的高程值可靠性越高。因此可以通过距离平方反比加权[21]的方式得到相对可靠的数字高程矩阵B。从而得出矩阵BAB所有元素在加权中对应的权重参数矩阵wAB：

式中：dij为网格中第i行j列交叉点与扫掠起始线的欧氏距离。

同理，可分别求出BBC、BCD、BAD所有元素在加权中对应的权重参数矩阵wBC、wCD、wAD。设w1ij、w2ij、w3ij、w4ij分别为wAB、wBC、wCD、wAD中第i行j列元素，则可得到数字高程矩阵B中第i行j列元素Bij：

(9)

式中，BABij、BBCij、BCDij、BADij分别为数字高程矩阵BAB、BBC、BCD、BAD中第i行j列元素。

从而可以得到与采区煤层顶板面DEM所对应的数字高程矩阵T，进而得到与整个煤层DEM相对应的数字高程矩阵H：

(10)

3 动态修正煤层DEM

3.1 引入采区内分散煤层地质数据修正煤层DEM

上述构造煤层DEM的方法充分利用了采区边缘煤层地质信息，以下提出一种引入采区内分散煤层地质数据修正煤层DEM的方法，其数据通常来源于自上而下的、少量的、分散的钻孔勘测。钻孔的位置则需结合实际勘测条件在初始煤层DEM所包含的范围内选取，位置应尽量均匀分散，经试验发现选在煤层可能存在“谷”或“脊”的地方钻孔能够对模型的宏观走势修正起到更好的效果。

将DEM中某一节点的实际高程值与理论高程值之差作为ΔZ值，再结合该节点XY坐标，便可得到一系列含有DEM残差信息的点(Xi，Yi，ΔZi)。再将一部分位于采区边缘含有DEM残差信息的点(Xi，Yi，0)引入，便可得出所有参与残差曲面拟合的点。

根据采区内煤层变化的复杂程度，选取3～5次的多项式函数作为拟合函数(初始煤层模型越复杂则选择次数越高的多项式函数)，将所有残差点进行拟合得到具体的残差曲面函数：

(11)

式中，γ为拟合次数，3≤γ≤5；Ci为待定系数。

图6 边缘收敛后的残差曲面函数图像Fig.6 Image of residual surface function after edge convergence

(12)

rij=R(xij,yij)kij2

(13)

(14)

式中，xij为网格中第i行j列交叉点所对应的X坐标值；yij为网格中第i行j列交叉点所对应的Y坐标值。

最后，利用残差矩阵R修正数字高程矩阵H：

(15)

式中，RT为煤层顶板所对应的残差矩阵；RB为煤层底板所对应的残差矩阵。

3.2 引入工作面上连续煤层地质数据修正煤层DEM

随着煤层开采的进行，工作面逐步向采区深处推进，从暴露的煤层断面上能够更加容易获取到高质量的煤层地质数据，采用与之前处理巷道中煤层地质信息类似的方法，用3次样条差插值的方法可以得到工作面处完整的煤层顶底板曲线，利用已知煤层地质信息的工作面作为分界面可以把采区煤层分为多个部分。对每个部分采用双轨扫掠，采用加权融合的方法对模型进行重构，再对原有模型进行局部覆盖便可以实现煤层DEM动态精细修正。

4 试验结果及分析

为验证本文提出的煤层DEM的构建及修正方法的正确性和有效性，以某矿综放工作面为例，进行了煤层DEM模型的构建及动态修正试验。

4.1 综采工作面煤层概况

试验针对某矿综放工作面采区进行模拟建模，该工作面设计可采走向长度691.2 m，倾向长度180 m，煤层平均厚度6.67 m，直接顶为泥岩，基本顶为K2石灰岩，底板为铝质泥岩，煤层倾角2°～10°，平均倾角6°。当前该工作面采区已开采完成且保存了相对完整的采煤机工作记录(包括截割轨迹与采煤机位姿等)与钻探勘测记录。

在进行模拟建模试验之前，笔者依据采煤机工作记录与钻探勘测记录对煤层地质体进行了复原，以得到模拟地质体。在模拟地质体上进行虚拟勘测与建模，通过比较建模结果与模拟地质体的差异便可以实现对DEM建模方式的检验，图7为通过数据反演得到的模拟地质体(煤层顶板与煤层底板面模型)。

图7 模拟地质体Fig.7 Simulated geologic body

4.2 构建与修正煤层DEM

以采区边缘地质数据为已知条件，用本文第2节中提出的方法对煤层进行DEM建模，得到沿不同方向扫掠构建的煤层DEM，最终经加权融合后所构建的模型如图8所示。

图8 加权融合模型Fig.8 Weighted fusion model

选取垂直方向上的模型节点RMSE与MAE作为模型评价标准，对沿不同方向扫掠构建的煤层DEM进行误差分析，最终得到结果见表1。

表1中RMSE与MAE的计算方法如下：

表1 初始煤层DEM误差

(16)

(17)

由此可见，沿不同巷道方向进行扫掠建模将对应得到不同精度的煤层DEM，将多个单方向扫掠模型加权融合后，所得模型的精度优于大多数单方向扫掠模型。

按照第3节中所提出的方法，引入图9中红点所示位置处一系列采区内分散煤层地质数据对初始煤层DEM进行修正后，模型整体RMSE由0.239 m下降至0.135 m，MAE由0.180 m下降至0.114 m，由此可见，此种修正方法在本次试验中起到了很好的修正效果。

图9 钻孔勘测位置选择示意Fig.9 Schematic of selecting locations for drilling holes

采用本文第3节中所提出的方法，通过引入工作面上两处连续煤层地质数据，将原有模型分为3个区域分别建模并覆盖修正原有煤层DEM后(图10)，不同区域内煤层DEM误差变化见表2。

图10 分区重构法修正模型示意Fig.10 Schematic of modified model by zonal reconstruction method

表2 引入连续煤层地质数据修正模型前后煤层DEM误差

由此可见，该方法能够大幅降低数据源附近范围内模型误差，能够有效提高煤层DEM精度，通过反复引入逐步深入的工作面地质数据，反复分区重构模型便可实现煤层DEM的动态精细修正。

5 结 论

1)提出基于采区边缘煤层地质数据，结合双轨扫掠与加权融合构建初始煤层DEM方法。在模拟建模试验中，模型均方根误差与平均绝对误差分别达到0.25 m与0.2 m以下。

2)提出引入采区内离散煤层地质数据，在原有模型的基础上计算煤层DEM中相应节点处高程误差值，并基于这些误差值通过拟合的方法得到残差曲面，最后利用残差曲面对DEM中节点进行整体偏移的煤层DEM修正方法。在模拟建模试验中，通过该方法对煤层DEM进行修正后，模型均方根误差与平均绝对误差由0.25 m与0.2 m分别降至0.14 m与0.12 m。

3)提出引入工作面上连续煤层地质数据，在该工作面处产生新的模型分段面，将分段面两侧分为不同区域，在各自区域中重新构建区域煤层DEM，再用新构建的区域DEM覆盖原有模型的煤层DEM修正方法。在模拟建模试验中，通过该方法对煤层DEM进行修正后，模型局部均方根误差与平均绝对误差均达到0.12 m以下。最终模型精度水平能够为无人工作面的实现提供可靠的地质信息保障基础。