CN102044077B - 多分辨率多区域变分水平集图像分割方法 - Google Patents

Abstract

一种图像处理技术领域多分辨率多区域变分水平集图像分割方法。包括：设置分辨率级数以及分割区域的数目，将原始图像按照空间分辨率在每一维进行连续下采样生成分辨率为2^L的图像；利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集演化曲线方程；以2ⁱ为系数得到的该演化曲线作为下一分辨率构建初始化演化曲线，然后采用多分辨率水平集方法，进行曲线演化得到当前分辨率下N-1个零水平集演化曲线方程；最后，演化过程不断重复，直至达到原始分辨率图像，得到分割结果。本发明避免分割区域的重叠和漏分，降低了噪声的干扰、减小了搜索的空间。

Description

多分辨率多区域变分水平集图像分割方法

技术领域

本发明涉及的是一种图像处理技术领域图像分割的方法，具体是一种多分辨率多区域变分水平集(Multiresolution and Multiregion Level Set，MR-MRSL)图像分割方法。

技术背景

图像分割是图像特征提取和分类的重要环节，图像分割的目的就是将图像中的灰度同质区域分离开，并通过各个同质区域的边界来表达。近年来，水平集分割方法凭借其自由拓扑性及多信息共融性，被广泛应用于计算机视觉，例如图像分割、运动跟踪、三维重建。基于C-V模型的水平集图像分割方法具有几个优点：1、所涉及的图像函数的定义域是整个图像，具有全局特性。因此，该模型图像分割方法具有全局优化的特点，仅用一条初始闭合轮廓线就可以进行图像的分割；2、初始曲线的位置无关性，轮廓线经过数次循环可以正确地分割出目标和背景；3、与图像中的边缘信息的无关性，即使图像中的边缘呈模糊或离散状，仍然可以获得理想的分割效果。但是，对具有噪声的遥感图像、医学图像以及自然图像而言，该方法具有几个缺陷：1、该C-V模型对同质区域的划分时仅考虑灰度，对多通道图像处理无能为力；2、该模型每次更新后，需要重新初始化符号距离函数，而高分辨率数据量丰富的图像计算量非常大；3、该模型对于带较厚空洞和三合点的目标，不能稳定地检测内部区域。

经对现有技术文献的检索发现，王爽等提出基于水平集和分水蛉相结合的图像轮廓检测方法(专利号：CN101567084)以及基于邻域概率密度函数特征的水平集图像分割方法(专利号：CN101571951)；曹宗杰提出了基于概率论模型的水平集方法(专利号：CN101221239)；杨杰、周娟等提出的基于空间矩的水平集图像分割方法(专利号：CN101493942A)。这些方法均是基于两区域水平集方法进行改进来提高分割结果，但是针对多区域的图像分割非常困难，主要原因是多个闭合曲线分割会导致分割的图像区域产生重叠现象，而且对于数据量大的在循环次数多时非常耗时。

发明内容

本发明的目的在于克服现有技术中的不足，提供一种多分辨率多区域变分水平集图像分割方法。本发明基于变分水平集C-V模型，用N-1个水平集函数Φ＝{φ₁，φ₂，…，φ_n-1}将图像分割成N(N＞1)个区域，每个水平集函数表达一个区域，从而避兔分割区域的重叠和漏分，通过多分辨率技术执行分割图像的演化曲线来解决初始化水平集能量函数陷入局部能量最小值，降低了噪声的干扰、减小了搜索的空间。

本发明是通过以下的技术方案实现的：

本发明包括如下步骤：

首先，设置分辨率级数以及分割区域的数目，将原始图像按照空间分辨率在每一维进行连续下采样生成分辨率为2^L的图像，其中2^L为低分辨率级数；

然后，为当前分辨率图像建立能量模型，利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集演化曲线方程；

进一步以2ⁱ(i＝2，…L)为系数采用双线性插值方法上采样演化曲线，得到的该演化曲线作为下一分辨率构建初始化演化曲线，然后构建该分辨率图像总能量模型，利用变分水平集最小化能量模型，采用多分辨率水平集方法，进行曲线演化得到当前分辨率下N-1个零水平集演化曲线方程；

最后，演化过程不断重复，直至达到原始分辨率图像，得到分割结果。

