具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

如图1所示，本发明提供了一种基于SEM图像的煤储层孔隙结构参数定量分析方法，需要说明的是，本实施方式定量分析的煤储层孔隙结构参数包括孔隙半径分布图、平均孔隙半径、孔隙配位数分布图以及平均孔隙配位数，本方法能够精确计算孔隙半径较小的孔隙单元的结构参数，孔隙半径的计算精确度高，且能够通过圆周曲线交叉算法计算普通砂岩储层样品的孔隙配位数，增设孔隙半径小对孔隙配位数的计算影响，使得计算方式简单且准确，相比于CT法，处理非常便捷同时消耗的人力物力较低，同时实现对孔隙配位数情况以图像和直方图的直观方式进行定量输出。

具体包括以下步骤，如图2所示：

步骤100、准备不同孔隙类型的煤储层样品，获取每种煤储层样品在电镜扫描下的SEM图像，且对每种煤储层样品的SEM图像进行二值化处理，将SEM图像转化为二值图像。

在步骤100中，煤储层样品的SEM图像中，SEM图像中非孔隙范围内的像素点的像素值大于孔隙范围内的像素点的像素值，且采样抽取多个非孔隙范围内的像素点的像素值以及孔隙范围内的像素点的像素值，且确定采样的非孔隙范围内的像素点的最小像素值以及孔隙范围内的像素点的最大像素值；

将无孔隙范围内超过最小像素值的像素点设置为255且将孔隙范围内低于最大像素值的像素点设置为0以实现对SEM图像的二值化处理。

SEM图像仅能获取关于煤储层样品的灰度图，图像的二值化是将图像上的像素点的灰度值设置为0或255，从256个灰度级强制减少到只有0和255两个个灰度级，从而将感兴趣目标区域和背景分离，在本实施方式中，非孔隙范围内像素点的像素值为0，作为白色背景，对应孔隙范围内像素点的像素值为255，从而将不同孔隙类型的煤储层样品在扫描电镜下获取的SEM图像转化为富含孔隙信息的二值图像。

步骤200、建立关于二值图像的二维坐标系，按照矩阵遍历方式识别每种煤储层样品对应二值图像中所有像素点对应的像素值，将像素值突变的像素点作为孔隙边缘，并且基于孔隙边缘划分二值图像中的孔隙空间个体。

在步骤200中，识别每种煤储层样品对应二值图像中的孔隙边缘轮廓是为了更加精准的获取孔隙参数，在本实施方式中，由于二值图像的像素点仅包含两个像素值，分别为0和255，按照行列式遍历二值图像每个像素点时，可得到每个像素点的行列坐标以及对应每个像素点的像素值，将每次像素点的像素值与上一个像素点的像素值进行对比，当遍历到像素点的像素值由255突变到0时，则认为该像素点对应孔隙单元的边缘位置，当遍历到像素点的像素值由0突变到255时，则将像素值为0的像素点作为孔隙单元的边缘位置。

识别每种煤储层样品对应二值图像中的孔隙边缘和孔隙空间个体的实现方法为：

步骤201、以二值图像的垂直交叉边缘为坐标轴建立二维坐标系，且按照先行后列的方式遍历二值图像以确定每个像素点的像素值，选定行列遍历时的遍历间距均为一个像素。

步骤202、将每个像素点的像素值与同一列的前一个像素点的像素值进行对比，确定像素值突变的像素点且重新定义该像素点的RGB值。

步骤203、将所有重新定义RGB值的像素点首尾连接构成的曲线设定为孔隙边缘，且将所有重新定义RGB值的像素点包围的区域设定为孔隙空间个体。

即使每种煤储层的孔隙单元的形状不同，但是行列方式遍历每种煤储层对应的二值图像时，每个孔隙空间个体均存在仅有单个像素值突变的像素点，在本实施方式中，当遍历一行中仅存在单个像素值突变的像素点，则将该像素点作为孔隙空间个体的起点，另一个存在单个像素值突变的像素点作为一个孔隙空间个体的终点，因此将行列遍历得到的单个像素值突变的像素点作为一个孔隙空间个体的起点，并将像素点的行列坐标保存在一个集合内。

将下一行遍历得到的单个像素值突变的像素点与集合内的像素点行列坐标进行对比，以将行坐标与集合内的像素点的行坐标差距为1的像素点导入同一个集合内。

将行遍历仅得到一个像素值突变的像素点作为一个孔隙空间个体的终点，这样便得到了每个孔隙空间个体对应的孔隙边缘轮廓对应的像素点在行列遍历时的行列坐标。

在步骤202中，确定遍历过程中的每一行的第一个像素值突变的像素点，按照遍历顺序划分每一行的像素值突变的像素点的优先级，并且将不同行的像素值突变的像素点的优先级设定为与遍历顺序一一对应，并且将同一行的像素值突变的像素点设置为同一个优先级。

将仅有一个像素值突变的像素点按照优先级顺序分别设置为一个孔隙空间个体的起点和终点，按照从起点到终点的顺序分为两条测绘曲线，两条测绘曲线顺次将不同优先级的像素值突变的像素点依次连接形成一个孔隙空间个体，单个孔隙空间个体即为一个孔隙单元，从而实现对孔隙单元的分割，便于查看识别，并且通过将分割的孔隙单元尺寸与求出的每个孔隙单元的长轴短轴尺寸对比，也可以将每个孔隙单元的孔隙配位数与程序编程求出的孔隙配位数进行人工对比，以实现对煤储层孔隙定量分析的人工验证工作。

步骤300、确定每个像素点在二维坐标系内的坐标值，基于坐标值计算二值图像中的每个孔隙空间个体的孔隙半径，并获取二值图像中每个孔隙空间个体的孔隙配位数。

