CN113702257A - 一种基于ct三维数据体的砾岩孔隙结构表征方法 - Google Patents

Abstract

本发明公开了一种基于CT三维数据体的砾岩孔隙结构表征方法，克服了常规孔隙结构研究的缺点，在砾岩样品不被破坏的条件下能够实现岩石微观结构的三维可视化，可以真实的还原砾岩内部的结构；能较好的将强非均质性砾岩中微裂缝及弱结构的孔隙结构提取出来，可得到较为精确的砾岩孔隙结构，更加真实的反应砾岩内部的孔隙结构。

Description

一种基于CT三维数据体的砾岩孔隙结构表征方法

技术领域

本发明涉及石油与天然气勘探与开发领域，尤其涉及一种基于 CT三维数据体的砾岩孔隙结构表征方法。

背景技术

砾岩中砾石颗粒大小混杂、裂缝发育、非均质性强，导致砾岩储层孔隙结构复杂、表征困难。孔隙结构是表征砾岩储层孔、渗物性的决定因素，科学的表征砾岩储层的复杂孔隙结构有利于储层识别、评价及分类，对油气资源高效的勘探开发具有重要的意义。目前，对砾岩的微观结构评价方法有很多种，如压汞法测量强非均质性砾岩孔隙时，由于砾岩颗粒大小混杂、孔隙度小及渗透性差，汞很难进入孔渗条件较低、物性较差的强非均质性砾岩中，导致压汞法测量的效果不佳且不能直观的描述孔隙空间的形态和空间分布特征。铸体薄片法在做孔隙结构分析时，可能会受到砾石与基质的影响，砾石颗粒过大、基质较少可能会影响测试分析的结果，且研究结果主要从二维角度分析孔隙结构信息，未对岩石的三维孔隙结构进行描述分析。常规阈值分割提取表征孔隙结构的方法适用于均质性的砂岩等，由于砾岩有较强的非均质性导致CT扫描的灰度图像灰度不均匀，使得阈值分割方法对弱结构孔隙提取效果较差。因此，急需要一种测量方便、可靠性高的基于CT三维数据体的砾岩孔隙结构表征方法。

发明内容

针对上述技术中存在的不足之处，本发明提供一种基于CT三维数据体的砾岩孔隙结构表征方法，克服了常规孔隙结构研究的缺点，在砾岩样品不被破坏的条件下能够实现岩石微观结构的三维可视化，可以真实的还原砾岩内部的结构；能较好的将强非均质性砾岩中微裂缝及弱结构的孔隙结构提取出来，可得到较为精确的砾岩孔隙结构，更加真实的反应砾岩内部的孔隙结构。

为实现上述目的，本发明提供一种基于CT三维数据体的砾岩孔隙结构表征方法，包括以下步骤：

S1：选取直径为8mm-100mm的砾岩样品，并对被选取的砾岩进行三维扫描，获取砾岩样品的二维灰度图像及三维数字岩心数据体；

S2：对二维灰度图像进行图像处理，并消除图像噪声和增强图像中孔隙细节特征；

S3：采用交互阈值分割方法对图像中的砂砾颗粒间孔隙结构进行第一次提取；

S4：采用顶底帽变换相结合的方法，对灰度图像的前景和背景灰度进行拉伸，凸显砾岩灰度图中弱结构孔隙的细节信息，对强非均质性砾岩灰度图像中孔隙的二次提取，实现微裂缝及弱结构孔隙的提取。

S5：将交互阈值分割与顶低帽变换两种算法分别提取的孔隙结构进行有机融合，获得砾岩二维灰度图像中的孔隙结构，并对融合后的二值化图进行精细化处理；

S6：对二值化切片数据中孔隙结构进行叠加，计算并构建砾岩孔隙网络模型；

