基于OPENCV分水岭的球团分割

       因为项目需求，需要对某铁矿厂的球团进行粒径检测。采集系统就不详细说了，主要是颗粒运动很快，粒径在8-12mm，范围1米左右，所以既要高像素的相机，又要曝光时间很短的相机，前期拍出来的光照不足，图像很暗。

   

图1 现场图                                                  图2. 球团图

        查了下资料，相似的文献比如水果分拣啊，硬币分拣啊，都是用的分水岭。我也尝试这分割了下。刚开始效果不行，对光滑的表面效果挺好，但是对这种小颗粒，表面不是特光滑的，效果不理想。主要是噪声多，精确度差。后来通过对背景和前景的调整，达到如图4的效果，基本上满足了生产需求。

 

图3. 改进前                                                                                     图4. 改进后

下面简单介绍下过程。

 一. 图像预处理

        因为原图是短曝光+提高增益的来的，比较黑，噪声大。需要对图像进行增强处理，本文用了顶帽变换+中值滤波的方式。

cv::Mat kernel2 = cv::getStructuringElement(cv::MORPH_RECT, cv::Size(12, 12), cv::Point(-1, -1));
cv::morphologyEx(src, src, CV_MOP_TOPHAT, kernel2, cv::Point(-1, -1), 1);//顶帽变换
cv::medianBlur(src, src, 3);

图5. 增强后的球团图 

二. 分水岭算法

      分水岭的基本理论主要是漫水筑坝的过程，可以参考文章1和文章2。这里主要谈谈实现过程。opencv的分水岭算法主要依靠marks来标记待分割的边界。一方面减少噪声，一方面提高巡边速度。因此，marks的设置就很重要。

2.1. 构建marks

2.1.1  构造前景

      这部分主要是对球团进行标记。因为球团是层叠的，我们主要把完整的球团标记出来作为前景。这部分我借鉴了网络上其他人的，然后调整了下参数。主要的过程：

（1）二值化，分割目标

（2）开运算，减少噪声

（3）距离变换+二值化

（4）膨胀运算，精简前景

cv::threshold(src, src, 0, 255, cv::THRESH_BINARY|cv::THRESH_OTSU);

cv::Mat old_src = src.clone(); //保存处理前的图，以方便求背景
// 开操作
cv::morphologyEx(src, src, cv::MORPH_OPEN, kernel, cv::Point(-1, -1), 2);
//---------------距离变换-------------------
cv::Mat dist_img;
cv::distanceTransform(src, dist_img, CV_DIST_L2, 3/*distance type*/);
cv::normalize(dist_img, dist_img, 0, 1., cv::NORM_MINMAX);
cv::threshold(dist_img, dist_img, 0.1, 1., CV_THRESH_BINARY);
// 膨胀
cv::Mat kernel1 = cv::Mat::ones(3, 3, CV_8UC1);
cv::morphologyEx(dist_img, dist_img, cv::MORPH_DILATE, kernel1, cv::Point(-1, -1), 2);//
dist_img.convertTo(src, CV_8U);

这部分的可操作性蛮大，主要表现在：

（1）是用开运算还是其他的，我也尝试过先两次膨胀，再1次腐蚀，效果也有差别，在某些图像上效果更好。还有，形态学运算的次数，也会影响结果；

（2）距离变换后的二值化阈值范围：我这里选的是0.1，起点越小，颗粒越小。所以噪声也会多。所以，这里要结合先验条件设置。

图6 二值化                                 图7 距离变换                                   图8 前景

2.1.2  构造背景

        这里的背景，是指暗区无法分割出来的球团。如果不指定背景，在分水岭过程中，算法认为每个前景是挨着的，会把暗区也分到前景去，使边缘过分割。也就是说，设计marks时，把背景作为标识0（也就是待分区），最后的结果不会分割出背景来。因此，有必要把背景指定。

      这里，只需要把二值化后的前景进行膨胀，剩余部分作为背景。操作如下：

cv::Mat bk;//找背景，扩大目标区，压缩背景，使背景更精确（背景为黑色）
cv::morphologyEx(src_old, bk, CV_MOP_DILATE, kernel, cv::Point(-1, -1), 6);