所述的将原始图像按照空间分辨率在每一维进行连续下采样生成分辨率为2^L的图像：先设置分辨率的级数L，生成分辨率级数为L低分辨率的图像，作为初始分辨率的图像。

图像的级数不能太高(本发明采用的级数L＝3，4)，否则会导致图像信息的丢失，影响分割结果。

所述的变分水平集，其多区域分割方法，包括分割区域表示和图像能量模型，在分割区域表示中，Chan和Vese提出了多相水平集分割方法，用N个水平集函数将图像分割成为2^N个区域，会产生交叉区域。

本发明采用N-1个水平集函数函数Φ＝{φ₁，φ₂，…，φ_n-1}将图像分割成N(N＞1)个区域，每个水平集函数表达一个区域，从而避免造成重叠和漏分。

所述的变分水平集，其函数曲线围成的区域可表示为：

$R_1=R_{{\mathrm{\γ}}_1}$

$\begin{array}{c}{R}_2=R_{{\mathrm{\γ}}_1}^c\∩R_{{\mathrm{\γ}}_2}\\ ......\end{array}$

$R_k=R_{{\mathrm{\γ}}_1}^c\∩R_{{\mathrm{\γ}}_2}^c\∩\·\·\·\∩R_{{\mathrm{\γ}}_k}$

$R_N=R_{{\mathrm{\γ}}_1}^c\∩R_{{\mathrm{\γ}}_2}^c\∩\·\·\·\∩R_{{\mathrm{\γ}}_{N-1}}$

在图像能量模型中，能量模型中图像的多区域信息、边界信息以及边缘演化模型，其能量泛函定义为：

$E(\mathrm{\φ},{\left\{c_i\right\}}_{i=1}^{N-1})=E_R(\mathrm{\φ},{\left\{c_i\right\}}_{i=1}^{N-1})+E_E\left(\mathrm{\φ}\right)+E_C\left(\mathrm{\φ}\right)$

图像的多区域能量定义为：假设给定的图像

为图像域，d为图像的维数，x为图像像素对应的向量。能量函数如下：

$E^R({\left\{{\stackrel{\→}{\mathrm{\γ}}}_i\right\}}_{i=1}^{N-1},{\left\{c_i\right\}}_{i=1}^N)={\mathrm{\λ}}_1\underset{R_1}{\∫}{\mathrm{\χ}}_{R_1}{|I\left(x\right)-c_1|}^2\mathrm{dx}+{\mathrm{\λ}}_2\underset{R_2}{\∫}{\mathrm{\χ}}_{R_2}{|I\left(x\right)-c_2|}^2\mathrm{dx}+\·\·\·$

$+{\mathrm{\λ}}_k\underset{R_k}{\∫}{\mathrm{\χ}}_{R_k}{|I\left(x\right)-c_k|}^2\mathrm{dx}+\·\·\·+{\mathrm{\λ}}_N\underset{R_N}{\∫}{\mathrm{\χ}}_{R_N}{|I\left(x\right)-c_N|}^2\mathrm{dx}$

$=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\underset{R_i}{\∫}{{\mathrm{\χ}}_{R_i}|I\left(x\right)-c_i|}^2\mathrm{dx}$

其中λ_i＞0(i＝1，…N)是各个能量项权重系数，c_i(i＝1，…N)表示曲线像素所在区域的像素均值。

(i＝1，…N)是第i区域的特征函数，并且满足假设H是Heaviside函数，

定义如下：

