CN104537673A - 基于多阈值和自适应模糊聚类的红外图像分割方法 - Google Patents

Abstract

本发明公开了一种基于多阈值和自适应模糊聚类的红外图像分割方法，主要解决现有多阈值分割方法在分割过程中存在伪峰干扰，分割结果不理想的问题。其实现步骤是：(1)输入原始红外图像；(2)利用引入控制因子的多阈值算法对原始红外图像进行粗分割；(3)对粗分割后图像进行形态学平滑处理；(4)随机选取细分割图像所需的聚类中心并根据自适应函数确定聚类数目；(5)对平滑后图像的像素点进行模糊聚类，得到最终的分割结果图像。本发明能在保证分割效率的同时提高分割准确度，具有分割结果轮廓清晰、目标完整的优点，可有效应用于红外精确制导、目标识别与跟踪。

Description

基于多阈值和自适应模糊聚类的红外图像分割方法

技术领域

本发明属于图像信息处理领域，涉及红外图像分割方法，可应用于红外目标探测和跟踪系统。

背景技术

图像分割，是指将一幅图像分解成各个有意义的部分或物体，它是计算机视觉领域和图像信息处理领域中最底层的处理技术。图像分割在图像分析和模式识别方面起着重要的作用，是图像目标特征提取、识别、跟踪及分类的基础。其中，红外图像分割在目标对象自动识别中起着的特殊作用。近年来，中外学者在红外图像分割的技术探索方面做出了很多贡献，并提出了许多方法，如边缘检测法、阈值分割法、区域生长法等。其中，阈值分割法因其实现简单、计算量小、自适应性强、性能较稳定且不许需要先验知识而成为图像分割中最基本和应用最为广泛的分割技术。

1975年Otsu在一篇关于灰度直方图的阈值选择技术的文章中提出了最为经典、应用最为广泛的最大类间方差法。1980年Pun将信息论中熵是平均信息的一种标准这一理论首先引入到图像分割技术中来，将后验熵最大的灰度值作为最佳阈值。但上述两种方法均是二元分类方法，1989年Lim和Lee在一篇基于阈值和模糊C均值技术的彩色图像分割算法的文章中提出了通过求导的方式直接对高斯平滑后的灰度统计直方图求峰谷值，并由对应的波谷值位置作为分割图像的多个阈值，实现了多元分类。

传统的多阈值算法分割图像的实现过程为：首先，计算输入图像的灰度统计直方图；其次，选择合适的高斯平滑因子对灰度统计直方图曲线进行平滑，以滤除灰度统计直方图中的无意义波峰；最后，计算平滑后灰度统计直方图的波谷位置，并将其作为分割图像的多个阈值，完成图像的分割。由于传统的多阈值分割方法不需要任何先验知识，完全靠空间中的样本点进行驱动，并且处理速度快，因此该方法被广泛地应用到图像的分割领域中。但是，研究人员发现由于多阈值算法本身的特性，使得该算法存在严重的伪峰干扰，从而使得用多阈值方法对图像进行分割时无法将目标和背景很好地分割开，造成分割结果不准确。

李建萍等人在论文“一种复杂背景下的目标分割算法”(《弹箭与制导学报》，文章编号：1673-9728(2010)04-0197-04)中提出了一种基于多阈值的复杂背景下的目标分割算法。该算法的实现过程为：首先对复杂背景下的图像进行灰度形态学滤波和增强，增强目标和背景的灰度对比度；然后利用图像直方图独立峰搜索方法确定图像阈值，并通过引入独立峰的峰谷比、峰面积、峰宽三个特征量对图像阈值进行优化，从而实现对图像的分割。该方法的不足之处在于，虽然引入了独立峰的判断准则，能滤除一定的伪峰干扰，但其特征量的选取无法自适应得到，不利于实时处理。

