段志鑫，吴 侃，宋月婷，白志辉，李 亮，周大伟

(1.中国矿业大学 环境与测绘学院，江苏 徐州 221116；2.武汉大学 遥感信息工程学院，湖北 武汉 430079； 3.冀中能源峰峰集团有限公司，河北 邯郸 056000)

相似材料模拟实验是研究岩体内部压力与岩层运移规律的重要手段，由于其具有易于实现、观测周期短、实验结果直观等优点，被广泛地应用于岩土工程和矿山开采等领域[1]。在矿山开采沉陷研究中，相似材料模拟实验的观测内容主要包括模拟岩层的位移和裂隙分布情况。目前，相似材料模型岩层位移监测的方法比较成熟多样，主要包括百分表法[2]、灯泡透镜法[3]、全站仪法[4]、近景摄影测量法[5-6]、光纤光栅传感器监测法[7]等。但是对于相似材料模型在模拟推进过程中裂隙发育和破裂形态信息，例如裂隙的发育高度、离层尺寸及位置等，由于观测手段和数据处理方法的限制，研究较少，目前此类信息的获取方法仍以手工量测为主，不仅效率低，且精度较差。许家林等对基于图像的相似材料模型采动裂隙提取与量化方法进行了大量研究，并开发了对应的模型裂隙处理程序FIMAGE，该方法可以较好地实现对背景单一的相似模型裂隙检测，但是未考虑图像畸变的影响[8-9]；张云鹏利用白光扫描仪获取模型的点云并通过邻域法提取模型破裂边界，该方法数据处理效率较低，无法满足实际应用需求[10]；陈朋等将相似模型表面的三维点云数据转换成深度图像，并通过边缘检测算子提取模型破裂边界，该方法观测成本高，且对细小裂缝的提取状况不佳[11]。

视觉测量技术是近些年工业检测中最常用的手段之一，具有简单、高效、成本低廉等优点[12]，在道路病害监测[13]、工件缺损检测[14]、结构体测量[15]等领域得到了广泛地应用。针对目前相似材料模型裂隙检测及定位方法自动化程度低、精度较差等问题，本文拟从二维数码图像出发，对基于视觉测量技术的相似材料模型裂缝提取和定位方法进行研究，为研究相似模拟实验中模型覆岩层的破坏和“三带”发育规律提供技术支持。

1 模型裂隙检测

原始模型图像仍然为三通道彩色图像，直接基于该图像进行裂缝提取过程比较复杂，因此需要先将模型图像转换成灰度图像，并进行图像预处理和滤波等一系列的处理流程，得到完整的模型裂隙区域信息。

1.1 图像预处理

为了增强图像裂缝区的显著性，采用分段线性灰度变换的方法对图像进行增强处理。该方法先对图像灰度进行分段，然后在各个灰度区间内，分别建立灰度线性映射关系，在增强兴趣区的同时达到抑制背景和噪声的目的。选取模型局部图像进行实验，原始灰度图如图1(a)所示，灰度变换后结果如图1(b)所示。

图1 图像预处理

灰度图像中冗余数据较多，不便于进行图像分析，因此需要将灰度图像二值化。采用Otsu算法进行阈值分割，该算法具有较好的适用性和分割效果，其核心是通过不断对阈值迭代调整，使图像前景与背景之间的灰度方差达到最大。设G(x，y)为二值图像，g(x，y)为灰度图，T0为Otsu法分割阈值，则分割数学过程如式(1)所示，二值图像如图1(c)所示。

(1)

1.2 多级连通域滤波

模型的二值图像中，仍含有大量形态多样，分布无规律的噪声，传统的数学滤波模型无法有效滤除此类噪声，因而，本文根据噪声的分布和形状特点，设计了一种多级连通域滤波方法。主要的步骤如下：

(1)搜索二值图像上每个非0像素的4-邻域，得到图像的若干个4-邻域连通域，建立离散像素之间的拓扑关系。

(2)基于连通域面积滤波。模型二值图像上存在着大量面积较小的斑块噪声，可通过比较统计，设置合理的面积阈值，滤除面积较小噪声区域。

(3)基于连通域外观比滤波。模型裂隙一般大体呈水平向分布，而连通域的外观比能反映其在水平方向上的细长程度，因此可利用该指数对纵向分布的噪声进行滤除。设gi(x，y)为二值图像的连通域，n为连通域个数，其目标外接矩形的长宽分别为Li，Wi，该连通域外观比为ri，设置滤波阈值为Tr，则滤波过程如式(2)，(3)所示。

ri=Li/Wi

(2)

(3)

(4)基于连通域形状因子滤波。由于网格墨迹及材料反光，二值图像上会存在边缘复杂的噪声，而裂隙一般分布紧凑且边缘相对规则，因此可以利用连通域的形状因子进行滤波。首先获取每个连通域gi(x，y)的周长Bi和面积Si，并根据公式(4)计算形状因子Fi，设定合理的阈值Tf，则滤波的数学过程如式(5)所示。

(4)

(5)