$\left\{\begin{array}{c}\begin{array}{cc}{\mathrm{\χ}}_{R_1}={\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1}}& =[1-H\left({\stackrel{\→}{\mathrm{\γ}}}_0\right)]H\left({\stackrel{\→}{\mathrm{\γ}}}_1\right)H\left({\stackrel{\→}{\mathrm{\γ}}}_0\right)\≡0\end{array}\\ \begin{array}{cc}{\mathrm{\χ}}_{R_2}={\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2}}& =[1-H\left({\stackrel{\→}{\mathrm{\γ}}}_0\right)][1-H\left({\stackrel{\→}{\mathrm{\γ}}}_1\right)]H\left({\stackrel{\→}{\mathrm{\γ}}}_2\right)\end{array}\\ .\\ .\\ .\\ {\mathrm{\χ}}_{R_k}={\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2^c}}\·\·\·{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{k-1}^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_k}}=\underset{i=0}{\overset{k-1}{\mathrm{\Π}}}[1-H\left({\stackrel{\→}{\mathrm{\γ}}}_i\right)]H\left({\stackrel{\→}{\mathrm{\γ}}}_k\right)\\ .\\ .\\ .\\ {\mathrm{\χ}}_{R_{N-1}}={\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2^c}}\·\·\·{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-2}^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-1}}}=\underset{i=0}{\overset{N-2}{\mathrm{\Π}}}[1-H\left({\stackrel{\→}{\mathrm{\γ}}}_i\right)]H\left({\stackrel{\→}{\mathrm{\γ}}}_{N-1}\right)\\ {\mathrm{\χ}}_{R_N}={\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2^c}}\·\·\·{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-2}^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-1}^c}}=\underset{i=0}{\overset{N-1}{\mathrm{\Π}}}[1-H\left({\stackrel{\→}{\mathrm{\γ}}}_i\right)]\end{array}\right.$

其中， $H_{\mathrm{\ϵ}}\left(x\right)=\frac{1}{2}[1+\frac{2}{\mathrm{\π}}\arctan \left(\frac{x}{\mathrm{\ϵ}}\right)]$

从上式可以看出，

表示所有独立的区域的交集为整个图像，即

${\mathrm{\χ}}_{R_1}+{\mathrm{\χ}}_{R_2}+\·\·\·{\mathrm{\χ}}_{R_k}+\·\·\·{\mathrm{\χ}}_{R_{N-1}}+{\mathrm{\χ}}_{R_N}={\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1}}+{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2}}+\·\·\·$

$+{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2^c}}\·\·\·{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{k-1}^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_k}}+\·\·\·+{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2^c}}\·\·\·{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-2}^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-1}}}+{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_1^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_2^c}}\·\·\·{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-2}^c}}{\mathrm{\χ}}_{R_{{\stackrel{\→}{\mathrm{\γ}}}_{N-1}^c}}=1.$

图像的边界信息能量中，由于水平集函数作为符号距离函数仅在水平集函数初始化时有定义，而在水平集函数的演化方程中未体现该约束，为了确保水平集近似符号距离函数(SignedDistance Function，SDF)。

本发明采用了Li Chunming(李纯明)等人在CVPR 2005(计算机视觉和模式识别国际会议2005)上发表的Level set evolution without reinitialization：A new variationalformulation文章中定义的边缘模型：

E_C(φ)＝vP(φ)

其中， $P\left(\mathrm{\φ}\right)=\frac{1}{2}\underset{\mathrm{\Ω}}{\∫}{\left(\right|\▿\mathrm{\φ}|-1)}^2\mathrm{dxdy}$

图像的边缘模型能量中，应为所有水平集曲线能量之和，表示为：

$E_E\left(\mathrm{\φ}\right)=\mathrm{\μL}\left(\mathrm{\φ}\right)=\mathrm{\μ}\underset{j=2}{\overset{N-1}{\mathrm{\Σ}}}\underset{R_j}{\∫}|\▿H|\mathrm{dxdy}=\mathrm{\μ}\underset{j=2}{\overset{N-1}{\mathrm{\Σ}}}\underset{R_j}{\∫}\mathrm{\δ}\left(\mathrm{\φ}\right)|\▿\mathrm{\φ}|\mathrm{dxdy}.$

所述的采用多区域变分水平集最小化能量模型，进行图像分割，具体为：将曲线表示成水平集的形式，利用水平集的曲线演化方法进行演化分割，可以得到变分泛函图像总能量的水平集函数的演化方程为：

$\frac{d{\stackrel{\→}{\mathrm{\γ}}}_1}{\mathrm{dt}}=-[{\mathrm{\λ}}_1{|I\left(x\right)-c_1|}^2-{\mathrm{\Φ}}_1\left({\stackrel{\→}{\mathrm{\γ}}}_1\right)+\mathrm{\μ}{k}_1]{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\stackrel{\→}{\mathrm{\γ}}}_1\right)+v(\▿{\stackrel{\→}{\mathrm{\γ}}}_1^2-k_1)$