南京理工大学王琼等人申请的专利“基于灰度直方图的多阈值分割法”(专利申请号：CN201210257574.8，公开号：CN103578095A)公开了一种基于灰度直方图的多阈值分割方法。该方法的实现过程为：首先根据类间方差最大法，找到图像的第一个阈值，然后以二分法为基础，对分割的图像的两部分分别求其类间的最大方差，再把两方差进行比较，找到最大的方差对应的阈值作为下一次分割的阈值点，以此类推，直至找到给定的阈值数为止；然后，根据平滑后的直方图，找到所有的波谷点；最后，把得到的所有阈值与波谷点进行比较，找到最靠近阈值的波谷点作为最终的阈值。该方法的不足之处在于，虽然以类间方差最大的思想寻找图像的最优阈值，改善了错分割或过分割的现象，但是该算法是以局部最优思想为基础的，没有考虑到图像的全局信息，且阈值数目得预先设定，在保证了算法运算效率的情况下，无法保证算法的实时性和分割结果的准确性。

发明内容

本发明方法的目的在于针对上述已有技术的不足，提出一种基于多阈值和自适应模糊聚类的红外图像分割方法，以在保证分割效率的同时提高分割的准确性。

实现本发明目的的技术方案是：利用引入控制因子的多阈值方法对待分割红外图像进行粗分割并形态学处理；随机选取模糊聚类的聚类中心并利用自适应函数确定聚类数目，对图像进行模糊聚类，实现对红外图像的最终分割，其具体步骤包括如下：

(1)输入原始红外图像I，计算其所有像素点数Co；

(2)利用原始红外图像I的一维灰度统计直方图函数H(l)，计算其灰度级H；

(3)对原始红外图像I进行粗分割：

3a)对一维灰度统计直方图函数H(l)进行平滑处理，得到平滑后的灰度统计直方图函数H₁(l)，并求取H₁(l)的差分D；

3b)寻找统计直方图函数H₁(l)的谷值点，即差分D由负到正的过零点位置，并将所有的谷值点连接起来，按从小到大排序，得到大小为1×p的谷值点矩阵M，其中p为谷值点个数；

3c)确定阈值点和阈值矩阵Th：

3c1)在多阈值算法中引入峰面积控制因子ψ和峰宽控制因子μ，通过峰面积控制因子ψ和峰宽控制因子μ分别计算得到灰度统计直方图H₁(l)中的相邻谷值点的最小峰面积MA和最小峰宽MW；

3c2)根据谷值点矩阵M计算相邻谷值点之间灰度统计直方图之和CT(n)及相邻谷值点的距离WT(n)，其中n＝1,2,...,p-1；

3c3)判断CT(n)＞MA且WT(n)＞MW是否满足，若满足，则将前一个谷值点标记为阈值点；若不满足，则合并两个谷值点，消去该波段，直至所有的谷值点全部标记完毕，得到阈值矩阵Th；

3d)对原始红外图像I中的像素点进行遍历，将像素点灰度值在阈值矩阵Th中相邻值域之间的像素点合并成一类，标记每一个像素点的类别，并分别计算总类别数q、每一类的像素点总个数N_t以及每一类的像素点的灰度值总和S_t，其中t＝1,2,...,q；

3e)计算合并后每个类别中的全部像素点的值域均值并用g_t代替原像素点的值，得到粗分割后图像I_mt及其像素点的值为g_j，其中j＝1,2,...,Co；

(4)对粗分割后图像I_mt进行形态学平滑处理，得到平滑后图像I_ms；

(5)对平滑后图像I_ms进行模糊聚类细分割：

5a)设定迭代停止条件ε，其中0＜ε＜1，聚类数目的初始值c＝2，自适应函数初始值L(1)＝0，迭代次数的初始值k＝0，模糊加权指数m＝2；

5b)随机初始化聚类中心矩阵其中i＝1,2,...,c；

5c)求取平滑后图像I_ms中第j个像素点与第i个聚类中心的第k次迭代的距离并计算第k次迭代的隶属度矩阵U^(k)和第k次迭代后的聚类中心矩阵V^(k+1)；

5d)对V^(k)和V^(k+1)进行比较，如果||V^(k+1)-V^(k)||＜ε，则迭代停止，执行步骤5e)，否则令k＝k+1，返回步骤5c)；