(5)特殊噪声滤除。相似模型表面布设有大量圆形摄影测量标志点，经过相机透镜中心投影后变形为椭圆，当这些标志点和裂隙连在一起时，会难以直接通过上述三步滤波方法滤除。在尽量保证正直拍摄的情况下，投影椭圆的扁率接近0，因而，可以通过近似圆形的形状因子指标，即Fi≈1，寻找出独立的标志点连通域，然后以此为连通域滤波模板，基于相似度判断，将混杂在裂隙区的标志点识别并滤除。

2 标志点识别与定位

2.1 标志点分类

本次试验中采用的标志点包括12位编码点和非编码点。其中12位编码点主要由白色圆形标志点和12等分的同心圆环编码块组成，白色圆点域和编码圆环内外径之比为3∶7∶12，每一个30°编码块取黑色或白色，分别对应0和1，如图2(a)所示，因此12个编码块按照一定顺序排列便可得到一个12位二进制数，将该二进制数进行移位运算，以其对应的最小十进制数值作为编码值，12位编码标志点通过编码块组合可以包含351个编码方案。非编码点由中心白色圆域和外围黑色同心圆环构成，其内外径之比为1∶2，如图2(b)所示。

图2 三维视觉测量标志点类型

2.2 标志点中心定位

为了确定裂缝的空间位置，需要对标志点中心进行精确定位。本文根据标志点的几何特点，设计了一种半径约束下最小二乘椭圆拟合的中心定位方法，该方法的主要思路为：

(1)首先基于连通域面积滤波的方法滤除二值图像中面积较大或者较小的噪声。

(2)根据实验中圆形标志点变形的特点，推导确定标志点连通域外观形状因子的范围为0.955

(3)利用形心法计算每个标志点连通域的中心坐标(u0，v0)，并计算连通域边界点到中心的平均距离r。

(4)利用canny算子对模型灰度图像进行边缘检测，并根据圆形标志点内外径之比为1∶2，确定约束半径为3/2r，滤除连通域中心约束半径以外的点，以此得到每个标志点的内圆边缘点(xi，yi)。

(5)利用公式(6)对步骤(4)中得到的边缘点进行最小二乘法椭圆拟合，得到椭圆的圆心坐标(xc，yc)，即为标志点的中心坐标，如图3(a)所示。

图3 标志点识别与定位

(6)

式中，f为目标函数；A，B，C，D，E为椭圆方程的5个系数；N为边缘点个数。

2.3 编码标志点识别

编码标志点是进行坐标空间转换的关键，因此对其进行准确识别具有重要的意义。根据编码标志点的几何和变形特点，设计了一种有效的编码点识别方法，该方法的主要步骤包括：

(1)在1.2节连通域面积滤波的基础上，根据编码标志点编码块的几何特点及变形特征，根据式(4)计算了编码块的外观形状因子系数变化范围为1.186

(2)根据2.2节中所求取各个标志点白色圆域的圆心坐标及椭圆方程，以一定的扩充系数σ对每个拟合椭圆进行扩充，并获取扩充椭圆的边缘点坐标及其灰度值，根据编码点的设计比例，确定椭圆扩充系数为σ=3。

(3)为了消除投影变形的影响，根据式(7)对扩充椭圆进行仿射变换，将其变换成圆，如图3(c)所示。

(7)

式中，(x，y)和(x’，y’)分别为变换前后边缘点坐标；a，b分别为扩充椭圆的长短半轴；α为椭圆倾角。

(4)按照顺时针顺序对圆边缘进行采样，寻找灰度分别为0和1的相邻点，取其中点作为分界点，并计算前后相邻两个分界点之间的圆心角及该角内包含的12位编码块数量。

(5)根据圆心角两侧灰度变化规律，确定夹角内编码块的灰度为0或1，据此得到一个12位的二进制编码序列，由于采样起点不同，可能得到不同的二进制序列，因而需要对非0二进制序列进行移位运算，取其对应的最小十进制数值作为该编码标志点的编码序号。由于本实验采用的编码规则与XJTUDP近景摄影测量系统不同，因此编码值也不相同，但仍是一一对应关系。点位识别效果如图3(d)所示。

3 模型裂隙定位

3.1 畸变校正

普通数码相机在成像的过程中都会存在一定程度的几何失真，即相机畸变，畸变的存在会使图像上裂隙轮廓点及标志点的真实位置产生误差，因此在进行坐标空间转换之前需要进行畸变校正，建立像点的坐标改正模型，如式(8)所示[16]。

(8)

式中，(Δx，Δy)为像点坐标改正值；(x，y)为像点实际量测坐标；(x0，y0)为像主点坐标；k1，k2，k3为径向畸变系数；p1，p2为偏心畸变系数；r为像点的畸变半径。相机畸变系数可通过相机标定得到。

3.2 直接线性变换

经过本文方法检测获取的裂隙轮廓坐标和标志点坐标均为像方坐标，无法直接用于分析研究，需要将像方坐标转换为物方坐标。为了确定裂隙在模型表面的实际位置，首先对模型裂隙图像上的编码标志点进行识别和中心定位，然后利用近景摄影测量系统测出标志点的物方坐标，进而建立裂隙坐标空间的直接线性变换模型。