S7：对砾岩孔隙结构特征参数提取及分析。克服了常规孔隙结构研究的缺点，在砾岩样品不被破坏的条件下能够实现岩石微观结构的三维可视化，可以真实的还原砾岩内部的结构；能较好的将强非均质性砾岩中微裂缝及弱结构的孔隙结构提取出来，可得到较为精确的砾岩孔隙结构，更加真实的反应砾岩内部的孔隙结构。

常规阈值分割的方法对灰度差异不明显及不同目标有重叠灰度值的灰度图像分割不明显，因此对强非均质性砾岩灰度图像中微裂缝及弱结构孔隙的提取效果差，阈值较低时部分裂缝及弱结构的孔隙无法被提取出来，阈值提高又会使得孔隙与基质的错误识别，导致孔隙体积的溢出影响分割的效果，对砾岩孔隙结构的实验结果产生一定影响。鉴于此，在提取砾岩样品内部孔隙时，采用交互阈值分割与顶低帽算法相结合的方式，可得到较为精确的岩石孔隙结构。

优选的，在步骤S1中，观察砾岩样品表观结构，选取具有代表性且满足分析需要的砾岩样品。基于样品尺寸大于其中最大砾石直径的原则，可选取直径为8mm-100mm的砾岩样品。

优选的，在步骤S1中，还包括利用CT扫描系统对选取的砾岩样品进行三维扫描，扫描分辨率的选取应满足最小空隙空间的识别需求，获取砾岩样品的二维灰度图像及三维数据体。

优选的，在步骤S2中，针对砾岩样品的灰度图像在其扫描环境下存在的各类噪声，对灰度图像进行中值滤波处理，消除图像中存在的噪声，使灰度图像变得更加清晰，同时对图像的细节特征及边缘信息有较好的保存。

优选的，在步骤S3中，利用交互阈值分割方法进行第一次孔隙结构的提取，从二维灰度图像的灰度级出发，将所有像素的灰度级与阈值进行比较进而选取一个合适的阈值，将大部分易提取的孔隙结构提取出来，在提取时尽量避免孔隙体积的溢出。

优选的，在步骤S4中，利用顶低帽变换相结合的方法，使得灰度图像前景和背景灰度进一步拉伸，突显出砾岩灰度图中弱结构孔隙的细节信息，起到对灰度图像增强的作用。

优选的，在步骤S4包括以下步骤：

S41:顶帽算法针对暗背景上的亮目标，具有一定的高通滤波器的特性，对不均匀光照的灰度图进行校证，可突出强非均质性砾岩灰度图像中弱结构孔隙的细节信息，增强图像的边缘信息，顶帽算法定义为T_hat(f)＝f-(f·b)；

S42：低帽变换算法用于亮背景上的暗目标，可以求出灰度图像中的谷值，突显出比原图轮廓周围更暗的区域，进而突出相互连接目标间的界限，分离出比邻近点更暗的一些弱结构孔隙，低帽算法定义为B_hat(f)＝(f·b)-f；

S43：基于顶低帽变换相结合的方法，突显出砾岩灰度图中弱结构孔隙的细节信息，可实现对强非均质性砾岩灰度图像中孔隙的二次提取，进而将微裂缝及弱结构孔隙提取出来。

优选的，在步骤S5中，对融合后的二值化图进行精细化处理，运用开、闭运算进行二次降噪，消除二值化图中的孤立的点以及填平小孔，得到高质量的岩石孔隙结构。

优选的，在步骤S6中包括，

S61:基于精细化处理的二值化图像，分析岩样的孔隙度随分析二值化图像切片尺寸的变化规律，获取砾岩孔隙结构的表征单元体尺寸大小；

S62:将图像处理后的二值化切片数据中孔隙结构叠加成一个三维数据场，建立三维孔隙空间网络模型；

S63:基于拓扑理论分水岭算法将三维孔隙空间网络模型中不同的孔隙进行提取；

S64:基于“最大球法”对分水岭算法后的三维孔隙结构进行孔喉的提取；

S65:将孔隙近似看作球体，喉道看作管道，建立一系列不同尺寸的球体与管道相互交叠形成的“球管”模型，实现砾岩内部孔隙结构的表征。

