溶蚀孔洞是碳酸盐岩储层中主要的储集空间类型，其成因复杂，分布非均质性强，预测难度很大(Loucks，1999；钟建华等，2012)，往往因古岩溶作用形成，通常和裂缝伴生，一定程度上控制着油气的分布(牛玉静等，2011；石书缘等，2014).研究孔洞连通域定量标记及参数拾取对有效评价碳酸盐岩储层具有重要的意义.全井眼微电阻率成像测井(FMI)图像分辨率高、信息量大(闫建平等，2011)，是一种重要的反映井筒剖面的高分辨率信息，它在地质学和岩石物理学中有着广泛的应用(Assous et al., 2013)，因此是研究非均质碳酸盐岩储层溶蚀孔洞、裂缝的重要资料，前人通过FMI图像对裂缝识别(陈钢花等，1999；Chitale，2005；柴华等，2012)、裂缝信息拾取(闫建平等，2009)、溶洞识别(钟广法等, 2009, 2004)、溶洞自动检测(田金文等，1999)、礁滩储层沉积相及岩性识别(Chai et al., 2009)、缝洞面孔率计算(Xiao et al., 2010)等方面进行了较广泛的研究，也有学者对岩石CT、薄片、电镜等多种微观图像进行处理，从而研究图像中的矿物特征(叶润青等，2011)以及结合分形维数计算研究孔隙结构(Ge et al., 2015；彭瑞东等，2011).但对于成像测井(FMI)图像上溶蚀孔洞连通域标记及其信息拾取方面的研究工作较少，尤其是对井筒剖面非均质溶蚀孔洞信息连续定量拾取方面(如连通域定量标记，连通域半径、圆度、孔洞面孔率及孔洞密度等相关定量参数)还不够深入.

连通域标记属于图像处理中的算法，已有的连通域标记算法可分为三类：基于像素标记的两遍扫描算法(Jiang et al., 2009)；基于线标记的两遍扫描算法(Tu and Lu, 2010)；以及基于模块的标记算法(马江林等，2008).已有一些学者在早期标记算法的基础上采用多种手段来提高连通域的标记速度(牛连强等，2015；王凯等，2014；马益杭等，2013)，并将其应用在医疗、交通、测绘等领域，目前还很少有将这一类图像处理方法应用于成像测井图像处理方面.

笔者以塔里木盆地奥陶系溶蚀孔洞型碳酸盐岩地层的FMI高分辨率彩色图像为研究对象，通过灰度、滤波图像处理后，结合岩心物性特征刻度图像，选取合适的阈值进行分割得到能够反映井壁孔洞特征的二值图像，溶蚀孔洞表现为黑色暗斑.由于溶蚀孔洞发育的FMI图像存在连通域数量多、形态多样的特征，同时鉴于基于等价对的图像连通域标记算法具有快速、不重复标记的优点(左敏等，2011)，因此利用该算法，对FMI二值化图像定量标记其溶蚀孔洞连通域，进而可对每个连通域进行目标信息拾取，包括孔洞半径尺寸、面孔率、圆度等，进一步计算图像面孔率曲线并利用其对溶蚀孔洞发育层段进行分层，即可提取每一层段溶蚀孔洞相关结构参数值分布的非均质信息，为定量快速地评价非均质性强、溶蚀孔洞发育的碳酸盐岩储层提供了技术支撑，也是利用成像测井定量拾取井筒岩石孔洞结构信息的有益探索.

岩溶型储层广泛存在于碳酸盐岩地层剖面中，发育的类型有裂缝型、裂缝-孔洞型、孔洞型(图 1)(Zhang et al., 2007).裂缝型储层中基岩孔隙度极低，倾斜裂缝面与井筒的交线在平面图上呈现为单周期的正弦曲线，因此，利用哈夫变换的点线对偶性可定量拾取裂缝的参数信息.裂缝-孔洞型、溶洞型储层中，孔洞发育非均质性强，有时和裂缝伴生，其参数信息定量拾取一直没有很好的方法来解决.

图 1

Fig. 1

图 1 碳酸盐岩岩溶型储层类型(Zhang et al., 2007) Fig. 1 Types of carbonate karst reservoirs (Zhang et al., 2007)

