CN105160662A - 基于局部高斯和玻尔兹曼的水平集图像分割方法 - Google Patents

Abstract

本发明涉及一种基于局部高斯和玻尔兹曼的水平集图像分割方法，主要解决现有技术对图像分割效率低及耗时的问题，其实现步骤是：1.输入待分割图像I到局部高斯模型中，得到图像I的能量函数E(I)；2.根据能量函数E(I)得到对应的演化方程Y(φ)；3.将演化曲线初始化为初始演化曲线φ(X,n)；4.求解初始演化曲线φ(X,n)的演化方程Y(φ)，得到新的演化曲线φ(X,n+1)；5.计算停止函数S，并将其余设置阈值T比较得出结果演化曲线φ_r；6.对结果演化曲线φ_r进行形态学开操作，得到最终分割结果φ_z。本发明提高了分割效率，加快了分割速度，可用于对自然图像、医学图像的快速分割。

Description

基于局部高斯和玻尔兹曼的水平集图像分割方法

技术领域

本发明属于图像处理技术领域，更进一步涉及一种基于水平集的图像分割方法，可用于自然图像、医学图像中目标的快速分割。

背景技术

在人们的实际生活中，绝大多数的信息都来源于图像。因为图像具有直观易懂，信息量大等特点，所以它是我们认知客观世界最重要的方式之一。

图像技术可分为图像前处理、图像分析、图像理解三部分。图像前处理主要是采集图像信息，并对其进行预处理的过程，这一过程主要包括图像采样、图像滤波、图像增强等内容；图像分析主要是提取感兴趣的目标和图像特征的过程，这一过程主要包括图像分割和特征提取；图像理解主要是对提取到的目标和特征进行研究的过程，这一过程主要包括目标跟踪和语义识别。图像分割在图像技术中占有十分重要的地位，一方面，图像分割是进行图像理解的基础，分割的结果将直接影响后续图像理解的处理效果；另一方面，图像分割出感兴趣的目标之后，把图像变为更为抽象的形式，大大减少了图像理解的复杂度及处理量，这使得进一步的高层分析、图像理解以及人工智能成为可能。图像分割技术已经广泛运用于实际生活的方方面面，如医学影像、人脸识别、交通控制等。

基于水平集的图像分割技术是近年来分割领域中的一个研究热点，它可以分割图像中具有复杂形状的目标对象，并能得到很好的结果。现有水平集分割模型有如下几种方法：

1.基于图像区域的水平集分割方法。这种方法通过图像的区域信息来构造能量函数。例如Chan等人[ChanTF,VeseLA.Activecontourswithoutedges[J],IEEETransactionsonImageProcessing,2001.pp:266-277]提出的无边界活动轮廓法，可用于分割灰度均匀的图像，但由于构造函数时默认图像灰度分布均匀，因此对灰度不均匀的自然图像分割效果不理想。

2.基于图像边界的水平集分割方法。这种方法通过图像的梯度信息来构造能量函数。例如Li等人[LiC,XuC,GuiC,etal.Levelsetevolutionwithoutre-initialization:anewvariationalformulation[C],IEEEConferenceonComputerVisionandPatternRecognition,2005.pp:430-436.]提出的无需重新初始化的变分法，可用于分割边界明显及梯度有意义的图像，但由于定义能量函数时使用的是图像梯度信息，因此对边界模糊及梯度无意义的医学图像图像分割效果不理想。

3.基于图像形状先验的水平集分割方法。这种方法通过图像的形状先验信息来构造能量函数。例如Chan等人[ChanT,ZhuW.Levelsetbasedshapepriorsegmentation[C].IEEEConferenceonComputerVisionandPatternRecognition,2005.pp:1164-1170.]提出的基于形状先验的水平集分割方法，可用于分割规则形状的目标，对形状不规则的目标分割效果不理想，且处理形状先验时计算量很大，使得分割效率不高。

发明内容

