【理论算法】基于数字图像处理的堆积体结构特征提取与建模

        堆积体作为一种特殊成因的第四系堆积物，在我国西南地区大型水利水电工程建设极为常见。土体的结构是指土体内部颗粒间的排列与连接方式，常见的结构类型有单粒结构、蜂窝结构和絮状结构三种类型。图1所示为澜沧江古水争岗冰水堆积体。

图1 澜沧江古水争岗冰水堆积体

       冰水堆积体是由冰水沉积形成的一类第四系堆积物，与一般的土体不同，其在结构上具有自身显著的特点：

（1）结构的二元性

整体上颗粒组成中巨粒和粗粒构成了冰水堆积体内部主体“骨架”，而细粒物质则多以“基质”的形式充填于“骨架”之间，从而整体上形成了“骨架”＋“基质”的结构特征。

（2）结构的随机性与无序性

堆积体在形成过程中巨粒的漂石和碎石、粗粒的碎石和角砾、细粒的粉土和粘土在毫无分选的条件下快速混杂堆积的产物，宏观上表现为无分选、无定向、无磨圆、无层理。

（3）细颗粒物质的充填作用和胶结作用

堆积体内部细颗粒物质含量较低，大多高度分散填充于巨颗粒、粗颗粒的空隙之间，在长期的固结作用下，常常与巨颗粒、粗颗粒一起构成了堆积体的紧密结构特征，并且充填于巨颗粒、粗颗粒之间的细颗粒物质具有一定的“粘结”作用^[1]。

堆积体的结构特性与节理岩体不同，内部块体的几何和发布特征是关注的重点。为了较好地获取堆积体的细观结构特性，本文提出了一种基于连通域标记算法的数字图像处理技术^[2,3,4]。

（1） 数字图像预处理

由于受外界各种因素的影响，利用数码相机拍摄的冰水堆积体图像往往存在大量噪声，下图2所示为本文采用中值滤波的去噪处理。

图2堆积体图像滤波去噪：

（a）去噪前图像，（b）去噪后图像

图像二值化处理在数字图像处理中占有非常重要的地位，通过二值化处理可以将灰度图像中感兴趣的目标凸显出来，便于对目标做进一步的处理。在计算机中，灰度图像中所有像素点的灰度值是由一个二维矩阵g(x,y)来存储的。灰度图像二值化是将图像中灰度值小于某个阀值的像素点的灰度值重置为0，把大于等这个阀值的像素点的灰度值重置为1，从而实现灰度值二值化，可用公式表述为：

式中：              表示二值图像中像素在      位置的灰度值，x和y分别表示像素点在图像像素矩阵中的行、列号，              为某个临界阀值，通常也被称为二值化阀值。      

（2） 连通域标记

连通区域标记是二值图像分析最重要的方法，它通过对二值图像中白色像素的标记，让每个单独的连通区域形成一个被标识的块，进一步的我们就可以获取这些块的轮廓、外接矩形、质心、不变矩等几何参数。

图像在计算机中是以像素为基本单元存储的，每个像素周围均有若干个邻接像素，通常将这些邻接像素所组成的集 合称为该像素的邻域。如图3a所示，对于一个像素来说，其常见的邻域类型有两种：4-邻域和8-邻域。

图3 常见邻域类型：（a）4-邻域；（b）8-邻域

连通域，在数学上顾名思义是指一个闭合连通的区域。一般而言，按照不同邻域类型，则标记出的连通域不同。对于图4a所示的一个二值图像，如果考虑4-邻域，则该图像包含了3个连通区域；如果考虑8-邻域，则该图像包含了2个连通区域。对于冰水堆积体来说，考虑4-邻域的连通域标记是比较合理的。

图4 邻域标记算法图例：（a） 原始二值图像；

（b）4-邻域识别连通域；（c）8-邻域识别连通域

（3） 连通域边界提取

对冰水堆积体二值图像经过连通域标记后，可将图像中每个石块所对应的连通标识出来，所标识出的每个连通域均为一个像素集 合，通过对像素集 合进行处理可以提取连通域的边界。文中获取连通域的边界的方法可以分为三个主要步骤：

下面以图5为例，来介绍连通域边界提取步骤。图5a所示为由3个像素组成的一个连通域。如图5b所示，连通域边界提取的第一步是将三个像素进行分解，得到12条边，其中P2P5及P5P6为重复的边；第二步是将P2P5及P5P6重复的边删除；第三步是按照边的拓扑关系将剩余边连接起来得到了连通域的边界P1-P4-P5-P7-P8-P6-P3-P2，如图5c所示。

图5 连通域边界提取：（a）连通域；

（b）像素分解；（c）删除重复边

（4） 边界光滑处理

根据连通域标记及边界提取获得了冰水堆积体细观结构模型。由于连通域边界是由一系列像素的边构成，因此提取到的石块边界呈一个锯齿状的形状。为了得到形状更为“光滑”的石块，本文采用了岳中琦等人^[5]所提出的边界光滑处理方法对石块边界进行光滑处理，其具体步骤如下：

[1] 选定一个适当的阀值T；

[2] 查找边界中距离最大的两个点，连接这两个点将石块分割为两半部分，任意选择其中一部分进行后续操作；

[3] 搜索该所选部分中所有点中到分割线距离最大的点，并记录该点的编号及该点到分割线的最大距离        ；

[4] 若        小于T，则将分割线作为该部分的边界，开始考虑另一部分；

[5] 若        大于T，则将步骤[3]中所搜索出的点与步骤[2]中所查找的两个点分别连续，并将该部分边界再分为两部分，并重复步骤[3]、[4]；

[6] 重复执行步骤[3]-[5]，直到所有部分的边界点到其相应分割线的最大距离均小于T为止，此时生成的新的闭合分界线就作为石块的新边界。

图6  连通域边界光滑处理技术

基上述理论开发了DigiSim程序（见图7）并在Advances in Engineering Software（JCR Q1，IF=4.255）上发表”DigiSim — An Open Source Software Package for Heterogeneous Material Modeling Based on Digital Image Processing” (doi.org/10.1016/j.advengsoft.2020.102836)。

此外相关代码已公布在github，地址为https://github.com/GeoGroup/DigiSim.git。

图7  DigSim程序用户图形界面