溶蚀孔洞按其形态及产状可以分为斑点、斑块状溶洞、层状大洞及片状和网状缝洞等类型(钟广法等，2009).①斑点-斑块状溶洞：此类溶洞外形呈等轴状，但尺寸较小，单个溶洞的直径一般从2 mm~20 cm不等，其中规模较小者(2 mm~2 cm)在岩心或成像测井图像上表现为小斑点；规模较大者(数厘米至20 cm)则呈斑块状(图 2).②层状大洞：指洞径大于井筒直径的等轴状大溶洞或厚度大于10 cm且在各极板图像上连续分布的大型溶洞，因其在岩心或成像测井图像上似层状而得名，其顶、底界面多不规则，通常切割围岩层理.③片状或网状缝洞：单个溶洞的形态呈片状，在岩心表面或展开的二维成像测井图像(FMI)上呈线性，与裂缝类似，但其形状多不规则，宽窄不一，局部膨大，解释为被溶蚀作用改造的裂缝或层理.按产状可进一步细分为水平缝洞、倾斜缝洞和垂直缝洞等类型.正是由于孔洞、缝洞的随机、非均质性，其连通情况非常复杂，那么数字图像处理中的连通域标记算法是拾取孔洞连通域的有效方法.

图 2

Fig. 2

图 2 斑点-斑块状溶洞 Fig. 2 Karst carve of spot-patch shape

(1)连通域

在图像自动目标识别和跟踪过程中，需要对图像进行灰度图转换并将目标进行阈值分割提取后得到二值(或多值)图像，图像中通常包含多个连通区域.为了提取不同区域的特征，就需要对连通区域进行标记.

所谓连通域，是指由若干个像素组成的集合，该集合中的像素具有以下特性：①所有像素的灰度级别均小于或等于连通域级别；②同一个连通域中的像素两两相通，即在任意两个像素之间存在一条完全由这个集合的元素构成的通路(左敏等，2011).而二值图像中的连通域就是互相连在一起的1像素或0像素区域(王凯等，2014).对于FMI经二值化后的图像，感兴趣的是像素值为0(黑色)的连通成分即孔洞连通域.

研究图像的连通域，其实就是研究像素间的连通性问题.在二维图像中假设目标像素点与周围某相邻的像素点其像素值相同，则称这两个像素点连通.在研究连通性时需要首先确定邻域连通规则为4邻域连通还是8邻域连通.4邻域连通关注的是目标像素点的上、下、左、右4个位置点.8邻域连通则选取目标像素在二维空间中3×3矩阵中所有的相邻像素点，即除了上、下、左、右点外，还包括左上、右上、左下、右下4个位置点(图 3).实际上4邻域连通包含于8邻域连通，但是由于两种连通规则关注的像素点的数目与位置的差异，4邻域连通会影响FMI二值图像实际标记操作的准确性.所以本文的算法采用了8邻域连通规则.

图 3

Fig. 3

图 3 4邻域和8邻域连通性 Fig. 3 Four-and eight-neighborhood connectivities

(2)图像标记

如果想从井筒剖面上连续、定量地评价溶蚀孔洞的发育程度，需要对各个连通区域进行分别检测判断.连通域标记算法可通过对FMI二值化后的图像进行分析处理，依次对图像上的溶蚀孔洞进行编号标记.该方法在碳酸盐岩储层洞穴连通性二维表征(Izgec et al., 2008)、孔径分布网络模拟(Rodríguez et al., 2011)方面已有一定的应用.目前，对连通域标记方法已有的改进和研究主要集中在如何降低等价区域合并过程的时间复杂度和空间复杂度；或者是通过某种特殊手段，免去等价区域合并过程，一次性完成连通域标记(王晶等，2005).

二值图像中可能存在一个连通域也可能存在多个连通域，FMI图像经二值化后，往往存在多个连通域，因此需要按次序将图像中符合某种连通规则的(8邻域连通)像素用相同的标号标记出来.而只有经过准确标记处理的二值图像才能更好地进行图像信息的定量提取.

图像的标记方法有很多种，常见的方法大致可分为像素标记法、线性标记法和区域增长法，其中像素标记法使用最为广泛.迷宫算法和种子算法是连通域检测比较成熟的算法，两个算法均采用像素标记法.笔者考虑到以上算法的优缺点，参考像素标记的连通域标识方法，对FMI二值图像进行一次扫描，对每个像素进行标识，同时处理等价对解决相同连通域不同标识的问题，最终标识出所有的连通域.

(3)基于等价对处理的连通域标记算法原理、思路