5e)计算聚类数目为c时的自适应函数L(c)，当2＜c＜q时，若存在满足L(c-1)＞L(c-2)且L(c-1)＞L(c)的点，则取聚类数目c＝c-1时的聚类结果作为图像细分割的结果，否则令c＝c+1，返回步骤5b)；当c＞q时，则取聚类数目时的聚类结果作为图像细分割的结果，得到最终的分割结果图像。

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

第一、本发明利用引入控制因子的多阈值方法对红外图像进行粗分割后，采用模糊聚类方法对粗分割后的图像进行细分割，有效地改善了多阈值方法中的伪峰干扰现象，克服了单一使用多阈值方法对红外图像进行分割时分割结果不理想的问题。

第二、本发明对粗分割后的图像进行细分割时采用的改进的模糊聚类方法能自适应获得待分割图像的聚类数目，有效地避免了传统模糊聚类方法需要对聚类结果的有效性进行验证，从而使得本发明在提高了图像分割准确度的同时，保证了图像分割的实时性。

附图说明

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

图2为用本发明对射击人红外图像的分割结果图；

图3为用本发明对丛林人红外图像的分割结果图；

图4为用本发明对持枪人红外图像的分割结果图；

图5为用本发明对树林单人红外图像的分割结果图；

以下结合附图对本发明具体实施方式和效果做进一步详细描述。

具体实施方式

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

步骤1.输入原始红外图像I，并计算其所有像素点数Co。

1a)统计原始红外图像I的行数作为其长度L，其列数作为其宽度W；

1b)计算原始红外图像I的所有像素点数Co＝L*W。

步骤2.利用原始红外图像I的一维灰度统计直方图函数H(l)，计算其灰度级H。

将原始红外图像I从像素空间映射到其灰度直方图特征空间得到其一维灰度统计直方图函数H(l)，并将该直方图函数H(l)的横坐标长度作为其灰度级H。

步骤3.对原始红外图像I进行粗分割。

3a)使用Witkin的高斯卷积平滑运算对一维灰度统计直方图函数H(l)进行平滑处理，得到平滑后的灰度统计直方图函数 $H_1\left(l\right)=H\left(l\right)*g(l,\mathrm{\τ})=H\left(l\right)*\frac{1}{\sqrt{2\mathrm{\π}}\mathrm{\τ}}\exp \left(\frac{{-x}^2}{{2\mathrm{\τ}}^2}\right),$ 其中，l为原始红外图像I中像素点映射到灰度统计直方图特征空间的灰度值，τ为高斯平滑滤波器的标准差；

3b)求取直方图函数H₁(l)的差分D：

D＝[H₁(2)-H₁(1)H₁(3)-H₁(2)...H₁(H)-H₁(H-1)]，

并寻找直方图函数H₁(l)的谷值点，即差分D由负到正的过零点位置，将所有的谷值点连接起来，按从小到大排序，得到大小为1×p的谷值点矩阵M，其中p为谷值点个数；

3c)确定阈值点和阈值矩阵Th；

3c1)统计相邻谷值点间灰度统计直方图之和CT(n)，计算相邻谷值点的距离WT(n)＝M(n+1)-M(n)，其中n＝1,2,...,p-1；

3c2)计算灰度统计直方图函数H₁(l)中相邻谷值点的最小峰面积MA＝ψ*Co和最小峰宽MW＝μ*H，其中，ψ为峰面积控制因子，μ为峰宽控制因子；

3c3)判断是否满足CT(n)＞MA且WT(n)＞MW这两个条件，若满足，则将前一个谷值点标记为阈值点；若不满足，则合并两个谷值点，消去该波段，直至所有的谷值点全部标记完毕，得到阈值矩阵Th；

3d)对原始红外图像I中的像素点进行遍历，将像素点灰度值在阈值矩阵Th中相邻值域之间的像素点合并成一类，标记每一个像素点的类别，并分别计算总类别数q、每一类的像素点总个数N_t以及每一类的像素点的灰度值总和S_t，其中t＝1,2,...,q；

3e)计算合并后每个类别中的全部像素点的值域均值并用g_t代替原像素点的值，得到粗分割后的图像I_mt及其像素点的值为g_j，其中j＝1,2,...,Co。

步骤4.对粗分割后图像I_mt进行形态学平滑处理。

