CN105046689A - 一种基于多层次图结构的交互式立体图像快速分割方法 - Google Patents

Abstract

一种基于多层次图结构的交互式立体图像快速分割方法，首先输入一组立体图像，通过立体图像匹配算法得到视差图。在原始图像左右任意一图中指定部分前、背景。根据指定部分应用CUDA并行计算的方法建立前、背景的颜色以及视差分布的先验统计模型。通过对原始图像进行高斯滤波、下采样得到粗糙尺度较小的图像，然后将粗糙图像与原始图像一起构成多层次图结构。鉴于目前立体图像分割存在分割模型复杂，计算效率低的问题。本发明在基于视差图的立体图像同步分割的理论框架下，探索新的分割方法。力图简化模型的复杂度，并行处理计算密集型的任务，提高立体图像分割速度，实现实时分割常见尺寸立体图像的目的。

Description

一种基于多层次图结构的交互式立体图像快速分割方法

技术领域

本发明属于图像处理、计算机图形学和计算机视觉等交叉领域，涉及一种基于多层次图结构的交互式立体图像快速分割方法。

背景技术

近年来3D技术不断发展，从3D立体电视到3D立体电影，对3D内容的创作以及3D编辑工具的开发提出了迫切的需求。交互式立体图像分割是其中一项重要工作，它是许多应用最重要处理的环节，像物体识别、追踪，图像分类，图像编辑以及图像重建等。目前立体图像分割已经应用于医疗图像中器官的分割与分析，物体的追踪，场景的理解等实际生活中。因此，立体图像分割效率成为重要的研究方向。

相比单幅图像的分割，交互式立体图像的智能分割起步较晚。目前图像分割方法主要存在两方面的挑战：计算准确率和计算速度。这是一对矛盾的问题，很难在两者之间达到较好的平衡。在提高计算准确率方面，人们做了很多的努力。Price等人在2011年的ICCV上发表的“StereoCut：ConsistentInteractiveObjectSelectioninStereoImagePairs”中，利用立体图像对间的视差信息来提高立体图像分割的准确率度。其将图像中每个像素的颜色、梯度、视差等信息融入传统的图割理论中，通过求解最大流来得到立体图像边界优化的结果。这种方法虽然分割精度较高，但是构建的分割模型边与节点的数目庞大，计算复杂，效率低下。目前分割算法多通过改变graphcut算法的具体实施过程来提高分割速度。对于立体图像像素数目多，边结构复杂的问题，仅改变graphcut算法的实施过程无法从根本上解决。同时，在立体图像分割过程中，存在很多单指令流多数据流计算密集型的任务。传统方法没有很好的利用这种任务可并行执行的特点，串行处理，使效率低，消耗大量的时间，从而使分割效率低下。

发明内容

鉴于目前立体图像分割存在分割模型复杂，计算效率低的问题。本发明在基于视差图的立体图像同步分割的理论框架下，探索新的分割方法。力图简化模型的复杂度，并行处理计算密集型的任务，提高立体图像分割速度，实现实时分割常见尺寸立体图像的目的。

为实现这个目标，本发明的技术方案为：首先输入一组立体图像，通过立体图像匹配算法得到视差图。在原始图像左右任意一图中指定部分前、背景。根据指定部分应用CUDA并行计算的方法建立前、背景的颜色以及视差分布的先验统计模型。通过对原始图像进行高斯滤波、下采样得到粗糙尺度较小的图像，然后将粗糙图像与原始图像一起构成多层次图结构。以此为基础，在图割理论框架下形式化多层次图结构中的颜色、梯度以及视差等约束，构造能量函数。为了提高效率，应用CUDA并行计算的方法处理建图过程。采用图的最大流/最小割算法求解多层次图的全局最优化结果。然后统计边界处误差较大的像素点，采用传统的图割理论，对统计的边界像素点进行局部优化。将全局处理与局部优化的结果融合在一起，构成最后的分割结果。若用户没有得到理想的效果，还可以继续对图中错误区域进行勾画，直到得到理想结果。

与现有技术相比，本发明具有以下优点：本发明通过构架基于多层次图结构的立体图像分割模型，简化了边的复杂度，显著提高了处理的速度。同时，将一些计算密集型的单指令流多数据流的任务用CUDA技术并行处理，节省大量时间。实验证明：相比现有方法，在同等交互量的前提下，本发明所述方法在分割准确率以及一致性变化不大的情况下，可以显著提高分割速度。

附图说明

图1为本发明所涉及方法的流程图；

图2为本发明应用实例实验结果：(a)、(b)为输入的左、右图像，(c)、(d)是采用Price等人在2011年的ICCV上发表的“StereoCut：ConsistentInteractiveObjectSelectioninStereoImagePairs”中的方法分割的结果；(e)、(f)为本发明的分割结果；两种方法所用的用户输入在(c)、(e)图中显示，其中第一线条标识前景，第二线条标识背景。同时给出了两种方法分割的准确率以及分割的时间。本实施例测试所用笔记本电脑配置为：CPU处理器Intel(R)Pentium(R)CPUB9502.10GHz2.10GHz；Gpu处理器NVIDIAGeForceGT540M。

具体实施方式

下面结合附图和具体实施方式对本发明做进一步说明。