$\begin{array}{c}\frac{d{\stackrel{\→}{\mathrm{\γ}}}_2}{\mathrm{dt}}=-[{\mathrm{\λ}}_2{|I\left(x\right)-c_2|}^2-{\mathrm{\Φ}}_2\left({\stackrel{\→}{\mathrm{\γ}}}_2\right)+\mathrm{\μ}{k}_2]{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\stackrel{\→}{\mathrm{\γ}}}_2\right)+v(\▿{\stackrel{\→}{\mathrm{\γ}}}_2^2-k_2)\\ .\\ .\\ .\end{array}$

$\begin{array}{c}\frac{d{\stackrel{\→}{\mathrm{\γ}}}_k}{\mathrm{dt}}=-[{\mathrm{\λ}}_k{|I\left(x\right)-c_k|}^2-{\mathrm{\Φ}}_k\left({\stackrel{\→}{\mathrm{\γ}}}_k\right)+\mathrm{\μ}{k}_k]{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\stackrel{\→}{\mathrm{\γ}}}_k\right)+v(\▿{\stackrel{\→}{\mathrm{\γ}}}_k^2-k_k)\\ .\\ .\\ .\end{array}$

$\frac{d{\stackrel{\→}{\mathrm{\γ}}}_{N-1}}{\mathrm{dt}}=-[{\mathrm{\λ}}_{N-1}{|I\left(x\right)-c_{N-1}|}^2-{\mathrm{\Φ}}_{N-1}\left({\stackrel{\→}{\mathrm{\γ}}}_{N-1}\right)+\mathrm{\μ}{k}_{N-1}]{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\stackrel{\→}{\mathrm{\γ}}}_{N-1}\right)+v(\▿{\stackrel{\→}{\mathrm{\γ}}}_{N-1}^2-k_{N-1})$

其中，

表示第i条水平集演化曲线，δ_ε(x)为平滑函数H_ε(x)的导数，

${\mathrm{\δ}}_{\mathrm{\ϵ}}\left(x\right)=H_{\mathrm{\ϵ}}^{\′}\left(x\right)=\frac{1}{\mathrm{\π}}\frac{\mathrm{\ϵ}}{{\mathrm{\ϵ}}^2+x^2},$

${\mathrm{\&Phi;}}_i\left(x\right)={\mathrm{\&lambda;}}_{i+1}{|I\left(x\right)-c_{i+1}|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)>0}+{\mathrm{\&lambda;}}_{i+2}{|I\left(x\right)-c_{i+2}|}^2{\mathrm{\&chi;}}_{u_{j+1}(x,t)<0}{\mathrm{\&chi;}}_{u_{i+2}(x,t)>0}+\&CenterDot;\&CenterDot;\&CenterDot;$

$+{\mathrm{\&lambda;}}_{N-1}{|I\left(x\right)-c_{N-1}|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)>0},$

$+{\mathrm{\&lambda;}}_N{|I\left(x\right)-c_N|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}$

为水平集函数的曲率， $k=\stackrel{\&RightArrow;}{\&dtri;}\&CenterDot;\frac{\stackrel{\&RightArrow;}{\&dtri;}{u}_i}{\left|\right|\stackrel{\&RightArrow;}{\&dtri;}{u}_i\left|\right|}=\frac{u_{\mathrm{xx}}{u}_y^2-2u_x{u}_y{u}_{\mathrm{xy}}+u_{\mathrm{yy}}{u}_x^2}{{(u_x^2+u_y^2)}^{3/2}},$

$\left\{\begin{array}{cc}{\mathrm{\&chi;}}_{u_i(x,t)>0}=H\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_i\right)& u_i(x,t)>0\\ {\mathrm{\&chi;}}_{u_i(x,t)\&le;0}=1-H\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_i\right)& u_i(x,t)\&le;0\end{array}\right.,$

为各局部区域的像素均值，