本发明的目的在于针对上述已有技术的不足，提出一种基于局部高斯和玻尔兹曼的水平集图像分割方法，以提高分割精度和分割效率。

为实现上述目的，本发明的技术方案包括如下步骤：

(1)输入待分割图像I到局部高斯模型中，得到关于图像I的能量函数E(I)；

(2)通过Euler-Lagrange方程得到能量函数E(I)对应的演化方程Y(φ)，其中，φ为演化曲线；

(3)将演化曲线初始化为初始演化曲线φ(X,n)，其中X为图像上的点，状态n＝0；

(4)求解演化曲线为初始演化曲线φ(X,n)的演化方程Y(φ)，得到新的演化曲线φ(X,n+1)；

(5)设置阈值T，计算停止函数S，比较阈值T和停止函数S的大小：

(5.1)根据新的演化曲线φ(X,n+1)和初始演化曲线φ(X,n)计算停止函数S：

$S=\frac{m_1-m_3}{\sqrt{{\mathrm{\σ}}_1-{\mathrm{\σ}}_3}}+\frac{m_2-m_4}{\sqrt{{\mathrm{\σ}}_2-{\mathrm{\σ}}_4}}$

其中，m₁,m₂,σ₁,σ₂分别为新的演化曲线φ(X,n+1)内外均值和内外方差，m₃,m₄,σ₃,σ₄分别为初始演化曲线φ(X,n)的内外均值和内外方差；

(5.2)设置阈值T，比较阈值T和停止函数S的大小：

当T＜S时，n＝n+1，继续迭代，返回步骤(4)；

当T≥S时，迭代终止，得到演化曲线φ_r＝φ(X,n+1)；

(6)对结果演化曲线φ_r进行形态学开操作，得到最终分割结果φ_z。

本发明与现有技术相比具有以下优点：

1)本发明由于采用了局部高斯模型，提高了灰度不均匀图像及边界模糊图像的分割效果；

2)本发明由于采用了玻尔兹曼技术，加快了分割速度，提高了分割效率。

附图说明

图1为本发明的实现流程图；

图2为用本发明和LGD模型对自然图像进行分割的对比图；

图3为用本发明和LGD模型对医学图像进行分割的对比图；

图4为用本发明和LGD模型对遥感图像进行分割的对比图。

具体实施方式

下面结合附图对本发明的实施率及效果作进一步详细描述。

参照图1，本发明的实施步骤如下：

步骤1：输入待分割图像I到局部高斯模型中，得到关于图像I的能量函数E(I)。

能量函数可以由如下模型得到：CV模型，图像边界模型，形状先验模型，局部高斯模型等。

本实例采用局部高斯模型，局部高斯模型的能量函数为：

$\begin{array}{c}E\left(I_0\right)=\mathrm{\α}\underset{i=1}{\overset{2}{\Σ}}\underset{\mathrm{\Ω}}{\∫}w(X-Y){\mathrm{logP}}_{i,X}\left(I_0\right(Y\left)\right)M_i\left(\mathrm{\φ}\right(Y\left)\right)dYdX\\ +\mathrm{\β}\underset{\mathrm{\Ω}}{\∫}|\▿H\left(\mathrm{\φ}\right(X\left)\right)|dX+\mathrm{\γ}Q\left(\mathrm{\φ}\right(X\left)\right)\end{array}$

其中，Ω为图像域，X为图像上的点，Y∈|X-Y|≤ρ,I₀(Y)为点Y的像素值，ρ为局部半径，α为局部高斯权重，β为长度参数，γ为惩罚项参数，w(X-Y)为高斯核函数，P_i,X(I₀(Y)),i＝1,2为局部高斯函数，φ为演化曲线，M₁＝H(φ(X))，M₂＝1-H(φ(X))，H为阶跃函数，为惩罚项，▽为拉普拉斯算子；

(1.1)把待分割图像I带入到高斯模型的能量函数表达式中，得到关于图像I的能量函数E(I)：

