Demo Video:
注:图片来自[2]
作者:@Mapupcal
前言:
智能剪刀算法(Intelligent Scissors)是由Eric N. Mortensen和William A. Barrett于1995年提出来的一种图像分割的交互算法,可以用于2D图像分割,论文见[1][2]。
该算法可以用于辅助用户精确地勾勒出感兴趣的区域。
其可以在运行时快速地定位到该图像区域的边缘上,通过与用户交互,更精确地完成整个感兴趣区域的勾画、分割。(如上图所示)
算法理论基础:
①将图像转换成一个庞大的加权有向图。
②在运行时,通过特定的寻路算法,在出发点[a]开始,寻找一条最优路径,该路径通往用户指定的目标点[b]。
[a] 论文中为种子点,SeedPoint,为了和论文相对应,以下均使用SeedPoint,下同。
[b] 论文中为*点,FreePoint.
注:智能剪刀算法的核心在于①,论文[1][2]着重讨论了如何确定加权有向图的边权重。
注:关于②,论文也提供了相应的算法实现,论文当中的算法就是Dijkstra最短路径搜索算法的一种变种。
将图像变成一个庞大的加权有向图:
① 图像中每个像素,在加权有向图中都将成为该图中的一个结点Node。
② 相邻的像素,互相之间存在一条加权有向边。
③ 确定每条加权有向边的权重。
注:在内存表示加权有向图的方式有很多,有每个节点基于顺序容器存放有向边的方式,也有基于矩阵的方式。在实现中,为了数据的局部性,本文采用了基于矩阵的方式去实现可提高性能。
对于某一个非边界像素p,其有八个邻居像素,如:
当转换成加权有向图之后,就变成如下这个样子:
注:上面有向边上面的数字为边的编号,并不是边的实际权重。
根据论文给出的计算边权重的计算,实际影响的因子有以下六个:
拉普拉斯交叉零点值 Fz(q)
像素点梯度值 Fg(q)
像素点梯度方向 Fd(p,q)
注:以上三个为静态因子,当不涉及Training时,有向边权重只由这三个公式确定[1][2]。
像素值 Fp(q)
“内部”像素值 Fi(q)
“外部”像素值 Fo(q)
注:以上三个为动态因子,Traning部分核心,这些因子的函数值会随用户切割图像的时候发生改变[2],本文只会简单提及Training部分。
p表示上两图中的中心结点,q指与p相邻的结点。
[1]确定边权重的表达式如下:
其中Wz = 0.43 Wg=0.43 Wd=0.14
[2]确定边权重的表达式如下:
其中Wz = 0.3 Wg=0.3 Wd=0.1 Wp=0.1 Wi=0.1 Wo=0.1
拉普拉斯交叉零点
拉普拉斯算子可以作为图像锐化的算子,最简单的三阶3x3拉普拉斯算子如下:
用这个算子去处理图像的像素后,得到每个像素对应的拉普拉斯算子处理后的值。
对于每个结果值,我们使用如下函数去处理每一个结果值,得到拉普拉斯交叉零点的成本映射(论文中称CostMap)。
疑惑:边界像素点经过拉普拉斯算子处理后,处于边界的像素的结果值应该比较高的,而处于边界像素的有向边权重应该比较低才对。而不知为何,根据论文的公式,恰恰与这个思想相反。
注:可能是我没有找到正确的拉普拉斯算子去处理图像,我觉得论文的这个公式是反直觉的。
注:在实现中,我使用了符合我自己直觉的方法,效果还算是满意。
注:对于具备RGB彩色的图像,其像素的Fz应该为RGB三个通道中Fz的最大值。
像素的梯度值
“沿着梯度方向,函数变化得最快!”
高数的知识大多都已经丢给回高数老师了,但这句话一直都还记得。基于此,那么是不是可以认为,梯度值大的像素,处于边缘的几率就比梯度小的像素大得多呢?论文[1][2]确实也将梯度值作为权重因子考虑进去了。
梯度值计算方式个根据以下公式来计算的:
其中Ix为该像素在X方向的偏导数,Iy为该像素在Y方向的偏导数。
那么在二维图像上,如何对一个像素求其偏导数呢。
以下两个算子可用于求Ix和Iy:
注:假设横轴为x轴,竖轴为y轴,根据上面你也一定知道哪个是求x偏导,哪个是求y偏导数。(第一个为x,第二个为y)
到此为止,我们应该能够计算所有像素的对应的梯度值了。但是,这些梯度值现在还不能作为Fg的CostMap。因为对于梯度值大(有可能很大,因为求偏导的卷积运算结果是处于255*255*255~0区间)的像素点,我们要降低其梯度值对起有向边权重的影响,而对于梯度值小的,我们也要提高之。
于是我们有以下公式去达到目的:
其中
最后,让切割线段走的路径更加倾向于经过横竖像素点比斜像素点。
根据相邻结点与中心结点p的欧几里得距离,我们还要加上欧几里得距离的权重,权重如下:
注:根据欧几里得距离加权,简单来说,就是横竖边的权重应该要比斜边的权重低。(简单几何学推导就能得出该加权公式)
好了,至此,我们已完成梯度值Fg的CostMap计算。
注意:对于RBG彩色图片,请取RBG三个通道中Fg的最大值。
梯度方向
这里的梯度方向,其实也可以理解为向量PQ(对应点p、q,起点为p,终点为q)与梯度向量之间的量化权重。与梯度方向接近的边PQ权重应该越小。
注:以上均为个人猜测和理解,尚且不清楚论文的公式是否为这种意思。
见[2]中的量化公式(其中[1]此处的公式不尽相同,详情见[1]):
其中
注:D(p)为单位向量。实现的时候记得将之归一化。
注:||p-q||表示向量的模,即向量的平方再开方。
注:对于彩色图片,在实现上,我简单使用了三通道的组合和来代替Ix,Iy。
FP、Fi、Fo
关于Fp、Fi、Fo权重因子,是和training有关的。
关于什么是training,简单点来讲,就是当用户选定了一部分分割边后(这条路径是通过智能剪刀算法来提取的),就可以用这些分割边的像素点作为Samples进行有关像素值因子的学习,从而继续影响每个像素的Fp、Fi、Fo三个CostMap的值(不仅仅是像素值相关的因子,连梯度值因子也可以做成动态的。),这样一来,每条边的权重将进行一次更新,而这次更新后的边权重,执行寻路算法后,将能更精确地选取用户想要的边界。因为有向边的权重会根据用户的交互而进行动态更新,所以说这些因子是动态的。
算法根据用户选定了的分割边的像素进行学习,这就是training。
注:动态更新每条边的权重,会增加运算操作,这些运算的成本将在每一次SeedPoint更新的时候体现。
当有向图的权重确定下来之后,就可以利用加权有向图的寻路算法来进行寻路了。
算法的思想很简单:只要每次设置SeedPoint(某个开始点),就重新计算一次该SeedPoint到其它结点的最短路径。这个时候,所有路径的集合将构成一棵由SeedPoint扩散出去的最短路径树,有了这个最短路径树(每个结点结构中,都有一个指向路径前一个结点的指针,以便快速回溯这条最优路径),交互式的用户响应就呼之欲出了。
因为我们设置有向边权重的时候,区域边界的像素结点的边总是比其它边的权重要低,所以最短路径总是倾向于沿着区域的分割边界产生,这就得到了一把只会切割区域边界的“剪刀”。
注:注意,我用词是“倾向”而不是“总是”。具体原因你可以根据最短路径的知识来推断。
寻路算法
论文[2]中给出了寻路算法:
注:这种算法其实也算是Dijkstra最短路径搜索算法的一种变种。
注:在实现中,我使用了[3]中优化的Dijkstra算法来实现寻路算法。
快速Dijkstra算法的具体思想,就是使用一个“侵入式”的索引优先队列来从当论文[2]算法中的list,以快速获取最小成本的结点,并且,该优先队列可以非常快速地移除其底层堆中具备某个索引值的元素。
详情可以见[4]或者任何一本基础算法书当中有关最快的Dijkstra算法的实现。实现代码中使用的索引优先队列结构在[4]的211页中的习题中有所提及,个人实现是使用了[3]中提供的skeleton代码当中的优先队列和有向图结点数据结构。
实现效果
注:红色的X表示在该点设置了SeedPoint.
注:可以看到效果还是挺好的,最后来抠图。
注:关于抠图算法,在我的代码实现上,并没有对凹多边形做凸多边形分解,只是是简(偷)单(懒)地使用两点之间“扫描线”算法的实现。所以凹多边形的抠图会出现如下问题。
代码实现:智能剪刀——交互式图像分割 实现
注:项目代码基于Opencv3.2,主要用在图像卷积运算,CostMap的数据结构表示,以及图形界面的显示操作。
注:代码注重在功能的实现上,仍有较大的优化空间。
注:本代码并未实现有关training部分,如果你想知道详情,可以参考论文[2]自行实现。
Reference:
[1] Mortensen E N, Barrett W A. Intelligent scissors for image composition[C]//Proceedings of the 22nd annual conference on Computer graphics and interactive techniques. ACM, 1995: 191-198.
[2] Mortensen E N, Barrett W A. Interactive segmentation with intelligent scissors[J]. Graphical models and image processing, 1998, 60(5): 349-384.
[3] CS 4670/5670, Project 1: Image Scissors
[4] Robert Sedgewick Kevin Wayne《算法(中文第四版)》