​一、分水岭理论概述

任何灰度图像都可以看作是一个地形表面，其中高强度表示峰和丘陵，而低强度表示山谷。您开始用不同颜色的水（标签）填充每个孤立的山谷（局部最小值）。随着水位的上升，根据附近的山峰（梯度），来自不同山谷的水，显然具有不同的颜色将开始合并。为避免这种情况，您可以在水合并的位置建立障碍。你继续注水和建造障碍物，直到所有的山峰都被淹没。然后，您创建的障碍会为您提供分割结果。这就是分水岭背后的“哲学”。

但是由于图像中的噪声或任何其他不规则性，这种方法会给您带来过度分割的结果。因此OpenCV实现了一个基于标记的分水岭算法，您可以在其中指定要合并哪些谷点，哪些不合并。它是一种交互式图像分割。我们所做的是为我们所知道的对象赋予不同的标签。用一种颜色（或强度）标记我们确定为前景或对象的区域，用另一种颜色标记我们确定为背景或非对象的区域，最后标记我们不确定的区域，用 0 标记它。那是我们的标记。然后应用分水岭算法。然后我们的标记将使用我们提供的标签进行更新，并且对象的边界将具有 -1 的值。

在将图像传递给函数之前，您必须用正 (>0) 索引粗略地勾勒出图像标记中所需的区域。因此，每个区域都表示为一个或多个具有像素值 1、2、3 等的连通分量。可以使用findContours和drawContours从二进制掩码中检索此类标记（请参阅 watershed.cpp 演示）。标记是未来图像区域的“种子”。标记中的所有其他像素，其与轮廓区域的关系未知，应由算法定义，应设置为 0。在函数输出中，标记中的每个像素都设置为“种子”组件的值，或者在区域之间的边界处设置为 -1。

任何两个相邻的连通分量不一定被分水岭边界（-1 的像素）分开；例如，它们可以在传递给函数的初始标记图像中相互接触。

分水岭算法基于Fernand Meyer. Color image segmentation. In Image Processing and its Applications, 1992., International Conference on, pages 303–306. IET, 1992.

基于分水岭的算法一般使用方式如下：

（1）使用诸如 OpenCV 函数cv::filter2D执行一些拉普拉斯滤波以进行图像锐化。

（2）使用 OpenCV 函数cv::distanceTransform以获得二值图像的派生表示，其中每个像素的值被其到最近的背景像素的距离替换。

（3）用 OpenCV 函数cv::watershed将图像中的对象与背景隔离开来。

二、watershed函数

1、函数原型

void cv::watershed(InputArray image,InputOutputArray markers )

2、参数详解

三、OpenCV源码

1、源码路径

opencv\modules\imgproc\src\segmentation.cpp

2、源码代码