$\begin{array}{c}E\left(I\right)=\mathrm{\α}\underset{i=1}{\overset{2}{\Σ}}\underset{\mathrm{\Ω}}{\∫}w(X-Y){\mathrm{logP}}_{i,X}\left(I\right(Y\left)\right)M_i\left(\mathrm{\φ}\right(Y\left)\right)dYdX\\ +\mathrm{\β}\underset{\mathrm{\Ω}}{\∫}|\▿H\left(\mathrm{\φ}\right(X\left)\right)|dX+\mathrm{\γ}Q\left(\mathrm{\φ}\right(X\left)\right)\end{array}$

步骤2：通过Euler-Lagrange方程得到能量函数E(I)对应的演化方程Y(φ)。

Euler-Lagrange方程表达式为：

L为被积函数，g为被导函数，g'为被导函数g的一阶导，q为被导函数g中的变量。

根据Euler-Lagrange方程得到能量函数E(I)对应的演化方程Y(φ)的实施步骤如下：

(2.1)对能量函数E(I)求导，得到被积函数G，其表达式为：

$\begin{array}{c}G=\mathrm{\α}\underset{i=1}{\overset{2}{\Σ}}\underset{\mathrm{\Ω}}{\∫}w(X-Y){\mathrm{logP}}_{i,X}\left(\right(Y\left)\right)M_i\left(\mathrm{\φ}\right(Y\left)\right)dY\\ +\mathrm{\β}|\▿H\left(\mathrm{\φ}\right(X\left)\right)|+\frac{1}{2}\mathrm{\γ}{\left(\right|\▿\mathrm{\φ}\left(X\right)|-1)}^2\end{array};$

(2.2)当能量函数E(I)达到极小时，被积函数G对应的Euler-Lagrange方程为：

$\frac{\∂G}{\∂\mathrm{\φ}}-\left(\frac{\∂}{\∂x}\right(\frac{\∂G}{\∂{\mathrm{\φ}}_x})+\frac{\∂}{\∂y}\left(\frac{\∂G}{\∂{\mathrm{\φ}}_y}\right))=0,$ (x,y)为图像上的点。

(2.3)引入时间变量t，得到能量函数E(I)对应的演化方程Y(φ)：

$\begin{array}{c}Y\left(\mathrm{\φ}\right)=\frac{\∂\mathrm{\φ}}{\∂t}=\frac{\∂G}{\∂\mathrm{\φ}}-\left(\frac{\∂}{\∂x}\right(\frac{\∂G}{\∂{\mathrm{\φ}}_x})+\frac{\∂}{\∂y}\left(\frac{\∂G}{\∂{\mathrm{\φ}}_y}\right))\\ =\mathrm{\α}\mathrm{\δ}\left(\mathrm{\φ}\right)(e_1-e_2)+\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right)div\left(\frac{\▿\mathrm{\φ}}{|\▿\mathrm{\φ}|}\right))+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div(\frac{\▿\mathrm{\φ}}{|\▿\mathrm{\φ}|}))\end{array}$

其中， $e_j\left(X\right)=\underset{\mathrm{\Ω}}{\∫}w(X-Y)\[\log \left({\mathrm{\σ}}_j\right(Y\left)\right)+\frac{{\left({\mathrm{\μ}}_j\right(X)-I(Y\left)\right)}^2}{2{\mathrm{\σ}}_j^2\left(X\right)}\]dY,j=1,2.$

μ_j(X),σ_j(X)分别为点X的高斯均值和高斯标准差，δ(φ)是冲击函数。

步骤3：将演化曲线初始化为初始演化曲线φ(X,n)，其中X为图像上的点，设状态n＝0。

步骤4：求解演化曲线为初始演化曲线φ(X,n)的演化方程Y(φ)，得到新的演化曲线φ(X,n+1)。