本发明的流程如图1所示，具体包括如下步骤：

步骤一，匹配立体图像。

读入一对立体图像I＝{I^l，I^r}，I^l与I^r分别表示左、右图像。通过立体匹配算法计算得到左、右图像对应的视差图，分别用D^l与D^r表示。立体匹配算法采用的是Felzenszwalb等人在CVPR04上发表的论文“EfficientBeliefPropagationforEarlyVision”中提出的算法。

步骤二，添加前、背景线索。

用户通过所设计的界面在其中任意一张图像中指定部分前、背景。本发明实施采用类似于Price等人在2011年的ICCV上发表的“StereoCut：ConsistentInteractiveObjectSelectioninStereoImagePairs”中用到的方法，利用鼠标、触摸屏或者手写笔等输入设备，通过在图像上勾画不同颜色的线条指定部分前、背景像素。如图2(e)所示，第一线条覆盖的像素属于前景，第二线条覆盖的像素属于背景。本发明的后续步骤对于该步骤中所用的前、背景像素指定方式并无限制，其它方式亦可使用。

步骤三，建立前、背景的颜色、视差先验模型。

用F表示用户指定的前景像素集合，B表示用户指定的背景像素集合；前、背景的颜色、视差的先验模型采用GMM、直方图以及多个类簇的形式表达。本发明采用的是多类簇形式，通过统计对应像素集合的颜色、视差得到类簇。为了提高处理速度，采用基于CUDA并行的Kmeans算法，对F与B中的像素对应的颜色值、视差值分别进行聚类。处理颜色模型的具体过程如下：每个线程处理一个像素，计算每个像素到所有前景、背景类簇的距离，选择最近的距离，将像素聚类到对应的类簇中。得到N_c个前景颜色类簇M_c个背景颜色类簇上述颜色类簇分别表示前景、背景的颜色分布统计模型；同时，用同样的方法，对F和B中的像素对应的视差值分别进行聚类，得到N_d个前景视差类簇M_d个背景视差类簇上述视差类簇分别表示前景、背景的视差分布统计模型；在本实施例中，N_c＝M_c＝64；N_d＝M_d＝16。

步骤四，基于多层次图结构的全局优化；

由于图像中前景、背景各自的分布比较聚集即前、背景内部像素差异较小，边界处像素差异较大。利用这一特性，用区域具有代表性的像素来表示邻域所有像素。本方法采用高斯滤波、下采样的方式，得到代表性的像素点。进而得到粗糙的尺度较小的图像。将粗糙图像与原始图像融合，构成多层次图结构。对多层次图结构的模型进行全局处理。将原始立体图像对表示为I＝{I^l,I^r}，粗糙的立体图像对表示为I^τ＝{I^l,τ,I^r,τ}，I^l、I^l,τ与I^r、I^r,τ分别表示左、右图像。将原始立体图像与粗糙立体图像共同表示成一个无向图G＝<ν,ε>；其中，ν为无向图G中的节点集合，ε为边的集合；无向图G中的每个顶点对应立体图像I与I^τ中的一个像素；交互式立体图像快速分割是在输入笔画的约束下，为原始立体图像对中的每个像素p_i赋予一个标签x_i；x_i∈{1,0}，分别表示前、背景；无向图G中的边包含每个像素与源点、汇点的连接边，图像内相邻像素的连接边，以及视差图决定的立体图像对应点之间的连接边；同时还包含粗糙层与原始图像的父子节点之间的连接边。令为粗糙层图像像素点。由于粗糙层是对原始层进行下采样得到的，所以一个代表采样前的I图像中N_l*N_l的区域内的像素，在本实施例中N_l＝3。

把求解上述基于多层次图结构的立体图像快速分割问题定义为以下目标能量函数的最优化问题：

$\begin{array}{c}E\left(X\right)=\\ w_{unary}\underset{p_i^{\mathrm{\&tau;}}\&Element;I^{\mathrm{\&tau;}}}{\mathrm{\&Sigma;}}{E}_{unary}\left(p_i^{\mathrm{\&tau;}}\right)+w_{\mathrm{int}ra}\underset{(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})\&Element;N_{\mathrm{int}ra}}{\mathrm{\&Sigma;}}{E}_{\mathrm{int}ra}(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})+w_{\mathrm{int}er}\underset{(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})\&Element;N_{\mathrm{int}er}}{\mathrm{\&Sigma;}}{E}_{\mathrm{int}er}(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})\\ +w_{paternity}\underset{(p_i^{\mathrm{\&tau;}},p_{i,j})\&Element;N_{paternity}}{\mathrm{\&Sigma;}}{E}_{paternity}(p_i^{\mathrm{\&tau;}},p_{ij})\end{array}---\left(1\right)$

其中是一元项，表示粗糙层像素的颜色、视差与前、背景颜色和视差统计模型的相似性，也叫做数据项；相似性越高，值越大；是粗糙层图像内二元项，反映了粗糙层图像所有像素与四邻域之间的差异，Ν_intra表示包含左右粗糙层图中所有像素点的邻接关系的集合；差异越大，则该项越小；根据图割算法的原理，此时邻域像素之间倾向于取不同的标签；是粗糙图像间的二元项，定义了对应点的匹配的结果，匹配度越高，则该项越大；Ν_inter表示含有左右粗糙层像素点对应关系的集合。