优选的，在步骤S7中包括以下步骤：

S71:通过相应位置的最大球多簇体积，确定“球管”模型中每个孔隙的体积。每个孔隙体积确定后，根据公式计算出所分析砾岩样品的绝对孔隙度与连通孔隙度；

S72:由分析得到孔隙体积，利用公式确定砾岩三维孔隙模型中每个孔隙与喉道的等效直径：

式中：Vpore为单个孔隙体积，单位pix³；Vvoxel为总体素的体积，单位pix³；

S73:基于建立的“球管”模型，通过对每个孔隙所连通的喉道进行遍历搜索并统计出每个孔隙连通喉道的个数，实现所分析砾岩样品配位数的计算分析；

S74:基于三维连通孔隙网络模型，针对所分析两个空间位置点，计算其直线距离，然后将其除以沿分析方向的连通喉道路径长度，进而计算得到喉道的迂曲度。

本发明的有益效果是：与现有技术相比，本发明提供一种基于CT三维数据体的砾岩孔隙结构表征方法，克服了常规孔隙结构研究的缺点，在砾岩样品不被破坏的条件下能够实现岩石微观结构的三维可视化，可以真实的还原砾岩内部的结构；能较好的将强非均质性砾岩中微裂缝及弱结构的孔隙结构提取出来，可得到较为精确的砾岩孔隙结构，更加真实的反应砾岩内部的孔隙结构。

附图说明

图1为本发明的砾岩孔隙结构表征流程图；

图2为本发明的CT扫描原始二维灰度图像；

图3为本发明的CT扫描原始三维数据体；

图4为本发明的CT扫描中值滤波后的灰度图；

图5为本发明的孔隙结构的提取体积溢出灰度图；

图6为本发明的孔隙结构的提取图；

图7为本发明的孔隙结构的提取弱结构及微裂缝图；

图8为本发明的孔隙结构的提取融合后的二值化图；

图9为本发明的孔隙结构的提取精细化处理后的二值化图；

图10为本发明的表征单元体二维切片分析示意图；

图11为本发明的表征单元体孔隙度与REV尺寸大小的关系图；

图12为本发明的岩石内部三维孔隙网络模型图；

图13为本发明的图12的孔隙分割后模型图；

图14为本发明的三维球管模型图；

图15为本发明的孔隙体积的孔隙半径分布图；

图16为本发明的孔隙体积的喉道半径分布图；

图17为本发明的孔隙体积的孔隙体积分布图。

具体实施方式

为了更清楚地表述本发明，下面结合附图对本发明作进一步地描述。

实施例1：请参阅图1至图17，本实施例的一种基于CT三维数据体的砾岩孔隙结构表征方法，包括以下步骤：S1：选取直径为8 mm-100mm的砾岩样品，并对被选取的砾岩进行三维扫描，获取砾岩样品的二维灰度图像及三维数字岩心数据体；S2：对二维灰度图像进行图像处理，并消除图像噪声和增强图像中孔隙细节特征；S3：采用交互阈值分割方法对图像中的砂砾颗粒间孔隙结构进行第一次提取； S4：采用顶底帽变换相结合的方法，对灰度图像的前景和背景灰度进行拉伸，凸显砾岩灰度图中弱结构孔隙的细节信息，对强非均质性砾岩灰度图像中孔隙的二次提取，实现微裂缝及弱结构孔隙的提取。S5：将交互阈值分割与顶低帽变换两种算法分别提取的孔隙结构进行有机融合，获得砾岩二维灰度图像中的孔隙结构，并对融合后的二值化图进行精细化处理；S6：对二值化切片数据中孔隙结构进行叠加，计算并构建砾岩孔隙网络模型；S7：对砾岩孔隙结构特征参数提取及分析。克服了常规孔隙结构研究的缺点，在砾岩样品不被破坏的条件下能够实现岩石微观结构的三维可视化，可以真实的还原砾岩内部的结构；能较好的将强非均质性砾岩中微裂缝及弱结构的孔隙结构提取出来，可得到较为精确的砾岩孔隙结构，更加真实的反应砾岩内部的孔隙结构。