新的演化曲线φ(X,n+1)可以通过数值迭代方程组和玻尔兹曼方程组得到，本实例采用但不限于玻尔兹曼方程求解新的演化曲线φ(X,n+1)。

玻尔兹曼方程组表达式包括如下主、次两个方程：

主方程为： $f_k(X,n+1)=f_k(X,n)+z_n\[f_k^{eq}-f_k(X,n)\]+{\mathrm{\ΔtF}}_n,k=1,2,3,4,5,$

次方程为： $\mathrm{\φ}(X,n+1)=\underset{k=0}{\overset{4}{\Σ}}{f}_k(X,n+1),$

其中f_k(X,n+1)为新的分布函数，f_k(X,n)为初始分布函数，为均衡分布函数，z_n为玻尔兹曼参数、F_n为外力，Δt为时间步长；

根据玻尔兹曼方程组得到新的演化曲线φ(X,n+1)的实施步骤如下：

(4.1)令|▽φ|＝1，将演化方程Y(φ)化简为：

$Y\left(\mathrm{\φ}\right)=div\left(\frac{2}{5\mathrm{\Δ}t}\right(\frac{1+5\mathrm{\Δ}t\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right)}{2}-\frac{1}{2})\·\▿\mathrm{\φ})+\frac{25}{2}\·\frac{2}{25}\left(\mathrm{\δ}\right(\mathrm{\φ}\left)\mathrm{\α}\right(e_1-e_2)+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div\left(\▿\mathrm{\φ}\right)\left)\right);$

(4.2)将扩散方程C_E与化简后的演化方程Y(φ)相对应，得出玻尔兹曼参数z，外力F。

Chapman-Enskog扩散方程C_E为：

$C_E=\frac{\∂\mathrm{\φ}}{\∂t}=div\left(\frac{2}{5\mathrm{\Δ}t}\right(\frac{1}{z}-\frac{1}{2})\·\▿\mathrm{\φ})+\frac{25}{2}F,$

玻尔兹曼参数z，外力F为：

$z=\frac{2}{1+5\mathrm{\Δ}t\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right)},F=\frac{2}{25}\left(\mathrm{\δ}\right(\mathrm{\φ}\left)\mathrm{\α}\right(e_1-e_2)+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div\left(\▿\mathrm{\φ}\right)\left)\right);$

(4.3)把初始演化曲线φ(X,n)代入到(4.2)的波尔兹曼参数z和外力F等式中，得到在状态n时的玻尔兹曼参数z_n和外力F_n：

$z_n=\frac{2}{1+5\mathrm{\Δ}t\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right(X,n))},$

$F_n=\frac{2}{25}\left(\mathrm{\δ}\right(\mathrm{\φ}\left(X,n\right)\left)\mathrm{\α}\right(e_1-e_2)+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div\left(\▿\mathrm{\φ}\right)\left)\right);$

(4.4)令初始分布函数 $f_k(X,0)=f_k^{eq}=\frac{1}{5}\mathrm{\φ}(X,0),k=1,2,...,5,$ 为均衡分布函数，把玻尔兹曼参数z_n、外力F_n、初始分布函数f_k(X,n)和均衡分布函数代入玻尔兹曼的主方程中，计算新的分布函数f_k(X,n+1)：

$f_k(X,n+1)=f_k(X,n)+z_n\[f_k^{eq}-f_k(X,n)\]+{\mathrm{\ΔtF}}_n,k=1,2,3,4,5;$

(4.5)把新的分布函数f_k(X,n+1)代入玻尔兹曼的次方程中，计算新的演化曲线φ(X,n+1)：

$\mathrm{\φ}(X,n+1)=\underset{k=0}{\overset{4}{\Σ}}{f}_k(X,n+1).$

步骤5：计算停止函数S，设置阈值T，比较阈值T和停止函数S的大小。

(5.1)根据新的演化曲线φ(X,n+1)和初始演化曲线φ(X,n)计算停止函数S：