孔隙半径为单个孔隙空间个体的长轴与短轴的平均值的一半，根据上述，按照行列遍历二值图像时，将像素值由255突变到0的像素点作为孔隙边缘的像素点，并将此像素点的行列坐标保存在单个孔隙空间个体对应的集合内，因此当计算每个孔隙空间个体的孔隙半径时，先查找出每个集合内的行坐标最大值、行坐标最小值、列坐标最大值和列坐标最小值，然后将行坐标最大值、行坐标最小值、列坐标最大值和列坐标最小值转化为在二维坐标系的最大X坐标，最小X坐标，最大Y坐标，最小Y坐标，求取单个孔隙空间个体的长轴与短轴，具体的实现步骤为：

选择每个孔隙空间个体对应集合内像素值突变的像素点的横坐标最大值对应的像素点、横坐标最小值对应的像素点，纵坐标最大值的像素点和坐标表最小值的像素点，每个孔隙空间个体的长轴半径为Max[(Xmax-Xmin)，(Ymax-Ymin)]，而每个孔隙空间个体的短轴半径为Min[(Xmax-Xmin)，(Ymax-Ymin)]，将通过公式[(Xmax-Xmin)+(Ymax-Ymin)]/2计算每个孔隙空间个体的孔隙半径。

而孔隙配位数是指一个中心孔隙与其周边相邻的孔隙数的比值，如一个中心孔隙周边相邻有5个孔隙，则该孔隙的配位数为5。

在步骤300中，获取二值图像中每个孔隙空间个体的孔隙配位数的实现方法为：①将所述孔隙空间个体的孔隙半径中心位置作为所述孔隙空间个体的质点，以所述孔隙空间个体的质点为圆心并获取不同半径的圆周曲线；

②依次统计孔隙空间个体与不同半径的圆周曲线之间的交叉数量，确定每次与圆周曲线交叉时的孔隙空间个体，通过累加不同半径的圆周曲线对应的不同孔隙空间个体总数作为孔隙配位数；

③顺次对所有孔隙空间个体依次计算对应的孔隙配位数，且统计相同孔隙配位数对应的孔隙空间个体数量。

由于孔隙空间个体的像素点与不同半径的圆周曲线的交叉数量不同，且设定圆周曲线的半径选定范围，仅统计半径选定范围内的圆周曲线与像素点的交叉个数，且识别每次与所述圆周曲线交叉时时所述孔隙空间个体，将不同半径的圆周曲线对应的不同孔隙空间个体的数量叠加作为孔隙配位数。

其中，确定与不同圆周曲线交叉的孔隙空间个体是否为同一个的具体操作方式为：将与圆周曲线的交叉的像素点与统计的像素值突变的像素点所在集合进行对比，确定属于同一个孔隙空间个体的像素点，并将不同的孔隙空间个体作为孔隙配位数。

通过上述方法通过改变圆周曲线的半径，可以获得对应一个孔隙单元的多个孔隙配位数，且将圆周曲线与孔隙空间个体交叉数量的叠加结果作为孔隙配位数，将孔径半径较小的孔隙空间个体也增加到孔隙配送数的计算过程中，提高了计算精度，另外，本实施方式还可以设定圆周曲线的选定范围，避免将圆周曲线的选定方式过大而造成计算错误的问题，从而得到的孔隙配位数更加精准，且可通过程序调控，实现更加快捷简单。

步骤400、输出每种煤储层样品的孔隙半径分布直方图以及孔隙配位数直方图，并计算平均孔隙半径和平均孔隙配位数。

在步骤400中，统计每种煤储层样品的孔隙配位数进行分类，确定对应孔隙配位数相同的孔隙空间个体的数量，并计算每种煤储层样品的平均孔喉配位数；

统计每种煤储层样品的不同孔隙半径对应的孔隙空间个体的数量，确定孔隙半径相同的孔隙空间个体的数量，并计算每种煤储层样品的孔隙半径平均值。

为了验证上述对于煤储层孔隙结构参数定量分析方法，本发明以沁水盆地山西组3号煤层为研究对象，如图3所所示，基于SEM图像对孔隙分布、孔隙半径大小、和孔隙配位数进行定量计算。

如图4所示，通过扫描电镜显微镜下观察，对煤样品进行观察及捕捉富含孔隙特征的SEM图像。

根据定量分析流程图，首先对SEM图像进行二值化处理获取相应的二值图像，如图5所示。

基于孔隙空间个体划分结果之上，计算图4中每个孔隙空间个体的孔隙半径，获取的孔隙半径分布直方图及平均孔隙大小(18.0697微米)，如图6所示。

在获取孔隙半径相关参数及填充后的二值图像后，基于对图6中孔隙的进行分析，获取每个单元孔隙的配位分配中每个彩色单元空间进行计算，联通的图像为图7所示，而每个孔隙空间个体及附近的配位情况分布为图8所示。

最终如图9所示，输出孔隙配位数直方图以及平均孔喉配位数，煤岩样品的平均孔隙配位数为1.1009个。

因此本方法能够精确计算孔隙半径较小的孔隙单元的结构参数，孔隙半径的计算精确度高，且能够通过圆周曲线交叉算法计算普通砂岩储层样品的孔隙配位数，排除孔隙半径小对孔隙配位数的计算影响，使得计算方式简单且准确，相比于CT法，处理非常便捷同时消耗的人力物力较低，同时实现对孔隙配位数情况以图像和直方图的直观方式进行定量输出。

以上实施例仅为本申请的示例性实施例，不用于限制本申请，本申请的保护范围由权利要求书限定。本领域技术人员可以在本申请的实质和保护范围内，对本申请做出各种修改或等同替换，这种修改或等同替换也应视为落在本申请的保护范围内。