摘要
裂隙参数是围岩稳定性分析的参考因素，为了实现围岩裂隙的可视化分析，该文提出一种基于可视化工具包 VTK 的围岩裂隙三维重建方法。首先将采集的隧道围岩图像使用 OpenCV 进行图像灰度化、降噪的预处理；其次运用 Otsu 法和形态学中的开运算进行图像分割处理，得到较为完整的围岩裂隙图；然后将 VTK 中传统的 Marching Cubes 算法采用中值法代替线性插值法计算等值点，并运用改进后的 Marching Cubes 算法对处理后的裂隙图像进行三维重建；最后将改进后的算法与原算法的三维模型重建速率进行对比。研究结果表明：① 将 VTK 应用于隧道围岩裂隙的三维重建，可以得到较完整的裂隙三维模型，实现局部围岩裂隙的三维重建；② 两种算法重建的裂隙三维模型效果无明显差别，但改进算法的三维模型重建速率提高了 18%。
Abstract
Fracture parameters are important reference factors for surrounding rock stability analysis. To realize the visualization analysis of surrounding rock fractures， this paper proposed a three-dimensional （3D） reconstruction method of surrounding rock fractures based on the visualization toolkit （VTK）. Firstly， the collected tunnel surrounding rock images were preprocessed by OpenCV for image grayscale and noise reduction. Secondly， the Otsu method and the open operation in morphology were used to segment the image， and the complete surrounding rock fracture map was obtained. Then， the median method was used to replace the linear interpolation method with the traditional Marching Cubes algorithm in VTK to calculate the equivalent point， and the improved Marching Cubes algorithm was used to carry out the 3D reconstruction of the processed fracture image. Finally， the 3D model reconstruction rates of the improved algorithm were compared with that of the original algorithm. The results show that ① by applying VTK to 3D reconstruction of tunnel surrounding rock fractures， a complete 3D model of fractures can be obtained to realize 3D reconstruction of local surrounding rock fractures；② there are no significant differences between the two algorithms in the 3D model of fractures， but the 3D model reconstruction rate of the improved algorithm is increased by 18%.
Keywords
0 引言
随着社会和科学技术的进步，中国交通运输领域得到了前所未有的发展。隧道是中国交通运输网的重要组成部分，然而在隧道施工过程中，常常由于隧道围岩受到扰动，导致隧道围岩稳定性降低，从而引发围岩体坍塌事故，造成了巨大的经济损失和负面影响。近些年越来越多的专家学者对隧道围岩稳定性进行深入研究，在隧道围岩稳定性的诸多影响因素中，围岩裂隙是重要的因素之一，因此对围岩裂隙进行可视化研究具有重要的意义[1-3]。
随着不同学科之间的交叉融合，计算机技术也越来越多地运用在土木工程领域。计算机视觉一直是计算机技术中比较热门的方向，其中的三维重建技术在近几年取得了巨大的发展，国内外的专家学者不断地尝试将三维重建技术与隧道工程联系起来，进行了多方面研究。彭守建等[4] 提出了一种基于 CGAL（计算几何算法库）的裂隙面三维重构方法，为裂隙面的几何特征量化提供了基础；谢壮等[5] 采用机器视觉的方法获取隧道掌子面的三维点云模型，以此来得到隧道掌子面的节理特征信息，对节理面进行量化描述；李赤谋等[6] 将深度学习与隧道工程相结合，提出了在 Unet 网络上的隧道掌子面三维重建，并以此来实现隧道围岩的快速分级； García-Luna 等[7] 提出了一种在 SFM（运动恢复结构）摄影测量技术下的隧道掌子面三维点云模型获取方法；Xue 等[8]将 SFM-DLT 运用在地铁隧道上，通过获得三维点云重建隧道来监测隧道的变形等定量特征；Zhao 等[9] 利用深度卷积生成对抗网络完成了对致密砂岩的三维数字重构，所得的三维模型与真实砂岩的孔隙大小，几何结构大致相同；赵红华等[10]通过自主编写的体绘制三维重构程序，建立了三维位移场，用于测量透明土在静压试验下的位移变形；李丽萍[11]基于双目立体视觉原理，构建了一个三维重建系统来监测隧道围岩变形，以减少隧道事故的发生。综上所述，三维重建技术已经在隧道工程领域得到了广泛的应用。
本文将三维重建技术运用于隧道围岩裂隙的可视化研究中，将采集的隧道围岩裂隙图片经过图像处理和特征提取，得到裂隙的主体部分。然后使用 VTK 可视化工具包，在改进的 MC 算法基础上，将经过处理后的围岩裂隙图像进行三维重建，得到围岩裂隙的三维模型，为隧道围岩稳定性的分析提供了另一种思路。
1 裂隙图像处理
1.1 灰度化处理
图像上的最小单元通常称为像素，而图像上的每一个像素的颜色数值常被称为灰度，灰度的数值范围一般为 0~255，其中 0 指黑色，颜色最深；255 指白色，颜色最浅[12]。
彩色图像通常是由 R、G、B 3 个通道组成，可用函数 f（ x，y）来描述，R、G、B 分别表示每个像素的 3 个分量。所以图像灰度化指图像上的 R、G、B 值相等，且图像的颜色只有黑、白、灰 3 种灰度图像。在进行图像处理时，通常都会对目标彩色图像进行灰度化，以得到它的灰度图像。图像灰度化通常采用以下几种方法：
（1）最大值灰度处理法。其原理是使 R、G、B 数值等于 3 个数值中最大的一个，公式表达为：
(1)
（2）分量值法。其原理是分别采用 R、G、B 3 种不同的灰度值作为图像的灰度，从而得到 3 张不同灰度的图像，公式表达为：
(2)
（3）平均值法。其原理是把图像的 R、G、B 3 个分量值相加，然后取平均值，公式表达为：
(3)
（4）加权平均法。其原理是分别在 3 个分量前面加上一个系数求和，其公式表达为：
(4)
式（4）为 OpenCV 开源库所采用的灰度权值。
由于本文在 python 的基础上采用 OpenCV 对目标图像进行处理，所以使用加权平均法对目标图像灰度化，经过处理的裂隙灰度图像和裂隙灰度直方图如图1、2所示。
图1 裂隙灰度示意图
Figure1 Fracture grayscale
由裂隙灰度直方图可知：该围岩裂隙图像中像素灰度值存在明显的峰值区域，且在峰值区域两侧分布较平缓，因此说明灰度图像中的特征区域相对背景区域更加明显，在后续的图像识别过程中效果更好。
图2 裂隙灰度直方图
Figure2 Fracture grayscale histogram
1.2 图像降噪处理
在采集和处理图像的过程中，由于设备因素、外界环境等影响，使采集到的图片受到噪声的干扰，导致图像质量下降，这也为接下来的研究带来了不利[13-14]，所以将采集的图像进行降噪处理也是图像处理中非常关键的步骤。在图像处理中常用的图像降噪方法有均值滤波、中值滤波、高斯滤波等方法。
均值滤波是一种典型的线性滤波算法，其主要原理是领域平均法，即用一部分图像区域的各个像素的平均值来替换整个图像的各个像素值。
高斯滤波是另一种线性滤波算法，其原理是将一个扫描窗口中确定区域内像素进行加权平均，再将得到的加权平均灰度值替换窗口中心的像素值。
中值滤波是一种非线性滤波算法，其原理是将当前像素点和周围领域的像素点进行排序，然后再取中间的那个像素点作为当前像素点的像素值。
经过对 3 种滤波的比较分析发现，中值滤波在图像降噪的过程中，几乎不会影响图像的原来形状，能够很好地达到去除噪声的效果，所以本文决定使用中值滤波对目标裂隙图像进行处理。经过处理后的裂隙灰度图像如图3所示。
2 裂隙图像特征提取
2.1 裂隙图像分割
在图像处理的过程中，需要提取出图像中的特征目标，由于图像中的特征目标与背景具有相对分明的灰度值。因此，得到图像的灰度直方图后，可以根据灰度直方图中灰度的分布情况，选择适当的方法对不同的灰度集合进行分割，这样就能分辨出图像中的特征目标区域与背景区域[15]。
图3 中值滤波处理后的裂隙图
Figure3 Fracture after median filtering
目前在计算机视觉领域最常用的图像分割方法为：区域分割方法、阈值分割方法、边缘检测分割方法等[16]。本文考虑图像灰度分布情况，选择了一种基于阈值的图像分割方法——最大类间方差法（Otsu 法）[17]。
最大类间方差法因其简单、快速、分割明确、适用性广等优点而被广泛采用。其基本思想为：假设图像的像素为 N，图像灰度范围是[0，L-1]，灰度级为 i∈[0，L-1]的像素值为 ni。概率 pi为：
(5)
并且：
(6)
C1 类的灰度值 ∈[0，T]，C2 类的灰度值 ∈[T+1， L-1]，因此图像的总均值、C1 类、C2 类的均值为 μ1、 μ2 和 μ3，则：
(7)
可知：
(8)
类间方差定义：
(9)
将式（8）代入式（9），得到：
(10)
上述 θ、μ3、μ2 都为阈值 T 的因变量，所以使用遍历 T 的方法就能求得最大类间方差的阈值 T。
在使用 Otsu 法进行图像分割时，需要确定一个阈值，并对图像灰度直方图中灰度的分布情况进行分析，得到大致的阈值范围，图4为不同阈值情况下的处理图像。
图4 不同阈值的裂隙分割示意图
Figure4 Fracture segmentation with different thresholds
经过对上述不同阈值图像比较观察之后，当阈值 T=50 时，处理后的图像效果较好，分割出的裂隙主体部分较完整，能较好地保留裂隙的形状信息，并且干扰因素少。所以本文决定使用阈值为 50 的分割图像作为后续的处理图像。
2.2 形态学处理
数学形态学是一门包含了几何积分、集合代数、拓扑论的数学理论学科。在计算机视觉领域，最常用的形态学运算有膨胀、腐蚀、开运算、闭运算等[18]。
腐蚀是采用滤波对图像内部进行处理，而膨胀是利用结构元素对图像的外部进行处理。
膨胀：假设输入图像为集合 A，结构元素为集合 B，则集合 A 被集合 B 膨胀可以表示为：
(11)
腐蚀：假设输入图像为集合 A，结构元素为集合 B，则集合 A 被集合 B 腐蚀可以表示为：
(12)
在形态学处理中，除了腐蚀和膨胀之外，还有两个互为对偶运算的运算，开运算和闭运算。开运算和闭运算在原理上即是膨胀和腐蚀的组合运算。
开运算：假设输入图像为集合 A，结构元素为集合 B，则 B 对 A 的开运算可以表示为：
(13)
闭运算：假设输入图像为集合 A，结构元素为集合 B，则 B 对 A 的闭运算可以表示为：
(14)
在对裂隙图像进行图像处理时，由于灰度的分布不均匀，导致处理后的裂隙图像不连续，有断开的痕迹。所以针对此问题，本文采用了数学形态学中的开运算，来对分割的裂隙图像进行处理，在 python 基础上调用 OpenCV 函数库，得到经过开运算处理后的图像，如图5所示。
图5 开运算处理后的裂隙示意图
Figure5 Fracture after opening operation
由图5可知：开运算在原来裂隙图像基础上将中间断裂部分连接起来，并且将多余干扰部分去除，使得裂隙图像更加完整，丰富了裂隙连接处细节，有效地提升了后续模型重建精度。
3 裂隙三维重建
VTK 可视化工具包是最近几年在计算机可视化领域中，得到科研学者广泛应用的一种科研和开发工具，里面包含了超过 500 多种图形图像和可视化算法，使用者可以在函数库中找到对应的目标函数进行二次开发和利用[19]。
3.1 传统的面绘制 MC 算法
在使用 VTK 绘制三维模型时，最常用的方法就是 Marching Cubes（MC）算法，它是通过对二维图像进行识别和处理，区分开建模目标和背景，从而进行三维重建[18]。MC 算法是一种基于体素级的面绘制重建方法，其主要的核心是“体元”。在图像中相邻的 8 个体素点构成的一个立方体就叫做体元。其结构如图6所示。
图6 体元示意图
Figure6 Voxel
在采用 MC 算法进行重建时，重建过程中所生成的三维曲面网格被称为等值面，面绘制也因此被称为等值面提取。等值面与每个体元棱边的交点，可以连接成三角片，三角片的连接类型有 28 种[20]。面绘制重构出的三维模型则是由若干个三角片连接而成的。其步骤框图如图7所示。
图7 MC 算法步骤示意图
Figure7 MC algorithm step
但采用传统的 MC 算法进行三维重建时，依然存在诸多问题。在重建的过程中，由于 VTK 中自带的 MC 算法，需要对输入目标图像的所有体素进行遍历，但是由于上文提到的 MC 算法中的三角片类型可以简化为 15 种拓扑构型[21]，所以存在许多无效运算。而且传统的面绘制 MC 算法在计算三角片顶点和法向量时采用线性插值计算方法，造成计算量大、过程复杂等问题，最终导致三维重建效率低、精度不高。
3.2 改进的面绘制 MC 算法
本文采用中值法代替线性插值方法，其原理简单，采用中值法计算等值面交点和法向量时，误差值小于交点所在棱边边长的一半，并且随着照片分辨率的不断提高，采用中值法所重建的三维模型与原本采用线性插值方法所获得的三维模型在视觉上并没有明显差异。
假设 x1、x2 为体元棱边上的两顶点，n1、n2 为顶点的法向量，x0 为棱边上的中点，n0 为中点法向量。则用中值法可表示为：
(15)
(16)
通过上式可知：中值法能够有效地减少线性插值法所带来的计算量，缩短等值面生成的时间，改善传统 MC 算法的效率问题。
3.3 三维重建试验
本文所采用的三维重建设备为：AMD-RYZEN 5 笔记本电脑、8G 运行内存。图像数据来自那丘、六月田隧道。重建环境为 win 10-64 位操作系统、python 3.8 以及 VTK 9.1.0 版本。
试验分别使用传统的 MC 算法和改进的 MC 算法建立三维模型，比较两种方法模型构建的耗时。
采用传统 MC 算法和改进的 MC 算法进行三维重建，模型分别如图8、9所示。
图8 传统 MC 算法的裂隙三维重建示意图
Figure8 3D reconstruction of fracture based on traditional MC algorithm
图9 改进 MC 算法的裂隙三维重建示意图
Figure9 3D reconstruction of fracture based on improved MC algorithm
重建完成后，分别对两种重建方法消耗时间进行对比，结果如表1所示。
表1 传统 MC 算法与改进 MC 算法速率对比
Table1 Rate comparison between traditional MC algorithm and improved MC algorithm
由表1可知：改进后的 MC 算法在重建速度上得到了显著提高，三维重建效率提高了 18%。虽然使用了改进算法重构裂隙三维模型，但由于改变了等值面插值计算方法，导致三维模型在连接处缺失细节。因此，在整体视觉效果上，与传统 MC 算法并无明显区别。
4 结论
（1）本文将计算机视觉中的三维建模技术运用在隧道围岩裂隙上，并在 python 基础上对采集的裂隙图像进行图像处理，然后使用可视化工具包 VTK 中的 Marching Cubes 算法对经过图像处理后的裂隙图像进行三维重建，发现得到的裂隙三维模型与原始图像中的裂隙形状并无太大差别，实现了从二维图像到三维模型的转化，为隧道围岩裂隙的研究提供了一种新思路。
（2）本文对 VTK 中的传统 Marching Cubes 算法进行了改进，使用中值法替代线性插值法对等值面交点和法向量进行计算，在一定程度上减少了三维模型生成的时间，使三维重建效率提高了 18%。但是改进的 Marching Cubes 算法重建出来的三维模型在细节上还有一些不足，也是以后改进的方向。