是粗糙层图像与原始图像间的二元约束关系，表示父子节点的相似性，父子节点差异越小，该值越大，边界经过两者的可能性越小。Ν_paternity表示父子对应关系的集合。w_unary，w_intra，w_inter，w_paternity调节各能量项之间的权值；w_unary＝1，w_intra＝4000，w_inter＝8000，w_paternity＝1000000。

(1)定义一元约束项

一元约束项包括颜色一元项和视差一元项两部分，定义如下：

$E_{unary}\left(p_i^{\mathrm{\&tau;}}\right)=w_c(1-P_c(x_i^{\mathrm{\&tau;}}\left|c_i^{\mathrm{\&tau;}}\right))+w_d(1-P_d(x_i^{\mathrm{\&tau;}}\left|d_i^{\mathrm{\&tau;}}\right))---\left(2\right)$

其中，表示给定像素的颜色取前景或背景标签的概率值；因为概率越大，能量函数应当越小，所以取1-P_c表示颜色一元项；同样地，表示给定像素的视差值取前景或背景标签的概率值；取1-P_d表示视差一元项；w_c、w_d分别代表颜色与视差的影响权值，w_c+w_d＝1；

本方法以类簇形式表示前、背景的颜色和视差模型，包括N_c个前景颜色类簇M_c个背景颜色类簇N_d个前景视差类簇M_d个背景视差类簇给出一元项的计算方法；

颜色一元项的计算方式如下：本方法采用基于CUDA并行方法来计算。将CPU端的所有像素的颜色值传到GPU端。在GPU中，并行处理所有像素。每个线程表示一个未标记像素。线程相互独立，所有线程同时计算像素颜色到前景、背景颜色模型的类簇中心的距离，找到其中最小的距离；用这个最小的距离描述像素颜色与前、背景颜色的相似度；离前景或背景颜色距离越小，则颜色越相近，根据图割理论，该像素越倾向于选择前景或背景标签；待所有线程结束，将GPU端每个像素的求解结果传到CPU端，在CPU端进行详细的建图过程。颜色一元项的数学形式描述为：

$1-P_c\left(x_i^{\mathrm{\&tau;}}\right|c_i^{\mathrm{\&tau;}})=\left\{\begin{array}{c}\frac{s_i^{\min }}{s_i^{\min }+t_i^{\min }},x=1\\ \frac{t_i^{\min }}{s_i^{\min }+t_i^{\min }},x=0\end{array}\right.---\left(3\right)$

其中，分别表示像素的颜色到前景和背景颜色的各类簇中心的最小距离，其表达式分别为：

$s_i^{\min }=\min \left({||c_i^{\mathrm{\&tau;}}-C_n^F||}^2\right),n=1,...,N_c$

$t_i^{\min }=\min \left({||c_i^{\mathrm{\&tau;}}-C_m^B||}^2\right),m=1,...,M_c$

视差一元项与颜色一元项的计算过程相同；

(2)定义图像内二元约束项

图像内二元约束项包含两项，分别描述像素点周围颜色变化和视差变化，即颜色梯度和视差梯度，定义如下：

$E_{\mathrm{int}ra}(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})=f_c(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})f_d(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})|x_i^{\mathrm{\&tau;}}-x_j^{\mathrm{\&tau;}}|---\left(4\right)$

其中，表示相邻像素间颜色的相似性，颜色越相近其值越大，根据图割算法的原理，边界穿过二者的几率就较小；表示像素相对于邻接像素点视差的相似性；二者视差越相近，其值越大，根据图割算法的原理，二者取不同标签的几率就较小；为了减少视差产生的误差，视差项中的视差，本步骤采用的是经过高斯滤波以及下采样得到的粗糙层的视差信息。两项的定义形式如下：

$f_c(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})=\frac{1}{{||c_i^{\mathrm{\&tau;}}-c_j^{\mathrm{\&tau;}}||}^2+1},(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})\&Element;N_{\mathrm{int}ra}---\left(5\right)$

$f_d(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})=\frac{1}{{||d_i^{\mathrm{\&tau;}}-d_j^{\mathrm{\&tau;}}||}^2+1},(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})\&Element;N_{\mathrm{int}ra}---\left(6\right)$

(3)定义图像间二元约束项

图像间二元项约束图像间对应像素取相同标签，定义如下：

$E_{\mathrm{int}er}(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})=\frac{C(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})+C(p_i^{r,\mathrm{\&tau;}},p_i^{l,\mathrm{\&tau;}})}{2}|x_i^{l,\mathrm{\&tau;}}-x_j^{r,\mathrm{\&tau;}}|---\left(7\right)$

其中，C表示立体图像中之间作为对应点的可能性，是非对称函数：

$C(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})=P\left(x_i^{l,\mathrm{\&tau;}}\right|M\left(p_i^{l,\mathrm{\&tau;}}\right)=p_j^{r,\mathrm{\&tau;}},p_j^{l,\mathrm{\&tau;}})P\left(M\right(p_i^{l,\mathrm{\&tau;}})=p_j^{r,\mathrm{\&tau;}})---\left(8\right)$

