图像分割就是利用相关方法将图像中的目标和背景区分开来, 以便进一步的研究和分析. 在图像分割的研究和应用中, 水平集方法由于其鲁棒性和效率高等特点受到越来越多的研究 者重视. 水平集方法的基本思想是用连续曲线来表达目标边缘, 并定义一个能量泛函使得其 自变量包括边缘曲线, 因此分割过程就转变为求解能量泛函的最小值过程, 一般可以通过求 解函数对应的欧拉方程来实现, 能量达到最小时的曲线位置就是目标的轮廓所在. 目前应用 较广的水平集算法是 Chan 和 Vese 在 2000 年提出的无边缘活动轮廓模型（C-V 模型） ^[1，2] ,但该模型仅仅利用了图像全局的区域信息, 没有利用图像的边缘信息. 文献 [4]及文献 [5] 在C-V 模型的基础上添加了图像的边缘信息, 提出了一种改进的模型, 该模型被证明对提高曲线的演化速度和图像的边缘检测有一定的作用. 文献 [7] 在 C-V 模型中通过定义惩罚项, 避免了在演化过程中水平集函数的重新初始化过程, 提高了算法效率. 本文针对 SAR 图像一般多是灰度不均匀, 且多有凹陷区域等特点, 提出了一种改进的水平集模型.

水平集方法的基本思想是将图像空间中的二维曲线嵌入到三维曲面中作为其零水平集. 当曲面演变时, 嵌入到三维曲面中的零水平集曲线也随之改变, 最终只要确定零水平集即可 确定移动曲面演化的结果. C-V 模型 [1] 是由 Chan 和 Vese 在 2001 年提出的一种区域性的水平集方法, 它利用曲线的内外灰度均值来促进水平集的演化. C-V 模型的能量泛函如下所示:

$E({{c}_{1}},{{c}_{2}},C)=\alpha \int\limits_{inside(C)}{\left| uo(x,y)-{{c}_{1}}\left| ^{2}dxdy+\beta \right. \right.}\int\limits_{outside(C)}{\left| uo(x,y)-{{c}_{2}}\left| ^{2}dxdy+\mu .Length(C), \right. \right.}$

其中 C 表示光滑闭合的轮廓曲线, u₀ 表示源图像, c₁ 和 c₂ 分别是演化曲线 C 内部和外部的图像灰度均值, 一般情形下取$\alpha =\beta =1$,$\mu $为正常数. 引入水平集函数$\varphi (x,y)$来代替演化曲线 C, 并且规定点 (x,y)在 C 的内部时,$\varphi (x,y)>0$；点(x,y)在C的外部时,$\varphi (x,y)﹤0$；点(x,y)在C上时,$\varphi (x,y)=0$.从而上述能量泛函可写为如下形式:

$E({{c}_{1}},{{c}_{2}},C)=\mu \int\limits_{\Omega }{{{\delta }_{\varepsilon }}(\varphi (x,y))\left| \nabla \varphi (x,y)\left| dxdy+\alpha \right. \right.}\int\limits_{\Omega }{\left| uo(x,y)-{{c}_{1}}\left| ^{2}{{H}_{\varepsilon }}(\varphi (x,y))dxdy \right. \right.} \\ +\beta \int\limits_{\Omega }{\left| {{u}_{0}}(x,y)-{{c}_{2}}\left| ^{2} \right.(1-H\varepsilon (\varphi (x,y))dxdy, \right.}$

其中${{H}_{\varepsilon }}(z)$和${{\delta }_{\varepsilon }}(z)$分别是海氏函数 $H(z)=\left\{ _{0,z\le 0}^{1,z\ge 0,} \right.$和狄拉克函数$\delta (z)=\frac{d}{dz}H(z)$的正则化形式,$\Omega $表示图像域. 上述能量泛函最小化问题可以通过求解能量泛函对应的欧拉方程来实现, 从而得到如下的水平集演化方程

$\frac{\partial \varphi }{\partial t}={{\delta }_{\varepsilon }}(\varphi )\left[ \mu .div(\frac{\nabla \varphi }{\left| \nabla \varphi \right|})-\alpha {{({{u}_{0}}-{{c}_{1}})}^{2}}+\beta {{({{u}_{0}}-{{c}_{2}})}^{2}} \right],$

其中$\varphi (0,x.y)={{\varphi }_{0}}(x,y)\in \Omega $. 上式中的灰度均值 c₁ 和 c₂ 可分别在每次迭代中采用如下的方式进行更新

${{c}_{1}}(\varphi )=\frac{\int\limits_{\Omega }{{{u}_{0}}(x,y){{H}_{\varepsilon }}(\varphi (x,y))dxdy}}{\int\limits_{\Omega }{{{H}_{\varepsilon }}(\varphi (x,y))dxdy}},{{c}_{2}}(\varphi )=\frac{\int\limits_{\Omega }{{{u}_{0}}(x,y)(1-{{H}_{\varepsilon }}(\varphi (x,y)))dxdy}}{\int\limits_{\Omega }{(1-{{H}_{\varepsilon }}(\varphi (x,y)))dxdy}}.$

(2.1)

C-V模型对噪声不敏感, 并且该模型对边缘模糊的图像以及含有内部轮廓的图像很有效^[2]. 但C-V模型也存在一些不足, 首先C-V模型没有利用图像的边缘信息, 从而在分割一些灰度不均匀的图像时效果不好^[3]; 其次C-V模型对初始轮廓的设置有一定的依赖性, 不同的初始轮廓经常会得到不同的分解结果 ^[4]最后在迭代过程中需要不断地重新初始化水平集函数, 计算量很大. 基于此, 本文在C-V模型的基础上, 提出了一种加入边缘信息、并且无需重新初始化的改进模型.

1) 由于水平集函数$\varphi $通常是取由初始轮廓线生成的符号距离函数, 而水平集函数在迭代时往往会发生退化, 使其不再保持为符号距离函数, 所以在进行若干次迭代之后, 必须对水平集函数$\varphi $重新初始化, 使其恢复为对当前零水平集的符号距离函数. 但重新初始化的计算量很大, 从而导致分割效率较低. 为了避免重新初始化, 很多文献^[5]提出了一种改进方案, 即在关于水平集函数$\varphi $的能量泛函中增加一项内部约束能量 $E(\varphi )=\displaystyle\int\limits_\Omega {\dfrac{1}{2}} (\left| {\nabla \varphi } \right|-1)dxdy$ 它的梯度下降流为

$\frac{\partial \varphi }{\partial t}=\Delta \varphi -{\rm div}(\frac{\nabla \varphi }{\left| {\nabla \varphi } \right|}).$

由于水平集函数$\varphi $为符号距离函数必须满足$\left| {\nabla \varphi } \right|=1$, 而上面的方程是一个非线性的热传导方程, 当$\left| {\nabla \varphi } \right| > 1$时, 热量正向传递, 使$\left| {\nabla \varphi } \right|$减小; 当$\left| {\nabla \varphi } \right|<1$时, 热量负向传递, 使$\left| {\nabla \varphi } \right|$增大. 这样对于偏离$\left| {\nabla \varphi } \right|=1$的水平集函数$\varphi $, 在后续的迭代过程中都将被纠正过来, 使其保持为符号距离函数.

2) 由于传统的C-V模型只考虑了图像的全局区域信息, 而没有考虑图像的边缘梯度信息, 从而水平集的演化速度和边缘信息检测会受到影响, 从而导致在演化时往往容易陷入到局部极小值, 于是考虑引入边缘检测函数^[6]

$g(\left| {\nabla G_\sigma \ast u} \right|)=\dfrac{a}{1+\left| {\nabla G_\sigma \ast u} \right|^2}$

来替换上式中的$\delta _\varepsilon (\varphi )$, 其中$G_\sigma $为均值为0, 方差为$\sigma $的高斯函数, $\nabla $为梯度算子, $\ast $为二维图像卷积算子, 通过对原始图像做卷积运算, 可以降低原始图像的噪声, 从而有效抑制了噪声对图像边缘提取的影响; 而参数$a$是用来控制曲线的演化速度, 在图像灰度变化较平缓的区域, $\left| {\nabla G\ast u} \right|$的值趋于0, 于是取较大的$a$值, 从而加快曲线的演化速度; 在靠近图像的边缘处, 图像灰度变化较大, $\left| {\nabla G\ast u} \right|$的值也很大, 这时相应地取较小的$a$值, 来保证演化的稳定性, 引导演化曲线向目标图像边缘靠近, 避免分割陷入局部极小值. 结合1), 2)得到改进后的C-V水平集演化方程如下:

$\begin{equation} \label{eq2} \frac{\partial \varphi }{\partial t}=g(\left| {\nabla G_\sigma \ast u} \right|)[\mu \cdot {\rm div}(\frac{\nabla \varphi }{\left| {\nabla \varphi } \right|})-\alpha (u_0 -c_1 )^2+\beta (u_0 -c_2 )^2]+\Delta \varphi -{\rm div}(\frac{\nabla \varphi }{\left| {\nabla \varphi } \right|}). \end{equation}$