连通域标识算法的目的是通过对图像的扫描和分析得到图像中的连通域划分和连通域的数目，进而识别出图像中单个或多个活动目标.为了便于处理，本文算法的输入为FMI图像经灰度、滤波、二值化并取反后的图像(Er)数组Count[][]，其中像素值非0即1，该数组在图像显示时分别代表了黑色(背景颜色)和白色(目标颜色)两种颜色，即白色区域表示溶洞.输出的是根据Count[][]做了连通域标记后的新结构数组NCount[][]以及由NCount[][]还原出的图像，该结果数组中不同的值表示属于不同的连通域.图 4是以对一个像素点进行连通域检测为例的检测流程.

图 4

Fig. 4

图 4 FMI二值图像连通域标记算法流程 Fig. 4 Flow chart of algorithm for connected domain marks in FMI binary image

当确定8邻域连通规则之后，对图像的二值数值进行扫描，顺序为从左到右、从上至下.在处理时，每一个目标像素点只需判断该像素点和周围已确定连通性的像素点之间的关系就可确定自己的连通性.然而不同的位置点所需关注的邻域点不同，图 5和图 6分别为所有可能的位置点及各自确定连通性时所需关注的邻域点.

图 5

Fig. 5

图 5 二值图像数组不同位置标识 Fig. 5 Identifications of binary image array at different locations

图 6

Fig. 6

图 6 不同位置点连通域检测时关注邻域点 Fig. 6 Focused neighborhood points in connected domain detection of different points

但由于溶洞形态以及扫描顺序的影响，导致开始认为是两个不连通的区域，随着标记扫描过程的深入，发现它们实际上属于同一个连通区域，这一现象可以将其称为等价对现象，图 6中，(b)、(c)、(d)三种类型都有可能产生等价对.我们必须在程序中对等价对进行处理，否则这将严重影响标记效果，从而也会影响到后期图像中溶洞信息的准确拾取.

图 7可以看出，在对图 7a进行扫描标记时，虚线方框区域会产生等价对，如图 7b所示，图中白色区域实际为同一个连通区域，但是由于扫描顺序和标记算法的影响，这一区域被标记为‘8’和‘14’两个连通域，这一结果不仅影响溶洞、溶孔个数的准确统计，也影响溶蚀孔洞表征信息的准确提取.所以在图像标记过程中，当出现等价对时需要对其进行及时的处理，以免在之后的信息拾取工作中遇到问题.

图 7

Fig. 7

图 7 等价对处理前后效果 Fig. 7 Effect comparison of before and after equivalent processing

在此，以图 6c为例说明等价对的处理过程.将已完成标记的四个像素点的标记值放在新的数组B[]中，假设B[]中有两个非0且不相等的值，那么它们就是一对等价对.此时，先令NCount[x][y](4号位置处的像素点标记值)等于B[]中第一个非0值，然后找到NCount[][]中所有标记值等于B[]中第二个非0值的点并记录它们的位置，最后将NCount[x][y]的值逐一赋给它们.至此，这一等价对被处理完成.至于B[]中非0值大于两个的情况，均可按照这种方法进行处理.

图 7c显示的是虚线方框区域进行了等价对处理的结果，该区域被统一标记为‘7’，排除了等价对的影响.在对全图进行图像标记和等价对处理之后，连通域标记的最大数值就可以认为是图中溶蚀孔洞连通域的个数.

在进行FMI成像测井图像溶蚀孔洞、裂缝信息拾取之前，FMI的原始彩色图像需经过图像预处理将其他干扰信息去除，从而更有利于检测孔洞、裂缝信息.显然，这是一个很复杂的过程，需要经过一系列相关的图像处理步骤来完成.

(1)标记、拾取信息步骤

①图像预处理：首先对原始图像(图 8a)进行灰度化处理，得到灰度图像(图 8b)；灰度化处理后，为了消除图像噪声并同时保护边缘信息，再对图像进行中值滤波处理实现对图像的非线性平滑处理(图 8c)；接着对滤波后的图像进行阈值分割(阈值选择时需要尽可能结合岩心孔隙度进行限制和修正，从而选取最能反映真实溶孔的阈值)，得到二值图像(图 8d)，为孔洞连通域标记、信息参数拾取奠定了基础.

图 8

Fig. 8

图 8 成像测井图像预处理及孔洞标记、信息拾取 Fig. 8 Preprocessing of imaging logging image and hole labeling and information picking up