选取半径为3的圆盘结构元素，先对粗分割后图像I_mt进行闭运算，再进行开运算，得到平滑后图像I_ms。

步骤5.设定模糊聚类细分割初始值。

5a)设定迭代停止条件ε，其中0＜ε＜1，聚类数目的初始值c＝2，自适应函数的初始值L(1)＝0，迭代次数的初始值k＝0，模糊加权指数m＝2；

5b)随机初始化聚类中心 $V^{\left(0\right)}=\left[\begin{array}{ccccc}{v}_1^{\left(0\right)}& ...& v_i^{\left(0\right)}& ...& v_c^{\left(0\right)}\end{array}\right],$ 其中i＝1,2,...,c。

步骤6.对平滑后图像I_ms进行模糊聚类细分割。

6a)计算平滑后图像I_ms中第j个像素点与第i个聚类中心第k次迭代的距离：

$d(g_j,v_i^{\left(k\right)})=\mathrm{In}(1+\mathrm{\β}||g_j-v_i^{\left(k\right)}|\left|\right),$

其中，β表示调节曲线曲率的自由度参数，||·||表示计算欧氏距离；

6b)根据第k次迭代的距离计算第k次迭代的隶属度矩阵U^(k)的第i行第j列的元素其中，j'表示平滑后图像I_ms聚类数目中第j'类；

6c)根据得到第k次迭代的隶属度矩阵U^(k)：

6d)根据计算第k次迭代后的聚类中心矩阵V^(k+1)的第i个元素

$v_i^{(k+1)}=\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m{g}_j/\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m,$

由得到第k次迭代后的聚类中心矩阵 $V^{(k+1)}=\left[\begin{array}{ccccc}{v}_1^{(k+1)}& ...& v_i^{(k+1)}& ...& v_c^{(k+1)}\end{array}\right];$

6e)对V^(k)和V^(k+1)进行比较，如果||V^(k+1)-V^(k)||＜ε，则迭代停止，执行步骤7，否则令k＝k+1，返回步骤6a)。

步骤7.计算聚类数目为c时的自适应函数L(c)：

$L\left(c\right)=\frac{\underset{i=1}{\overset{c}{\mathrm{\Σ}}}\left(\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m\right)d^2(v_i^{\left(k\right)},X)/(c-1)}{\underset{i=1}{\overset{c}{\mathrm{\Σ}}}\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m{d}^2(g_j,v_i^{\left(k\right)})/(\mathrm{Co}-c)},$

其中，表示平滑后图像I_ms中所有像素点的中心向量。

步骤8.确定聚类数目c和聚类分割结果，输出最终的分割结果图像。

当2＜c＜q时，判断是否存在满足L(c-1)＞L(c-2)且L(c-1)＞L(c)这两个条件的点，若存在，则取聚类数目c＝c-1时的聚类结果作为图像细分割的结果，否则，令c＝c+1，返回步骤5b)；当c＞q时，取聚类数目时的聚类结果作为图像细分割的结果，得到最终的分割结果图像。

以下通过仿真对本发明的效果做进一步的描述。

1.仿真条件：

选用射击人红外图像、丛林人红外图像、持枪人红外图像和树林单人红外图像四组数据进行实验，实验在CPU为Core(TM)i3-41303.40GHz、内存为4.00GB的WINDOWS 7系统上使用MatlabR2012b进行仿真。

实验中所使用的原始数据如图2(a)、图3(a)、图4(a)和图5(a)所示，其中图2(a)所示的射击人红外图像和图4(a)所示的持枪人红外图像来源于Terravic Weapon红外数据库；图3(a)所示的丛林人红外数据和图5(a)所示的树林单人红外图像来源于TerravicMotion红外数据库。

实验参数的选取：下列所有实验中，高斯平滑滤波器的标准差τ＝0.5，峰面积控制因子ψ＝0.0001，峰宽控制因子μ＝0.15，自由度参数β＝0.1，模糊聚类迭代停止阈值ε＝0.00001，FCM聚类数目设定与本发明算法自适应得到的聚类数目一致。

2.仿真内容：