void cv::watershed( InputArray _src, InputOutputArray _markers ){CV_INSTRUMENT_REGION();// Labels for pixelsconst int IN_QUEUE = -2; // Pixel visitedconst int WSHED = -1; // Pixel belongs to watershed// possible bit values = 2^8const int NQ = 256;Mat src = _src.getMat(), dst = _markers.getMat();Size size = src.size();// Vector of every created nodestd::vector<WSNode> storage;int free_node = 0, node;// Priority queue of queues of nodes// from high priority (0) to low priority (255)WSQueue q[NQ];// Non-empty queue with highest priorityint active_queue;int i, j;// Color differencesint db, dg, dr;int subs_tab[513];// MAX(a,b) = b + MAX(a-b,0)#define ws_max(a,b) ((b) + subs_tab[(a)-(b)+NQ])// MIN(a,b) = a - MAX(a-b,0)#define ws_min(a,b) ((a) - subs_tab[(a)-(b)+NQ])// Create a new node with offsets mofs and iofs in queue idx#define ws_push(idx,mofs,iofs)\{ \if( !free_node )\free_node = allocWSNodes( storage );\node = free_node; \free_node = storage[free_node].next;\storage[node].next = 0; \storage[node].mask_ofs = mofs;\storage[node].img_ofs = iofs; \if( q[idx].last ) \storage[q[idx].last].next=node; \else \q[idx].first = node; \q[idx].last = node; \}// Get next node from queue idx#define ws_pop(idx,mofs,iofs) \{ \node = q[idx].first;\q[idx].first = storage[node].next; \if( !storage[node].next ) \q[idx].last = 0;\storage[node].next = free_node;\free_node = node; \mofs = storage[node].mask_ofs;\iofs = storage[node].img_ofs; \}// Get highest absolute channel difference in diff#define c_diff(ptr1,ptr2,diff) \{\db = std::abs((ptr1)[0] - (ptr2)[0]);\dg = std::abs((ptr1)[1] - (ptr2)[1]);\dr = std::abs((ptr1)[2] - (ptr2)[2]);\diff = ws_max(db,dg);\diff = ws_max(diff,dr); \assert( 0 <= diff && diff <= 255 ); \}CV_Assert( src.type() == CV_8UC3 && dst.type() == CV_32SC1 );CV_Assert( src.size() == dst.size() );// Current pixel in input imageconst uchar* img = src.ptr();// Step size to next row in input imageint istep = int(src.step/sizeof(img[0]));// Current pixel in mask imageint* mask = dst.ptr<int>();// Step size to next row in mask imageint mstep = int(dst.step / sizeof(mask[0]));for( i = 0; i < 256; i++ )subs_tab[i] = 0;for( i = 256; i <= 512; i++ )subs_tab[i] = i - 256;// draw a pixel-wide border of dummy "watershed" (i.e. boundary) pixelsfor( j = 0; j < size.width; j++ )mask[j] = mask[j + mstep*(size.height-1)] = WSHED;// initial phase: put all the neighbor pixels of each marker to the ordered queue -// determine the initial boundaries of the basinsfor( i = 1; i < size.height-1; i++ ){img += istep; mask += mstep;mask[0] = mask[size.width-1] = WSHED; // boundary pixelsfor( j = 1; j < size.width-1; j++ ){int* m = mask + j;if( m[0] < 0 ) m[0] = 0;if( m[0] == 0 && (m[-1] > 0 || m[1] > 0 || m[-mstep] > 0 || m[mstep] > 0) ){// Find smallest difference to adjacent markersconst uchar* ptr = img + j*3;int idx = 256, t;if( m[-1] > 0 )c_diff( ptr, ptr - 3, idx );if( m[1] > 0 ){c_diff( ptr, ptr + 3, t );idx = ws_min( idx, t );}if( m[-mstep] > 0 ){c_diff( ptr, ptr - istep, t );idx = ws_min( idx, t );}if( m[mstep] > 0 ){c_diff( ptr, ptr + istep, t );idx = ws_min( idx, t );}// Add to according queueassert( 0 <= idx && idx <= 255 );ws_push( idx, i*mstep + j, i*istep + j*3 );m[0] = IN_QUEUE;}}}// find the first non-empty queuefor( i = 0; i < NQ; i++ )if( q[i].first )break;// if there is no markers, exit immediatelyif( i == NQ )return;active_queue = i;img = src.ptr();mask = dst.ptr<int>();// recursively fill the basinsfor(;;){int mofs, iofs;int lab = 0, t;int* m;const uchar* ptr;// Get non-empty queue with highest priority// Exit condition: empty priority queueif( q[active_queue].first == 0 ){for( i = active_queue+1; i < NQ; i++ )if( q[i].first )break;if( i == NQ )break;active_queue = i;}// Get next nodews_pop( active_queue, mofs, iofs );// Calculate pointer to current pixel in input and marker imagem = mask + mofs;ptr = img + iofs;// Check surrounding pixels for labels// to determine label for current pixelt = m[-1]; // Leftif( t > 0 ) lab = t;t = m[1]; // Rightif( t > 0 ){if( lab == 0 ) lab = t;else if( t != lab ) lab = WSHED;}t = m[-mstep]; // Topif( t > 0 ){if( lab == 0 ) lab = t;else if( t != lab ) lab = WSHED;}t = m[mstep]; // Bottomif( t > 0 ){if( lab == 0 ) lab = t;else if( t != lab ) lab = WSHED;}// Set label to current pixel in marker imageassert( lab != 0 );m[0] = lab;if( lab == WSHED )continue;// Add adjacent, unlabeled pixels to corresponding queueif( m[-1] == 0 ){c_diff( ptr, ptr - 3, t );ws_push( t, mofs - 1, iofs - 3 );active_queue = ws_min( active_queue, t );m[-1] = IN_QUEUE;}if( m[1] == 0 ){c_diff( ptr, ptr + 3, t );ws_push( t, mofs + 1, iofs + 3 );active_queue = ws_min( active_queue, t );m[1] = IN_QUEUE;}if( m[-mstep] == 0 ){c_diff( ptr, ptr - istep, t );ws_push( t, mofs - mstep, iofs - istep );active_queue = ws_min( active_queue, t );m[-mstep] = IN_QUEUE;}if( m[mstep] == 0 ){c_diff( ptr, ptr + istep, t );ws_push( t, mofs + mstep, iofs + istep );active_queue = ws_min( active_queue, t );m[mstep] = IN_QUEUE;}}}

四、官方参考代码1

import numpy as npimport cv2 as cvfrom matplotlib import pyplot as pltimg = cv.imread('C:\\Users\\xiaomao\\Desktop\\water_coins.jpg')gray = cv.cvtColor(img,cv.COLOR_BGR2GRAY)ret, thresh = cv.threshold(gray,0,255,cv.THRESH_BINARY_INV+cv.THRESH_OTSU)# noise removalkernel = np.ones((3,3),np.uint8)opening = cv.morphologyEx(thresh,cv.MORPH_OPEN,kernel, iterations = 2)# sure background areasure_bg = cv.dilate(opening,kernel,iterations=3)# Finding sure foreground areadist_transform = cv.distanceTransform(opening,cv.DIST_L2,5)ret, sure_fg = cv.threshold(dist_transform,0.7*dist_transform.max(),255,0)# Finding unknown regionsure_fg = np.uint8(sure_fg)unknown = cv.subtract(sure_bg,sure_fg)# Marker labellingret, markers = cv.connectedComponents(sure_fg)# Add one to all labels so that sure background is not 0, but 1markers = markers+1# Now, mark the region of unknown with zeromarkers[unknown==255] = 0markers = cv.watershed(img,markers)img[markers == -1] = [255,0,0]cv.imshow("1", img)cv.waitKey(0)

五、官方参考代码2

from __future__ import print_functionimport cv2 as cvimport numpy as npimport random as rngrng.seed(12345)src = cv.imread('C:\\Users\\zyh\\Desktop\\card.png')# Show source imagecv.imshow('Source Image', src)src[np.all(src == 255, axis=2)] = 0# Show output imagecv.imshow('Black Background Image', src)kernel = np.array([[1, 1, 1], [1, -8, 1], [1, 1, 1]], dtype=np.float32)# do the laplacian filtering as it is# well, we need to convert everything in something more deeper then CV_8U# because the kernel has some negative values,# and we can expect in general to have a Laplacian image with negative values# BUT a 8bits unsigned int (the one we are working with) can contain values from 0 to 255# so the possible negative number will be truncatedimgLaplacian = cv.filter2D(src, cv.CV_32F, kernel)sharp = np.float32(src)imgResult = sharp - imgLaplacian# convert back to 8bits gray scaleimgResult = np.clip(imgResult, 0, 255)imgResult = imgResult.astype('uint8')imgLaplacian = np.clip(imgLaplacian, 0, 255)imgLaplacian = np.uint8(imgLaplacian)#cv.imshow('Laplace Filtered Image', imgLaplacian)cv.imshow('New Sharped Image', imgResult)bw = cv.cvtColor(imgResult, cv.COLOR_BGR2GRAY)_, bw = cv.threshold(bw, 40, 255, cv.THRESH_BINARY | cv.THRESH_OTSU)cv.imshow('Binary Image', bw)dist = cv.distanceTransform(bw, cv.DIST_L2, 3)# Normalize the distance image for range = {0.0, 1.0}# so we can visualize and threshold itcv.normalize(dist, dist, 0, 1.0, cv.NORM_MINMAX)cv.imshow('Distance Transform Image', dist)_, dist = cv.threshold(dist, 0.4, 1.0, cv.THRESH_BINARY)# Dilate a bit the dist imagekernel1 = np.ones((3,3), dtype=np.uint8)dist = cv.dilate(dist, kernel1)cv.imshow('Peaks', dist)dist_8u = dist.astype('uint8')# Find total markerscontours, _ = cv.findContours(dist_8u, cv.RETR_EXTERNAL, cv.CHAIN_APPROX_SIMPLE)# Create the marker image for the watershed algorithmmarkers = np.zeros(dist.shape, dtype=np.int32)# Draw the foreground markersfor i in range(len(contours)):cv.drawContours(markers, contours, i, (i+1), -1)# Draw the background markercv.circle(markers, (5,5), 3, (255,255,255), -1)#cv.imshow('Markers', markers*10000)cv.watershed(imgResult, markers)#mark = np.zeros(markers.shape, dtype=np.uint8)mark = markers.astype('uint8')mark = cv.bitwise_not(mark)# uncomment this if you want to see how the mark# image looks like at that point#cv.imshow('Markers_v2', mark)# Generate random colorscolors = []for contour in contours:colors.append((rng.randint(0,256), rng.randint(0,256), rng.randint(0,256)))# Create the result imagedst = np.zeros((markers.shape[0], markers.shape[1], 3), dtype=np.uint8)# Fill labeled objects with random colorsfor i in range(markers.shape[0]):for j in range(markers.shape[1]):index = markers[i,j]if index > 0 and index <= len(contours):dst[i,j,:] = colors[index-1]# Visualize the final imagecv.imshow('Final Result', dst)cv.waitKey()

六、其它参考代码3

# import the necessary packagesfrom skimage.feature import peak_local_maxfrom skimage.morphology import watershedfrom scipy import ndimageimport numpy as npimport imutilsimport cv2# construct the argument parse and parse the arguments# load the image and perform pyramid mean shift filtering# to aid the thresholding stepimage = cv2.imread('C:\\Users\\xiaomao\\Desktop\\234.png')shifted = cv2.pyrMeanShiftFiltering(image, 21, 51)cv2.imshow("Input", image)# convert the mean shift image to grayscale, then apply# Otsu's thresholdinggray = cv2.cvtColor(shifted, cv2.COLOR_BGR2GRAY)thresh = cv2.threshold(gray, 0, 255,cv2.THRESH_BINARY | cv2.THRESH_OTSU)[1]cv2.imshow("Thresh", thresh)# compute the exact Euclidean distance from every binary# pixel to the nearest zero pixel, then find peaks in this# distance mapD = ndimage.distance_transform_edt(thresh)localMax = peak_local_max(D, indices=False, min_distance=20,labels=thresh)# perform a connected component analysis on the local peaks,# using 8-connectivity, then appy the Watershed algorithmmarkers = ndimage.label(localMax, structure=np.ones((3, 3)))[0]labels = watershed(-D, markers, mask=thresh)print("[INFO] {} unique segments found".format(len(np.unique(labels)) - 1))# loop over the unique labels returned by the Watershed# algorithmfor label in np.unique(labels):# if the label is zero, we are examining the 'background'# so simply ignore itif label == 0:continue# otherwise, allocate memory for the label region and draw# it on the maskmask = np.zeros(gray.shape, dtype="uint8")mask[labels == label] = 255# detect contours in the mask and grab the largest onecnts = cv2.findContours(mask.copy(), cv2.RETR_EXTERNAL,cv2.CHAIN_APPROX_SIMPLE)cnts = imutils.grab_contours(cnts)c = max(cnts, key=cv2.contourArea)# draw a circle enclosing the object((x, y), r) = cv2.minEnclosingCircle(c)cv2.circle(image, (int(x), int(y)), int(r), (0, 255, 0), 2)cv2.putText(image, "#{}".format(label), (int(x) - 10, int(y)),cv2.FONT_HERSHEY_SIMPLEX, 0.6, (0, 0, 255), 2)# show the output imagecv2.imshow("Output", image)cv2.waitKey(0)