$c_i=\left\{\begin{array}{cc}\frac{\underset{R_1}{\&Integral;}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)>0}\mathrm{dx}}{\underset{R_1}{\&Integral;}{\mathrm{\&chi;}}_{u_1(x,t)>0}\mathrm{dx}}& i=1\\ \frac{\underset{R_i}{\&Integral;}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{i-1}(x,t)<0}{\mathrm{\&chi;}}_{u_i(x,t)>0}\mathrm{dx}}{\underset{R_i}{\&Integral;}{\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{i-1}(x,t)<0}{\mathrm{\&chi;}}_{u_i(x,t)>0}\mathrm{dx}}& i\&Element;[2,N-1]\\ \frac{\underset{R_N}{\&Integral;}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}\mathrm{dx}}{\underset{R_N}{\&Integral;}{\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}\mathrm{dx}}& i=N\end{array}\right..$

所述的多分辨率水平集方法，包括两个过程：首先，对该图像的每一维进行下采样生成分辨率为2ⁱ(i＝2，…，L)的低分辨率图像。这样的好处在于，经过空间下采样得到不同分辨率的图像在低分辨率图像进行曲线演化过程中降低了噪声的干扰、减小搜索的空间以及局部最大值数。其次，相对应地，以2ⁱ(i＝2，…L)为系数上采样该分辨率下的演化曲线，然后构建该分辨率图像总能量模型，利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集分割图像。演化过程不断重复直至达到原始分辨率图像，得到分割结果。

在图像处理过程通常采用原始分辨率的图像，这会导致给定的初始化水平集能量函数陷入局部能量最小值，从而导致分割的效果。特别对于还有大量噪声的图像而言，这种情况发生的概率非常大。本发明采用了多分辨率技术来解决上述问题。

本发明基于变分水平集C-V模型，采用N-1个水平集函数将图像分割成N(N＞1)个区域，每个水平集函数表达一个区域，从而避免分割区域的重叠和漏分。为了防止给定的初始化水平集能量函数陷入局部能量最小值，降低噪声的干扰，减小搜索的空间，提出了多分辨率技术来执行分割图像的演化曲线来解决这个问题。本发明弥补了传统两区域的水平集方法容易产生冗余轮廓缺陷，具有重要的实用价值。

附图说明

图1表示多区域图像表示，采用4个水平集函数划分5个独立区域。

图2为C-V模型和多分辨率多区域水平集方法的分割效果：

其中：(a)初始化图像，两个圆表示初始化的轮廓线；(b)采用c-v模型分割的结果：(c)和(d)为本发明方法分割结果的两目标区域图像。

图3为单分辨率和多分辨率多区域水平集方法的分割结果：

其中：(a)中三幅图像分别为初始图像、本发明方法分割的轮廓图以及分割结果图；(b)中两幅图像分别表示单分辨率下多区域水平集方法分割的轮廓图和分割结果；(c)中四幅图像分别表示分辨率的级数为L＝4，3，2，1的分割结果；(d)中三幅图像表示各目标区域分割的结果。

具体实施方式

下面结合附图对本发明的实施例作详细说明：本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

实施例

本实施例基于变分水平集C-V模型，获得N个区域的区域信息、边界信息以及边缘演化模型为图像建立总的能量模型，利用变分水平集方法最小化能量模型，用N-1个水平集函数将图像分割成N(N＞1)个区域，每个水平集函数表达一个区域，得到各区域分割的结果。为了防止给定的初始化水平集能量函数陷入局部能量最小值，降低噪声的干扰，减小搜索的空间，采用多分辨率技术得到比单一分辨率多区域水平集方法更好的分割结果。

本实施例包括如下步骤：

1)设置分辨率级数设置为L以及分割区域数N(N≥2)，演化曲线数初始值m＝N-1，本实施例设置L＝3或4。将分割图像按空间分辨率在每一维进行连续下采样生成分辨率为2^L的图像。

2)利用图像的全局信息，为当前分辨率图像建立总的总能量模型：

$E(\mathrm{\&phi;},{\left\{c_i\right\}}_{i=1}^{N-1})=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&lambda;}}_i\underset{R_i}{\&Integral;}{\mathrm{\&chi;}}_{R_i}{|I\left(x\right)-c_i|}^2\mathrm{dx}+\mathrm{\&mu;}\underset{j=2}{\overset{N-1}{\mathrm{\&Sigma;}}}\underset{R_j}{\&Integral;}\mathrm{\&delta;}\left(\mathrm{\&phi;}\right)|\&dtri;\mathrm{\&phi;}|\mathrm{dxdy}+\mathrm{vP}\left(\mathrm{\&phi;}\right)$

