基于分水岭算法的孔洞拾取与参数计算

2022-04-23石新罗谋兵朱萌林强徐苗

测井技术 2022年1期

石新,罗谋兵,朱萌,林强,徐苗

(1.中国石油川庆钻探工程有限公司地质勘探开发研究院,四川成都610051;2.中国石油川庆钻探工程有限公司,四川成都610051)

0 引言

孔洞型碳酸盐岩储层通常发育大量孔洞,为油气储层提供了良好的储集条件,具有较大的开发潜力,一直是油气开发的重点研究对象。对于碳酸盐岩储层,应用常规测井资料无法准确地识别储层类型以及储层级别,而电成像测井资料较好地弥补了这一缺陷,特别是对于孔洞型碳酸盐岩储层,通过对溶蚀孔洞的统计与计算能够更准确快速地识别储层段,并对储层做出定量评价,从而围绕溶蚀孔洞的发育状况进行孔洞型储层评价。

多年来,国内外围绕孔洞的发育状况对孔洞型储层开展了大量研究,孔洞型油气储层的电成像测井图相关研究也逐步增多。Meger 和Becher提出标准层控制下的图像处理分割技术,开创了利用图像处理技术拾取孔洞的先河[1],国外这方面的技术较国内更加成熟。国内很多学者相继研究了储层孔洞空间的拾取。张丽莉等[2]将分割算法应用于地层微电阻率扫描成像测井资料的数据处理,在一定程度上分割出了有意义的简单缝洞信息,取得了良好的处理效果;秦巍等[3]提出将腐蚀膨胀的算法原理以及图像拾取的方法应用于电成像测井资料,实现了对地层信息的自动拾取与获取;李雪英等[4]提出了经典的轮廓跟踪法对孔洞信息进行提取,并对孔洞参数进行计算,但由于边缘点坐标信息存储计算量较大,导致计算速度较慢;赵军等[5]通过比对大量图像分割算法,发现利用奇异多阈值分割算法、基于过渡区的分割算法、霍普菲尔德神经网络法、基于图像间模糊散度阈值化算法能够较好地将地层中的有用信息从背景中提取出来,对后期缝洞的定量计算提供了可靠保证,同时他们利用边缘标识法以及边缘填充法对缝洞参数进行了定量计算,但由于对内存的频繁操作导致计算速度较慢;李茂兵等[6]采用弗里曼链码来记录目标边缘,并在此基础上计算孔洞面积、孔洞密度、孔洞平均大小及面孔率等参数,有效地减少了计算所占用的存储空间并提高了计算速度,但其孔洞轮廓提取精度有一定的局限性;张林等[7]提出了利用梯度方向直方图统计量和熵来提取裂缝、溶蚀孔洞、燧石结核3种典型地质现象的特征量,用BP神经网络对这些地质现象进行自动拾取,但由于研究样品较少且地层背景复杂,造成该方法应用局限性较大。

以上方法存在孔洞拾取精度低、处理速度慢以及应用局限性大等问题,本文基于孔洞型碳酸盐岩储层电成像测井资料,提出利用优化的分水岭算法对孔洞进行拾取,提高对该类储层电成像测井资料的处理效率与精度。

1 电成像测井孔洞响应特征

电成像测井主要通过探测井壁附近的地层电阻率变化,经过标准化处理校正后得到一系列阵列式井壁电阻率数据,将标准化后的电阻率数据刻度后转换成一定色谱下不同的色标数据,以电成像测井图的形式输出[8-9]。在电成像图上,颜色越亮表明电阻率越高,颜色越暗表明电阻率越低[10-11]。由于碳酸盐岩岩石骨架呈现高电阻率特征,而孔洞中充满钻井液呈现低电阻率特征,且碳酸盐岩储层的孔洞由次生溶蚀作用形成,孔径变化范围大,因此,在电成像测井图上孔洞一般表现为大小不一的暗色斑点状图形(图1)。基于电成像测井图的孔洞响应特征,准确拾取电成像测井图中的暗色区域是准确计算面孔率的重要前提,本文以此为基础展开对于孔洞的研究。

图1 孔洞型碳酸盐岩储层的电成像测井响应特征图

2 基于分水岭算法的孔洞拾取

2.1 分水岭算法原理

