CN109658378B - 基于土壤ct图像的孔隙辨识方法及系统 - Google Patents

Abstract

本发明实施例提供一种基于土壤CT图像的孔隙辨识方法及系统，该方法包括：建立土壤CT图像的灰度‑梯度直方图，并通过筛除、投影、卷积的方式获得灰度直方图，基于灰度直方图的极点数量，确定土壤图像的聚类数目；根据聚类数目对灰度‑梯度直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；基于邻域空间信息构造目标函数，并通过迭代优化完成对土壤CT图像的模糊辨识；对模糊辨识结果进行去模糊化，获得土壤孔隙的辨识结果。本发明实施例通过确定聚类数目、构造初始隶属度矩阵、构造加入空间信息的目标函数和去模糊化处理，对土壤CT图像具有普适性，解决了初始值制约辨识精确度的问题，保证了土壤孔隙辨识的执行效率。

Description

基于土壤CT图像的孔隙辨识方法及系统

技术领域

本发明实施例涉及图像处理领域，更具体地，涉及一种基于土壤CT图像的孔隙辨识方法及系统。

背景技术

目前，基于土壤CT图像研究者大多采用传统FCM法、快速FCM法来辨识孔隙结构。但由于在传统FCM算法中，所有像素点均参与计算，每次迭代需消耗大量时间，并且其稳定效果和速度极易受初始条件的影响，因此，降低了辨识方法的运算效率。而快速FCM法采用灰度级参与迭代运算，虽极大地提高运算速度，但忽略了细节信息，仍无法准确描述较为复杂的孔隙模糊边界图像(土壤CT图像)的特征。因此，研究一种能在保证孔隙辨识精度的前提下具有较高的执行效率的辨识方法成为亟待解决的问题。

发明内容

为了解决上述问题，本发明实施例提供一种克服上述问题或者至少部分地解决上述问题的基于土壤CT图像的孔隙辨识方法及系统。

根据本发明实施例的第一方面，提供一种基于土壤CT图像的孔隙辨识方法，该方法包括：建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目；根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识；对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果。

根据本发明实施例第二方面，提供了一种基于土壤CT图像的孔隙辨识系统，该系统包括：确定模块，用于建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于灰度直方图的极点数量，确定土壤图像的聚类数目；划分模块，用于根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；构造模块，用于基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识；辨识模块，用于对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果。

根据本发明实施例的第三方面，提供了一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，处理器执行程序时实现如第一方面的各种可能的实现方式中任一种可能的实现方式所提供的基于土壤CT图像的孔隙辨识方法。

根据本发明实施例的第四方面，提供了一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现如第一方面的各种可能的实现方式中任一种可能的实现方式所提供的基于土壤CT图像的孔隙辨识方法。

本发明实施例提供的基于土壤CT图像的孔隙辨识方法及系统，通过建立灰度-梯度二维直方图、确定聚类数目、构造初始隶属度矩阵、构造加入空间信息的目标函数和去模糊化处理，从而相比于现有技术对土壤CT图像具有普适性，解决了初始值制约辨识精确度的问题，在保证了土壤孔隙辨识精度的前提下具有较高的执行效率。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍。显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些图获得其他的附图。

图1为本发明实施例提供的基于土壤CT图像的孔隙辨识方法的流程示意图；

图2为本发明实施例提供的灰度-梯度二维直方图的示意图；

图3为本发明实施例提供的去模糊化处理的效果示意图；

图4为本发明实施例提供的基于土壤CT图像的孔隙辨识系统的结构示意图；

图5为本发明实施例提供的电子设备的实体结构示意图。

具体实施方式

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

现有技术中，关于土壤CT图像的孔隙辨识有传统FCM法和快速FCM法；其中：

传统FCM法遵循像素点以不同概率值属于不同集群的原则，通过迭代优化完成像素点的模糊划分。当目标函数通过迭代优化趋于极小值时，则认为所有像素点都趋于某个聚类中心并远离其他聚类中心，聚类结果达到理想状态；