其中， $P\left(\mathrm{\&phi;}\right)=\frac{1}{2}\underset{\mathrm{\&Omega;}}{\&Integral;}{\left(\right|\&dtri;\mathrm{\&phi;}|-1)}^2\mathrm{dxdy}$

3)利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集演化曲线方程。

演化方程为：

$\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1}{\mathrm{dt}}=-[{\mathrm{\&lambda;}}_1{|I\left(x\right)-c_1|}^2-{\mathrm{\&Phi;}}_1\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1\right)+\mathrm{\&mu;}{k}_1]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1^2-k_1)$

$\begin{array}{c}\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2}{\mathrm{dt}}=-[{\mathrm{\&lambda;}}_2{|I\left(x\right)-c_2|}^2-{\mathrm{\&Phi;}}_2\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2\right)+\mathrm{\&mu;}{k}_2]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2^2-k_2)\\ .\\ .\\ .\end{array}$

$\begin{array}{c}\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_k}{\mathrm{dt}}=-[{\mathrm{\&lambda;}}_k{|I\left(x\right)-c_k|}^2-{\mathrm{\&Phi;}}_k\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_k\right)+\mathrm{\&mu;}{k}_k]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_k\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_k^2-k_k)\\ .\\ .\\ .\end{array}$

$\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}}{\mathrm{dt}}=-[{\mathrm{\&lambda;}}_{N-1}{|I\left(x\right)-c_{N-1}|}^2-{\mathrm{\&Phi;}}_{N-1}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}\right)+\mathrm{\&mu;}{k}_{N-1}]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}^2-k_{N-1})$

${\mathrm{\&Phi;}}_i\left(x\right)={\mathrm{\&lambda;}}_{i+1}{|I\left(x\right)-c_{i+1}|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)>0}+{\mathrm{\&lambda;}}_{i+2}{|I\left(x\right)-c_{i+2}|}^2{\mathrm{\&chi;}}_{u_{j+1}(x,t)<0}{\mathrm{\&chi;}}_{u_{i+2}(x,t)>0}+\&CenterDot;\&CenterDot;\&CenterDot;$

$+{\mathrm{\&lambda;}}_{N-1}{|I\left(x\right)-c_{N-1}|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)>0}$

$+{\mathrm{\&lambda;}}_N{|I\left(x\right)-c_N|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}$

其中，

表示第i条水平集演化曲线，δ_ε(x)为平滑函数H_ε(x)的导数，

为水平集函数的曲率，

$k=\stackrel{\&RightArrow;}{\&dtri;}\&CenterDot;\frac{\stackrel{\&RightArrow;}{\&dtri;}{u}_i}{\left|\right|\stackrel{\&RightArrow;}{\&dtri;}{u}_i\left|\right|}=\frac{u_{\mathrm{xx}}{u}_y^2-2u_x{u}_y{u}_{\mathrm{xy}}+u_{\mathrm{yy}}{u}_x^2}{{(u_x^2+u_y^2)}^{3/2}},$

$\left\{\begin{array}{cc}{\mathrm{\&chi;}}_{u_i(x,t)>0}=H\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_i\right)& u_i(x,t)>0\\ {\mathrm{\&chi;}}_{u_i(x,t)\&le;0}=1-H\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_i\right)& u_i(x,t)\&le;0\end{array}\right.,$

为各局部区域的像素均值，

$c_i=\left\{\begin{array}{cc}\frac{\underset{R_1}{\&Integral;}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)>0}\mathrm{dx}}{\underset{R_1}{\&Integral;}{\mathrm{\&chi;}}_{u_1(x,t)>0}\mathrm{dx}}& i=1\\ \frac{\underset{R_i}{\&Integral;}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{i-1}(x,t)<0}{\mathrm{\&chi;}}_{u_i(x,t)>0}\mathrm{dx}}{\underset{R_i}{\&Integral;}{\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{i-1}(x,t)<0}{\mathrm{\&chi;}}_{u_i(x,t)>0}\mathrm{dx}}& i\&Element;[2,N-1]\\ \frac{\underset{R_N}{\&Integral;}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}\mathrm{dx}}{\underset{R_N}{\&Integral;}{\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}\mathrm{dx}}& i=N\end{array}\right.$