常规阈值分割的方法对灰度差异不明显及不同目标有重叠灰度值的灰度图像分割不明显，因此对强非均质性砾岩灰度图像中微裂缝及弱结构孔隙的提取效果差，阈值较低时部分裂缝及弱结构的孔隙无法被提取出来，阈值提高又会使得孔隙与基质的错误识别，导致孔隙体积的溢出影响分割的效果，对砾岩孔隙结构的实验结果产生一定影响。鉴于此，在提取砾岩样品内部孔隙时，采用交互阈值分割与顶低帽算法相结合的方式，可得到较为精确的岩石孔隙结构。

实施例2：本实施例在步骤S1中，观察砾岩样品表观结构，选取具有代表性且满足分析需要的砾岩样品。基于样品尺寸大于其中最大砾石直径的原则，可选取直径为8mm-100mm的砾岩样品。本实施例在步骤S1中，还包括利用CT扫描系统对选取的砾岩样品进行三维扫描，扫描分辨率的选取应满足最小空隙空间的识别需求，获取砾岩样品的二维灰度图像及三维数据体。本实施例在步骤S2中，针对砾岩样品的灰度图像在其扫描环境下存在的各类噪声，对灰度图像进行中值滤波处理，消除图像中存在的噪声，使灰度图像变得更加清晰，同时对图像的细节特征及边缘信息有较好的保存。

实施例3：本实施例在步骤S3中，利用交互阈值分割方法进行第一次孔隙结构的提取，从二维灰度图像的灰度级出发，将所有像素的灰度级与阈值进行比较进而选取一个合适的阈值，将大部分易提取的孔隙结构提取出来，在提取时尽量避免孔隙体积的溢出。本实施例在步骤S4中，利用顶低帽变换相结合的方法，使得灰度图像前景和背景灰度进一步拉伸，突显出砾岩灰度图中弱结构孔隙的细节信息，起到对灰度图像增强的作用。

实施例4：本实施例在步骤S4包括以下步骤：S41:顶帽算法针对暗背景上的亮目标，具有一定的高通滤波器的特性，对不均匀光照的灰度图进行校证，可突出强非均质性砾岩灰度图像中弱结构孔隙的细节信息，增强图像的边缘信息，顶帽算法定义为 T_hat(f)＝f-(f·b)；S42：低帽变换算法用于亮背景上的暗目标，可以求出灰度图像中的谷值，突显出比原图轮廓周围更暗的区域，进而突出相互连接目标间的界限，分离出比邻近点更暗的一些弱结构孔隙，低帽算法定义为B_hat(f)＝(f·b)-f；S43：基于顶低帽变换相结合的方法，突显出砾岩灰度图中弱结构孔隙的细节信息，可实现对强非均质性砾岩灰度图像中孔隙的二次提取，进而将微裂缝及弱结构孔隙提取出来。本实施例的在步骤S5中，对融合后的二值化图进行精细化处理，运用开、闭运算进行二次降噪，消除二值化图中的孤立的点以及填平小孔，得到高质量的岩石孔隙结构。

实施例5：本实施例在步骤S6中包括，S61:基于精细化处理的二值化图像，分析岩样的孔隙度随分析二值化图像切片尺寸的变化规律，获取砾岩孔隙结构的表征单元体尺寸大小；S62:将图像处理后的二值化切片数据中孔隙结构叠加成一个三维数据场，建立三维孔隙空间网络模型；S63:基于拓扑理论分水岭算法将三维孔隙空间网络模型中不同的孔隙进行提取；S64:基于“最大球法”对分水岭算法后的三维孔隙结构进行孔喉的提取；S65:将孔隙近似看作球体，喉道看作管道，建立一系列不同尺寸的球体与管道相互交叠形成的“球管”模型，实现砾岩内部孔隙结构的表征。