$S=\frac{m_1-m_3}{\sqrt{{\mathrm{\σ}}_1-{\mathrm{\σ}}_3}}+\frac{m_2-m_4}{\sqrt{{\mathrm{\σ}}_2-{\mathrm{\σ}}_4}}$

其中，m₁为新的演化曲线φ(X,n+1)内均值，m₂为新的演化曲线φ(X,n+1)外均值，σ₁为新的演化曲线φ(X,n+1)内方差，σ₂为新的演化曲线φ(X,n+1)外方差，m₃为初始演化曲线φ(X,n)的内均值，m₄为初始演化曲线φ(X,n)的外均值，σ₃为初始演化曲线φ(X,n)的内方差，σ₄为初始演化曲线φ(X,n)的外方差；

(5.2)设置阈值T，比较阈值T和停止函数S的大小：

当T＜S时，n＝n+1，继续迭代，返回步骤(4)；

当T≥S时，迭代终止，得到迭代演化曲线φ_r＝φ(X,n+1)。

步骤6：对迭代演化曲线φ_r进行形态学开操作，得到最终分割结果φ_z。

(6.1)对迭代演化曲线φ_r进行形态学腐蚀，去除演化曲线φ_r中的孤立点，缩小其目标区域，得到临时演化曲线φ_f；

(6.2)对临时演化曲线φ_f进行形态学膨胀，恢复目标区域大小，得到最终分割结果φ_z。

本发明的效果可以通过使用以下仿真实验进一步说明

1、仿真条件

本发明是在中央处理器为Intel(R)Core(TM)i32.93GHZ、内存2G、WINDOWS7操作系统上，运用MATLAB软件进行的仿真。

2、仿真内容

仿真1，用本发明和现有LGD模型对自然图像进行分割，结果如图2所示，其中：

图2(a)为初始自然图像，

图2(b)为使用LGD模型对自然图像的初始化，

图2(c)为使用LGD模型对自然图像进行分割的结果，

图2(d)为使用本发明对自然图像进行分割的结果，

对比图2(c)和图2(d)，可以看出本发明对自然图像分割效果好。

仿真2，用本发明和现有LGD模型对医学图像进行分割，结果如图3所示，其中：

图3(a)为初始医学图像，

图3(b)为使用LGD模型对医学图像的初始化，

图3(c)为使用LGD模型对医学图像进行分割的结果，

图3(d)为使用本发明对医学图像进行分割的结果，

对比图3(c)和图3(d)，可以看出本发明对医学图像分割效果好。

仿真3，用本发明和现有LGD模型对遥感图像进行分割，结果如图4所示，其中：

图4(a)为初始遥感图像，

图4(b)为使用LGD模型对遥感图像的初始化，

图4(c)为使用LGD模型对遥感图像进行分割的结果，

图4(d)为使用本发明对遥感图像进行分割的结果，

对比图4(c)和图4(d)，可以看出本发明对遥感图像分割效果好,可以取得全局最优分割结果。

仿真1、仿真2和仿真3所有的迭代次数和时间如表1所示：

表1

由表1可知，与LGD模型相比，本发明只需要迭代几次便可得到结果，图像分割所有时间远远小于LGD模型，图像分割效率较高。

Claims (4)

1.一种基于局部高斯和玻尔兹曼的水平集图像分割方法，包括如下步骤：

(1)输入待分割图像I到局部高斯模型中，得到关于图像I的能量函数E(I)；

(2)通过Euler-Lagrange方程得到能量函数E(I)对应的演化方程Y(φ)，其中φ为演化曲线；

(3)将演化曲线初始化为初始演化曲线φ(X,n)，其中X为图像上的点，状态n＝0；

(4)求解演化曲线为初始演化曲线φ(X,n)的演化方程Y(φ)，得到新的演化曲线φ(X,n+1)；

(5)计算停止函数S，设置阈值T，比较阈值T和停止函数S的大小：