快速FCM法用灰度直方图特征空间代替像素空间参与迭代计算，计算图像中灰度相同的像素点数量和灰度等级，在优化过程中，由灰度等级与相应数量的乘积代替像素灰度，在此基础上，通过不同等级灰度对应的隶属度向量，可以反求出原始图像中所有像素点的隶属度。

但是，由上面的传统FCM法对孔隙结构进行辨识，由于其大量的数据集和迭代次数使得执行效率偏低，不适用于大批量的土壤CT图像；快速FCM法虽执行效率提高，但像素的细节信息丢失，孔隙辨识精度降低。因此，研究一种能精确辨析孔隙，适用于孔隙结构微小、灰度值均一化的土壤CT图像的辨识方法成为亟待解决的问题。

基于此，本发明实施例提供一种基于土壤CT图像的孔隙辨识方法。参见图1，该方法包括：

101、建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目。

其中，土壤孔隙是指土壤颗粒之间、团聚体之间或团聚体内部的孔隙。孔隙结构包括孔隙数目、大小等几何形态。CT(Computed Tomography)，即电子计算机断层扫描，它是利用精确准直的X线束、γ射线、超声波等，与灵敏度极高的探测器一同围绕人体的某一部位作一个接一个的断面扫描，具有扫描时间快，图像清晰等特点。CT扫描技术是孔隙量化研究的有效手段，在步骤101之前，可以利用CT扫描技术获得土壤CT图像。然后建立土壤CT图像对应的灰度-梯度二维直方图。参见图2，灰度-梯度二维直方图的横坐标为灰度，纵坐标为梯度。并基于灰度-梯度二维直方图确定聚类的数目。

将灰度-梯度二维直方图向灰度轴作投影，获得一维灰度直方图f(x)。之后，采用高斯模板对一维灰度直方图f(x)进行卷积，并基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目。

102、根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3。

具体地，根据聚类数据对灰度-梯度二维直方图进行区域划分可采用如下划分规则：

分区数＝梯度块数(m)*灰度块数(n)

分区数＝4*4聚类数目≥16

分区数＝3*4聚类数目≥11

分区数＝3*3聚类数目>7

分区数＝2*3聚类数目≤7

如图2所示，例如当聚类数目为9时，基于上述划分规则，将灰度-梯度二维直方图按照3*3的进行区域划分，获得9个区域，每个区域代表一个类。其中，T和G分布表示图像的梯度和灰度。T和G的取值规则为：

在划分完成后，可以根据划分得到的类别构造初始隶属度矩阵。

103、基于初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据目标函数对土壤CT图像进行模糊辨识。

其中，应当说明的是，传统的目标函数是不包括像素点邻域的空间信息。而本发明实施例中，由于CT机器的部分容积效应，土壤CT图像中的某一个像素点受其邻域的影响，以灰度值的均值状态存在。因此，对某一像素点所属类别进行判断时，应该考虑该像素点的邻域像素点对其产生的影响。即，在传统的目标函数基础上，将其邻域内的空间信息作为约束条件，创建了新的目标函数，能够更加准确的判断像素点的类别。在构造目标函数后，基于目标函数可以重新进行隶属度矩阵和聚类中心的计算，获得目标隶属度矩阵和目标聚类中心。从而能够根据每个像素点的隶属度对像素点进行分类，即可完成模糊辨识。

104、对模糊辨识结果进行去模糊化处理，获得土壤CT图像中的孔隙结构辨识结果。

具体地，参见图3，在对像素点进行分类后，去模糊化处理的目的是提取出土壤CT图像中的孔隙结构。由于土壤CT图像主要分类白色、浅灰色、深灰色和黑色四类，而孔隙结构的像素点相对于其他部分的颜色较深，因此可以将灰度值最小的一类像素点提取出来作为孔隙结构辨识结果，从而完成基于土壤CT图像的孔隙辨识。