仿真实验1：分别运用多阈值法、FCM法和本发明方法对射击人红外图像进行分割，实验结果如图2所示，其中图2(a)是原始射击人红外图像、图2(b)是用多阈值法的分割结果图像、图2(c)是用FCM法的分割结果图像、图2(d)是用本发明方法的分割结果图像。

仿真实验2：分别运用多阈值法、FCM法和本发明方法对丛林人红外图像进行分割，实验结果如图3所示，其中图3(a)是原始丛林人红外图像、图3(b)是用多阈值法的分割结果图像、图3(c)是用FCM法的分割结果图像、图3(d)是用本发明方法的分割结果图像。

仿真实验3：分别运用多阈值法、FCM法和本发明方法对持枪人红外图像进行分割，实验结果如图4所示，其中图4(a)是持枪人原始红外图像、图4(b)是用多阈值法的分割结果图像、图4(c)是用FCM法的分割结果图像、图4(d)是用本发明方法的分割结果图像。

仿真实验4：分别运用多阈值法、FCM法和本发明方法对树林单人红外图像进行分割，实验结果如图5所示，其中图5(a)是原始树林单人红外图像、图5(b)是用多阈值法的分割结果图像、图5(c)是用FCM法的分割结果图像、图5(d)是用本发明方法的分割结果图像。

3.仿真实验结果分析：

从图2(b)、图3(b)、图4(b)、图5(b)可以看出，采用传统的多阈值方法对红外图像进行分割时，不能将图像目标正确分割出来，存在严重的伪峰干扰，图像目标会被继续分割成许多小块，分割结果不准确；

从图2(c)、图3(c)、图4(c)、图5(c)可以看出，采用FCM法对红外图像进行分割时，图像目标虽然不会被分割成很多小块，但图像目标与图像背景之间的区域连通，产生了明显的误分现象；

从图2(d)、图3(d)、图4(d)、图5(d)可以看出，本发明方法能够较好的去除背景干扰，准确地将图像目标分割出来，分割出的图像目标轮廓清晰、目标完整，分割结果无错分现象。

Claims (5)

1.基于多阈值和自适应模糊聚类的红外图像分割方法，包括如下步骤：

(1)输入原始红外图像I，计算其所有像素点数Co；

(2)利用原始红外图像I的一维灰度统计直方图函数H(l)，计算其灰度级H；

(3)对原始红外图像I进行粗分割：

3a)对一维灰度统计直方图函数H(l)进行平滑处理，得到平滑后的灰度统计直方图函数H₁(l)，并求取H₁(l)的差分D；

3b)寻找直方图函数H₁(l)的谷值点，即差分D由负到正的过零点位置，并将所有的谷值点连接起来，按从小到大排序，得到大小为1×p的谷值点矩阵M，其中p为谷值点个数；

3c)确定阈值点和阈值矩阵Th：

3c1)在多阈值算法中引入峰面积控制因子ψ和峰宽控制因子μ，通过峰面积控制因子ψ和峰宽控制因子μ分别计算得到灰度统计直方图H₁(l)中相邻谷值点的最小峰面积MA和最小峰宽MW；

3c2)根据谷值点矩阵M计算相邻谷值点之间灰度统计直方图之和CT(n)及相邻谷值点的距离WT(n)，其中n＝1,2,…,p-1；

3c3)判断CT(n)＞MA且WT(n)＞MW是否满足，若满足，则将前一个谷值点标记为阈值点；若不满足，则合并两个谷值点，消去该波段，直至所有的谷值点全部标记完毕，得到阈值矩阵Th；

3d)对原始红外图像I中的像素点进行遍历，将像素点灰度值在阈值矩阵Th中相邻值域之间的像素点合并成一类，标记每一个像素点的类别，并分别计算总类别数q、每一类的像素点总个数N_t以及每一类的像素点的灰度值总和S_t，其中t＝1,2，…，q；

3e)计算合并后每个类别中的全部像素点的值域均值并用g_t代替原像素点的值，得到粗分割后图像I_mt及其像素点的值为g_j，其中j＝1,2，…，Co；

(4)对粗分割后图像I_mt进行形态学平滑处理，得到平滑后图像I_ms；