是基于视差图确定的之间作为对应点的概率分布函数；函数表示是左粗糙层像素在右粗糙层上的对应点，对应关系根据原始视差图决定；采用一致的Delta函数，定义方式如下；

$P(M\left(p_i^{l,\mathrm{\&tau;}}\right)=p_j^{t,\mathrm{\&tau;}})=\{\begin{array}{c}1,\left|p_i^{l,\mathrm{\&tau;}}-p_j^{r,\mathrm{\&tau;}}\right|=d_i^{l}and\left|p_j^{r,\mathrm{\&tau;}}-p_i^{l,\mathrm{\&tau;}}\right|=p_j^r\\ 0,other\end{array}---\left(9\right)$

其中，为左粗糙层中像素与右图中对应点的视差值；为右粗糙层中像素与左图对应点的视差；为了更好的确定左右图像素的对应关系，在此采用的是未经处理的原始视差图的视差。

式(8)中表示与之间的颜色相似的概率，在视差完全准确的情况下，但目前的视差计算方法存在误差，为了更好的确定左右图的对应关系，摒弃了视差项。仅利用颜色项，采取如下形式：

$P\left(x_i^{l,\mathrm{\&tau;}}\right|M\left(p_i^{l,\mathrm{\&tau;}}\right)=p_j^{r,\mathrm{\&tau;}},x_j^{r,\mathrm{\&tau;}})=\frac{1}{{||c_i^{l,\mathrm{\&tau;}}-c_j^{r,\mathrm{\&tau;}}||}^2+1}---\left(10\right)$

其中，为左粗糙层图像素的颜色值，是在右粗糙层对应点的值；

(4)定义上下层间的父子约束关系

图像分割最终的结果应在像素层中表示出来。为了将粗糙层的处理结果传递到像素层，同时保持上下层图像间的父子像素的一致性，将上下层间的父子约束关系定义为：

$E_{paternity}(p_i^{\mathrm{\&tau;}},p_{i,j})=\mathrm{\&infin;},(p_i^{\mathrm{\&tau;}},p_{i,j})\&Element;N_{paternity}---\left(11\right)$

表示上下层父子像素间的相似性。由于粗糙层的像素代表原始像素层N_l*N_l区域的所有像素，粗糙层像素的标签即代表像素层对应区域的所有像素标签，因此将父子像素间的边权定义为无穷大。非父子节点像素的边不再考虑。

(5)求解能量函数最小值

对于上下层间的父子约束关系，本发明中定义为无穷大，因此父子之间的边永不会被分割，父节点的标签会直接传递到子节点。由于计算父子节点的边会消耗大量的内存，同时增加计算的时间。在具体优化求解过程中，不再详细计算父子节点间的边。采用图割算法，例如YuriBoykov等人于2004年在《IEEETransactiononPAMI》上发表的论文“AnExperimentalComparisonofMin-Cut/Max-FlowAlgorithmsforEnergyMinimizationinVision”中所提出的最大流/最小割算法，通过最优化本发明所定义的能量函数(式(1))，得到最优的标记结果，即粗糙层分割结果。然后根据粗糙层像素的标签，直接确定像素层对应的区域像素标签。通过这种方法在准确率不变的情况下，可以显著提高分割的速度。由于直接将粗糙层的标签传递到像素层，对于边界处邻域像素差异较大的像素点存在较大的误差。为了提高分割的准确率，统计边界处误差较大的点，进行局部优化。

步骤五，基于原始图像的边界处局部优化

经过步骤四的全局优化，得到粗糙的分割边界。由于粗糙层像素对应原始像素层的N_l*N_l区域内像素的集合，将的标签直接传递到像素层N_l*N_l的区域。在本实施例中N_l＝3。对于边界处，邻域像素的差异大，直接把粗糙层像素的标签赋给区域的所有像素，会存在较大的误差。因此，对边界处进行单独的局部优化。

进行局部优化前，先***部边界信息。首先将得到的粗糙的分割边界分为上、下边界与左、右的边界两部分。然后将上、下边界向边界线的上面与下面分别扩充N_l个像素，将左、右边界分别向边界线的左面与右面扩充N_l个像素，在本实施例中N_l＝3。对统计的边界像素，采用传统图割理论进行局部优化。局部优化是在像素层上进行的，由于视差计算存在误差，在局部优化时放弃了视差信息。在全局处理时，保证了立体图像分割的一致性，而且局部优化是对局部像素点进行的处理。因此，在局部优化时，同时在左右两幅图像上独立进行。若I^e为统计的局部待处理图。定义局部的能量函数为：

$E^e\left(X\right)=w_{unary}\underset{p_i\&Element;I^e}{\&Sigma;}{E}_{unary}^e\left(p_i\right)+w_{\mathrm{int}ra}\underset{(p_i,p_j)\&Element;N_{\mathrm{int}ra}^e}{\&Sigma;}{E}_{\mathrm{int}ra}^e(p_i,p_j)---\left(12\right)$

是一元项即数据项，表示边界处的像素与前、背景颜色模型的相似性，相似性越大，值越大。是二元项即平滑项，表示邻域像素的相似性，二者越相似，值越小。边界经过二者的可能性就越小。代表边界图中所有邻接关系的结合。其中，w_unary+w_intra＝1

一元项具体定义如下：

$E_{unary}^e\left(p_i\right)=P\left(x_i\right|c_i)=\frac{p\left(c_i\right|x_i)}{P\left(c_i\right|x_i=1)+P\left(c_i\right|x_i=0)}---\left(13\right)$

边界处的优化是局部的精确的优化，应尽可能减少误差，因此，一元项仅采用颜色项。一元项的具体计算同全局优化中一元项颜色的计算。

二元项为了减少误差，也仅采用颜色项。具体定义如下所示：

$E_{\mathrm{int}ra}^e(p_i,p_j)=\frac{1}{{||c_i-c_j||}^2+1}|x_i-x_j|,(p_i,p_j)\&Element;N_{\mathrm{int}ra}^e---\left(14\right)$

局部能量函数定义好后，采用步骤四提到的最大流/最小割优化算法，最优化局部能量函数即式(12)，得到最优的标记结果，即分割结果；同步骤四分割的结果相融合，构成整个图像对的分割结果。

步骤六，交互

如对分割结果不满意，返回步骤二，继续添加前、背景线索；每添加一笔，将触发一次完整的分割过程。在已分割的基础上，进行进一步的分割，直到得到满意的结果。

以Price等人在2011年的ICCV上发表的“StereoCut：ConsistentInteractiveObjectSelectioninStereoImagePairs”中的方法为对比对象，说明本发明方法的有效性。两种方法均采用一致的Delta函数(式(9))作为对应点之间的概率分布函数。图2给出了效果对比。图2(a)、(b)为输入的左、右图像。(c)、(d)是采用StereoCut方法分割的结果；图2(e)、(f)为本发明的分割结果；下面两列给出了两种方法分割的准确率以及分割的总时间。准确率(用A表示)的具体定义如下：

$A=\frac{1}{2}(\frac{{\mathrm{\&Sigma;}}_{i=1}^{N_L}{f}_A(c_i^l-c_i^{\lg })}{N_L}+\frac{{\mathrm{\&Sigma;}}_{j=1}^{N_r}{f}_A(c_j^r-c_j^{rg})}{N_r}),---\left(15\right)$

其中 $f_A\left(\mathrm{\&Delta;}c\right)=\left\{\begin{array}{cc}1,& \mathrm{\&Delta;}c=0\\ 0,& \mathrm{\&Delta;}c\&NotEqual;0\end{array}\right.$

其中，N_L和N_r分别表示左图和右图图像的像素总数，为分割后左图中第i个像素的标签(0或1),相应的表示分割后右图第j个像素的标签。分别代表左、右图真值，则反映了左图某一像素的标签与真值之间的差异。函数f_A是关于差异的函数，差异为0时，函数为1，否则记为0。从公式(15)可看出，单幅图像中与真值的无差异总数与图像大小的比值即为分割的准确率，立体图像的分割准确性就是左右两图准确率的平均值。

两种方法所用的用户输入分别在图(c)、(e)中显示，目标物内部的第一线条的线标记前景，目标物外部的第二线条的线标记背景。对比图(c)、(d)和图(e)、(f)，以及所给出的两种方法的计算时间和准确率值，可看出：本方法在同等交互量的前提下，在分割准确率变化不大的情况下，可以显著提高图像分割的速度。

Claims (1)

1.一种基于多层次图结构的交互式立体图像快速分割方法，其特征在于：该方法首先输入一组立体图像，通过立体图像匹配算法得到视差图；在原始图像左右任意一图中指定部分前、背景；根据指定部分应用CUDA并行计算的方法建立前、背景的颜色以及视差分布的先验统计模型；通过对原始图像进行高斯滤波、下采样得到粗糙尺度较小的图像，然后将粗糙图像与原始图像一起构成多层次图结构；以此为基础，在图割理论框架下形式化多层次图结构中的颜色、梯度以及视差等约束，构造能量函数；为了提高效率，应用CUDA并行计算的方法处理建图过程；采用图的最大流/最小割算法求解多层次图的全局最优化结果；然后统计边界处误差较大的像素点，采用传统的图割理论，对统计的边界像素点进行局部优化；将全局处理与局部优化的结果融合在一起，构成最后的分割结果；若用户没有得到理想的效果，继续对图中错误区域进行勾画，直到得到理想结果；

其特征在于：该方法具体包括如下步骤：

步骤一，匹配立体图像；

读入一对立体图像I＝{I^l，I^r}，I^l与I^r分别表示左、右图像；通过立体匹配算法计算得到左、右图像对应的视差图，分别用D^l与D^r表示；

步骤二，添加前、背景线索；

用户通过所设计的界面在其中任意一张图像中指定部分前、背景；利用鼠标、触摸屏或者手写笔等输入设备，通过在图像上勾画不同颜色的线条指定部分前、背景像素；第一线条覆盖的像素属于前景，第二线条覆盖的像素属于背景；本方法的后续步骤对于该步骤中所用的前、背景像素指定方式并无限制，其它方式亦可使用；

步骤三，建立前、背景的颜色、视差先验模型；

用F表示用户指定的前景像素集合，B表示用户指定的背景像素集合；前、背景的颜色、视差的先验模型采用GMM、直方图以及多个类簇的形式表达；本方法采用的是多类簇形式，通过统计对应像素集合的颜色、视差得到类簇；为了提高处理速度，采用基于CUDA并行的Kmeans算法，对F与B中的像素对应的颜色值、视差值分别进行聚类；处理颜色模型的具体过程如下：每个线程处理一个像素，计算每个像素到所有前景、背景类簇的距离，选择最近的距离，将像素聚类到对应的类簇中；得到N_c个前景颜色类簇M_c个背景颜色类簇上述颜色类簇分别表示前景、背景的颜色分布统计模型；同时，用同样的方法，对F和B中的像素对应的视差值分别进行聚类，得到N_d个前景视差类簇M_d个背景视差类簇上述视差类簇分别表示前景、背景的视差分布统计模型；在本实施例中，N_c＝M_c＝64；N_d＝M_d＝16；

步骤四，基于多层次图结构的全局优化；

由于图像中前景、背景各自的分布比较聚集即前、背景内部像素差异较小，边界处像素差异较大；利用这一特性，用区域具有代表性的像素来表示邻域所有像素；本方法采用高斯滤波、下采样的方式，得到代表性的像素点；进而得到粗糙的尺度较小的图像；将粗糙图像与原始图像融合，构成多层次图结构；对多层次图结构的模型进行全局处理；将原始立体图像对表示为I＝{I^l,I^r}，粗糙的立体图像对表示为I^τ＝{I^l,τ,I^r,τ}，I^l、I^l,τ与I^r、I^r,τ分别表示左、右图像；将原始立体图像与粗糙立体图像共同表示成一个无向图G＝<ν,ε>；其中，ν为无向图G中的节点集合，ε为边的集合；无向图G中的每个顶点对应立体图像I与I^τ中的一个像素；交互式立体图像快速分割是在输入笔画的约束下，为原始立体图像对中的每个像素p_i赋予一个标签x_i；x_i∈{1,0}，分别表示前、背景；无向图G中的边包含每个像素与源点、汇点的连接边，图像内相邻像素的连接边，以及视差图决定的立体图像对应点之间的连接边；同时还包含粗糙层与原始图像的父子节点之间的连接边；令为粗糙层图像像素点；由于粗糙层是对原始层进行下采样得到的，所以一个代表采样前的I图像中N_l*N_l的区域内的像素，在本实施例中N_l＝3；

把求解上述基于多层次图结构的立体图像快速分割问题定义为以下目标能量函数的最优化问题：

$\begin{array}{c}E\left(X\right)=\\ w_{unary}\underset{p_i^{\mathrm{\&tau;}}\&Element;I^{\mathrm{\&tau;}}}{\&Sigma;}{E}_{unary}\left(p_i^{\mathrm{\&tau;}}\right)+w_{\mathrm{int}ra}\underset{(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})\&Element;N_{\mathrm{int}ra}}{\&Sigma;}{E}_{\mathrm{int}ra}(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})+w_{\mathrm{int}er}\underset{(p_i^{l,\mathrm{\&tau;}},p_j^{r,\mathrm{\&tau;}})\&Element;N_{\mathrm{inter}}}{\&Sigma;}{E}_{\mathrm{inter}}(p_i^{l.\mathrm{\&tau;}},p_j^{r,\mathrm{\&tau;}})\\ +w_{paternity}\underset{(p_i^{\mathrm{\&tau;}},p_{i,j})\&Element;N_{paternity}}{\&Sigma;}{E}_{paternity}(p_i^{\mathrm{\&tau;}},p_{i,j})\end{array}---\left(1\right)$

