图像分割—分水岭算法及Opencv实现

分水岭算法-图像分割

1.原理

在分水岭算法中，一幅图像中灰度值高的区域被看作山峰，灰度值低的区域被看作山谷。 然后从山谷的最低点灌水，水会慢慢在不同的地方汇合，而这些汇合的地方就是需要对图像分割的地方。 分水岭算法的核心思想就是建立堤坝阻止不同盆地的水汇合。 在一般分水岭算法中，通常是把一副彩色图像灰度化，然后再求梯度图， 最后在梯度图的基础上进行分水岭算法，求得分段图像的边缘线。如下所示是一个“地形”的剖面示意图。

以绿色虚线为界，左边表示地物1，右边表示地物2。A为地物1在当前范围内的最小值点， E为地物2在当前范围内的最小值点，两个盆地的交汇点为D。在地物1中C为小范围内的极小值点。

首先在两个盆地的最小值点A、E开始向盆地中注水，水会缓慢上升。 在两个盆地的水汇集的时刻，在交接的边缘线上(D点所在位置，也即分水岭线)， 建一个堤坝(图中黑色线段)，来阻止两个盆地的水汇集成一片水域。 这样图像就被分成2个像素集，一个是注水盆地像素集，一个是分水岭线像素集。

但仔细观察就会发现问题，传统的基于图像梯度的分水岭算法由于存在太多极小区域而产生许多小的集水盆地， 带来的结果就是图像过分割。 如图C点所在的极小值区域会形成一个小盆地，从而让地物1被分成两部分， 当C盆地的水和A盆地的水要汇合时，会在B点建立个水坝(图中灰色虚线)， 这显然不是我们想要的结果。 所以必须对分割相似的结果进行合并。 举个例子如一个桌面的图片，由于光照、纹理等因素，桌面会有很多明暗变化，反映在梯度图上就是一个个圈， 此时利用分水岭算法就会出现很多小盆地，从而分割出很多小区域。但这显而易见是不符合常识的。 因为桌面是一个整体，应该属于同一类，而不是因为纹理而分成不同的部分。

因此需要对分水岭算法进行改进。在OpenCV中采用的是基于标记的分水岭算法。 水淹过程从预先定义好的标记图像（像素）开始， 这样可以减少很多极小值点盆地产生的影响。 较好的克服了过度分割的不足。 本质上讲，基于标记点的改进算法是利用先验知识来帮助分割的一种方法。 对比如下图所示。

2.Opencv实现

从调用API的角度而言，在整个过程中并没有对原始影像求梯度，而是直接进行二值化， 进行像素标记。然后将标记图像与原图一起传入函数。 具体求梯度的操作封装在了函数中，用户无需关心。
在OpenCV中实现分水岭算法可以使用cv2.watershed()函数实现，主要有以下步骤：
步骤1.输入图像，对图像进行二值化
步骤2.对二值化后的图像进行噪声去除
步骤3.通过腐蚀、膨胀运算对图像进行前景、背景标注
步骤4.运用分水岭算法对图像进行分割

具体实例如下，下面是待分割的影像。影像中有很多彼此连接的硬币， 我们需要将这些硬币彼此分开。

2.1 图像读取，二值化操作

import numpy as np
import cv2
from matplotlib import pyplot as plt
"步骤1：二值化图像"
img = cv2.imread("E:/ruanjianDM/jupyternoerbookDM/Opencv3/data/coins.jpg")
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
# cv2.THRESH_BINARY_INV 背景和物体灰度都比较浅，只用它，会反色，全部为黑色
#cv2.THRESH_OTSU：背景变白色，物体黑色
#两者结合，先otsu，再反色
ret, thresh = cv2.threshold(gray, 0, 255, cv2.THRESH_BINARY_INV+cv2.THRESH_OTSU)
plt.subplot(121),plt.imshow(255-thresh,cmap='binary')
plt.axis("off")
plt.title("二值化图像")
"""
步骤2:
仔细观察二值化后的图像，会发现有一些噪声点。如果这些噪声点不去除，则在后续步骤中会被标记成前景或背景， 从而影响分割。
此处只有背景黑色中有小白点，可以使用开运算(先腐蚀后膨胀)去除，效果很好。
"""
kernel = np.ones((3, 3), np.uint8)
open = cv2.morphologyEx(thresh, cv2.MORPH_OPEN, kernel, iterations=2)
plt.subplot(122),plt.imshow(255-open,cmap='binary')
plt.axis("off")
plt.title("去燥之后的二值化图像")