②连通域标记、信息拾取：首先利用上述基于等价对的图像连通域标记算法对二值图像(图 8d)进行连通域标记处理，得到标记后图像(图 8e)，同时对标记图像中的每一连通域拾取其长、宽、内切圆半径、外接圆半径、圆度等参数信息(其详细的定义描述见下文)；连通域标记图像实际上仍是二值图像，但不同的是对连通域进行了连续标记，可通过程序对每个连通域进行信息定量化拾取，其放大后图像即图 8f，可清楚看到连通域标记情况；接着再对连通域标记图像(或二值化图像图 8d)按深度域像素拾取其面孔率(图 8g)、间隔20个像素(20 pixel)的面孔率直方图(图 8h)，为在深度域了解溶蚀孔洞发育程度及实际处理过程中对图像进行分层处理奠定了基础.

(2)标记、拾取信息参数

①连通域长、宽、内切圆半径、外接圆半径：Ⅰ.连通域长，指的是连通域最左边像素到最右边像素之间的像素个数(包括左右边像素)；Ⅱ.连通域宽，指的是连通域最顶端像素到最底端像素之间的像素个数(包括顶底端像素)；Ⅲ.连通域内切圆半径，连通域最左右两边像素和最顶低端像素可绘制成一个矩形，内切圆半径有两个值，一个指的是这个矩形短边的二分之一像素(R_内S)，另一个指的是这个矩形长边的二分之一像素(R_内L)；Ⅳ.连通域外接圆半径，连通域最左右两边像素和最顶低端像素可绘制成一个矩形，外接圆半径指的是这个矩形对角线的二分之一像素(R_外).由此定义，可通过程序设计输出图像(图 8e)中每个连通域的这些参数，同时可统计一段图像中连通域长和宽分布、内切圆半径分布(图 9a)、外接圆半径分布(图 9b)、连通域圆度分布(图 9c)等信息，这些参数的分布能够反映井筒地层孔洞发育程度(如果知道井眼直径，可利用图像宽度与井眼周长的关系，进一步得到上述参数像素值与实际尺寸(cm、in)之间的关系).

图 9

Fig. 9

图 9 连通域不同参数信息分布 Fig. 9 Different parameter information distributions of connected domain

②连通域圆度：圆度本来是沉积岩石学中描述岩石碎屑颗粒的原始棱角被磨圆的程度，它是碎屑颗粒的重要结构特征.它与颗粒的形状无关，只是棱角尖锐程度的函数.圆度在几何上反映了颗粒最大投影面影像中的隅角曲率(朱筱敏，2008)，韦尔德在1932年提出下列圆度计算公式：

(1)

式中，r为隅角的内切圆半径，n为隅角数，R为颗粒的最大内切圆半径.

成像测井FMI连通域标识图像中，想描述溶蚀孔洞连通域跟圆形的接近程度，于是借鉴碎屑颗粒圆度的定义及公式，但实际上该定义和想描述的连通域圆度在概念上是不同的，且该公式涉及的隅角参数也不便于图像处理过程中程序设计，因此，本文将连通域最左右两边像素和最顶低端像素组成的矩形与正方形的接近程度，定义成“连通域圆度”，具体由(2)式来计算：

(2)

式中，R_内S为短的内切圆半径，R_外为外接圆半径.分母意义：在正方形中内切圆半径与外接圆半径的比即为可以将其视为圆形的标准.在得到所有的连通域矩形内切圆和外接圆半径之比后，将其与标准值作比较可以更直观地看出连通域形状与圆的接近程度.

③连通域分选系数：分选系数原本是表示沉积物分选程度的参数，它反映颗粒大小的均匀程度，或者说是表现沉积物围绕集中趋势的离差(朱筱敏，2008)，给出的分选系数公式为：S_O=P₂₅/P₇₅，式中P₂₅和P₇₅分别代表累积曲线上颗粒累积含量25%和75%处所对应的颗粒直径.笔者想借鉴该公式来表征图像中溶蚀孔洞连通域尺度的分选性，因此将上述公式修改为：

(3)

式中，PC₂₅/PC₇₅表示在连通域内切圆半径、外接圆半径累积曲线上累计频率25%对应的半径值与累计频率75%对应的半径值之比.通过计算，图 8e图像中内切圆半径(R_内S)分选系数为2.25，外接圆半径分选系数为2.6.当然实际图像处理过程中可能并不会出现PC₂₅或PC₇₅，一定程度上用概率统计参数“方差”来表征分选性也是可行的.结合沉积物分选系数的概念，根据S_CO的大小也可以划分孔洞连通域的分选等级：S_CO=1~2.5，分选好；S_CO=2.5~4.0，分选中等；S_CO>4.0，分选差.

④溶洞密度：为了更好地表征溶洞发育程度，引入溶洞密度的概念，即单位长度上的溶洞个数，公式为：

(4)