(5)对平滑后图像I_ms进行模糊聚类细分割：

5a)设定迭代停止条件ε，其中0＜ε＜1，聚类数目的初始值c＝2，自适应函数初始值L(1)＝0，迭代次数的初始值k＝0，模糊加权指数m＝2；

5b)随机初始化聚类中心矩阵： $V^{\left(0\right)}=\left[\begin{array}{ccccc}{v}_1^{\left(0\right)}& ...& v_i^{\left(0\right)}& ...& v_c^{\left(0\right)}\end{array}\right],$ 其中i＝1,2，…，c；

5c)求取平滑后图像I_ms中第j个像素点与第i个聚类中心的第k次迭代的距离并计算第k次迭代的隶属度矩阵U^(k)和第k次迭代后的聚类中心矩阵V^(k+1)；

5d)对V^(k)和V^(k+1)进行比较，如果||V^(k+1)-V^(k)||＜ε，则迭代停止，执行步骤5e)，否则令k＝k+1，返回步骤5c)；

5e)计算聚类数目为c时的自适应函数L(c)，当2＜c＜q时，若存在满足L(c-1)＞L(c-2)且L(c-1)＞L(c)的点，则取聚类数目c＝c-1时的聚类结果作为图像细分割的结果，否则令c＝c+1，返回步骤5b)；当c＞q时，则取聚类数目时的聚类结果作为图像细分割的结果，得到最终的分割结果图像。

2.根据权利要求1所述的基于多阈值和自适应模糊聚类的红外图像分割方法,其特征在于，步骤3c)所述的最小峰面积MA和最小峰宽MW，分别按下式计算：

MA＝ψ*Co，

MW＝μ*H，

其中，ψ表示峰面积控制因子，μ表示峰宽控制因子。

3.根据权利要求1所述的基于多阈值和自适应模糊聚类的红外图像分割方法,其特征在于，所述步骤5c)的第j个像素点与第i个聚类中心的第k次迭代的距离通过下式计算：

$d\left(g_j{v}_i^{\left(k\right)}\right)=\mathrm{In}(1+\mathrm{\β}||g_j-v_i^{\left(k\right)}|\left|\right),$

其中，β表示调节曲线曲率的自由度参数，||·||表示计算欧氏距离。

4.根据权利要求1所述的基于多阈值和自适应模糊聚类的红外图像分割方法,其特征在于，所述步骤5c)中第k次迭代的隶属度矩阵U^(k)和第k次迭代后的聚类中心矩阵V^(k+1)，分别表示为：

$V^{(k+1)}=\left[\begin{array}{ccccc}{v}_1^{(k+1)}& ...& v_i^{(k+1)}& ...& v_c^{(k+1)}\end{array}\right],$

其中，U^(k)的第i行第j列的元素和V^(k+1)第i个元素分别按下式计算：

$u_{\mathrm{ij}}^{\left(k\right)}=1/\underset{j^{\′}=1}{\overset{c}{\mathrm{\Σ}}}{\left(\frac{d(g_j,v_i^{\left(k\right)})}{d(g_j,v_{j^{\′}}^{\left(k\right)})}\right)}^{2/(m-1)},$

$v_i^{(k+1)}=\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m{g}_j/\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m,$

其中，j'表示平滑后图像I_ms聚类数目中第j'类。

5.根据权利要求1所述的基于多阈值和自适应模糊聚类的红外图像分割方法，其特征在于，步骤5e)所述的自适应函数L(c)，通过下式计算：

$L\left(c\right)=\frac{\underset{i=1}{\overset{c}{\mathrm{\Σ}}}\left(\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m{d}^2(v_i^{\left(k\right)},X)\right)/(c-1)}{\underset{i=1}{\overset{c}{\mathrm{\Σ}}}\underset{j=1}{\overset{\mathrm{Co}}{\mathrm{\Σ}}}{\left(u_{\mathrm{ij}}^{\left(k\right)}\right)}^m{d}^2\left(g_j{v}_i^{\left(k\right)}\right)/(\mathrm{Co}-c)}$

其中，表示平滑后图像I_ms中所有像素点的中心向量。