实施例6：本实施例在步骤S7中包括以下步骤：S71:通过相应位置的最大球多簇体积，确定“球管”模型中每个孔隙的体积。每个孔隙体积确定后，根据公式计算出所分析砾岩样品的绝对孔隙度与连通孔隙度；

S72:由分析得到孔隙体积，利用公式确定砾岩三维孔隙模型中每个孔隙与喉道的等效直径：

式中：Vpore为单个孔隙体积，单位pix³；Vvoxel为总体素的体积，单位 pix³；S73:基于建立的“球管”模型，通过对每个孔隙所连通的喉道进行遍历搜索并统计出每个孔隙连通喉道的个数，实现所分析砾岩样品配位数的计算分析；S74:基于三维连通孔隙网络模型，针对所分析两个空间位置点，计算其直线距离，然后将其除以沿分析方向的连通喉道路径长度，进而计算得到喉道的迂曲度。

实施例7：第一步：CT扫描及图像预处理

1)观察砾岩样品表观结构，选取具有代表性且满足分析需要的砾岩样品。基于样品尺寸大于其中最大砾石直径的原则，可选取直径为 8mm～100mm的砾岩样品；2)利用CT扫描系统对选取的砾岩样品进行三维扫描，获取砾岩样品的二维灰度图像(图2)及三维数据体(图 3)；3)针对砾岩样品的灰度图像在其扫描环境下存在的各类噪声，对灰度图像进行中值滤波处理，消除图像中存在的噪声，使灰度图像变得更加清晰，同时对图像的细节特征及边缘信息有较好的保存(图 4)。

第二步：微裂缝及弱结构的提取

常规阈值分割的方法对灰度差异不明显及不同目标有重叠灰度值的灰度图像分割不明显，对强非均质性砾岩灰度图像中微裂缝及弱结构孔隙的提取效果差，阈值较低时部分裂缝及弱结构的孔隙无法被提取出来，阈值提高又会使得孔隙与基质的错误识别(图5)，导致孔隙体积的溢出影响分割的效果，对砾岩孔隙结构的实验结果产生一定影响。鉴于此，在提取砾岩样品内部孔隙时，采用交互阈值分割与顶低帽算法相结合的方式，可得到较为精确的岩石孔隙结构。

1)利用交互阈值分割方法进行第一次孔隙结构的提取，从二维灰度图像的灰度级出发，将所有像素的灰度级与阈值进行比较进而选取一个合适的阈值，将大部分易提取的孔隙结构提取出来，如图6所示，在提取时尽量避免孔隙体积的溢出；

2)利用顶低帽变换相结合的方法，使得灰度图像前景和背景灰度进一步拉伸，突显出砾岩灰度图中弱结构孔隙的细节信息，起到对灰度图像增强的作用。原始灰度图像上加上顶低帽变换的结果，可有效地提高灰度图像的对比度，有利于对砾岩灰度图中微裂缝及弱结构孔隙的提取。

a.顶帽算法针对暗背景上的亮目标，具有一定的高通滤波器的特性，对不均匀光照的灰度图进行校证，可突出强非均质性砾岩灰度图像中弱结构孔隙的细节信息，增强图像的边缘信息，顶帽算法定义为

b.低帽变换算法用于亮背景上的暗目标，可以求出灰度图像中的谷值，突显出比原图轮廓周围更暗的区域，进而突出相互连接目标间的界限，分离出比邻近点更暗的一些弱结构孔隙，低帽算法定义为

B_{h at}(f)＝(f·b)-f (2)

c.基于顶低帽变换相结合的方法，突显出砾岩灰度图中弱结构孔隙的细节信息，可实现对强非均质性砾岩灰度图像中孔隙的二次提取(图7)，进而将微裂缝及弱结构孔隙提取出来。