式中，N为发育溶洞的个数，L为图像长度.由于溶蚀孔洞发育的大小均匀程度存在差异，单独使用这一参数评价溶洞发育程度容易产生偏差，所以溶洞密度需要结合连通域分选系数才能比较全面地评价地层的溶洞发育情况.

按照上述介绍的详细步骤，通过实例来说明连通域标识法在FMI成像测井图像中标记、拾取溶蚀孔洞信息的应用.

首先，进行图像预处理，对较长的一段原始图像(图 10a)依次进行灰度化处理、中值滤波和阈值分割，得到与之对应的灰度图像(图 10b)、中值滤波图像(图 10c)和二值化图像(图 10d).

图 10

Fig. 10

图 10 图像连通域标记实例 Fig. 10 Example of connected domain image labeling

然后，进行连通域标记，利用连通域标记算法对二值图像(图 10d)进行标记处理，得到标记后图像(图 10e)；接着再对连通域标记图像(或二值化图像图 10d)按深度域像素拾取其面孔率(图 10f)、间隔20个像素的面孔率直方图(图 10g)，为在图像深度域了解溶蚀孔洞发育程度及图像分层奠定了基础.

最后，进行图像分层、连通域信息拾取，利用面孔率信息(图 10f、10g)反映的溶蚀孔洞发育程度对图像进行分层(图 11a)；同时提取每一连通域的长、宽、内切圆半径、外接圆半径、圆度等参数信息；接着分别统计得到每一层中连通域内切圆半径分布(图 11b)、外接圆半径分布(图 11c)、圆度分布(图 11d)等相关信息(表 1)，从图 11及表 1中溶蚀孔洞连通域的定量参数信息，可以较好地定量反映出FMI图像深度域溶蚀孔洞发育程度及非均质性，有助于更精确地评价井筒地质特征.

图 11

Fig. 11

图 11 连通域标记后多参数信息拾取实例 Fig. 11 Example of picking up parameter after connected domain labeling

表 1(Table 1)

表 1 图像分层中连通域参数信息拾取数据 Table 1 Picking up connected domain parameters information in different layers 分层 连通域长 (pixel) 连通域宽 (pixel) 内切圆半径 R内S(pixel) 外接圆半径 R外(pixel) 圆度 溶洞密度 (n/cm) 1 最大 121 122 60.2 85.9 0.99 0.132 最小 1 1 0.5 0.7 0.06 平均 23.4 29.2 11.7 19.2 0.79 方差 1967.3 1933.2 491.8 957.3 0.11 2 最大 31 36 15.5 23.8 0.99 0.267 最小 1 2 0.5 1.1 0.56 平均 9.6 11.9 4.7 7.7 0.8 方差 109.7 158.7 27.4 66.1 0.02 3 最大 157 127 63.5 101.0 0.99 0.201 最小 1 2 0.5 1.1 0.34 平均 32 23.2 10.9 20.1 0.76 方差 1906.4 922.2 226.1 693.0 0.03 4 最大 54 17 8.5 28.3 0.99 1.077 最小 1 1 0.5 0.7 0.14 平均 7.6 4.8 2.1 4.7 0.64 方差 103.4 17.8 4.6 28.4 0.05

表 1 图像分层中连通域参数信息拾取数据 Table 1 Picking up connected domain parameters information in different layers

(1)全井眼微电阻率成像测井(FMI)经数据处理后可得到全井壁高分辨率的彩色图像，在图像灰度化、中值滤波处理后，通过阈值分割可得到能够反映井壁溶蚀孔洞特征的二值图像，孔洞表现为黑色暗斑，只能定性地反映不同层段溶蚀孔洞发育的变化情况以及形态分布特征.

(2)利用基于等价对处理的图像连通域标记算法，可准确从二值图像中标记溶蚀孔洞连通域，进而可对每个连通域进行目标信息定量拾取，包括连通域半径尺寸、连通域面积(面孔率)、连通域圆度等孔洞结构参数.

(3)将反映溶蚀孔洞发育程度的面孔率曲线进行分层，可拾取每一层段溶蚀孔洞的连通域分选系数、溶洞密度、溶洞数量等结构参数，能够较好地表征FMI图像深度域溶蚀孔洞的发育程度及非均质性，是FMI图像应用于井筒孔洞信息拾取的有益探索，今后进一步研究图像三维信息处理方法并实现软件商业化，将更有助于精细定量地评价井筒地质特征.

本文完成过程中，西南石油大学地球科学与技术学院张云峰老师在提供溶洞成像测井图像方面给予了很多帮助，在此表示感谢！