本发明实施例提供的方法，通过建立灰度-梯度二维直方图、确定聚类数目、构造初始隶属度矩阵、构造加入空间信息的目标函数和去模糊化处理，从而相比于现有技术对土壤CT图像具有普适性，解决了初始值制约辨识精确度的问题，在保证了土壤孔隙辨识精度的前提下具有较高的执行效率。

基于上述实施例的内容，作为一种可选实施例，提供一种建立土壤CT图像的灰度-梯度二维直方图的方法，包括但不限于：基于八邻域的拉普拉斯算子计算获得土壤CT图像的灰度-梯度二维直方图。

具体地，八邻域的拉普拉斯算子的运算模板为：

F(i,j)＝8f(i,j)-f(i+1,j)-f(i-1,j)-f(i,j+1)-f(i,j-1)-f(i-1,j-1)-f(i-1,j+1)-f(i+1,j)-f(i+1,j+1)

(0＜i＜w,0＜j＜l)

式中，w为图像宽度，l为图像长度，f为梯度。

相比于四邻域的拉普拉斯算子，八邻域的拉普拉斯算子能够提取出较为细节的边缘信息，能够很好地解决对于受部分容积效应影响造成的边界均值模糊的现象，从而可以准确提取出孔隙的轮廓梯度信息。通过融合土壤结构的灰度信息和梯度信息，可准确获取土壤图像中灰度和梯度的分布情况，从而为后续基于灰度-梯度特征的分类提供初始条件。

基于上述实施例的内容，作为一种可选实施例，提供一种通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目的方法，包括但不限于：

步骤1、剔除灰度-梯度二维直方图中梯度大于10的像素点，并将更新后的灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图；

在确定聚类的数目前，可先去除灰度-梯度二维直方图中梯度大于梯度阈值的像素点。

具体地，首先可计算土壤CT图像中每个像素点的梯度值，具体可采用如下公式：

设自身的灰度为i、梯度为j、像素点的个数为N_ij，则具有这样属性的点出现的频率为

阵列(H_ij)就为该图像的灰度-梯度二维直方图。

在获得每个像素点的梯度值后，可设定一个梯度阈值，例如10，去除二维直方图中大于梯度阈值的高梯度像素点，以排除噪声。

步骤2、采用高斯模板对灰度直方图进行卷积运算，获得每一像素点对应的卷积结果Φ(x)。

作为一种可选实施例，可采用如下高斯模板进行卷积运算，

式中，g(x)为高斯模板，Φ(x)为卷积结果。

基于上述可以获得每一像素点对应的卷积结果。

步骤3、筛选获得卷积结果满足设定条件的像素点，将满足设定条件的像素点的数目作为聚类数目，其中，设定条件为Φ′(x_i)＝0，Φ″(x_i)<0且Φ(x_i)>0.005*max(Φ(x))。

具体地，满足设定条件的像素点即为极点。找出所有的满足Φ′(x_i)＝0，Φ″(x_i)<0且Φ(x_i)>0.005*max(Φ(x))的x_i，另C为集合{x_i}中x_i的个数，则可选择C个聚类中心。

基于上述实施例的内容，作为一种可选实施例，根据划分结果构造初始隶属度矩阵，包括但不限于：

按照如下规则生成初始隶属度矩阵，

式中，u_j为第j类的隶属度，C_j为划分结果中的第j类，c为聚类数目，x_i为第i个像素点的灰度值。

其中，对于每个像素点，其隶属度矩阵需要满足：

即每个像素点对于每个类的隶属度之和为1。

基于上述实施例的内容，作为一种可选实施例，提供一种基于初始隶属度矩阵及像素点的空间信息构造目标函数，并根据目标函数对土壤CT图像进行模糊辨识的方法，包括但不限于：

构造目标函数如下(该目标函数中，加号的后一项即为空间信息)，

式中，J为目标函数，v_j为聚类中心灰度，||·||为欧式距离的向量，c为聚类数目，u_ij为概率向量，m为常数，n为像素点个数，m为控制模糊度的常量，(m的值取可以为2)，x_r为像素点邻域的灰度中值，α为中心像素点的灰度值占中心像素点的邻域内像素点灰度值总和的比重(α的范围为0到1之间)；