3)将常规阈值分割与顶低帽变换两种算法分别提取的孔隙结果进行有机融合，获得砾岩二维灰度图像中的孔隙结构。如图8所示，为提取的孔隙结构的二值化图。

4)对融合后的二值化图进行精细化处理，运用开、闭运算进行二次降噪，消除二值化图中的孤立的点以及填平小孔，得到较高质量的岩石孔隙结构(图9)。

第三步：孔隙网络模型的建立

1)考虑岩心尺寸会对计算机运行和储存有较高的要求，既要满足岩样的宏观性，又要考虑现有的计算水平，折衷的方法就是选取合理的表征单元体。基于精细化处理的二值化图像，分析岩样的孔隙度随二值化尺寸的变化规律，进而分析出表征单元体的尺寸大小为400 pix×400pix×400pix(图10和图11)。

2)将图像处理后的二值化切片数据中孔隙结构叠加成一个三维数据场，建立三维孔隙空间网络模型。基于拓扑理论分水岭算法将三维孔隙空间网络模型中不同的孔隙独立的分隔开，并对不同的孔隙进行颜色渲染，使得每个孔隙度有自己独立的颜色(图12)。

3)基于“最大球法”对分水岭算法后的三维孔隙结构进行孔喉的提取，将孔隙近似看作球体，喉道看作管道，使得砾岩内部孔隙结构由一系列不同尺寸的球体与管道相互交叠形成，建立“球管”模型 (图13)。

第四步：孔隙结构分析

最大球算法通过在三维数字岩心的孔隙中选点作为球心，增大球体半径向四周扩展，周围可重叠相邻的较小球体，形成一个充满整个空间的最大球多簇；若1个球体属于2个最大球多簇，则这个公共的最大球体被认为是喉道。

1)基于建立的“球管”模型中每个孔隙的体积都与相应位置假想的最大球多簇体积相等，每个孔隙体积确定后，根据公式(3)可计算出数字岩心绝对孔隙度与连通孔隙度，如表1所示；

2)三维孔隙模型中每个孔隙与喉道的等效直径可通过公式(4) 确定，如表1与图14至图17所示；

式中：V_pore为单个孔隙体积，单位pix³；V_voxel为总体素的体积，单位pix³。

3)基于“球管”模型，在三维空间上转变不同的坐标方向统计出每个孔隙连通喉道的个数，即孔隙的配位数；

4)基于三维连通孔隙网络模型计算每个平面的质心及质心的路径长度，然后将其除以沿z轴的平面数可以得到喉道的迂曲度。

式中：

为质心的路径长度，z_n为z轴的平面个数。

表1孔喉结构参数分布

以上公开的仅为本发明的几个具体实施例，但是本发明并非局限于此，任何本领域的技术人员能思之的变化都应落入本发明的保护范围。

Claims (10)

1.一种基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，包括以下步骤：

S1：选取直径为8mm-100mm的砾岩样品，并对被选取的砾岩进行三维扫描，获取砾岩样品的二维灰度图像及三维数字岩心数据体；

S2：对二维灰度图像进行图像处理，并消除图像噪声和增强图像中孔隙细节特征；

S3：采用交互阈值分割方法对图像中的砂砾颗粒间孔隙结构进行第一次提取；

S4：采用顶底帽变换相结合的方法，对灰度图像的前景和背景灰度进行拉伸，凸显砾岩灰度图中弱结构孔隙的细节信息，对强非均质性砾岩灰度图像中孔隙的二次提取，实现微裂缝及弱结构孔隙的提取。

S5：将交互阈值分割与顶低帽变换两种算法分别提取的孔隙结构进行有机融合，获得砾岩二维灰度图像中的孔隙结构，并对融合后的二值化图进行精细化处理；

S6：对二值化切片数据中孔隙结构进行叠加，计算并构建砾岩孔隙网络模型；

S7：对砾岩孔隙结构特征参数提取及分析。

2.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S1中，观察砾岩样品表观结构，选取具有代表性且满足分析需要的砾岩样品。基于样品尺寸大于其中最大砾石直径的原则，可选取直径为8mm-100mm的砾岩样品。