其中是一元项，表示粗糙层像素的颜色、视差与前、背景颜色和视差统计模型的相似性，也叫做数据项；相似性越高，值越大；是粗糙层图像内二元项，反映了粗糙层图像所有像素与四邻域之间的差异，Ν_intra表示包含左右粗糙层图中所有像素点的邻接关系的集合；差异越大，则该项越小；根据图割算法的原理，此时邻域像素之间倾向于取不同的标签；是粗糙图像间的二元项，定义了对应点的匹配的结果，匹配度越高，则该项越大；Ν_inter表示含有左右粗糙层像素点对应关系的集合；是粗糙层图像与原始图像间的二元约束关系，表示父子节点的相似性，父子节点差异越小，该值越大，边界经过两者的可能性越小；Ν_paternity表示父子对应关系的集合；w_unary，w_intra，w_inter，w_paternity调节各能量项之间的权值；在本发明中w_unary＝1，w_intra＝4000，w_inter＝8000，w_paternity＝1000000。

(1)定义一元约束项

一元约束项包括颜色一元项和视差一元项两部分，定义如下：

$E_{unary}\left(p_i^{\mathrm{\&tau;}}\right)=w_c(1-P_c(x_i^{\mathrm{\&tau;}}\left|c_i^{\mathrm{\&tau;}}\right))+w_d(1-P_d(x_i^{\mathrm{\&tau;}}\left|d_i^{\mathrm{\&tau;}}\right))---\left(2\right)$