5)以2ⁱ(i＝2，…L)为系数采用双线性插值方法上采样演化曲线，得到的该演化曲线作为下一分辨率构建初始化演化曲线，然后构建该分辨率图像总能量模型，利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集分割图像。演化过程不断重复直至达到原始分辨率图像，得到分割结果。本实施例采用N-1个水平集函数函数Φ＝{φ₁，φ₂，…，φ_n-1}将图像分割成N(N＞1)个区域。如图1所示，采用4个水平集函数将图像分成5个区域，每个水平集函数表达一个区域，从而避免造成重叠和漏分。

如图2、3所示，分别为两个实例的分割结果，其中设置分辨率的级数均为4，图2比较了G-V模型和多分辨率多区域水平集方法的分割效果，图3比较了单分辨率和多分辨率多区域水平集方法的分割效果。在图2中，其中：(a)初始化图像，两个圆表示初始化的轮廓线；(b)采用c-v模型分割的结果；(c)和(d)为本发明实施分割结果的两目标区域图像。在图3中，其中：(a)中三幅图像分别为初始图像、本发明实施分割的轮廓图以及分割结果图；(b)中两幅图像分别表示单分辨率下多区域水平集方法分割的轮廓图和分割结果；(c)中四幅图像分别表示分辨率的级数为L＝4，3，2，1的分割结果；(d)中三幅图像表示各目标区域分割的结果。

从图中可以看出，上述实施例采用了独立区域水平集方法，减少了冗余的轮廓以及噪声的干扰，得到较好的分割效果。

Claims (2)

1.一种多分辨率多区域变分水平集图像分割方法，其特征在于，包括如下步骤：

第一步：设置分辨率级数以及分割区域的数目，将原始图像按照空间分辨率在每一维进行连续下采样生成分辨率为2^L的图像,其中L为分辨率级数；

第二步：为当前分辨率图像建立能量模型，利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集演化曲线方程，进行图像分割；

所述的变分水平集最小化能量模型，具体为：将曲线表示成水平集的形式，利用水平集的曲线演化方法进行演化分割，得到变分泛函图像总能量的水平集函数的演化方程为：

$\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1}{\mathrm{dt}}=-\left[{\mathrm{\&lambda;}}_1\right|I\left(x\right)-c_1{|}^2-{\mathrm{\&Phi;}}_1\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1\right)+{\mathrm{\&mu;k}}_1]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_1^2-k_1)$

$\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2}{\mathrm{dt}}=-\left[{\mathrm{\&lambda;}}_2\right|I\left(x\right)-c_2{|}^2-{\mathrm{\&Phi;}}_2\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2\right)+\mathrm{\&mu;}{k}_2]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_2^2-k_2)$

.

.

.

$\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_j}{\mathrm{dt}}=-\left[{\text{\&lambda;}}_j\right|I\left(x\right)-c_j{|}^2-{\mathrm{\&Phi;}}_j\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_j\right)+\mathrm{\&mu;}{k}_j]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_j\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_j^2-k_j)$

.

.

.

$\frac{d{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}}{\mathrm{dt}}=-\left[{\text{\&lambda;}}_{N-1}\right|I\left(x\right)-c_{N-1}{|}^2-{\mathrm{\&Phi;}}_{N-1}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}\right)+\mathrm{\&mu;}{k}_{N-1}]{\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}\right)+v(\&dtri;{\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_{N-1}^2-k_{N-1})$

其中，

表示第i条水平集演化曲线，λ_i,i＝1,2,…,N为大于0的系数，δ_ε(x)为平滑函数H_ε(x)的导数， ${\mathrm{\&delta;}}_{\mathrm{\&epsiv;}}\left(x\right)=H_{\mathrm{\&epsiv;}}^{\&prime;}\left(x\right)=\frac{1}{\mathrm{\&pi;}}\frac{\mathrm{\&epsiv;}}{{\mathrm{\&epsiv;}}^2+x^2},$

${\mathrm{\&Phi;}}_i\left(x\right)={\mathrm{\&lambda;}}_{i+1}|I\left(x\right)-c_{i+1}{|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)>0}+{\mathrm{\&lambda;}}_{i+2}|I\left(x\right)-c_{i+2}{|}^2{\mathrm{\&chi;}}_{u_{j+1}(x,t)<0}{\mathrm{\&chi;}}_{u_{i+2}(x,t)>0}+\&CenterDot;\&CenterDot;\&CenterDot;$

