CN108734714B - 一种基于Matlab分析碳酸盐岩结构的方法 - Google Patents

一种基于Matlab分析碳酸盐岩结构的方法 Download PDF

Info

Publication number
CN108734714B
CN108734714B CN201810571924.5A CN201810571924A CN108734714B CN 108734714 B CN108734714 B CN 108734714B CN 201810571924 A CN201810571924 A CN 201810571924A CN 108734714 B CN108734714 B CN 108734714B
Authority
CN
China
Prior art keywords
carbonate rock
matlab
pores
image
function
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201810571924.5A
Other languages
English (en)
Other versions
CN108734714A (zh
Inventor
久博
黄文辉
何明倩
刘凯
贾瑜
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China University of Geosciences Beijing
Original Assignee
China University of Geosciences Beijing
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China University of Geosciences Beijing filed Critical China University of Geosciences Beijing
Priority to CN201810571924.5A priority Critical patent/CN108734714B/zh
Publication of CN108734714A publication Critical patent/CN108734714A/zh
Application granted granted Critical
Publication of CN108734714B publication Critical patent/CN108734714B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/11Region-based segmentation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/60Analysis of geometric attributes

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Geometry (AREA)
  • Analysing Materials By The Use Of Radiation (AREA)

Abstract

本发明公开了一种基于Matlab分析碳酸盐岩结构的方法,包括如下步骤:首先读取碳酸盐岩的铸体薄片;然后选择色彩模式,将碳酸盐岩的铸体薄片图像二值化;再识别二值图像的黑色像素点,计算孔隙度,同时识别二值图像的白色像素点,分辨颗粒大小,计算粒度;然后求取碳酸盐岩的粒度与孔径大小频率分布图,最后求取碳酸盐岩粒度与孔径的平均直径及方差,并输出二值图像和求取的参数,通过Matlab定量计算不同碳酸盐孔隙结构,平均误差率为+4.661%,其中粒间孔平均误差率为‑4.370%,铸模孔平均误差率为+6.734%,生屑孔误差率为+6.144%,微裂缝平均误差率为‑0.314%,铸模孔与生屑孔的孔隙结构最为适用。

Description

一种基于Matlab分析碳酸盐岩结构的方法
技术领域
本发明涉及致密砂岩储层孔隙表征分析技术领域,具体为一种基于Matlab分析碳酸盐岩结构的方法。
背景技术
致密砂岩储层通常为储层渗透率低的砂岩储层,致密低渗透储层是一个相对概念,世界上并无一固定的标准和界限,它是由不同国家不同时期的资源和技术经济条件和世界油气价格而决定的。
储集层是指具有连通孔隙、允许油气在其中储存和渗滤的岩层。世界上已发现的油气储量大多数来自沉积岩层,其中以砂岩和碳酸盐岩储集层最为重要,裂缝性泥岩和煤层也可作为储集层;火成岩和变质岩储集层中也有工业性油气发现。
储层(Reservoir)的储集能力是由储集层的岩石物理性质决定的,通常包括其孔隙性、渗透性;孔隙性决定了储层储存能力的大小,渗透性决定了储集物的渗流能力。
碳酸盐岩是主要由方解石和白云石等碳酸盐矿物组成的沉积岩,其岩石结构包括颗粒、胶结物、泥晶、孔隙几部分,颗粒可分为内碎屑、鲕粒、球粒、藻粒、团块这几部分,岩石微观孔隙结构是指岩石中孔隙和吼道的几何形状、大小、分布及其相互连通关系。国内外多名学者对碳酸盐岩孔隙结构、表征等进行研究,多种手段被应用于研究致密砂岩孔隙,以核磁共振NMR测井技术、流体注入实验法和图像观测法这三大类为主。基于岩石电阻率参数的核磁共振测井技术,以测定储层孔径分布为主。以压汞法为代表的流体注入实验法能够间接地通过压汞曲线获得储集空间分布、孔径大小等参数,但并不能像以扫描电镜SEM及纳米CT技术等图像观测法直接地观测到储集空间分布、孔喉特征等特征。然而以扫描电镜及纳米CT为实验指导,研究结果大多为对孔喉的几何形状及联通关系进行定性描述,对孔隙表征的定量研究仍有不足。
发明内容
为了克服现有技术方案的不足,本发明提供一种基于Matlab分析碳酸盐岩结构的方法,能有效的解决背景技术提出的问题。
本发明解决其技术问题所采用的技术方案是:
一种基于Matlab分析碳酸盐岩结构的方法,包括如下步骤:
步骤100、读取碳酸盐岩的铸体薄片;
步骤200、选择色彩模式,将碳酸盐岩的铸体薄片图像二值化;
步骤300、识别二值图像的黑色像素点,计算孔隙度,同时识别二值图像的白色像素点,分辨颗粒大小,计算粒度;
步骤400、求取碳酸盐岩的粒度与孔径大小频率分布图;
步骤500、求取碳酸盐岩粒度与孔径的平均直径及方差,并输出二值图像和求取的参数。
作为本发明一种优选的技术方案,所述步骤200中,铸体薄片图像二值化的方法为:将铸体薄片图像导入Matlab的内置Color Threshod应用中通过HSV成像模式进行图像二值化转换。
作为本发明一种优选的技术方案,不同的成像模式适用于不同图像色彩空间,最终目的均为最大限度地识别碳酸盐岩结构并将其转为二值图像,所述HSV模式因其可以调控图像色调H、饱和度S以及明度V。
作为本发明一种优选的技术方案,所述孔隙度的求取方法为:通过Matlab中的size函数对二值图像中的像素进行识别提取,通过函数可以分别获得图像的总像点数Na,再进一步调用能够识别0像素值的函数bwarea获取孔隙像素点数Np,即黑色像素点Np,孔隙度P即黑色像素点与总像素点的比值即:
Figure GDA0003844263060000021
作为本发明一种优选的技术方案,所述粒度的求取方法包括:
首先,通过Matlab内置函数bwlabel对相同像素值的区域进行识别统计,同时可以确定区域的边界,基于此原理确定孔隙区域,即黑色或白色像素值连通区域;
然后对碳酸盐岩的每一个颗粒以及孔洞调用bwlabel函数,求出第i个区域的X和Y轴的最大值以及最小值即ximax、ximin、yimax、yimin,所以对于第i个颗粒或孔洞的直径大小di为X和Y方向的差值的平均值,即直径大小di
Figure GDA0003844263060000031
作为本发明一种优选的技术方案,所述步骤400中,碳酸盐岩的粒度与孔径大小频率分布图通过调用内部函数imhist对所有颗粒及孔洞进行叠加进行绘制。
作为本发明一种优选的技术方案,所述碳酸盐岩粒度与孔径的平均直径及方差的求取方法为:
分别调用函数mean以及函数var求取平均孔隙直径
Figure GDA0003844263060000034
以及方差S:
Figure GDA0003844263060000032
Figure GDA0003844263060000033
作为本发明一种优选的技术方案,该方法适用的颗粒包括内碎屑颗粒、鲕粒颗粒、球粒灰岩及云质灰岩。
作为本发明一种优选的技术方案,所述二值图像成像精度的误差主要是:利用函数bwlabel对二值图像中孔隙与颗粒进行识别会产生误差,通过Matlab对60-100个样品的孔隙度进行测量,其中每个样品取三个铸体薄片进行均值计算,最终结果与高压压汞测试结果进行对比,分析误差。
作为本发明一种优选的技术方案,在应用Matlab内置函数bwlabel对粒径大小及孔径大小进行计算时,函数bwlabel会对所有的孔隙及图像中的噪点进行统计,并对该统计过程加以限定调节即最小识别范围,最小识别范围为2μm。
与现有技术相比,本发明的有益效果是:本发明通过Matlab定量计算不同碳酸盐孔隙结构,平均误差率为+4.661%,其中粒间孔平均误差率为-4.370%,铸模孔平均误差率为+6.734%,生屑孔误差率为+6.144%,微裂缝平均误差率为-0.314%;对于不同碳酸盐岩石孔隙而言,通过Matlab定量分析铸模孔与生屑孔的孔隙结构最为适用。
附图说明
图1为本发明Matlab中ColorThreshold色彩模式选择界面图;
图2为本发明Color Threshold功能界面图;
图3为本发明实施方式中岩石结构定量分析相关图像;
图4为本发明的流程图;
图5为本发明读取铸体薄片/SEM图像的方法;
图6为本发明实施方式中对不同类型颗粒的碳酸盐岩定量分析示意图;
图7为本发明实施方式中对不同碳酸盐岩孔隙类型定量分析示意图;
图8为本发明实施方式中不同孔隙类型孔隙度与高压压汞孔隙度示意图;
图9为本发明实施方式中孔隙度实测值与真实值差值比较示意图;
图10为本发明实施方式中不同孔隙类型误差率示意图;
图11为本发明实施方式中RGB、HSV、L*a*b*三种色彩空间对比示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
如图4所示,本发明提供了一种基于Matlab分析碳酸盐岩结构的方法,包括如下步骤:
步骤一、读取碳酸盐岩的铸体薄片;
步骤二、选择色彩模式,将碳酸盐岩的铸体薄片图像二值化;
步骤三、识别二值图像的黑色像素点,计算孔隙度,同时识别二值图像的白色像素点,分辨颗粒大小,计算粒度;
步骤四、求取碳酸盐岩的粒度与孔径大小频率分布图;
步骤五、求取碳酸盐岩粒度与孔径的平均直径及方差,并输出二值图像和求取的参数。
在本实施方式步骤一中,利用铸体薄片图像分析岩石学结构是常见的手段,但是应用软件Matlab对铸体薄片图像分析进行相比现有技术具有突出的特点,能定量计算不同碳酸盐孔隙结构,平均误差率为+4.661%,其中粒间孔平均误差率为-4.370%,铸模孔平均误差率为+6.734%,生屑孔误差率为+6.144%,微裂缝平均误差率为-0.314%,铸模孔与生屑孔的孔隙结构最为适用,粒间孔与微裂缝相对较差。
应用软件Matlab分析处理图像时,需将图像预处理为软件可识别的二值图像,图像二值化就是将图像上的像素点转化为0或255,二值图像中灰度只有两种,非黑即白,也就是将整个图像呈现出明显的黑白效果,以这种方式来操作图像可以更容易识别出图像的结构特征。
因此将铸体薄片图像导入Matlab的内置Color Threshod应用中进行图像二值化转换,共有成像模式,分别为RGB模式、HSV模式、YCbCr模式、L*a*b*模式,如图1所示。
在本实施方式中,不同成像模式适用于不同图像色彩空间,最终目的均为最大限度地识别岩石学结构并将其转为二值图像,HSV模式因其可以调控图像色调H、饱和度S以及明度V,对于呈现褐黄色的碳酸盐颗粒及蓝色充铸孔隙的色彩空间而言,HSV模式能够准确地识别碳酸盐岩图像岩石学结构,所以,步骤二应用HSV模式将图像转化为二值图像,如图2所示。
图2中,左侧图像中黄色为岩石孔隙空间,右侧H为色彩空间,控制图像整个范围,将呈现呈蓝色的孔隙空间从色度盘中去掉,S与V分别为饱和度与明度,稍作调整,最终调整后的整个图像的色彩空间在右下方呈现,空间中呈现孔隙的蓝色基本消失,通过ShowBinary将其转为二值图像输出至Matlab。
通过预处理获取二值图像(如图3c)后,进一步计算图像孔隙度、粒度、孔喉分布、平均孔隙直径等参数。本实施方式以孔隙度、粒度、平均孔隙直径为例。
(1)孔隙度的参数求取方法为:
通过Matlab中的size函数对二值图像中的像素进行识别提取,通过函数可以分别获得图像的总像点数Na(所有黑白像素点总和),再进一步调用能够识别0像素值的函数bwarea获取孔隙像素点数Np,即黑色像素点Np,如图3d,孔隙度P即黑色像素点与总像素点的比值即:
Figure GDA0003844263060000061
对于岩石粒度以及孔隙度来说,Matlab内置函数bwlabel可对相同像素值的区域进行识别统计,同时可以确定区域的边界,基于此原理确定孔隙区域(黑色或白色像素值连通区域)。
(2)粒度的求取方法:
首先,通过Matlab内置函数bwlabel对相同像素值的区域进行识别统计,同时可以确定区域的边界,基于此原理确定孔隙区域,即黑色或白色像素值连通区域;
然后对碳酸盐岩的每一个颗粒以及孔洞调用bwlabel函数,求出第i个区域的X和Y轴的最大值以及最小值即ximax、ximin、yimax、yimin,所以对于第i个颗粒或孔洞的直径大小di为X和Y方向的差值的平均值,即直径大小di
Figure GDA0003844263060000071
(3)平均孔隙直径的求取方法:
调用内部函数imhist对所有颗粒及孔洞进行叠加绘制粒度及孔径分布直方图,直方图可清晰地展示孔径直方图。最终分别调用函数mean以及函数var求取平均孔隙直径
Figure GDA0003844263060000074
以及方差S:
Figure GDA0003844263060000072
Figure GDA0003844263060000073
在本实施方式中,通过内置函数,内置函数imhist最终可以获得粒度与孔径分布直方图,进一步对孔隙平均值非均质性以及连通性进行分析,通过调用Matlab内置函数mean函数可对平均孔隙半径进行求取。在获取每个颗粒及孔洞的基础上调用内置函数var函数对所有孔洞的方差进行求取,平均孔隙半径
Figure GDA0003844263060000075
及其方差S。
对于地层非均质性来说,方差越大,表现整个地层孔隙分布越复杂,颗粒变化大,非均质性越强。对碳酸盐岩图像二值化预处理后,通过Matlab相关内置函数可以定量求取储层砂岩孔隙度、孔隙半径分布、平均孔隙度、非均质性等表征参数。具体流程如图4所示。
碳酸盐岩的结构组分主要由颗粒、孔隙、泥晶及亮晶胶结物组成,其中碳酸盐岩颗粒包括内碎屑、鲕粒、藻粒、球粒、团块、生物颗粒,本实施方式分别以内碎屑、鲕粒、球粒及生物颗粒这四种颗粒类型为例,定量分析孔径分布。同时,孔隙也由原生孔隙、次生孔隙两种,分别以原生粒间孔、晶间孔、生物骨架孔、铸模孔、收缩孔这五类孔隙为例,定量分析地层孔隙度、孔径分布及大小、平均孔径大小及方差。
本实施方式应用Matlab对内碎屑、鲕粒、球粒及生物颗粒这四种颗粒类型分析,对颗粒含量、颗粒分布、颗粒粒径大小及分布、方差进行定量求解,并研究其方法的适用性。
其中,内碎屑颗粒主要是沉积盆地中沉积不久的汐水流、风暴流、重力流等的作用,破碎半固结或固结的各种碳酸盐沉积物,受波浪、潮搬运、磨蚀、再沉积而成的。以卡沙干油田KE-4井为例,应用图4流程图,首先通过Threshold将铸体薄片转为二值图像,如图6a-2,通过size函数测得地层颗粒含量为81.631%,进一步利用bwlabel函数对颗粒进行识别计算,所识别孔隙如图6a-3,粒径频率分布直方图如图6a-4,可以看出粒径主要集中在50μm左右,颗粒分布相对集中,引用Matlab内置函数mean及var函数求解平均孔隙度及粒径方差。
结果如表1所示,内碎屑粒径平均值为70.891μm,方差为47.520,数值较低,颗粒分选性较强。铸体薄片中大部分颗粒能够被系统所识别,但部分颗粒会被分割为几小部分或一些粒径较小的颗粒被拼合成较大的颗粒,这些误差主要来源于函数bwlabel对颗粒边缘的识别,而控制颗粒边缘主要受所选择的色彩模式及相关参数的影响,不同的色彩模式及参数会对生成的二值图像产生不同的误差。对于内碎屑颗粒,局部颗粒边缘会产生噪点,被系统误识别为颗粒边缘进而分割为几部分,而同时一些点-线接触的颗粒未被识别出,被函数一同识别,最终造成颗粒的拼合。但总体对内碎屑颗粒而言,误差较小,适用性较高。
鲕粒颗粒是具有核心和同心层结构的球状颗粒,根据鲕粒的结构和形态特征,可把鲕粒划分为以下类型:正常鲕、表皮鲕、复鲕、椭形鲕、放射鲕及单晶鲕这几种类型。本实施方式选取塔里木盆地塔中4井亮晶鲕粒灰岩为例,对其颗粒含量、颗粒分布、颗粒粒径大小及分布、方差进行定量分析。对二值图像进行读取分析,最终颗粒分布图如图6b-3所示,与6b-1图像相比,大部分鲕粒被还原较好,但部分鲕粒被Matlab错误识别为两个以上的颗粒,这会导致粒径计算结果与实际值偏低,平均粒径大小偏低,方差增大,利用Matlab所计算粒径分布图结果如图6b-4,所测的颗粒含量、平均粒径及粒径方差如表1所示,由于鲕粒形状非常规则,系统识别度非常高。
球粒是较细粒的(粗粉砂级或砂级)、由灰泥组成的、不具特殊内部结构的、球形或卵形的、分选较好的颗粒。以塔里木盆地中4井亮晶-生屑球粒灰岩为例,应用图4程序流程图对球粒灰岩进行定量分析,二值图像见图6c-2,与图6c-1相比,部分海百合颗粒由于颜色与亮晶胶结物相似,被Matlab误认为颗粒以外的物质,导致平均粒径偏小,但由于海百合粒径与球粒粒径相差较多,最终粒径方差将比理论值偏低,分选性也会变高。球粒颗粒分布图见图6c-3,图中部分球粒颗粒被系统分割为几小部分,这也是误差来源之一,适用性较高。
生物颗粒以以鄂尔多斯盆地延628井奥陶系云质灰岩为例,地层晶间孔发育,颗粒属于自形-半自形结构,由Matlab内置函数bwlabel对二值图像进行分析,得到图6d-3颗粒分布图,图中左上方部分晶间孔不发育的区域,系统未识别出细小的晶间孔,一些颗粒相联通,颗粒粒径与原铸体薄片相比偏高,利用内置函数imhist得到孔径频率分布直方图,由mean函数计算得到平均粒径为88.752μm。
表1:图6中不同类型地层的颗粒比例/%、平均粒径/μm及方差
Figure GDA0003844263060000091
如图5所示,应用Matlab对不同颗粒结构定量分析,在整个流程中,会产生误差,包括选择色彩模式、铸体图像的二值转化、颗粒被分割或重组等,每一步地操作都会对下一步产生影响,误差在整个过程中会叠加,因此对于图像的预处理非常重要,包括色彩模式的选择以及二值图像的转化。
对于不同碳酸盐岩孔隙类型定量分析如下:
碳酸盐岩孔隙类型分为原生孔隙、次生孔隙、裂缝三部分,包括原生粒间孔、原生晶间孔、铸模孔、生物孔、微裂缝等几个部分,主要的碳酸盐岩的颗粒类型包括原生粒间孔,所以下主要讨论晶间孔、铸模孔、生屑孔、微裂缝的定量分析及程序适用性。
(1)晶间孔
以鄂尔多斯盆地东缘延1757井奥陶系白云岩晶间孔为例,白云岩晶体自形程度较好,利用图5方法处理过程对白云岩晶间孔进行定量分析,图7a-2为Matlab生成的二值图像,与铸体薄片相比,大部分孔隙被识别出,但同时一些附着在薄片的黑色脏污也被误认为孔隙,这是误差来源之一。此外,一些孔洞也被系统分割成几部分或拼合几部分,这也是另一部误差来源。但对晶间孔的识别非常准确,相关定量分析结果如表2所示。
(2)铸模孔
样品取自卡莎干油田KE-4井,颗粒灰岩早期的强烈的选择性溶蚀形成了大量的铸模孔,如图7b-1所示,铸模孔形状较为规则,利用Matlab进行识别并利用相关内置函数进行定量分析后,孔隙重合度较高,除极少数孔隙别分割为几部分。
(3)生屑骨架孔
铸体薄片泥晶颗粒灰岩取样于卡沙干油田KED04井,颗粒主要为鲕粒、生物碎屑、砂屑等,生物颗粒见有孔虫、红藻、棘皮类,有孔虫体腔内的软体发生腐烂之后形成生物体腔孔如图7c-1,利用Matlab对铸体薄片二值化,孔隙包括有孔虫体腔内孔、粒间孔都能够被识别出,但一些孔隙会被分解为几小部分,影响孔径分布及平均孔隙度,孔径分布相对中,孔隙度为15.37%,程序非常适用。
(4)微裂缝-缝合线
以鄂尔多斯盆地东缘延1758井奥陶系泥晶灰岩为例,Matlab能够很好识别出微裂缝,如图7d-2,利用size函数及公式
Figure GDA0003844263060000111
能够便捷的求出孔隙度,但利用函数bwlalel对孔隙识别时,泥晶中一些细小的微裂缝未被识别出,同时大量微裂缝被分为细小的颗粒,所以应用Matlab定量分析微裂缝-缝合线时,只适用于定量分析孔隙度,对孔径分布、大小并不适用。
表2:图7中不同地层孔隙类型的孔隙度/%、平均粒径/μm及方差
不同孔隙类型 晶间孔 铸模孔 生屑骨架孔 微裂缝+缝合线
孔隙度/% 9.13 15.37 8.1675 11.094
平均孔径/μm 33.7323 49.1569 24.1185 22.854
孔径方差 23.7345 36.3214 16.6671 16.1559
关于二值图像的成像精度,首先,在整个流程中,从铸体图像转化为二值图像,到利用函数bwlabel对二值图像中孔隙与颗粒进行识别均会产生一定的误差,本实施方式通过Matlab对60个样品的孔隙度进行测量,其中每个样品取三个铸体薄片进行均值,最终结果与高压压汞测试结果进行对比,分析误差,结果见图8。
图8中,红色代表利用Matlab的测量结果,蓝色代表高压压汞孔隙度,四种被圈定的孔隙类型分别为粒间孔、铸模孔、生屑孔及微裂缝。测量值与真实值存在一定偏差,不同类型的偏差值也不一样,偏差结果如图9所示。
由图9可以看出,铸模孔与生屑孔两种孔隙类型测量值大部分大于真实值,这主要是在图像二值化过程中,由于色彩空间模式的选择及参数调控所产生的噪点引起的,噪点被Matlab误认为孔隙,由此造成孔隙度比真实值有所增加。粒间孔与微裂缝两种孔隙类型的一些颗粒同铸模孔与生屑孔一样,实际测量值有所偏大,但多数测量值低于理论值,这主要归根于粒间孔与微裂缝自身特征,由于一些粒间孔与微裂缝孔隙较窄,半径大小为纳米级孔隙,并未被Matlab识别,因此形成该种现象,具体误差率见图10。
大部分孔隙误差集中在-5%~10%之间,平均误差率为+4.661%,其中粒间孔平均误差率为-4.370%,标准差为8.818,铸模孔平均误差率为+6.734%,标准差为6.913,生屑孔误差率为+6.144%,标准差为6.731,微裂缝平均误差率为-0.314%,标准差为+9.271。所以对于不同碳酸盐岩石孔隙而言,通过Matlab定量分析铸模孔与生屑孔的孔隙结构最为适用,粒间孔与微裂缝相对较差。
在本实施方式中,图像二值化的精准度直接影响颗粒与孔隙分布,对定量分析产生最直接的影响。而图像二值化包括色彩空间的选择及参数的调控,Color Threshold中共有四种色彩空间模式RGB、HSV、YCbCr、L*a*b*,不同色彩空间模式适用于不同的色域的铸体薄片,以最常见的“黄蓝色”铸体薄片为例(黄色为碳酸盐岩干涉色,蓝色为孔隙被染色),分别以RGB、HSV、L*a*b*三种常规色彩空间模式对图像进行二值化处理,见图11。
对于RGB色彩模式而言,r、g、b三个变量并不能很好的消除色彩空间中的蓝色,调节变量,将会对色彩空间进行线性变换,导致图像中杂基、暗色矿物产生大量噪点,L*a*b*模式中L表示亮度(Luminosity),a表示从洋红色至绿色的范围,b表示从黄色至蓝色的范围,虽然控制参数b可以很好地消除图像中的蓝色,但由于同样是线性变换,一些代表颗粒的黄色也会受到波及,导致颗粒中产生大量噪点。对于HSV模式,H控制的色度盘可以轻易地将蓝色去除,为空间中的非线性变换,并不会影响其它色彩,基本无噪点产生,因此也最适用于碳酸盐岩铸体薄片。
在本实施方式中,应用Matlab对碳酸盐岩石结构定量分析的整个流程中,图像噪点是对计算结果最大的影响,因此去除噪点是必要的一部,应用Matlab内置函数bwlabel对粒径大小及孔径大小进行计算时,函数会对所有的孔隙及图像中的噪点进行统计,因此对于这个统计过程加以限定调节即最小识别范围,最小识别范围是根据研究图像而定的,本实施方式中的限定界限为2μm,但在对噪点去除的同时,一些微裂缝及较小的粒间孔也会被剔除,因此,对于粒径较低的碳酸盐岩如泥晶灰岩及微裂缝发育的碳酸盐岩并不适用。
在本实施方式中,通过Matlab定量计算不同碳酸盐孔隙结构,平均误差率为+4.661%,其中粒间孔平均误差率为-4.370%,铸模孔平均误差率为+6.734%,生屑孔误差率为+6.144%,微裂缝平均误差率为-0.314%。
对于不同碳酸盐岩石孔隙而言,通过Matlab定量分析铸模孔与生屑孔的孔隙结构最为适用,粒间孔与微裂缝相对较差。对于粒径较低的碳酸盐岩如泥晶灰岩及微裂缝发育的碳酸盐岩并不适用。
对于本领域技术人员而言,显然本发明不限于上述示范性实施例的细节,而且在不背离本发明的精神或基本特征的情况下,能够以其他的具体形式实现本发明。因此,无论从哪一点来看,均应将实施例看作是示范性的,而且是非限制性的,本发明的范围由所附权利要求而不是上述说明限定,因此旨在将落在权利要求的等同要件的含义和范围内的所有变化囊括在本发明内。不应将权利要求中的任何附图标记视为限制所涉及的权利要求。

Claims (8)

1.一种基于Matlab分析碳酸盐岩结构的方法,其特征在于,包括如下步骤:
步骤100、读取碳酸盐岩的铸体薄片;
步骤200、选择色彩模式,将碳酸盐岩的铸体薄片图像二值化;
步骤300、识别二值图像的黑色像素点,计算孔隙度,同时识别二值图像的白色像素点,分辨颗粒大小,计算粒度;
步骤400、求取碳酸盐岩的粒度与孔径大小频率分布图;
步骤500、求取碳酸盐岩粒度与孔径的平均直径及方差,并输出二值图像和求取的参数;
所述二值图像成像精度的误差分析方法包括:
利用函数bwlabel对二值图像中孔隙与颗粒进行识别会产生误差,通过Matlab对60-100个样品的孔隙度进行测量,其中每个样品取三个铸体薄片进行均值计算,最终结果与高压压汞测试结果进行对比,分析误差;
在应用Matlab内置函数bwlabel对粒径大小及孔径大小进行计算时,函数bwlabel会对所有的孔隙及图像中的噪点进行统计,并对该统计过程加以限定调节即最小识别范围,最小识别范围为2μm。
2.根据权利要求1所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:所述步骤200中,铸体薄片图像二值化的方法为:将铸体薄片图像导入Matlab的内置ColorThreshod应用中通过HSV成像模式进行图像二值化转换。
3.根据权利要求2所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:不同的成像模式适用于不同图像色彩空间,最终目的均为最大限度地识别碳酸盐岩结构并将其转为二值图像,所述HSV成像模式因其可以调控图像色调H、饱和度S以及明度V。
4.根据权利要求1所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:所述孔隙度的求取方法为:通过Matlab中的size函数对二值图像中的像素进行识别提取,通过函数可以分别获得图像的总像点数Na,再进一步调用能够识别0像素值的函数bwarea获取孔隙像素点数Np,即黑色像素点Np,孔隙度P即黑色像素点与总像素点的比值即:
Figure DEST_PATH_IMAGE002
5.根据权利要求1所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:所述粒度的求取方法包括:
首先,通过Matlab内置函数bwlabel对相同像素值的区域进行识别统计,同时可以确定区域的边界,基于此原理确定孔隙区域,即黑色或白色像素值连通区域;
然后对碳酸盐岩的每一个颗粒以及孔洞调用bwlabel函数,求出第i个区域的X和Y轴的最大值以及最小值即
Figure DEST_PATH_IMAGE004
Figure DEST_PATH_IMAGE006
Figure DEST_PATH_IMAGE008
Figure DEST_PATH_IMAGE010
,所以对于第i个颗粒或孔洞的直径大小
Figure DEST_PATH_IMAGE012
为X和Y方向的差值的平均值,即直径大小
Figure DEST_PATH_IMAGE012A
Figure DEST_PATH_IMAGE014
6.根据权利要求1所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:所述步骤400中,碳酸盐岩的粒度与孔径大小频率分布图通过调用内部函数imhist对所有颗粒及孔洞进行叠加进行绘制。
7.根据权利要求5所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:所述碳酸盐岩粒度与孔径的平均直径及方差的求取方法为:
分别调用函数mean以及函数var求取平均孔隙直径
Figure DEST_PATH_IMAGE016
以及方差S:
Figure DEST_PATH_IMAGE018
Figure DEST_PATH_IMAGE020
8.根据权利要求1所述的一种基于Matlab分析碳酸盐岩结构的方法,其特征在于:该方法适用的颗粒包括内碎屑颗粒、鲕粒颗粒、球粒灰岩及云质灰岩。
CN201810571924.5A 2018-06-06 2018-06-06 一种基于Matlab分析碳酸盐岩结构的方法 Expired - Fee Related CN108734714B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810571924.5A CN108734714B (zh) 2018-06-06 2018-06-06 一种基于Matlab分析碳酸盐岩结构的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810571924.5A CN108734714B (zh) 2018-06-06 2018-06-06 一种基于Matlab分析碳酸盐岩结构的方法

Publications (2)

Publication Number Publication Date
CN108734714A CN108734714A (zh) 2018-11-02
CN108734714B true CN108734714B (zh) 2022-11-25

Family

ID=63932055

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810571924.5A Expired - Fee Related CN108734714B (zh) 2018-06-06 2018-06-06 一种基于Matlab分析碳酸盐岩结构的方法

Country Status (1)

Country Link
CN (1) CN108734714B (zh)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109856029B (zh) * 2019-02-01 2021-07-30 中海石油(中国)有限公司上海分公司 一种基于图像分析的孔隙度评价方法
CN111504875B (zh) * 2020-04-28 2021-06-04 中国地质大学(北京) 一种致密砂岩孔喉参数的提取及计算方法及提取装置
CN111505038A (zh) * 2020-04-28 2020-08-07 中国地质大学(北京) 一种基于阴极发光技术定量分析砂岩胶结的实现方法
CN111413284B (zh) * 2020-05-18 2022-03-22 攀钢集团攀枝花钢铁研究院有限公司 钒钛烧结矿中主要物相的定量检测方法
CN112433069A (zh) * 2020-12-01 2021-03-02 安徽理工大学 一种奥陶系灰岩顶部相对隔水层的微观结构特征判别方法
CN113670791A (zh) * 2021-08-03 2021-11-19 中国地质大学(北京) 一种对储层孔隙及颗粒表面积的定量分析方法
CN113870200B (zh) * 2021-09-16 2022-05-31 中国地质大学(北京) 一种定量分析砂岩与碳酸盐岩储层中孔隙配位数的方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102749046A (zh) * 2012-07-23 2012-10-24 中国地质大学(武汉) 岩体结构面直剪试验中剪切面积的测量方法
CN106896213A (zh) * 2017-02-22 2017-06-27 中国地质大学(武汉) 一种基于点云数据的岩体结构面智能识别与信息提取方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6937744B1 (en) * 2000-06-13 2005-08-30 Microsoft Corporation System and process for bootstrap initialization of nonparametric color models
US20160203591A1 (en) * 2015-01-09 2016-07-14 Umm Al-Qura University System and process for monitoring the quality of food in a refrigerator

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102749046A (zh) * 2012-07-23 2012-10-24 中国地质大学(武汉) 岩体结构面直剪试验中剪切面积的测量方法
CN106896213A (zh) * 2017-02-22 2017-06-27 中国地质大学(武汉) 一种基于点云数据的岩体结构面智能识别与信息提取方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
分析碳酸盐岩孔隙系统数字图像的新方法;赵永刚 等.;《计算机应用研究》;20061001;第169-171页 *
应用MATLAB提取纳米模板特征几何参数;雷惊雷;《功能材料》;20111020;第1893-1897页 *

Also Published As

Publication number Publication date
CN108734714A (zh) 2018-11-02

Similar Documents

Publication Publication Date Title
CN108734714B (zh) 一种基于Matlab分析碳酸盐岩结构的方法
CN110458119B (zh) 一种非接触式测量的混凝土骨料级配快速识别方法
CN107228816B (zh) 一种泥页岩中不同类型孔隙演化评价方法
CN107133630B (zh) 一种基于扫描图像判断碳酸盐岩孔隙类型的方法
CN110264466A (zh) 一种基于深度卷积神经网络的钢筋检测方法
CN108961246B (zh) 一种基于人工智能的扫描电镜图像孔隙识别方法
CN102494976B (zh) 一种超细晶粒钢晶粒的自动测量及其形态分类方法
Quan et al. The method of the road surface crack detection by the improved Otsu threshold
CN104237103A (zh) 一种孔隙连通性定量表征方法及装置
CN114067114B (zh) 一种基于藻类细胞形态的面积核分割计数方法
CN112150430B (zh) 一种利用岩石细观结构数字图像的数值分析方法
CN107314957B (zh) 一种岩石块度尺寸分布的测量方法
CN112051200B (zh) 一种致密砂岩储层孔隙结构的定量评价方法和装置
CN104268600B (zh) 一种基于Minkowski距离的矿物浮选泡沫图像纹理分析及工况识别方法
YANG et al. Bubble size estimation using interfacial morphological information for mineral flotation process monitoring
CN110443793A (zh) 一种沥青混合料空隙分布均匀性评价方法
CN104318564A (zh) 一种矿物颗粒分相的方法
CN103971367B (zh) 水文资料图像分割方法
CN108986077A (zh) 基于双树复小波域共生增广矩阵的浮选泡沫工况识别方法
CN107220946B (zh) 一种岩石运输带上不良块度图像的实时剔除方法
CN113570652B (zh) 基于sem图像的砂岩储层矿物晶间孔的定量分析方法
CN112164029B (zh) 一种用于加气混凝土孔结构的图像分析方法
CN107977968A (zh) 基于建筑物阴影信息挖掘的建筑物分层检测方法
CN115359003A (zh) 两步式隧道灰度图像的裂缝识别方法、系统、介质及设备
CN115640546A (zh) 一种图像与特征信息融合的岩性识别方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20221125