分水岭算法是一种由L.Vincent提出的数学形态学分割方法,它基于拓扑理论将大地测量学的拓扑地貌应用于图像形态学之中,地形的高程对应于图像中每个像素点的灰度值,每个集水盆地由对应区域局部极小值点及其相关的像素集构成,分水岭则对应形态学图像中的局部极大值像素集[12-13]。在分水岭算法中,将整个图像看作一个大的集水盆地,其中包含众多小的集水盆地,注水过程不仅能在大盆地与每个小盆地之间形成分水岭,也能在紧邻的小盆地之间构筑分水岭堤坝。分水岭形成的过程是模拟向所有聚水盆地同时注入水,水从各个极小值点慢慢淹没盆地,随着淹没程度的加深,为了使各个聚水盆地的水不溢到其他聚水盆地,必须在2个聚水盆地汇合处构筑堤坝,而构筑的堤坝即为分水岭,这也是分水岭算法分割图像得到的边缘[14]。

基于图像形态学实现分水岭算法,其中洪水淹没的过程在梯度图像上进行,梯度图像同时也是灰度图像,图像的灰度空间可以类比为地球的地理结构,灰度值较大的像素连接起来,就是分水岭,灰度值较小的像素也就是常说的山谷。二值化阈值作为水平面,低于阈值的会被淹没,水从孤立山谷(局部最小值)开始充填。当水平面上升到一定高度时,水会溢出当前山谷,将图像其他位置的像素与当前山谷灰度值大小作为衡量修筑分水岭大坝的指标,从而避免两个山谷的水进行汇集。构建大坝的不同山谷中的水,也就是在起始点(灰度值最小的点)一定距离内像素值相差不大的点的集合,称作连通域。每一个灰度像素点称作连通分量(见图2),图中h代表水平面阈值,水平面阈值实际为一个像素阈值,范围在0～225,下面具体描述分水岭算法的实现。

图2 分水岭算法结构示意图

确定梯度图像中的像素,按照灰度值的大小进行分类,并且设定一个水平面阈值,水平面阈值类似于构筑分水岭的水平基准;找到灰度值最小的像素点,默认标记为灰度值最低点,以其作为起始点,将水平面阈值从最小值开始增长;设置像素之间距离阈值作为构建分水岭盆地条件,常用的分水岭算法利用梯度重构算法进行距离的衡量,本文采用欧几里德距离变换算法进行距离衡量。在水平面阈值增长的过程中,会碰到起始点周围的邻近像素,通过计算这些邻近像素与起始点的距离,在每一个连通分量中找到极小值G1,G2,…,GR(总共R个连通分量),从灰度的最小值开始泛滥,被淹没的点的集合记作T[n]