(2.2)

上述的水平集演化方程可以使用有限差分法 [7] 来实现数值求解, 即将所有空间域 上的偏导数用中心差分来实现, 所有时域上的偏导数用向前差分来实现, 具体表示如下: $\frac{\partial \varphi }{\partial }=\frac{\varphi _{i,j}^{n+1}-\varphi _{i,j}^{n}}{\Delta t}$其中$\Delta t$是时间步长. 演化方程中的曲率$div(\frac{\nabla \varphi }{\left| \nabla \varphi \right|})$可以使用二阶中心差分来进行求解:

其中${{\varphi }_{x}},{{\varphi }_{y}},{{\varphi }_{xx}},{{\varphi }_{xy}},{{\varphi }_{yy}}$可用如下所给公式进行计算

$\begin{align} && {{\varphi }_{x}}=\frac{1}{2h}({{\varphi }_{x}}+1,j-{{\varphi }_{i}}-1,j),{{\varphi }_{y}}=\frac{1}{2h}(\varphi i,j+1-\varphi i,j-1),{{\varphi }_{xx}}=\frac{1}{{{h}^{2}}}(\varphi j+1,j+\varphi i-1,j-2\varphi i,j), \\ && {{\varphi }_{yy}}=\frac{1}{{{h}^{2}}}(\varphi i,j+1+\varphi i,j-1-2\varphi i,j),{{\varphi }_{xy}}=\frac{1}{{{h}^{2}}}(\varphi i+1,j+1-\varphi i-1,j+1-\varphi i+1,j-1+\varphi i-1,j-1 , \end{align}$

其中 h 是离散网格的间隔, 演化方程中的 c₁; c₂ 可按上式 (2.1) 求解. 水平集演化终止准则:当曲线长度变化的绝对值在给定的迭代次数 (n) 里一直小于某个事先给定的阈值$\eta $, 水平集演化终止^[8]. 在本文的实验里, 固定 n=8, $\eta=4$.

该算法流程如下:

1) 输入原始图像$u_0 $.

2) 在$u_0 $中设置初始曲线C, 并设置如下相应的参数值: 时间步长$\Delta t$, 网格间隔$h$, 以及$\alpha ,\beta ,\mu ,a$的值.

3) 按上式(2.2)演化水平集函数$\varphi $.

4) 从水平集函数$\varphi $中提取零水平集即演化曲线.

5) 判断水平集演化是否终止, 如果是, 算法结束, 否则转第三步.

下面的实验均使用Matlab7进行编程仿真, 实验所用的计算机配置为: Intel Core4处理器, 10G内存, Windows7操作系统. 模型的参数设置为: $\alpha =\beta =1$, $\mu =0.01\ast 255^2$, 时间步长$\Delta t=0.1$, 网格间隔$h=1$, 正则化参数$\varepsilon =1$, $\sigma =2.0$, 刚开始演化时参数$a$取2, 分割到目标轮廓附近时$a$取1.

实验１为传统的C-V模型与改进后的模型对一副雷达河流图像分割实验对比. 图 2为传统的C-V模型实验结果, 图 3为改进后模型的实验结果. 其中传统的C-V模型需要迭代25次, 用时2.634s完成分割; 而改进后的模型只需迭代16次, 用时1.221s即可收敛, 并且可以检测出微小的边缘, 其分割效率和分割精确度都比传统的C-V模型有很大的提高.

实验2是对一幅SAR卫星云图像的分割. 图 5、图6和图 7为传统的C-V模型分别在$n=20$、$n=150$以及$n=240$时的分割结果; 图 8、图 9和图 10为改进后的模型分别在$n=20$、$n=150$以及$n=180$时的分割结果 . 通过对比分别用传统的C-V模型和改进后的模型在$n=20$ 和$n=150$的分割结果可以发现改进后的模型分割效率更高, 最终传统的C-V模型迭代了240次, 用时8.524s收敛; 而改进后的模型只迭代了180次, 用时5.456s就得到了分割结果, 并且由于改进后的模型引入了边缘检测算子, 降低了图像噪声对分割过程的影响, 从而得出了更为准确的分割结果.

本文依据传统的 C-V 模型在图像分割时只考虑图像的全局信息的基础上, 提出了一种加 入图像梯度边缘信息且无需重新初始化的改进模型, 并通过实验验证, 该模型比传统的 C-V 模型有稳定性好、收敛速度快、分割精确度高等优点. 但是本文引入的边缘检测函数仅与图 像的梯度信息有关, 而与图像的先验信息等无关, 从而影响了分割效率, 这是下一步研究的方向.