其中，表示给定像素的颜色取前景或背景标签的概率值；因为概率越大，能量函数应当越小，所以取1-P_c表示颜色一元项；同样地，表示给定像素的视差值取前景或背景标签的概率值；取1-P_d表示视差一元项；w_c、w_d分别代表颜色与视差的影响权值，w_c+w_d＝1；

本方法以类簇形式表示前、背景的颜色和视差模型，包括N_c个前景颜色类簇M_c个背景颜色类簇N_d个前景视差类簇M_d个背景视差类簇给出一元项的计算方法；

颜色一元项的计算方式如下：本方法采用基于CUDA并行方法来计算；将CPU端的所有像素的颜色值传到GPU端；在GPU中，并行处理所有像素；每个线程表示一个未标记像素；线程相互独立，所有线程同时计算像素颜色到前景、背景颜色模型的类簇中心的距离，找到其中最小的距离；用这个最小的距离描述像素颜色与前、背景颜色的相似度；离前景或背景颜色距离越小，则颜色越相近，根据图割理论，该像素越倾向于选择前景或背景标签；待所有线程结束，将GPU端每个像素的求解结果传到CPU端，在CPU端进行详细的建图过程；颜色一元项的数学形式描述为：

$1-P_c\left(x_i^{\mathrm{\&tau;}}\right|c_i^{\mathrm{\&tau;}})=\left\{\begin{array}{cc}\frac{s_i^{\min }}{s_i^{\min }+t_i^{\min }},& x=1\\ \frac{t_i^{\min }}{s_i^{\min }+t_i^{\min }},& x=0\end{array}\right.---\left(3\right)$

其中，分别表示像素的颜色到前景和背景颜色的各类簇中心的最小距离，其表达式分别为：

$s_i^{\min }=min\left(\right||c_i^{\mathrm{\&tau;}}-C_n^F|{|}^2),n=1,...,N_c$

$t_i^{\min }=min\left(\right||c_i^{\mathrm{\&tau;}}-C_m^B|{|}^2),m=1,...,M_c$

视差一元项与颜色一元项的计算过程相同；

(2)定义图像内二元约束项

图像内二元约束项包含两项，分别描述像素点周围颜色变化和视差变化，即颜色梯度和视差梯度，定义如下：

$E_{\mathrm{int}ra}(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})=f_c(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})f_d(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})|x_i^{\mathrm{\&tau;}}-x_j^{\mathrm{\&tau;}}|---\left(4\right)$

