CN113570652B - 基于sem图像的砂岩储层矿物晶间孔的定量分析方法 - Google Patents
基于sem图像的砂岩储层矿物晶间孔的定量分析方法 Download PDFInfo
- Publication number
- CN113570652B CN113570652B CN202110760753.2A CN202110760753A CN113570652B CN 113570652 B CN113570652 B CN 113570652B CN 202110760753 A CN202110760753 A CN 202110760753A CN 113570652 B CN113570652 B CN 113570652B
- Authority
- CN
- China
- Prior art keywords
- pixel
- pore
- image
- value
- radius
- 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.)
- Active
Links
- 239000011148 porous material Substances 0.000 title claims abstract description 193
- 238000000034 method Methods 0.000 title claims abstract description 48
- 238000001878 scanning electron micrograph Methods 0.000 title claims abstract description 40
- 238000004445 quantitative analysis Methods 0.000 title claims abstract description 18
- 229910052500 inorganic mineral Inorganic materials 0.000 title claims abstract description 16
- 239000011707 mineral Substances 0.000 title claims abstract description 16
- 239000011159 matrix material Substances 0.000 claims abstract description 38
- 238000012545 processing Methods 0.000 claims abstract description 14
- 238000005520 cutting process Methods 0.000 claims abstract description 8
- 230000001939 inductive effect Effects 0.000 claims abstract description 6
- 230000008439 repair process Effects 0.000 claims abstract description 6
- 238000010606 normalization Methods 0.000 claims description 16
- 239000012634 fragment Substances 0.000 claims description 15
- 230000035772 mutation Effects 0.000 claims description 11
- 238000005070 sampling Methods 0.000 claims description 11
- 230000008859 change Effects 0.000 claims description 9
- 238000009826 distribution Methods 0.000 claims description 9
- 238000013507 mapping Methods 0.000 claims description 8
- 230000008569 process Effects 0.000 claims description 8
- 230000008030 elimination Effects 0.000 claims description 7
- 238000003379 elimination reaction Methods 0.000 claims description 7
- 239000011800 void material Substances 0.000 claims description 4
- 230000000877 morphologic effect Effects 0.000 claims description 3
- 239000002734 clay mineral Substances 0.000 abstract description 27
- 238000005259 measurement Methods 0.000 abstract description 3
- 238000004364 calculation method Methods 0.000 abstract description 2
- 238000001035 drying Methods 0.000 abstract 1
- 238000001000 micrograph Methods 0.000 description 5
- 229910001919 chlorite Inorganic materials 0.000 description 4
- 229910052619 chlorite group Inorganic materials 0.000 description 4
- QBWCMBCROVPCKQ-UHFFFAOYSA-N chlorous acid Chemical compound OCl=O QBWCMBCROVPCKQ-UHFFFAOYSA-N 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 4
- NLYAJNPCOHFWQQ-UHFFFAOYSA-N kaolin Chemical compound O.O.O=[Al]O[Si](=O)O[Si](=O)O[Al]=O NLYAJNPCOHFWQQ-UHFFFAOYSA-N 0.000 description 4
- 229910052622 kaolinite Inorganic materials 0.000 description 4
- 229910052901 montmorillonite Inorganic materials 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 239000012530 fluid Substances 0.000 description 3
- 238000011161 development Methods 0.000 description 2
- QSHDDOUJBYECFT-UHFFFAOYSA-N mercury Chemical compound [Hg] QSHDDOUJBYECFT-UHFFFAOYSA-N 0.000 description 2
- 229910052753 mercury Inorganic materials 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 239000011435 rock Substances 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000005266 casting Methods 0.000 description 1
- 239000003245 coal Substances 0.000 description 1
- 239000003086 colorant Substances 0.000 description 1
- 239000013078 crystal Substances 0.000 description 1
- GUJOJGAPFQRJSV-UHFFFAOYSA-N dialuminum;dioxosilane;oxygen(2-);hydrate Chemical compound O.[O-2].[O-2].[O-2].[Al+3].[Al+3].O=[Si]=O.O=[Si]=O.O=[Si]=O.O=[Si]=O GUJOJGAPFQRJSV-UHFFFAOYSA-N 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 229910052900 illite Inorganic materials 0.000 description 1
- 230000006698 induction Effects 0.000 description 1
- 238000000691 measurement method Methods 0.000 description 1
- VGIBGUSAECPPNB-UHFFFAOYSA-L nonaaluminum;magnesium;tripotassium;1,3-dioxido-2,4,5-trioxa-1,3-disilabicyclo[1.1.1]pentane;iron(2+);oxygen(2-);fluoride;hydroxide Chemical compound [OH-].[O-2].[O-2].[O-2].[O-2].[O-2].[F-].[Mg+2].[Al+3].[Al+3].[Al+3].[Al+3].[Al+3].[Al+3].[Al+3].[Al+3].[Al+3].[K+].[K+].[K+].[Fe+2].O1[Si]2([O-])O[Si]1([O-])O2.O1[Si]2([O-])O[Si]1([O-])O2.O1[Si]2([O-])O[Si]1([O-])O2.O1[Si]2([O-])O[Si]1([O-])O2.O1[Si]2([O-])O[Si]1([O-])O2.O1[Si]2([O-])O[Si]1([O-])O2.O1[Si]2([O-])O[Si]1([O-])O2 VGIBGUSAECPPNB-UHFFFAOYSA-L 0.000 description 1
- 239000003208 petroleum Substances 0.000 description 1
- 238000003825 pressing Methods 0.000 description 1
- 239000004575 stone Substances 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
- 238000012800 visualization Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/60—Analysis of geometric attributes
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume, or surface-area of porous materials
- G01N15/08—Investigating permeability, pore-volume, or surface area of porous materials
- G01N15/088—Investigating volume, surface area, size or distribution of pores; Porosimetry
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N23/00—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00
- G01N23/22—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by measuring secondary emission from the material
- G01N23/225—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by measuring secondary emission from the material using electron or ion
- G01N23/2251—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by measuring secondary emission from the material using electron or ion using incident electron beams, e.g. scanning electron microscopy [SEM]
-
- G06T5/70—
-
- G06T5/77—
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/136—Segmentation; Edge detection involving thresholding
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/181—Segmentation; Edge detection involving edge growing; involving edge linking
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10056—Microscopic image
- G06T2207/10061—Microscopic image from scanning electron microscope
Abstract
本发明公开了一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，获取每种砂岩储层样品在电镜扫描下的SEM图像，将SEM图像切割为多个图像矩阵，利用自适应阈值方式对每个图像矩阵分别进行二值化处理；按照矩阵遍历方式识别每种砂岩储层样品对应图像矩阵中所有像素点对应的像素值，统一归纳二值图像的矩形阵列的像素值以同步完成修复噪点和消除噪点；重新除燥之后的二值图像，将像素值突变的像素点作为孔隙边界划分二值图像中的孔隙单元；确定每个像素点在二维坐标系内的坐标值，基于坐标值计算二值图像中的每个孔隙单元的孔隙半径以及喉道半径；本发明实现对SEM图像中各类粘土矿物的孔隙与喉道参数进行定量计算，且测量精度小于1微米的尺度。
Description
技术领域
本发明涉及砂岩储层粘土矿物晶间孔技术领域，具体涉及一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法。
背景技术
自生粘土矿物是砂岩中重要的填隙物，也是砂岩成岩阶段划分的重要指示矿物。该类矿物在油气勘探开发过程中备受关注，随着石油工业的大发展，砂岩中自生粘土矿物的研究，目前已取得了长足进展，孔隙是流体赋存于岩石中的基本储集空间，而喉道则是控制流体在岩石中渗流的重要通道。流体在自然界复杂的孔隙系统中流动时，都要经历一系列的交替着的孔隙和喉道，目前而言，大量试验集中于粘土矿物含量及孔隙形态学的测定与分析，砂岩中粘土矿物包括伊利石、绿泥石、蒙脱石、高岭石等。对于粘土矿物晶体内孔隙与喉道的分布及半径大小研究较少，尤其是针对单独的粘土矿物。
现存主要方法包括间接测量法，直接观察法以及数字演心法。间接测量法以高压压汞法为主，它能够针对砂岩的粘土矿物孔隙与喉道半径做出精确测量，但是无法对各类型粘土矿物的孔隙与喉道进行区分和测定；以铸体薄片或者扫描电镜的直接观察法能定性的对粘土矿物的孔隙与孔喉进行研究，但是无法实现定量分析；以X-CT数字岩芯孔隙结构三维模型重构技术为基础的数字岩芯技术能对各类孔隙喉道半径进行定量分析，但该方法非常高贵，费时，目前整体成本较高。
发明内容
本发明的目的在于提供一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，以解决现有技术中无法对各类型粘土矿物的孔隙与喉道进行区分和测定，无法实现对孔隙或者吼道半径定量分析的技术问题。
为解决上述技术问题，本发明具体提供下述技术方案：
一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，包括以下步骤：
步骤100、准备不同孔隙类型的砂岩储层样品，获取每种砂岩储层样品在电镜扫描下的SEM图像，将所述SEM图像切割为多个图像矩阵，且利用自适应阈值方式对每个所述图像矩阵分别进行二值化处理；
步骤200、按照矩阵遍历方式识别每种所述砂岩储层样品对应所述图像矩阵中所有像素点对应的像素值，统一归纳所述二值图像的矩形阵列的像素值以同步完成修复噪点和消除噪点；
步骤300、重新按照行列遍历消除噪点和修复噪点之后的所述二值图像，将像素值突变的像素点作为孔隙边界，并且基于孔隙边缘划分所述二值图像中的孔隙单元；
步骤400、确定每个像素点在所述二维坐标系内的坐标值，基于坐标值计算所述二值图像中的每个所述孔隙单元的孔隙半径以及喉道半径，并输出每种所述砂岩储层样品的孔隙半径分布直方图、平均孔隙半径、喉道半径直方图与平均喉道半径。
作为本发明的一种优选方案，在步骤100中，将所述SEM图像转化二值图像的具体实现步骤为：
将所述SEM图像转化为矩阵图像，并对所述矩阵图像平均切割为多个图像碎片；
抽样采集每个所述图像碎片中多个像素点的像素值，且确定对每个所述图像碎片二值化处理对应的分界像素值；
将每个所述图像碎片中低于分界像素值的像素点的像素值设置为0，且将高于分界像素值的像素点的像素值设置为255；
将所有的图像碎片按照行列式对齐，形成二值化处理后的二值图像。
作为本发明的一种优选方案，确定分界像素值的实现方式为：
将抽样采集的列间距和行间距形成抽样网格，并获取所述抽样网格的每个交叉点对应的像素值；
求取所有像素点的平均像素值，且将所述平均像素值作为分界像素值。
作为本发明的一种优选方案，在步骤300中，消除噪点和修复噪点可同时实现，具体的方法为：
利用归纳方式将每个自适应二值化处理的图像矩阵重新切割为多个归一化矩阵，统计每个所述归一化矩阵中像素值为0的像素点个数以及像素值为255的像素点个数；
将个数超过所述归一化矩阵中所有像素点数量一半的像素点作为多数元素，并将该归一化矩阵的所有像素点均重置为多数元素对应的像素值，以同步实现消除噪点和修复噪点。
作为本发明的一种优选方案，将每个所述图像矩阵拼接形成消除噪点和修复噪点后的二值图像，且重新遍历所述二值图像，以所述二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历所述二值图像以确定每个所述像素点的像素值，选定行列遍历时的遍历间距均为一个像素，确定出每行中像素值突变的像素点，且将像素值突变的像素点保存在一个集合内，对两个所述集合进行多次方数字形态学运算，将两个所述集合之间的像素点的像素值重置为0以实现二次修复噪点。
作为本发明的一种优选方案，对两个所述集合进行多次方数字形态学运算具体为：计算两个集合中同行或者同列像素点之间的邻域像素，并在所述邻域像素小于设定值时，将像素值为255的像素点设置为0以将两个集合对应的孔隙连接。
作为本发明的一种优选方案，在步骤300中，识别二次修复噪点的所述二值图像的孔隙边缘和孔隙单元的实现方法为：
步骤301、以所述二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历所述二值图像以确定每个所述像素点的像素值，选定行列遍历时的遍历间距均为一个像素；
步骤302、将每个像素点的像素值与同一列的前一个像素点的像素值进行对比，确定像素值突变的所述像素点且重新定义该像素点的RGB值；
步骤303、将所有重新定义RGB值的像素点首尾连接构成的曲线设定为所述孔隙边缘，且将所有重新定义RGB值的像素点包围的区域设定为孔隙单元。
作为本发明的一种优选方案，在步骤302中，确定遍历过程中的每一行的第一个像素值突变的所述像素点，按照遍历顺序划分每一行的像素值突变的所述像素点的优先级，并且将不同行的像素值突变的所述像素点的优先级设定为与遍历顺序一一对应，并且将同一行的像素值突变的所述像素点设置为同一个优先级；
将仅有一个像素值突变的所述像素点按照优先级顺序分别设置为一个所述孔隙单元的起点和终点，按照从起点到终点的顺序分为两条测绘曲线，两条所述测绘曲线顺次将不同优先级的所述像素值突变的所述像素点依次连接形成一个所述孔隙单元。
作为本发明的一种优选方案，在步骤400中，所述孔隙半径为单个所述孔隙单元的长轴与短轴的平均值的一半，所述喉道半径为单个所述孔隙单元的短轴，所述二值图像中三次划定的每个所述孔隙单元的所述孔隙半径利用每个所述孔隙单元的像素点在所述二维坐标系的最大X坐标，最小X坐标，最大Y坐标，最小Y坐标进行求取，具体的实现步骤为：
选择每个所述孔隙单元内像素值突变的所述像素点的横坐标最大值对应的像素点、横坐标最小值对应的像素点，纵坐标最大值的像素点和坐标表最小值的像素点；
每个所述孔隙单元的长轴半径为Max[(Xmax-Xmin)，(Ymax-Ymin)]，而每个所述孔隙单元的短轴半径为Min[(Xmax-Xmin)，(Ymax-Ymin)]；
将通过公式[(Xmax-Xmin)+(Ymax-Ymin)]/2计算每个所述孔隙单元的孔隙半径；
将Min[(Xmax-Xmin)，(Ymax-Ymin)]筛选的数据作为每个所述孔隙单元的喉道半径。
作为本发明的一种优选方案，在步骤400中，统计每种砂岩储层中的每个孔隙单元对应的所述孔隙半径，并确定孔隙半径相同的所述孔隙单元的数量，创建关于孔隙半径和孔隙单元数量的孔隙半径直方图，并计算每种砂岩储层样品的孔隙半径平均值；
统计每种砂岩储层中的每个孔隙单元对应的所述喉道半径，并确定喉道半径相同的所述孔隙单元的数量，创建关于喉道半径和孔隙单元数量的孔隙半径直方图，并计算每种砂岩储层样品的喉道半径平均值。
本发明与现有技术相比较具有如下有益效果：
本发明提供的定量分析方法结合直接观察法(扫描电镜图像)的可视化观察与数字图像技术，实现对砂岩矿物晶间孔的定量分析，且能够详细地对SEM图像中各类粘土矿物的孔隙与喉道参数进行定量计算，测量精度达到小于1微米的尺度。
附图说明
为了更清楚地说明本发明的实施方式或现有技术中的技术方案，下面将对实施方式或现有技术描述中所需要使用的附图作简单地介绍。显而易见地，下面描述中的附图仅仅是示例性的，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据提供的附图引伸获得其它的实施附图。
图1为本发明实施例提供的砂岩储层粘土矿物晶间孔定量分析的流程示意图；
图2为本发明实施例提供的砂岩储层粘土矿物晶间孔的SEM图像；
图3为本发明实施例提供的砂岩储层粘土矿物晶间孔的未降噪二值图像；
图4为本发明实施例提供的砂岩储层粘土矿物晶间孔的降噪二值图像；
图5为本发明实施例提供的砂岩储层粘土矿物晶间孔的孔隙半径直方图；
图6为本发明实施例提供的砂岩储层粘土矿物晶间孔的喉道半径直方图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
如图1所示，本发明提供了一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，本实施方式将直接观察法与算法相结合，通过将MATLAB中各类函数与编程算法的结合实现对不同类型的扫描电镜图像(SEM图像)中的粘土矿物的孔隙与喉道分布和半径进行定量研究，建立数字图像-函数的分析方法。
包括以下步骤：
步骤100、准备不同孔隙类型的砂岩储层样品，获取每种砂岩储层样品在电镜扫描下的SEM图像，将所述SEM图像切割为多个图像矩阵，且利用自适应阈值方式对每个所述图像矩阵分别进行二值化处理。
在步骤100中，将所述SEM图像转化二值图像的具体实现步骤为：
先将所述SEM图像转化为矩阵图像，并对所述矩阵图像平均切割为多个图像碎片。
再抽样采集每个所述图像碎片中多个像素点的像素值，且确定对每个所述图像碎片二值化处理对应的分界像素值，其中确定分界像素值的实现方式为：其后将抽样采集的列间距和行间距形成抽样网格，并获取所述抽样网格的每个交叉点对应的像素值，求取所有像素点的平均像素值，且将所述平均像素值作为分界像素值。
最后将每个所述图像碎片中低于分界像素值的像素点的像素值设置为0，且将高于分界像素值的像素点的像素值设置为255，将所有的图像碎片按照行列式对齐，形成二值化处理后的二值图像。
图像二值化就是将SEM图像上所有像素点的灰度值设置为0或255，也就是将整个图像呈现出明显的黑白效果的过程。生成的二值图像能够被计算机或者MATLAB所识别，二值图像的像素分布最终以矩阵的形式分部。图像二值化的基本原理主要是以设定一个阈值为依托，即将256个亮度等级的灰度图像选取适当的阈值，从而获得反映图像整体和局部特征的二值化图像。
由于SEM图像是一种灰度图像，自身并不像其他彩色图像如rgb图像(包含红绿蓝三个色彩参数)，SEM图像转换为二值图像的常规方式比较简单，选择一个像素值，将大于该像素值的像素点统一设置为255，小于该像素值的像素点统一设置为0。而本方法先对SEM图像进行分割图像像素，然后对每个分割图像进行多尺度灰度等级划分并区分每个分割图像的等级赋值，从而实现对整个SEM图像的多尺度渐变赋值。
步骤200、按照矩阵遍历方式识别每种所述砂岩储层样品对应所述图像矩阵中所有像素点对应的像素值，统一归纳所述二值图像的矩形阵列的像素值以同步完成修复噪点和消除噪点。
在本步骤中，消除噪点和修复噪点可同时实现，具体的方法为：
①利用归纳方式将每个自适应二值化处理的图像矩阵重新切割为多个归一化矩阵，统计每个所述归一化矩阵中像素值为0的像素点个数以及像素值为255的像素点个数。
②将个数超过所述归一化矩阵中所有像素点数量一半的像素点作为多数元素，并将该归一化矩阵的所有像素点均重置为多数元素对应的像素值，以同步实现消除噪点和修复噪点，具体的实现方法为：
扫描电镜方法由于具有超强的放大倍数，能够观察至纳米级孔隙。在这个过程中，受高压影响，SEM图像在拍摄过程中会产生一定数量的噪点，这对于纳米级与微米级孔隙的影响是不可忽视的。因此需要人为对获取的二值图像进行降噪。降噪包括受噪点影响而导致孔隙发生分割的区域进行修复与对图像中噪点的去除。
其中需要特别说明的是，当一个序列中存在一个占大多数量的元素(即该元素在序列中的数量超过50％)，比如说在一个归一化矩阵的开始位置获得的像素点的像素值为255，则以该像素值作为对比元素，如果在3×3的归一化阵列中累加像素值为255的像素点，如果存在至少有5个像素点的像素值为255，则统一3×3的归一化阵列中所有像素点的像素值均为255；否则将所述归一化矩阵的所有像素点的像素值置于0。
将每个所述图像矩阵拼接形成消除噪点和修复噪点后的二值图像，且重新遍历所述二值图像，以所述二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历所述二值图像以确定每个所述像素点的像素值，选定行列遍历时的遍历间距均为一个像素，确定出每行中像素值突变的像素点，且将像素值突变的像素点保存在一个集合内，对两个所述集合进行多次方数字形态学运算，将两个所述集合之间的像素点的像素值重置为0以实现二次修复噪点，其中，将每个像素值突变的像素点与上一行像素值突变的像素点的行坐标进行对比，如果此时的像素点的行坐标有两个，且与上一行像素值突变的像素点的行坐标的差距为1，那么将此时遍历到的像素点与上一行的像素点保存在同一个集合内。
对两个所述集合进行多次方数字形态学运算具体为：计算两个集合中同行或者同列像素点之间的邻域像素，并在所述邻域像素小于设定值时，将像素值为255的像素点设置为0以将两个集合对应的孔隙连接。
步骤300、重新按照行列遍历消除噪点和修复噪点之后的所述二值图像，将像素值突变的像素点作为孔隙边界，并且基于孔隙边缘划分所述二值图像中的孔隙单元。
在步骤300中，识别二次修复噪点的所述二值图像的孔隙边缘和孔隙单元的实现方法为：
步骤301、以所述二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历所述二值图像以确定每个所述像素点的像素值，选定行列遍历时的遍历间距均为一个像素。
步骤302、将每个像素点的像素值与同一列的前一个像素点的像素值进行对比，确定像素值突变的所述像素点且重新定义该像素点的RGB值，在步骤302中，确定遍历过程中的每一行的第一个像素值突变的所述像素点，按照遍历顺序划分每一行的像素值突变的所述像素点的优先级，并且将不同行的像素值突变的所述像素点的优先级设定为与遍历顺序一一对应，并且将同一行的像素值突变的所述像素点设置为同一个优先级。
将仅有一个像素值突变的所述像素点按照优先级顺序分别设置为一个所述孔隙单元的起点和终点，按照从起点到终点的顺序分为两条测绘曲线，两条所述测绘曲线顺次将不同优先级的所述像素值突变的所述像素点依次连接形成一个所述孔隙单元。
步骤303、将所有重新定义RGB值的像素点首尾连接构成的曲线设定为所述孔隙边缘，且将所有重新定义RGB值的像素点包围的区域设定为孔隙单元。
即使每种砂岩储层的孔隙单元的形状不同，但是行列方式遍历每种煤储层对应的二值图像时，每个孔隙空间个体均存在仅有单个像素值突变的像素点，在本实施方式中，当遍历一行中仅存在单个像素值突变的像素点，则将该像素点作为孔隙空间个体的起点，另一个存在单个像素值突变的像素点作为一个孔隙空间个体的终点，因此将行列遍历得到的单个像素值突变的像素点作为一个孔隙空间个体的起点，并将像素点的行列坐标保存在一个集合内。
将下一行遍历得到的单个像素值突变的像素点与集合内的像素点行列坐标进行对比，以将行坐标与集合内的像素点的行坐标差距为1的像素点导入同一个集合内。
将行遍历仅得到一个像素值突变的像素点作为一个孔隙空间个体的终点，这样便得到了每个孔隙空间个体对应的孔隙边缘轮廓对应的像素点在行列遍历时的行列坐标。
将仅有一个像素值突变的像素点按照优先级顺序分别设置为一个孔隙空间个体的起点和终点，按照从起点到终点的顺序分为两条测绘曲线，两条测绘曲线顺次将不同优先级的像素值突变的像素点依次连接形成一个孔隙空间个体，单个孔隙空间个体即为一个孔隙单元，从而实现对孔隙单元的分割，便于查看识别，并且通过将分割的孔隙单元尺寸与求出的每个孔隙单元的长轴短轴尺寸对比，以实现对砂岩储层孔隙定量分析的人工验证工作。
步骤400、确定每个像素点在所述二维坐标系内的坐标值，基于坐标值计算所述二值图像中的每个所述孔隙单元的孔隙半径以及喉道半径，并输出每种所述砂岩储层样品的孔隙半径分布直方图、平均孔隙半径、喉道半径直方图与平均喉道半径。
所述孔隙半径为单个所述孔隙单元的长轴与短轴的平均值的一半，所述喉道半径为单个所述孔隙单元的短轴，所述二值图像中三次划定的每个所述孔隙单元的所述孔隙半径利用每个所述孔隙单元的像素点在所述二维坐标系的最大X坐标，最小X坐标，最大Y坐标，最小Y坐标进行求取，具体的实现步骤为：
1、选择每个所述孔隙单元内像素值突变的所述像素点的横坐标最大值对应的像素点、横坐标最小值对应的像素点，纵坐标最大值的像素点和坐标表最小值的像素点。
2、每个所述孔隙单元的长轴半径为Max[(Xmax-Xmin)，(Ymax-Ymin)]，而每个所述孔隙单元的短轴半径为Min[(Xmax-Xmin)，(Ymax-Ymin)]。
3、将通过公式[(Xmax-Xmin)+(Ymax-Ymin)]/2计算每个所述孔隙单元的孔隙半径。
4、将Min[(Xmax-Xmin)，(Ymax-Ymin)]筛选的数据作为每个所述孔隙单元的喉道半径。
统计每种砂岩储层中的每个孔隙单元对应的所述孔隙半径，并确定孔隙半径相同的所述孔隙单元的数量，创建关于孔隙半径和孔隙单元数量的孔隙半径直方图，并计算每种砂岩储层样品的孔隙半径平均值；统计每种砂岩储层中的每个孔隙单元对应的所述喉道半径，并确定喉道半径相同的所述孔隙单元的数量，创建关于喉道半径和孔隙单元数量的孔隙半径直方图，并计算每种砂岩储层样品的喉道半径平均值。
为了验证上述对砂岩储层粘土矿物晶间孔的定量分析方法，本发明以鄂尔多斯盆地上古生界石盒子组砂岩储层中高岭石晶间孔、绿泥石晶间孔、伊利石晶间孔、伊利石-蒙脱石混层晶间孔以及云母化基质孔的孔隙与喉道的分布和大小进行定量分析。
(a)对富含粘土矿物晶间孔的扫描电镜图像二值化；分别对图2中的高岭石晶间孔、绿泥石晶间孔、伊利石晶间孔、伊利石-蒙脱石混层晶间孔以及云母化基质孔的SEM图像进行图像二值化。将得到的扫描电镜图像进行拆分，并获得多尺度灰度赋值图像，基于不同灰度值，进行二值化赋值，即得到了如图3所示的富含各个粘土矿物晶间孔孔隙形态的二值图像。
(b)对获取的二值图像进行降噪；
基于之前的分析，扫描电镜受高压影响，SEM图像在拍摄过程中会产生一定数量的噪点，这对于纳米级与微米级孔隙的影响是不可忽视的。因此需要人为的对获取的二值图像进行降噪处理，得到如图4所示的已降噪二值图像。
(c)对进行了降噪的二值图像中孔隙单元的划分；对各类粘土矿物孔隙特征的二值图像进行孔隙单元的识别与划分，得到的孔隙单元图像，利用不同颜色的孔隙边缘划分表示不同孔隙单元。
(d)对步骤c中的孔隙单元进行粘土矿物的孔隙半径计算，根据每个孔隙单元的最大X值，最小X值，最大Y值和最小Y值计算孔隙半径，分别为X轴间距和Y轴间距的一半，如图5所示输出孔隙半径的分布直方图；
(e)对二值图像中的喉道进行识别，并计算喉道半径，选择X轴间距和Y轴间距的最小值作为喉道半径，如图6所示输出喉道半径的分布直方图。
综上的到了以下各类粘土矿物晶间孔的孔隙半径即喉道半径参数(表1)。
表1各类粘土矿物平均孔隙与喉道大小
因此本实施方式提供的定量分析方法结合直接观察法(扫描电镜图像)的可视化观察与数字图像技术，做到了定量分析，这是直接观察法所不能做到的，能够详细地对SEM图像中各类粘土矿物(如上述分析的绿泥石、高岭石、伊利石、伊利石-蒙脱石混层、云母化基质)的孔隙与喉道参数进行定量计算，这是压汞法或者数字岩芯法不能够达到的小于1微米的尺度。
以上实施例仅为本申请的示例性实施例，不用于限制本申请，本申请的保护范围由权利要求书限定。本领域技术人员可以在本申请的实质和保护范围内，对本申请做出各种修改或等同替换，这种修改或等同替换也应视为落在本申请的保护范围内。
Claims (6)
1.一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，其特征在于，包括以下步骤：
步骤100、准备不同孔隙类型的砂岩储层样品，获取每种砂岩储层样品在电镜扫描下的SEM图像，将所述SEM图像切割为多个图像矩阵，且利用自适应阈值方式对每个所述图像矩阵分别进行二值化处理；
步骤200、按照矩阵遍历方式识别每种所述砂岩储层样品对应所述图像矩阵中所有像素点对应的像素值，统一归纳二值图像的矩形阵列的像素值以同步完成修复噪点和消除噪点；
步骤300、重新按照行列遍历消除噪点和修复噪点之后的所述二值图像，将像素值突变的像素点作为孔隙边界，并且基于孔隙边缘划分所述二值图像中的孔隙单元；
将每个所述图像矩阵拼接形成消除噪点和修复噪点后的二值图像，且重新遍历所述二值图像，以所述二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历所述二值图像以确定每个所述像素点的像素值，选定行列遍历时的遍历间距均为一个像素，确定出每行中像素值突变的像素点，且将像素值突变的像素点保存在一个集合内，对两个所述集合进行多次方数字形态学运算，将两个所述集合之间的像素点的像素值重置为0以实现二次修复噪点；
对两个所述集合进行多次方数字形态学运算具体为：计算两个集合中同行或者同列像素点之间的邻域像素，并在所述邻域像素小于设定值时，将像素值为255的像素点设置为0以将两个集合对应的孔隙连接；
在步骤300中，识别二次修复噪点的所述二值图像的孔隙边缘和孔隙单元的实现方法为：
步骤301、以所述二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历所述二值图像以确定每个所述像素点的像素值，选定行列遍历时的遍历间距均为一个像素；
步骤302、将每个像素点的像素值与同一列的前一个像素点的像素值进行对比，确定像素值突变的所述像素点且重新定义该像素点的RGB值；
步骤303、将所有重新定义RGB值的像素点首尾连接构成的曲线设定为所述孔隙边缘，且将所有重新定义RGB值的像素点包围的区域设定为孔隙单元；
在步骤302中，确定遍历过程中的每一行的第一个像素值突变的所述像素点，按照遍历顺序划分每一行的像素值突变的所述像素点的优先级，并且将不同行的像素值突变的所述像素点的优先级设定为与遍历顺序一一对应，并且将同一行的像素值突变的所述像素点设置为同一个优先级；
将仅有一个像素值突变的所述像素点按照优先级顺序分别设置为一个所述孔隙单元的起点和终点，按照从起点到终点的顺序分为两条测绘曲线，两条所述测绘曲线顺次将不同优先级的所述像素值突变的所述像素点依次连接形成一个所述孔隙单元；
步骤400、确定每个像素点在所述二维坐标系内的坐标值，基于坐标值计算所述二值图像中的每个所述孔隙单元的孔隙半径以及喉道半径，并输出每种所述砂岩储层样品的孔隙半径分布直方图、平均孔隙半径、喉道半径直方图与平均喉道半径。
2.根据权利要求1所述的一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，其特征在于：在步骤100中，将所述SEM图像转化二值图像的具体实现步骤为：
将所述SEM图像转化为矩阵图像，并对所述矩阵图像平均切割为多个图像碎片；
抽样采集每个所述图像碎片中多个像素点的像素值，且确定对每个所述图像碎片二值化处理对应的分界像素值；
将每个所述图像碎片中低于分界像素值的像素点的像素值设置为0，且将高于分界像素值的像素点的像素值设置为255；
将所有的图像碎片按照行列式对齐，形成二值化处理后的二值图像。
3.根据权利要求2所述的一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，其特征在于，确定分界像素值的实现方式为：
将抽样采集的列间距和行间距形成抽样网格，并获取所述抽样网格的每个交叉点对应的像素值；
求取所有像素点的平均像素值，且将所述平均像素值作为分界像素值。
4.根据权利要求1所述的一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，其特征在于，在步骤300中，消除噪点和修复噪点可同时实现，具体的方法为：
利用归纳方式将每个自适应二值化处理的图像矩阵重新切割为多个归一化矩阵，统计每个所述归一化矩阵中像素值为0的像素点个数以及像素值为255的像素点个数；
将个数超过所述归一化矩阵中所有像素点数量一半的像素点作为多数元素，并将该归一化矩阵的所有像素点均重置为多数元素对应的像素值，以同步实现消除噪点和修复噪点。
5.根据权利要求1所述的一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，其特征在于：在步骤400中，所述孔隙半径为单个所述孔隙单元的长轴与短轴的平均值的一半，所述喉道半径为单个所述孔隙单元的短轴，所述二值图像中三次划定的每个所述孔隙单元的所述孔隙半径利用每个所述孔隙单元的像素点在所述二维坐标系的最大X坐标，最小X坐标，最大Y坐标，最小Y坐标进行求取，具体的实现步骤为：
选择每个所述孔隙单元内像素值突变的所述像素点的横坐标最大值对应的像素点、横坐标最小值对应的像素点，纵坐标最大值的像素点和坐标表最小值的像素点；
每个所述孔隙单元的长轴半径为Max[(Xmax-Xmin)，(Ymax-Ymin)]，而每个所述孔隙单元的短轴半径为Min[(Xmax-Xmin)，(Ymax-Ymin)]；
将通过公式[(Xmax-Xmin)+(Ymax-Ymin)]/2计算每个所述孔隙单元的孔隙半径；
将Min[(Xmax-Xmin)，(Ymax-Ymin)]筛选的数据作为每个所述孔隙单元的喉道半径。
6.根据权利要求5所述的一种基于SEM图像的砂岩储层矿物晶间孔的定量分析方法，其特征在于：在步骤400中，统计每种砂岩储层中的每个孔隙单元对应的所述孔隙半径，并确定孔隙半径相同的所述孔隙单元的数量，创建关于孔隙半径和孔隙单元数量的孔隙半径直方图，并计算每种砂岩储层样品的孔隙半径平均值；
统计每种砂岩储层中的每个孔隙单元对应的所述喉道半径，并确定喉道半径相同的所述孔隙单元的数量，创建关于喉道半径和孔隙单元数量的孔隙半径直方图，并计算每种砂岩储层样品的喉道半径平均值。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110760753.2A CN113570652B (zh) | 2021-07-06 | 2021-07-06 | 基于sem图像的砂岩储层矿物晶间孔的定量分析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110760753.2A CN113570652B (zh) | 2021-07-06 | 2021-07-06 | 基于sem图像的砂岩储层矿物晶间孔的定量分析方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113570652A CN113570652A (zh) | 2021-10-29 |
| CN113570652B true CN113570652B (zh) | 2023-03-10 |
Family
ID=78163750
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110760753.2A Active CN113570652B (zh) | 2021-07-06 | 2021-07-06 | 基于sem图像的砂岩储层矿物晶间孔的定量分析方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113570652B (zh) |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN114441405B (zh) * | 2021-12-22 | 2023-06-06 | 中国地质大学（北京） | 基于压实和胶结减孔趋势的次生增孔幅度定量评价方法 |
| CN117409408B (zh) * | 2023-12-15 | 2024-03-08 | 北京大学 | 层理缝参数获取方法、装置、设备及可读存储介质 |
Family Cites Families (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108763997A (zh) * | 2018-04-13 | 2018-11-06 | 广州中大微电子有限公司 | 一种基于交叉方向快速二值化的二维码降噪方法及系统 |
| CN108956416A (zh) * | 2018-06-06 | 2018-12-07 | 中国地质大学（北京） | 一种基于Matlab分析致密砂岩储层孔隙表征的方法 |
| CN112965485B (zh) * | 2021-02-03 | 2022-10-04 | 武汉科技大学 | 一种基于二次区域划分的机器人全覆盖路径规划方法 |
| CN112801846B (zh) * | 2021-02-09 | 2024-04-09 | 腾讯科技（深圳）有限公司 | 一种水印嵌入、提取方法、装置、计算机设备和存储介质 |
-
2021
- 2021-07-06 CN CN202110760753.2A patent/CN113570652B/zh active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN113570652A (zh) | 2021-10-29 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN113570652B (zh) | 基于sem图像的砂岩储层矿物晶间孔的定量分析方法 | |
| CN105445160B (zh) | 一种沥青混合料的空隙特征及其提取方法 | |
| CN113222992B (zh) | 基于多重分形谱的裂纹特征表征方法及系统 | |
| CN108734714B (zh) | 一种基于Matlab分析碳酸盐岩结构的方法 | |
| CN112150430B (zh) | 一种利用岩石细观结构数字图像的数值分析方法 | |
| CN112051200B (zh) | 一种致密砂岩储层孔隙结构的定量评价方法和装置 | |
| CN103471973A (zh) | 一种图像法测定织物孔隙特征的方法 | |
| CN104268830A (zh) | 基于数字图像确定不均匀岩土材料渗透系数的方法 | |
| CN108109157B (zh) | 一种基于数字式全景钻孔图像的岩体评估分析方法 | |
| CN114581764B (zh) | 基于深度学习算法的地下结构裂纹病害判别方法 | |
| CN111626358B (zh) | 一种基于bim图片识别的隧道围岩分级方法 | |
| CN109146888B (zh) | 一种基于差商分析的土裂隙分形维数计算方法 | |
| CN115290053B (zh) | 一种高层建筑施工建设异常快速检测方法 | |
| CN113570651B (zh) | 基于sem图像的碳酸盐岩储层孔隙半径分布定量方法 | |
| CN112465743A (zh) | 一种周期性结构质量检测方法 | |
| CN112070788A (zh) | 基于分块梯度分割的快速统计形变孪晶的图像处理方法 | |
| CN107657620A (zh) | 一种带纹理的金属凝固区域识别的方法及系统 | |
| CN115597494B (zh) | 一种基于点云的预制构件预留孔的精度检测方法、系统 | |
| CN110930423B (zh) | 一种物体边缘特征识别提取方法 | |
| CN113705001B (zh) | 输电线路故障检测方法、系统及嵌入式设备 | |
| CN112862816B (zh) | 一种hrtem图像中煤芳香烃晶格条纹的智能提取方法 | |
| CN111784059B (zh) | 一种煤层宏观裂隙优势发育方位预测方法 | |
| CN104504702A (zh) | 基于方格搜索法的水泥刻槽路面裂缝识别方法 | |
| CN108765381B (zh) | 一种水利工程裂纹提取方法及系统、信息处理终端 | |
| CN105241886A (zh) | 一种基于hvs与引导滤波器的土遗址裂缝检测方法 |
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 |