直接线性变换(DLT)不需要事先测定相机的内外方位元素，而可以直接建立点位的物方坐标与像方坐标之间的线性变换关系，因此在非量测相机摄影测量中应用较多[17]。考虑本文的监测目标为二维相似材料模型，模型裂隙与标志点均布设在同一平面，因此只需要测量不少于4个编码标志点的物方坐标，便可以求解出二维DLT单应变换矩阵及其变换参数，进而通过变换模型将裂隙轮廓坐标转换到物方空间，具体的变换模型如式(9)所示。

(9)

式中，(xu，yu)为编码标志点去畸变像方坐标；(X，Y)为对应标志点物方坐标；L1～L8为二维DLT变换的8个参数。

4 实验结果及分析

4.1 模型裂隙检测

本次试验对模拟西部某矿区开采的相似材料模型进行裂隙检测与定位研究，模型的几何尺寸为3000mm×1500mm×300mm，正面布设有大量编码标志点与非编码标志点，试验采用普通数码相机对模型裂隙进行观测。选取2期(11.09，11.18)数据进行试验验证，试验效果如图4所示，其中图4(a)为原始裂隙图像，图4(b)为裂隙二值图像，图4(c)为裂隙检测与定位结果。

图4 模型裂隙检测结果

由上述实验结果可以看出，利用本文方法可以较好地检测并提取出模型的裂隙，检测结果与原图中裂缝位置基本一致。但是，受图像分辨率以及光照条件的影响，对于一些细小裂隙，会存在一定的“漏检”现象。

4.2 模型裂隙定位结果分析

在视觉测量中，受摄影方式、图像分辨率、标志点像方定位精度、标志点实测精度、裂隙像方坐标精度等因素的影响，相似材料模型裂隙定位会存在一定的误差。受限于现有的技术手段，实验中难以直接对裂隙的位置精确实测并进行精度分析，而在模型裂隙提取准确的前提下，其定位精度等于周围非编码点的定位精度，因此本文通过对裂隙周围的非编码标志点进行定位误差分析，来反映裂隙的位置精度。由于XJTUDP近景摄影测量系统单点测量精度可达0.03～0.1mm[18]，完全可以满足矿山相似材料模型监测的要求，本文将XJTUDP系统测量的点位坐标值作为实测值，与本文方法的定位结果进行对比分析。实验以119，15，459和349共4个编码点作为控制点，并随机选取模型图像中50个点进行定位误差分析，得到统计曲线如图5所示。

图5 模型非编码点定位误差

综合统计数据，按照公式(10)计算得到检查点在X，Y方向上坐标值相对于XJTUDP系统的中误差分别为：mx=0.840mm，my=0.714mm。XJTUDP系统点位测量误差最大为0.1mm，考虑该系统自身误差，本文方法在X，Y方向上定位误差小于0.940mm，完全可以满足一般相似材料模拟实验中裂隙定位的要求。利用本文方法对11.18期模型图像中的裂隙发育最高点定位，结果如表1所示，其中裂隙编号C1，C2，C3等见图4(c)(11.18)。

(10)

4.3 裂隙检测结果应用

由4.1可知，通过本文的模型裂隙检测方法，可以得到只包含裂隙的二值图像，而每一条裂隙都对应二值图像上的一个连通域，对每个连通域进行边界遍历，便可得到每条裂隙所对应的轮廓点像方坐标集合；然后利用本文2.2和2.3中方法对模型图像上摄影测量编码点进行识别与中心坐标检测，并根据XJTUDP系统测得的编码点实际坐标，解算该张图像的直接线性变换矩阵；最后通过 DLT变换将裂隙轮廓点的像方坐标转换到真实物方坐标空间，此时裂隙的轮廓尺寸、位置等信息均为实际值。为了便于模型裂隙的观察与量测分析，可以将裂隙轮廓点和开切眼的真实坐标导入到CAD中，并通过闭合曲线对每一条裂隙的轮廓点连接，便可得到如图6所示的三期裂隙分布图。后续通过CAD查询功能可方便快捷地获取每条裂隙的尺寸、面积、分布角度等信息，为研究岩层破坏空间、“三带”发育规律等研究提供直观的依据。

图6 模型裂隙分布

5 结 论

(1)对基于视觉测量的相似材料模型裂隙检测与定位方法进行了研究，通过模型图像预处理和多级连通域滤波实现了对裂隙的准确检测与提取；设计了有效的编码点识别与中心定位方法，并基于DLT变换实现了对模型裂隙位置的精确定位。

(2)试验研究表明，在裂隙检测准确的前提下，本文方法在X，Y方向上的裂隙定位误差小于0.940mm，可以较为准确高效地确定模型的裂隙发育位置；同时，本文方法具有成本低廉、可实时连续观测等特点，能够为相似材料模型中岩层破坏和 “三带”发育规律研究提供数据支持，对于矿山开采沉陷等的物理模拟研究具有一定的实用价值。