其中，表示相邻像素间颜色的相似性，颜色越相近其值越大，根据图割算法的原理，边界穿过二者的几率就较小；表示像素相对于邻接像素点视差的相似性；二者视差越相近，其值越大，根据图割算法的原理，二者取不同标签的几率就较小；为了减少视差产生的误差，视差项中的视差，本步骤采用的是经过高斯滤波以及下采样得到的粗糙层的视差信息；两项的定义形式如下：

$f_c(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})=\frac{1}{\left|\right|c_i^{\mathrm{\&tau;}}-c_j^{\mathrm{\&tau;}}|{|}^2+1},(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})\&Element;N_{\mathrm{int}ra}---\left(5\right)$

$f_d(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})=\frac{1}{\left|\right|d_i^{\mathrm{\&tau;}}-d_j^{\mathrm{\&tau;}}|{|}^2+1},(p_i^{\mathrm{\&tau;}},p_j^{\mathrm{\&tau;}})\&Element;N_{\mathrm{int}ra}---\left(6\right)$

(3)定义图像间二元约束项

图像间二元项约束图像间对应像素取相同标签，定义如下：

$E_{\mathrm{int}er}(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})=\frac{C(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})+C(p_i^{r,\mathrm{\&tau;}},p_i^{l,\mathrm{\&tau;}})}{2}|x_i^{l,\mathrm{\&tau;}}-x_j^{r,\mathrm{\&tau;}}|---\left(7\right)$

其中，C表示立体图像中之间作为对应点的可能性，是非对称函数：

$C(p_i^{l,\mathrm{\&tau;}},p_i^{r,\mathrm{\&tau;}})=P\left(x_i^{l,\mathrm{\&tau;}}\right|M\left(p_i^{l,\mathrm{\&tau;}}\right)=p_j^{r,\mathrm{\&tau;}},x_j^{r,\mathrm{\&tau;}})P\left(M\right(p_i^{l,\mathrm{\&tau;}})=p_j^{r,\mathrm{\&tau;}})---\left(8\right)$

是基于视差图确定的之间作为对应点的概率分布函数；函数表示是左粗糙层像素在右粗糙层上的对应点，对应关系根据原始视差图决定；采用一致的Delta函数，定义方式如下；