// 对我们已经确定分类的区域（无论是前景还是背景）使用不同的正整数标记，
// 对我们不确定的区域（要进行分水岭查找边界的区域）使用 0 标记。
cv::Mat marks(src_old.rows, src_old.cols, CV_32S, cv::Scalar(0));

for (int i = 0; i < src.rows; i++)
{
	LPBYTE lpBk = bk.ptr<BYTE>(i);
	INT32 *lpMask = marks.ptr<INT32>(i);
	for (int j = 0; j < src.cols; j++)
	{
		if (*(lpBk + j) == 0)//标记背景为1
			*(lpMask + j) = 1;
	}
}

 2.1.3  构造完整的marks

         marks里面已经指定了背景，再融合进前景就行了。这里用的是drawContours函数，把前景轮廓画到marks里。

vector<vector<cv::Point>> contours;
//此时，src为二值图，白色区域为前景
cv::findContours(src, contours, CV_RETR_EXTERNAL, CV_CHAIN_APPROX_SIMPLE);
for (int i = 0; i < contours.size() && i<1000; i++)
{
	double area = cv::contourArea(contours[i]);
	if (area < 200)//剔除小区域
		continue;
    // 标记从2开始，1已经画了背景
	cv::drawContours(marks, contours, static_cast<int>(i), cv::Scalar::all(static_cast<int>(i) + 2), 1, 8);
}

  

图9 背景marks                                                                                             图10 最终的marks

---

       其中，灰色部分标记为1的背景，白色部分为前景（这里只是示意，实际上灰色为1，每个白色封闭区实际上标号大于1，各自也不一样），这两者为注水区，中间黑色部分（标记0）为未知区，也就是待搜索区。

 2.2. 分水岭分割

     下一步很简单了，就是直调用函数了。但是，操作对象可以是前面增强后的图，也可以是原图。但是，增强后的图对比度大了，但是噪点也多了。所以我用的还是原图。

cv::watershed(input_src, marks);

       

       图11 处理后的marks

三. 统计显示  

         这一步marks里就已经把边界用-1，标识出来了。所以marks一定是一个32位整型矩阵。再下一步就是遍历统计了。我这里是遍历每个标识面积（不要遍历标识为1的背景）和区域中心，然后简单的换算一下圆的半径：

                    r =  2 * sqrt( nArea / PI)

typedef struct _ContourInfo
{
	float x;
	float y;
	int nArea;
	_ContourInfo()
	{
		x = y = -1;
		nArea = 0;
	}
}ContourInfo, *pContourInfo;

std::map<int, ContourInfo> indexContour;
for (int i = 0; i < marks.rows; i++)
{
	for (int j = 0; j < marks.cols; j++)
	{
		int index = marks.at<int>(i, j);//区域标记
		if (index != -1 && index <= static_cast<int>(contours.size()))//-1为边界像素,1为背景
		{
			indexContour[index].nArea++;
			indexContour[index].x += j;
			indexContour[index].y += i;
		}
	}
}

std::map<int, ContourInfo>::iterator itr;
for (itr = indexContour.begin(); itr != indexContour.end(); itr++)
{
	if (itr->second.nArea > 0)
	{
		itr->second.x /= itr->second.nArea;
		itr->second.y /= itr->second.nArea;
	}
}
float fRadio = 0;

//----------show--------------
if (bManuInput)
{
	std::map<int, ContourInfo>::iterator itr;
	for (itr = indexContour.begin(); itr != indexContour.end(); itr++)
	{
		if (itr->second.nArea > 100 && itr->first != 1/*1为背景*/)
		{
			fRadio = sqrt(double(itr->second.nArea) / PI) + 2;
			cv::circle(oImage, cv::Point(itr->second.x, itr->second.y), int(fRadio), cv::Vec3b(0, 0, 255), 1);
		}
	}
	cv::imshow("result", oImage);
}

效果：

四.参考文献

轮廓提取试验_元气少女缘结神的博客-CSDN博客

OpenCV---分水岭算法 - 山上有风景 - 博客园

标记控制的分水岭分割- MATLAB & Simulink Example- MathWorks 中国

分水岭算法的详细介绍（附c代码）_一步一个脚印的屌丝的博客-CSDN博客