获得满足目标函数的目标隶属度矩阵和目标聚类中心，

其中，目标隶属度矩阵为，

其中，目标聚类中心为，

根据目标隶属度矩阵及目标聚类中心，将每一像素点划分至隶属度最高的类。

其中，在计算获得目标隶属度矩阵和目标聚类中心后，u_ij、v_j更新并记录数据，根据最大隶属度矩阵的原则将其划分为相应的类。像素i属于隶属度最高的类C，即：

C_i＝arg{max(u_ij)}j＝1,2,…,c

基于上述实施例的内容，作为一种可选实施例，邻域为3*3邻域，相应地，

式中，x₅为中心像素点的灰度值，x₁至x₄以及x₆至x₉均为中心像素点的相邻像素的灰度值。

具体地，换言之，此时将x_i作为下式中的x₅进行α的计算，x₁到x₄以及x₆到x₉均为与x_i或x₅相邻的像素点的灰度值。

基于梯度-灰度综合信息赋值的初始条件与理想条件较为接近，因此，该方法可以较快地达到收敛状态，以自适应完成孔隙结构的辨识。

基于上述实施例的内容，作为一种可选实施例，提供一种对模糊辨识结果进行去模糊化处理，获得土壤CT图像中的孔隙结构辨识结果的方法，包括但不限于：

设定孔隙结构的辨识准则如下，

式中，c为聚类数目，k为一个变量，且k＝c-2，η为0至1之间的小数；

根据辨识准则，对土壤CT图像进行如下二值化处理，

式中，I为经过聚类后土壤CT图像中各像素点的灰度值，I的值为(0-1)，x_i为第i个像素点的灰度值。

具体地，设定孔隙辨识标准则对聚类结果进行优化，以将模糊聚类结果转化为确定性分类，完成孔隙结果的辨识。

综上，本发明实施例提供的方法引入梯度-灰度二维矩阵表征土壤不同物质间的差异。然后，通过矩阵分区确定聚类数目和初始化隶属度函数，指导传统FCM法进行迭代优化，以确保目标函数较快趋于稳定状态。最后，通过建立孔隙辨识准则保证孔隙结构辨识的精确性。该方法解决了初始值制约辨识精确度的问题，对于孔隙结构的定量研究具有重要意义。

本发明实施例提供的方法，是基于灰度-梯度特征的空间模糊C均值孔隙辨识方法，依次为建立灰度-梯度二维特征矩阵、确定聚类数目、自动分区并构造初始隶属度矩阵、构造加入空间约束的目标函数、确定孔隙辨识准则、优化模糊聚类结果，从而相较于现有技术，对土壤CT图像具有普适性，解决了初始值制约辨识精确度的问题，在保证孔隙辨识精度的前提下具有较高的执行效率。

基于上述实施例的内容，本发明实施例提供了一种基于土壤CT图像的孔隙辨识系统，该基于土壤CT图像的孔隙辨识系统用于执行上述方法实施例中的基于土壤CT图像的孔隙辨识方法。参见图4，该系统包括：确定模块401、划分模块402、构造模块403和辨识模块404；其中，确定模块401，用于建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目；划分模块402，用于根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；构造模块403，用于基于初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据目标函数对土壤CT图像进行模糊辨识；辨识模块404，用于对模糊辨识结果进行去模糊化处理，获得土壤CT图像中的孔隙结构辨识结果。