$+{\mathrm{\&lambda;}}_{N-1}|I\left(x\right)-c_{N-1}{|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)>0},$

$+{\mathrm{\&lambda;}}_N|I\left(x\right)-c_N{|}^2{\mathrm{\&chi;}}_{u_{i+1}(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}$

其中，为表示i区域对应各区域的示性函数，

为各局部区域的像素均值，

$c_i=\left\{\begin{array}{cc}\frac{{\&Integral;}_{R_1}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)>0}\mathrm{dx}}{{\&Integral;}_{R_1}{\mathrm{\&chi;}}_{u_1(x,t)>0}\mathrm{dx}}& i=1\\ \frac{{\&Integral;}_{R_i}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,y)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{i-1}(x,y)<0}{\mathrm{\&chi;}}_{u_i(x,t)>0}\mathrm{dx}}{{\&Integral;}_{R_i}{\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{i-1}(x,t)<0}{\mathrm{\&chi;}}_{u_i(x,t)>0}\mathrm{dx}}i\&Element;[2,N-1],& \\ \frac{{\&Integral;}_{R_N}I\left(x\right){\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}\mathrm{dx}}{{\&Integral;}_{R_N}{\mathrm{\&chi;}}_{u_1(x,t)<0}\&CenterDot;\&CenterDot;\&CenterDot;{\mathrm{\&chi;}}_{u_{N-2}(x,t)<0}{\mathrm{\&chi;}}_{u_{N-1}(x,t)<0}\mathrm{dx}}& i=N\end{array}\right.$

其中，R_i,i＝1,2,…,N表示第i区域；

k为水平集函数的曲率， $k=\stackrel{\&RightArrow;}{\&dtri;}\&CenterDot;\frac{\stackrel{\&RightArrow;}{\&dtri;}{u}_i}{\left|\right|\stackrel{\&RightArrow;}{\&dtri;}{u}_i\left|\right|}=\frac{u_{\mathrm{xx}}{u}_y^2-2u_x{u}_y{u}_{\mathrm{xy}}+u_{\mathrm{yy}}{u}_x^2}{{(u_x^2+u_y^2)}^{3/2}},$

$\left\{\begin{array}{cc}{\mathrm{\&chi;}}_{u_i(x,t)>0}=H\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_i\right)& u_i(x,t)>0\\ {\mathrm{\&chi;}}_{u_i(x,t)\&le;0}=1-H\left({\stackrel{\&RightArrow;}{\mathrm{\&gamma;}}}_i\right)& u_i(x,t)\&le;0\end{array}\right.,$

第三步：以2ⁱ，i＝2,…L为系数采用双线性插值方法上采样演化曲线，得到的该演化曲线作为下一分辨率构建初始化演化曲线，然后构建该分辨率图像总能量模型，利用变分水平集最小化能量模型，采用多分辨率水平集方法，进行曲线演化得到当前分辨率下N-1个零水平集演化曲线方程；

所述的多分辨率水平集方法，包括两个过程：首先，对该图像的每一维进行下采样生成分辨率为2ⁱ,i=2,…,L的低分辨率图像；其次，相对应地，以2ⁱ,i=2,…,L为系数上采样该分辨率下的演化曲线，然后构建该分辨率图像总能量模型，利用变分水平集最小化能量模型，进行曲线演化得到N-1个零水平集分割图像；演化过程不断重复直至达到原始分辨率图像，得到分割结果；

第四步：判断当前图像分辨率级数是否达到原始分辨率，当达到原始分辨率，则输出分割结果；否则，返回第二步。

2.根据权利要求1所述的多分辨率多区域变分水平集图像分割方法，其特征是，所述的将原始图像按照空间分辨率在每一维进行连续下采样生成分辨率为2^L的图像：先设置分辨率的级数L，生成分辨率级数为L低分辨率的图像，作为初始分辨率的图像。

Images