$P(M\left(p_i^{l,\mathrm{\&tau;}}\right)=p_j^{r,\mathrm{\&tau;}})=\left\{\begin{array}{cc}1,& |p_i^{l,\mathrm{\&tau;}}-p_j^{r,\mathrm{\&tau;}}|=d_i^{l}and|p_j^{r,\mathrm{\&tau;}}-p_i^{l,\mathrm{\&tau;}}|=d_j^r\\ 0,& others\end{array}\right.---\left(9\right)$

其中，为左粗糙层中像素与右图中对应点的视差值；为右粗糙层中像素与左图对应点的视差；为了更好的确定左右图像素的对应关系，在此采用的是未经处理的原始视差图的视差；

式(8)中表示与之间的颜色相似的概率，在视差完全准确的情况下，但目前的视差计算方法存在误差，为了更好的确定左右图的对应关系，摒弃了视差项；仅利用颜色项，采取如下形式：

$P\left(x_i^{l,\mathrm{\&tau;}}\right|M\left(p_i^{l,\mathrm{\&tau;}}\right)=p_j^{r,\mathrm{\&tau;}},x_j^{r,\mathrm{\&tau;}})=\frac{1}{\left|\right|c_i^{l,\mathrm{\&tau;}}-c_j^{r,\mathrm{\&tau;}}|{|}^2+1}---\left(10\right)$

其中，为左粗糙层图像素的颜色值，是在右粗糙层对应点的值；(4)定义上下层间的父子约束关系

图像分割最终的结果应在像素层中表示出来；为了将粗糙层的处理结果传递到像素层，同时保持上下层图像间的父子像素的一致性，将上下层间的父子约束关系定义为：

$E_{paternity}(p_i^{\mathrm{\&tau;}},p_{i,j})=\mathrm{\&infin;},(p_i^{\mathrm{\&tau;}},p_{i,j})\&Element;N_{paternity}---\left(11\right)$

表示上下层父子像素间的相似性；由于粗糙层的像素代表原始像素层N_l*N_l区域的所有像素，粗糙层像素的标签即代表像素层对应区域的所有像素标签，因此将父子像素间的边权定义为无穷大；非父子节点像素的边不再考虑；

(5)求解能量函数最小值

对于上下层间的父子约束关系，本方法中定义为无穷大，因此父子之间的边永不会被分割，父节点的标签会直接传递到子节点；由于计算父子节点的边会消耗大量的内存，同时增加计算的时间；在具体优化求解过程中，不再详细计算父子节点间的边；采用图割算法，通过最优化本方法所定义的能量函数(式(1))，得到最优的标记结果，即粗糙层分割结果；然后根据粗糙层像素的标签，直接确定像素层对应的区域像素标签；通过这种方法在准确率不变的情况下，可以显著提高分割的速度；由于直接将粗糙层的标签传递到像素层，对于边界处邻域像素差异较大的像素点存在较大的误差；为了提高分割的准确率，统计边界处误差较大的点，进行局部优化；

步骤五，基于原始图像的边界处局部优化

经过步骤四的全局优化，得到粗糙的分割边界；由于粗糙层像素对应原始像素层的N_l*N_l区域内像素的集合，将的标签直接传递到像素层N_l*N_l的区域；对于边界处，邻域像素的差异大，直接把粗糙层像素的标签赋给区域的所有像素，会存在较大的误差；因此，对边界处进行单独的局部优化；

进行局部优化前，先***部边界信息；首先将得到的粗糙的分割边界分为上、下边界与左、右的边界两部分；然后将上、下边界向边界线的上面与下面分别扩充N_l个像素，将左、右边界分别向边界线的左面与右面扩充N_l个像素；在本发明中N_l＝3；对统计的边界像素，采用传统图割理论进行局部优化；局部优化是在像素层上进行的，由于视差计算存在误差，在局部优化时放弃了视差信息；在全局处理时，保证了立体图像分割的一致性，而且局部优化是对局部像素点进行的处理；因此，在局部优化时，同时在左右两幅图像上独立进行；若I^e为统计的局部待处理图；定义局部的能量函数为：

$E^e\left(X\right)=w_{unary}\underset{p_i\&Element;I^e}{\&Sigma;}{E}_{unary}^e\left(p_i\right)+w_{\mathrm{int}ra}\underset{(p_i,p_j)\&Element;N_{\mathrm{int}ra}^e}{\&Sigma;}{E}_{\mathrm{int}ra}^e(p_i,p_j)---\left(12\right)$

是一元项即数据项，表示边界处的像素与前、背景颜色模型的相似性，相似性越大，值越大；是二元项即平滑项，表示邻域像素的相似性，二者越相似，值越小；边界经过二者的可能性就越小；代表边界图中所有邻接关系的结合；一元项具体定义如下：

$E_{unary}^e\left(p_i\right)=P\left(x_i\right|c_i)=\frac{P\left(c_i\right|x_i)}{p\left(c_i\right|x_r=1)+p\left(c_i\right|x_i=0)}---\left(13\right)$

边界处的优化是局部的精确的优化，应尽可能减少误差，因此，一元项仅采用颜色项；一元项的具体计算同全局优化中一元项颜色的计算；

二元项为了减少误差，也仅采用颜色项；具体定义如下所示：

$E_{\mathrm{int}ra}^e(p_i,p_j)=\frac{1}{\left|\right|c_i-c_j|{|}^2+1}|x_i-x_j|,(p_i,p_j)\&Element;N_{\mathrm{int}ra}^e---\left(14\right)$

局部能量函数定义好后，采用步骤四提到的最大流/最小割优化算法，最优化局部能量函数即式(12)，得到最优的标记结果，即分割结果；同步骤四分割的结果相融合，构成整个图像对的分割结果；

步骤六，交互

如对分割结果不满意，返回步骤二，继续添加前、背景线索；每添加一笔，将触发一次完整的分割过程；在已分割的基础上，进行进一步的分割，直到得到满意的结果。