3.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S1中，还包括利用CT扫描系统对选取的砾岩样品进行三维扫描，扫描分辨率的选取应满足最小空隙空间的识别需求，获取砾岩样品的二维灰度图像及三维数据体。

4.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S2中，针对砾岩样品的灰度图像在其扫描环境下存在的各类噪声，对灰度图像进行中值滤波处理，消除图像中存在的噪声，使灰度图像变得更加清晰，同时对图像的细节特征及边缘信息有较好的保存。

5.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S3中，利用交互阈值分割方法进行第一次孔隙结构的提取，从二维灰度图像的灰度级出发，将所有像素的灰度级与阈值进行比较进而选取一个合适的阈值，将大部分易提取的孔隙结构提取出来，在提取时尽量避免孔隙体积的溢出。

6.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S4中，利用顶低帽变换相结合的方法，使得灰度图像前景和背景灰度进一步拉伸，突显出砾岩灰度图中弱结构孔隙的细节信息，起到对灰度图像增强的作用。

7.根据权利要求6所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S4包括以下步骤：

S41:顶帽算法针对暗背景上的亮目标，具有一定的高通滤波器的特性，对不均匀光照的灰度图进行校证，可突出强非均质性砾岩灰度图像中弱结构孔隙的细节信息，增强图像的边缘信息，顶帽算法定义为T_hat(f)＝f-(f·b)；

S42：低帽变换算法用于亮背景上的暗目标，可以求出灰度图像中的谷值，突显出比原图轮廓周围更暗的区域，进而突出相互连接目标间的界限，分离出比邻近点更暗的一些弱结构孔隙，低帽算法定义为B_hat(f)＝(f·b)-f；

S43：基于顶低帽变换相结合的方法，突显出砾岩灰度图中弱结构孔隙的细节信息，可实现对强非均质性砾岩灰度图像中孔隙的二次提取，进而将微裂缝及弱结构孔隙提取出来。

8.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S5中，对融合后的二值化图进行精细化处理，运用开、闭运算进行二次降噪，消除二值化图中的孤立的点以及填平小孔，得到高质量的岩石孔隙结构。

9.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S6中包括，

S61:基于精细化处理的二值化图像，分析岩样的孔隙度随分析二值化图像切片尺寸的变化规律，获取砾岩孔隙结构的表征单元体尺寸大小；

S62:将图像处理后的二值化切片数据中孔隙结构叠加成一个三维数据场，建立三维孔隙空间网络模型；

S63:基于拓扑理论分水岭算法将三维孔隙空间网络模型中不同的孔隙进行提取；

S64:基于“最大球法”对分水岭算法后的三维孔隙结构进行孔喉的提取；

S65:将孔隙近似看作球体，喉道看作管道，建立一系列不同尺寸的球体与管道相互交叠形成的“球管”模型，实现砾岩内部孔隙结构的表征。

10.根据权利要求1所述的基于CT三维数据体的砾岩孔隙结构表征方法，其特征在于，在步骤S7中包括以下步骤：

S71:通过相应位置的最大球多簇体积，确定“球管”模型中每个孔隙的体积。每个孔隙体积确定后，根据公式计算出所分析砾岩样品的绝对孔隙度与连通孔隙度；

S72:由分析得到孔隙体积，利用公式确定砾岩三维孔隙模型中每个孔隙与喉道的等效直径：

式中：Vpore为单个孔隙体积，单位pix³；Vvoxel为总体素的体积，单位pix³。

S73:基于建立的“球管”模型，通过对每个孔隙所连通的喉道进行遍历搜索并统计出每个孔隙连通喉道的个数，实现所分析砾岩样品配位数的计算分析；

S74:基于三维连通孔隙网络模型，针对所分析两个空间位置点，计算其直线距离，然后将其除以沿分析方向的连通喉道路径长度，进而计算得到喉道的迂曲度。

Images