(5.1)根据新的演化曲线φ(X,n+1)和初始演化曲线φ(X,n)计算停止函数S：

$S=\frac{m_1-m_3}{\sqrt{{\mathrm{\σ}}_1-{\mathrm{\σ}}_3}}+\frac{m_2-m_4}{\sqrt{{\mathrm{\σ}}_2-{\mathrm{\σ}}_4}}$

其中，m₁,m₂,σ₁,σ₂分别为新的演化曲线φ(X,n+1)内外均值和内外方差，m₃,m₄,σ₃,σ₄分别为初始演化曲线φ(X,n)的内外均值和内外方差；

(5.2)设置阈值T，比较阈值T和停止函数S的大小：

当T＜S时，n＝n+1，继续迭代，返回步骤(4)；

当T≥S时，迭代终止，得到演化曲线φ_r＝φ(X,n+1)；

(6)对演化曲线φ_r进行形态学开操作，得到最终分割结果φ_z。

2.根据权利要求1所述的基于局部高斯和玻尔兹曼的水平集图像分割方法，其中所述步骤(2)中通过Euler-Lagrange方程得到能量函数E(I)对应的演化方程Y(φ)，按如下步骤进行：

(2.1)对能量函数E(I)求导，得到被积函数G，其表达式为：

$\begin{array}{c}G=\mathrm{\α}\underset{i=1}{\overset{2}{\Σ}}\underset{\mathrm{\Ω}}{\∫}w(X-Y){\mathrm{logP}}_{i,X}\left(I\right(Y\left)\right)M_i\left(\mathrm{\φ}\right(Y\left)\right)dY\\ +\mathrm{\β}|\▿H\left(\mathrm{\φ}\right(X\left)\right)|+\frac{1}{2}\mathrm{\γ}{\left(\right|\▿\mathrm{\φ}\left(X\right)|-1)}^2\end{array}$

其中，Ω为图像域，X为图像上的点，Y∈|X-Y|≤ρ,I(Y)为点Y的像素值，ρ为局部半径，α为局部高斯权重，β为长度参数，γ为惩罚项参数，w(X-Y)为高斯核函数，M₁＝H(φ(X))，M₂＝1-H(φ(X))，H(φ)为阶跃函数，▽为拉普拉斯算子；

(2.2)当能量函数E(I)达到极小时，被积函数G对应的Euler-Lagrange方程为：

$\frac{\∂G}{\∂\mathrm{\φ}}-\left(\frac{\∂}{\∂x}\right(\frac{\∂G}{\∂{\mathrm{\φ}}_x})+\frac{\∂}{\∂y}\left(\frac{\∂G}{\∂{\mathrm{\φ}}_y}\right))=0,$ (x,y)为图像上的点。

(2.3)引入时间变量t，得到能量函数E(I)对应的演化方程Y(φ)：

$\begin{array}{c}Y\left(\mathrm{\φ}\right)=\frac{\∂\mathrm{\φ}}{\∂t}=\frac{\∂G}{\∂\mathrm{\φ}}-\left(\frac{\∂}{\∂x}\right(\frac{\∂G}{\∂{\mathrm{\φ}}_x})+\frac{\∂}{\∂y}\left(\frac{\∂G}{\∂{\mathrm{\φ}}_y}\right))\\ =\mathrm{\α}\mathrm{\δ}\left(\mathrm{\φ}\right)(e_1-e_2)+\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right)div\left(\frac{\▿\mathrm{\φ}}{|\▿\mathrm{\φ}|}\right))+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div(\frac{\▿\mathrm{\φ}}{|\▿\mathrm{\φ}|}))\end{array}$

其中， $e_j\left(X\right)=\underset{\mathrm{\Ω}}{\∫}w(X-Y)\[\log \left({\mathrm{\σ}}_j\right(Y\left)\right)+\frac{{\left({\mathrm{\μ}}_j\right(X)-I(Y\left)\right)}^2}{2{\mathrm{\σ}}_j^2\left(X\right)}\]dY,j=1,2.$

