引言
图像分割中，采用腐蚀、膨胀等形态学操作一般可以对图像进行一定效果的分割。腐蚀采用的kernel往往决定了分割的效果，核太小则难以分割开来，核太大则会将某些部分腐蚀掉同时对原本图像形状产生较大的影响。所以，换了一个思路想，从外边分割不如从里面扩张。所以分割问题的重点放在了图像质心的查找。
图像处理中，针对二值化分割图可以使用距离变幻得到一张灰度图，而亮度中心则是质心，查找质心既可以实现简单的分割。Matlab中查找局部峰值的函数有imregionalmax，而在python中模块skimage中有peak_local_max函数同样实现了对灰度图的峰值查找。但是，在C/C++中，笔者还没有找到相关的图像处理库有这一方面的功能，所以只能对着论文写了。
1、距离变换
距离变换在图像处理应用十分广泛，特别在图像分割方面。简单来说，距离变换是针对二值化图像，其实就是计算非零点到最近零点的距离。所以，图像经过距离变换后，分割图像中心则呈现出峰值。距离的定义一般有以下3种类型：
OpenCV中，则有7中类型的距离定义，除此之外用于还可以使用自己定义的距离。（所以，有无限制种吧）
/** Distance types for Distance Transform and M-estimators */ enum { CV_DIST_USER =-1, /**< User defined distance */ CV_DIST_L1 =1, /**< distance = |x1-x2| + |y1-y2| */ CV_DIST_L2 =2, /**< the simple euclidean distance */ CV_DIST_C =3, /**< distance = max(|x1-x2|,|y1-y2|) */ CV_DIST_L12 =4, /**< L1-L2 metric: distance = 2(sqrt(1+x*x/2) - 1)) */ CV_DIST_FAIR =5, /**< distance = c^2(|x|/c-log(1+|x|/c)), c = 1.3998 */ CV_DIST_WELSCH =6, /**< distance = c^2/2(1-exp(-(x/c)^2)), c = 2.9846 */ CV_DIST_HUBER =7 /**< distance = |x|<c ? x^2/2 : c(|x|-c/2), c=1.345 */ };
以下是针对两个圆形粘连情况下采用CV_DIST_L2（常用欧式距离）作为距离进行距离变换：
2、形态学灰度重建（Morphological Grayscale Reconstruction）
一般对于局部峰值极值点的查找很容易联想到使用窗口判断的方法。但是，窗口的大小决定着算法的好坏。更致命的是，窗口这种方法很难应对出现连续峰值的问题（峰值平坦？不太会描述）。所以这种方法很难解决峰值查找的问题。于是，笔者在各种论坛中查找各种峰值查找的方法，最终还是回到Matlab参考文献。
Matlab中的imregionalmax函数是使用形态学灰度重建的方法实现了局部极大值的查找，具体参考论文：Morphological Grayscale Reconstruction in Image Analysis: Applications and Effcient Algorithms。这里只进行简单描述。
从上图可知，对于图像I，用其偏置I-1进行grayscale reconstruction灰度重建得到重建后的图像。然后在使用原图与重建后的图像相减即可得到局部极大值部分。上图也可以解释Matlab中imregionalmax等于原图与其imreconstruct相减。接下来，则是如何使用算法实现这个过程呢？论文末尾处已经给出了许多算法实现的伪代码，其中包括并行算法、串行算法等等。本文使用队列对算法进行实现，该算法在论文中已经给出伪代码：
特别说明，笔者取J=I-1，R(J)是遍历每个点取八邻近最大值来初始化队列。NG(p)也是取得p点的八邻近。部分实现代码如下：
Mat findpeeks(Mat I,Mat J) { queue<Point> que; /*最大值滤波*/ J = regionalmaximaImg(J); /*对满足条件的点入队列，初始化队列*/ for(int i=0;i<I.rows;i++) for(int j=0;j<I.cols;j++) if(J.at<uchar>(i,j) && isexistNBzero(J,Point(j,i))) que.push(Point(j,i)); /*循环队列，直至队列为空*/ while(que.empty()) { Point p = que.front(); que.pop(); for(int i=0;i<NEGHTBOR;i++) { /*get the neighbor point*/ Point q; q.x = p.x+x_offset[i]; q.y = p.y+y_offset[i]; if(q.x>I.cols-1) q.x = I.cols-1; if(q.x<0) q.x = 0; if(q.y>I.rows-1) q.y = I.rows-1; if(q.y<0) q.y = 0; if(J.at<uchar>(p.y,p.x)<J.at<uchar>(q.y,q.x) && I.at<uchar>(q.y,q.x)!=I.at<uchar>(q.y,q.x)) { J.at<uchar>(p.y,p.x)=min(J.at<uchar>(p.y,p.x), I.at<uchar>(q.y,q.x)); que.push(p); } } } return (I-J)>0; }
实现效果如下：
由图像可见峰值点可以找出来，但是不一定为质心，因为连接处也可能存在峰值，这里需要后续的优化。需要注意的是峰值点不是单点，可能出现连续或间断点；笔者认为后续使用区域生长、主动轮廓模型可以完成实现分割。
结语
形态学灰度重建提供了一种峰值查找的解决方法，论文比较长，我也看得不太懂，由于前阶段毕业设计急着需要用就直接按照后面的伪代码实现了，后续有时间再将毕设后续的内容总结上来。
另外：该算法的并行算法实现详见我的简书：OpenCV CUDA 图像峰值查找，我总觉得搬过来不太好。
本文源码下载：01_findpeeks
参考文献：
Morphological Grayscale Reconstruction in Image Analysis: Applications and Effcient Algorithms