T[n]={(x,y)|g(x,y)

(1)

式中,g(x,y)为原图像(x,y)坐标位置的像素值,无量纲;x、y为某个像素点对应的行索引和列索引,无量纲;n为水位值,无量纲。

令Q为T[n]中的连通分量集合,对每个连通分量q∈Q[n],在每一次水淹后需要判别集合K

K=q∩M[n-1]

(2)

式中,K为判别水淹是否会溢出的指标,无量纲;[n-1]为水位为n-1时被水淹部分的所有连通分量的并集,无量纲。

第1种判别情况:若K为空集或K包含M[n-1]的1个连通分量,则将q并入M[n-1]得到D[n],在D[n]的基础上进行分水岭构筑;第2种判别情况:若K包含M[n-1]的2个及以上连通分量则意味M[n]中的q在上一次水淹结束时是多个连通分量,而在水位为n时汇聚成一个汇水盆地,即需要返回到水位为n-1时,将此时的图像以水位作为灰度阈值进行二值化,M[n-1]对应的二值图像记为S[n-1],在S[n-1]的基础上进行分水岭构筑[15]。

对K集合判断的过程,实际上是筛选出需要构筑分水岭的连通分量,利用连通分量再对目标区域进行分水岭的构筑,由此实现对目标轮廓的拾取及分割。

2.3 分水岭算法优化

分水岭算法具有计算速度快,对图像的边缘能够精准定位的优点。该算法应用在简单的图像上效果良好,但对于复杂且噪声较多的图像,其强大的边缘提取能力对噪声及无意义纹理的响应也同样敏感,导致拾取不必要的细节,造成较严重的过度分割。

为解决过分割问题,选取中值滤波器进行去噪,使用形态学消除细节,对目标区域进行标记。其优点:利用欧几里得距离变换算法取代梯度重构算法,建立地形差异更明显的集水盆地,在利用阈值算法重构二值图的过程中,减少了部分伪极小值及其相关像素构成的聚水盆地数量,有效抑制过度分割现象;并且利用目标区域标记(前景目标标记以及分水岭构筑区域标记)取代循环队列标记水淹的过程,减小了原算法水淹过程中不断判断K集合的大量计算,从而提高了算法的计算速度。

2.3.1图像平滑预处理

将电成像测井响应特征图(见图1)转化为原始灰度图[见图3(a)],对原始灰度图进行中值滤波处理,在实际应用的过程中,需要在尽量保留原有信息的情况下,过滤掉图像内部的噪声,这一过程会得到较平滑的图像。

图3 中值滤波平滑处理前后对比图

由于中值滤波不存在均值滤波、高斯滤波等所带来的细节模糊问题,能够在不影响原有图像的情况下较好地去除“椒盐”噪声,在滤波器的适配上,中值滤波器能较好地适应实际需求。图3(b)是对图3(a)进行中值滤波的结果,可以看到原始灰度图内含有的噪声信息被有效地过滤掉了,减小了电成像测井图中因电阻率的微小波动造成的图像干扰。

2.3.2形态学开运算

将中值滤波后的灰度图[见图3(b)]通过阈值算法转换为初始二值图(见图4),通过在初始二值图像上平移结构元素并进行特定的逻辑运算,实现膨胀与腐蚀操作[16-17]。开运算为先腐蚀后膨胀的过程,闭运算为先膨胀后腐蚀的过程[13],开运算一般用来消除小的前景目标、在细长的点处分离目标、在不明显改变目标面积的情况下平滑其边缘轮廓,这样就不会影响到图像中孔洞的拾取和面积大小。在形成初始二值图像的基础上,选择3×3大小的结构元对初始二值图进行开运算,其对比效果如图4、图5所示,图4的许多噪点在开运算后(见图5)已经被消除,并且图5的孔洞边缘更加平滑。

图4 初始二值图

图5 开运算后二值图

2.3.3目标区域标记

当孔洞之间没有较多连接时,直接利用腐蚀操作就能确定前景目标(孔洞),但在储层中往往发育较多的孔洞彼此相互连接,所以确定前景目标就变得困难。

本文利用简单欧几里得距离变换算法确定前景目标,基于欧几里得距离变换算法原理,计算开运算后的二值图像(见图5)中每个前景目标像素点与其最邻近背景像素点的欧氏距离,将此距离设置为对应前景目标像素点的灰度值,前景目标像素点灰度越高,距离背景越远。形成的基于欧几里得距离变换算法的灰度图像见图6(a),对图6(a)进行简单阈值操作得到确定前景目标图[见图6(b)]。

图6 前景目标区域标记过程

对初始二值图像(见图4)的前景目标进行膨胀操作得到图6(c)所示图像,将膨胀图与确定前景目标图相减得到孔洞边缘的未知区域(包含前景目标边缘与膨胀边缘之间的区域)[见图6(d)]。

对确定前景目标图进行标记,将确定背景(包含未知区域部分)的像素点全部标记为0,前景目标从1开始进行正整数标记,其中标记的数字以不同颜色呈现,以代表不同的目标区域,形成的标记结果图见图7。

图7 标记结果

2.3.4分水岭构筑

利用标记后的结果图(见图8),可以确定孔洞边缘区域与背景区域:对标记结果图创建标记矩阵,用于表示图中的不同区域;利用标记结果图对真实区域进行矩阵表示,将标记为0的区域视为不确定区域,将标记为1的区域视为背景区域,将标记大于1的区域视为想得到的前景目标,也就是孔洞边缘;将各个标记值相同的区域连接起来从而形成连接图,进而完成分水岭的构筑,即对孔洞进行边缘轮廓的拾取与分割。优化后的分水岭算法大大减少了无地质意义的极小值点的干扰,除了少部分噪点,大多数的暗色孔洞都被拾取分割了出来,为计算面孔率提供了较准确的边缘(见图8)。

图8 孔洞拾取结果

3 应用效果及分析

基于优化的分水岭算法对X1井3 219～3 224 m层段电成像测井图进行孔洞边缘拾取并分割,统计得到该井段的孔洞总面积为0.763 1 m2,计算该井段的面孔率为9.278%。为了验证拾取结果的准确性,对孔洞处理拾取并将计算后的面孔率于图9第4道中以连续曲线呈现,将面孔率与对应深度的岩心分析孔隙度对比发现,其对应关系较好,优化的分水岭算法得到的面孔率与对应岩心分析孔隙度的符合率达83.67%,说明经过优化后的分水岭算法应用于孔洞型碳酸盐岩储层电成像测井图的孔洞拾取效果较好。