其中，可以利用CT扫描技术获得土壤CT图像。然后确定模块401建立土壤CT图像对应的灰度-梯度二维直方图，并基于灰度-梯度二维直方图确定聚类的数目。划分模块402根据聚类数据对灰度-梯度二维直方图进行区域划分，在划分完成后，可以根据划分得到的类别构造初始隶属度矩阵。构造模块403在传统的目标函数基础上，将其邻域内的空间信息作为约束条件，创建了新的目标函数，能够更加准确的判断像素点的类别。构造模块403在构造目标函数后，基于目标函数可以重新进行隶属度矩阵和聚类中心的计算，获得目标隶属度矩阵和目标聚类中心。从而能够根据每个像素点的隶属度对像素点进行分类，即可完成模糊辨识。由于孔隙结构的像素点相对于土壤CT图像中土壤的其他部分的颜色较深，因此辨识模块404可以将灰度值最大的一类像素点提取出来作为孔隙结构辨识结果，从而完成基于土壤CT图像的孔隙辨识。

本发明实施例提供的系统，通过建立灰度-梯度二维直方图、确定聚类数目、构造初始隶属度矩阵、构造加入空间信息的目标函数和去模糊化处理，从而相比于现有技术对土壤CT图像具有普适性，解决了初始值制约辨识精确度的问题，在保证了土壤孔隙辨识精度的前提下具有较高的执行效率。

本发明实施例提供了一种电子设备，如图5所示，该设备包括：处理器(processor)501、通信接口(Communications Interface)502、存储器(memory)503和通信总线504，其中，处理器501，通信接口502，存储器503通过通信总线504完成相互间的通信。处理器501可以调用存储器503上并可在处理器501上运行的计算机程序，以执行上述各实施例提供的基于土壤CT图像的孔隙辨识方法，例如包括：建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目；根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识；对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果。

此外，上述的存储器503中的逻辑指令可以通过软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读取存储介质中。基于这样的理解，本发明的技术方案本质上或者说对现有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个实施例方法的全部或部分步骤。而前述的存储介质包括：U盘、移动硬盘、只读存储器(ROM，Read-Only Memory)、随机存取存储器(RAM，Random Access Memory)、磁碟或者光盘等各种可以存储程序代码的介质。

本发明实施例还提供一种非暂态计算机可读存储介质，其上存储有计算机程序，该计算机程序被处理器执行时实现以执行上述各实施例提供基于土壤CT图像的孔隙辨识方法，例如包括：建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目；根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识；对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果。

以上所描述的电子设备等实施例仅仅是示意性的，其中作为分离部件说明的单元可以是或者也可以不是物理上分开的，作为单元显示的部件可以是或者也可以不是物理单元，即可以位于一个地方，或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本实施例方案的目的。本领域普通技术人员在不付出创造性的劳动的情况下，即可以理解并实施。

通过以上的实施方式的描述，本领域的技术人员可以清楚地了解到各实施方式可借助软件加必需的通用硬件平台的方式来实现，当然也可以通过硬件。基于这样的理解，上述技术方案本质上或者说对现有技术做出贡献的部分可以以软件产品的形式体现出来，该计算机软件产品可以存储在计算机可读存储介质中，如ROM/RAM、磁碟、光盘等，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行各个实施例或者实施例的某些部分方法。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明各实施例技术方案的精神和范围。

Claims (8)

1.一种基于土壤CT图像的孔隙辨识方法，其特征在于，包括：

建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目；

根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；

基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识；

对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果；

所述基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识，包括：

构造目标函数如下，

式中，J为目标函数，v_j为聚类中心灰度，||·||为欧式距离的向量，u_ij为概率向量，m为常数，n为像素点个数，m为控制模糊度的常量，

为像素点邻域的灰度中值，α为中心像素点的灰度值占所述中心像素点的邻域内像素点灰度值总和的比重，x_i为第i个像素点的灰度值；

获得满足所述目标函数的目标隶属度矩阵和目标聚类中心，

其中，所述目标隶属度矩阵为，

其中，所述目标聚类中心为，

根据所述目标隶属度矩阵及所述目标聚类中心，将每一所述像素点划分至隶属度最高的类；

所述对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果，包括：

设定所述孔隙结构的辨识准则如下，

式中，c为聚类数目，k为一个变量，且k＝c-2，η为0至1之间的小数；