μ_j(X),σ_j(X)分别为点X的高斯均值和高斯标准差，δ(φ)是冲击函数。

3.根据权利要求1所述的基于局部高斯和玻尔兹曼的水平集图像分割方法，其中所述步骤(4)求解演化曲线为初始演化曲线φ(X,n)的演化方程Y(φ)，得到新的演化曲线φ(X,n+1)，按如下步骤进行：

(4.1)令|▽φ|＝1，将演化方程Y(φ)化简为：

$Y\left(\mathrm{\φ}\right)=div\left(\frac{2}{5\mathrm{\Δ}t}\right(\frac{1+5\mathrm{\Δ}t\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right)}{2}-\frac{1}{2})\·\▿\mathrm{\φ})+\frac{25}{2}\·\frac{2}{25}\left(\mathrm{\δ}\right(\mathrm{\φ}\left)\mathrm{\α}\right(e_1-e_2)+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div\left(\▿\mathrm{\φ}\right)\left)\right)$

(4.2)将Chapman-Enskog扩散方程C_E与化简后的演化方程Y(φ)相对应，得出玻尔兹曼参数z，外力F：

$z=\frac{2}{1+5\mathrm{\Δ}t\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right)},F=\frac{2}{25}\left(\mathrm{\δ}\right(\mathrm{\φ}\left)\mathrm{\α}\right(e_1-e_2)+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div\left(\▿\mathrm{\φ}\right)\left)\right)$

其中，Δt为时间步长， $C_E=\frac{\∂\mathrm{\φ}}{\∂t}=div\left(\frac{2}{5\mathrm{\Δ}t}\right(\frac{1}{z}-\frac{1}{2})\·\▿\mathrm{\φ})+\frac{25}{2}F;$

(4.3)把初始演化曲线φ(X,n)代入到(4.2)的波尔兹曼参数z和外力F等式中，得到在状态n时的玻尔兹曼参数z_n和外力F_n：

$z_n=\frac{2}{1+5\mathrm{\Δ}t\mathrm{\β}\mathrm{\δ}\left(\mathrm{\φ}\right(X,n))},F_n=\frac{2}{25}\left(\mathrm{\δ}\right(\mathrm{\φ}\left(X,n\right)\left)\mathrm{\α}\right(e_1-e_2)+\mathrm{\γ}({\▿}^2\mathrm{\φ}-div\left(\▿\mathrm{\φ}\right)\left)\right)$

(4.4)令初始分布函数 $f_k(X,0)=f_k^{eq}=\frac{1}{5}\mathrm{\φ}(X,0),k=1,2,\mathrm{...},5,$ 为均衡分布函数，把参数z_n、外力F_n、初始分布函数f_k(X,n)和均衡分布函数代入玻尔兹曼的主方程中，计算新的分布函数f_k(X,n+1)：

$f_k(X,n+1)=f_k(X,n)+z_n\[f_k^{eq}-f_k(X,n)\]+{\mathrm{\ΔtF}}_n,k=1,2,3,4,5$

(4.5)把新的分布函数f_k(X,n+1)代入玻尔兹曼的次方程中，计算新的演化曲线φ(X,n+1)：

$\mathrm{\φ}(X,n+1)=\underset{k=0}{\overset{4}{\Σ}}{f}_k(x,n+1).$

4.根据权利要求1所述的基于局部高斯和玻尔兹曼的水平集图像分割方法，其中所述步骤(6)中对演化曲线φ_r进行形态学开操作，按如下步骤进行：

(6.1)对演化曲线φ_r进行形态学腐蚀，去除演化曲线φ_r中的孤立点，缩小其目标区域，得到临时演化曲线φ_f；

(6.2)对临时演化曲线φ_f进行形态学膨胀，恢复目标区域的大小，得到最终分割结果φ_z。