2.2图像标记

上面完成了步骤1-2，现在开始对图像进行标记： 前景(物体)、背景以及未知区域。 我们现在知道靠近对象中心的区域肯定是前景，而远离对象中心的区域肯定是背景，不能确定的区域即是图像边界。 对于前景区域，可以在第二步得到的图像上进行多次腐蚀运算， 这样就可以得到肯定是前景的区域。同样，对该图像进行多次膨胀运算，可以得到比前景大的范围， 除去这些范围的部分肯定是背景。至于在两者之间的就是未知区域，也就是需要运用分水岭算法， 从而给出分水岭边界的地方。

前景标记

方案一：提取肯定是硬币的区域(前景)，可以使用腐蚀操作。腐蚀操作可以去除边缘像素，剩下就可以肯定是硬币了。 当硬币之间没有接触时，这种操作是有效的。但是由于硬币之间是相互接触的，效果并不是很好。

"若直接进行腐蚀操作，效果如下"
open2 = cv2.erode(open, kernel)  # 腐蚀
plt.imshow(255-open2,cmap="binary")
"显然这样的效果是不合理的，因为物体之间有连接，"

方案二：距离变换再加上合适的阈值。

"步骤3.1：前景标记"
# 第一个参数是二值化图像
# 第二个参数是类型distanceType
# 第三个参数是maskSize
# 返回的结果是一张灰度图像，但注意这个图像直接采用OpenCV的imshow显示是有问题的
# 所以采用Matplotlib的imshow显示或者对其进行归一化再用OpenCV显示
distance_transform = cv2.distanceTransform(open, cv2.DIST_L2, 5)
ret, sure_fg = cv2.threshold(distance_transform, 0.7 * distance_transform.max(), 255, cv2.THRESH_BINARY)
plt.subplot(121),plt.imshow(255-distance_transform,cmap='binary')
plt.axis("off")
plt.title('图像距离变换-骨架提取')
plt.subplot(122),plt.imshow(255-sure_fg ,cmap='binary')
plt.axis("off")
plt.title('获取前景')

背景标记、获取标签

在知道了哪些肯定是前景区域后，就可以给它们创建标签了(一个与原图像大小相同，数据类型为int32的数组)。 对我们已经确定分类的区域（无论是前景还是背景）使用不同的正整数标记，对不确定的区域使用0标记。 我们可以使用函数cv2.connectedComponents()来完成。 它会把将背景标记为0，其它对象使用从1开始的正整数标记。 但我们知道如果背景标记为0，那分水岭算法就会把它当成未知区域了。 因此我们还需要对返回的结果进行一些修改，如统一加1或其它数字。 但这里不能减。因为OpenCV会把边界标记成负数。因此这里都使用正数。 在把背景和前景标记完后，最后是将未知区域标记为0。

"步骤3.2：背景标记"#膨胀三次，黑色部分一定是背景
plt.rcParams["font.family"] = "SimHei"#直接修改配置字典，设置默认字体
sure_bg = cv2.dilate(open, kernel, iterations=3)
plt.subplot(141),plt.imshow(255-sure_fg,cmap="binary")
plt.title("获取前景")
plt.axis('off')
plt.subplot(142),plt.imshow(255-sure_bg,cmap="binary")
plt.title("获取背景")
plt.axis('off')
#距离变换后获得的图像数据类型是float32，因此需要转成uint8才能和sure_bg做运算
sure_fg = np.uint8(sure_fg)
unknown = cv2.subtract(sure_bg, sure_fg)
plt.subplot(143),plt.imshow(255-unknown ,cmap="binary")
plt.title("合成")
plt.axis('off')
"步骤3.3：创建标签"
ret, markers1 = cv2.connectedComponents(sure_fg)
markers = markers1 + 1#背景，前景的标记为非0
markers[unknown == 255] = 0#将未知区域标记为0
plt.subplot(144),plt.imshow(markers,cmap='jet')
plt.title("获取标签")
plt.axis('off')

2.3实现分水岭算法

"步骤4：实现分水岭算法"
markers3 = cv2.watershed(img, markers)
#opencv中将分水岭边界设为-1
img[markers3 == -1] = [255, 0, 0]
plt.title("分水岭图像")
plt.axis('off')
plt.subplot(122),plt.imshow(img)
plt.title("分割图像")
plt.axis('off')