根据所述辨识准则，对所述土壤CT图像进行如下二值化处理，

式中，I为经过聚类后所述土壤CT图像中各像素点的灰度值，x_i为第i个像素点的灰度值。

2.根据权利要求1所述的方法，其特征在于，所述建立土壤CT图像的灰度-梯度二维直方图，包括：

基于八邻域的拉普拉斯算子计算获得所述土壤CT图像的灰度-梯度二维直方图。

3.根据权利要求1所述的方法，其特征在于，所述通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目，包括：

剔除灰度-梯度二维直方图中梯度大于10的像素点，并将更新后的灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图；

采用高斯模板对所述灰度直方图进行卷积运算，获得每一所述像素点对应的卷积结果Φ(x)；

筛选获得所述卷积结果满足设定条件的所述像素点，将满足所述设定条件的像素点的数目作为所述聚类数目，其中，所述设定条件为Φ′(x_i)＝0，Φ″(x_i)<0且Φ(x_i)>0.005*max(Φ(x))。

4.根据权利要求1所述的方法，其特征在于，根据划分结果构造初始隶属度矩阵，包括：

按照如下规则生成所述初始隶属度矩阵，

式中，u_j为隶属度，C_j为划分结果中的第j类，c为聚类数目，x_i为第i个像素点的灰度值。

5.根据权利要求4所述的方法，其特征在于，所述邻域为3*3邻域，相应地，

式中，x₅为中心像素点的灰度值，x₁至x₄以及x₆至x₉均为中心像素点的相邻像素的灰度值。

6.一种基于土壤CT图像的孔隙辨识系统，其特征在于，包括：

确定模块，用于建立土壤CT图像的灰度-梯度二维直方图，通过将灰度-梯度二维直方图向灰度轴作投影获得一维灰度直方图，基于卷积后灰度直方图的极点数量，确定土壤图像的聚类数目；

划分模块，用于根据所述聚类数目对所述灰度-梯度二维直方图进行区域划分，并根据划分结果构造初始隶属度矩阵；其中，分区数＝梯度块数*灰度块数，聚类数目≥16时，分区数＝4*4；聚类数目≥11时，分区数＝3*4；聚类数目>7时，分区数＝3*3时；聚类数目≤7时，分区数＝2*3；

构造模块，用于基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识；

辨识模块，用于对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果；

所述基于所述初始隶属度矩阵及像素点邻域的空间信息构造目标函数，并根据所述目标函数对所述土壤CT图像进行模糊辨识，包括：

构造目标函数如下，

式中，J为目标函数，v_j为聚类中心灰度，||·||为欧式距离的向量，u_ij为概率向量，m为常数，n为像素点个数，m为控制模糊度的常量，

为像素点邻域的灰度中值，α为中心像素点的灰度值占所述中心像素点的邻域内像素点灰度值总和的比重，x_i为第i个像素点的灰度值；

获得满足所述目标函数的目标隶属度矩阵和目标聚类中心，

其中，所述目标隶属度矩阵为，

其中，所述目标聚类中心为，

根据所述目标隶属度矩阵及所述目标聚类中心，将每一所述像素点划分至隶属度最高的类；

所述对模糊辨识结果进行去模糊化处理，获得所述土壤CT图像中的孔隙结构辨识结果，包括：

设定所述孔隙结构的辨识准则如下，

式中，c为聚类数目，k为一个变量，且k＝c-2，η为0至1之间的小数；

根据所述辨识准则，对所述土壤CT图像进行如下二值化处理，

式中，I为经过聚类后所述土壤CT图像中各像素点的灰度值，x_i为第i个像素点的灰度值。

7.一种电子设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述程序时实现如权利要求1至5任一项所述基于土壤CT图像的孔隙辨识方法的步骤。

8.一种非暂态计算机可读存储介质，其上存储有计算机程序，其特征在于，该计算机程序被处理器执行时实现如权利要求1至5任一项所述基于土壤CT图像的孔隙辨识方法的步骤。

Images