数字图像处理综合练习——水瓶水位线合格检测

马上就要转到学习深度学习的主干线了，这也是大势所趋，但不能忘本，传统图像处理的知识也是非常重要的，特此记录一下之前学习时做过的小练习。

整个项目的资源放在：水瓶水位线合格检测

项目需求

题目来源于冈萨雷斯《数字图像处理》第11章练习题11.38，最终要解决的就是判断瓶中水量是否达到液位标准，液位标准就是瓶颈底部和肩部之间的中点。肩部是瓶子侧面与瓶子倾斜部分的交点。所以要解决的问题主要有以下几点：

• 要能正确识别到肩部和颈部关键点，以至于动态获得每个瓶子的液位标准线。
• 水瓶液位线并非直线，而是曲线，通过图像处理技术定位到液位区域提取曲线并求平均值来代替液位线。
• 正确处理部分瓶身的情况。
• 考虑系统的鲁棒性，尽可能的减少使用先验知识。

这个练习很简单，整个过程只使用VS2019+OpenCV4.20，有兴趣的小伙伴可以去尝试一下。

一、图像预处理

1. 滤波与阈值化

图像读入时是三通道的，需要将原图进行灰度化，然后进行高斯滤波，图像打光方式很明显是采用的背光，前景和背景灰度差很明显，只需要指定一个阈值简单的阈值化即可。

// 灰度化并高斯平滑cv::Mat GaussianSrc, GraySrc;cv::cvtColor(ColorSrc, GraySrc, cv::COLOR_BGR2GRAY);cv::GaussianBlur(GraySrc, GaussianSrc, cv::Size(7, 7), 0, 0);cv::Mat binarySrc;cv::threshold(GaussianSrc, binarySrc, 50, 255, cv::THRESH_BINARY);

2. 垂直投影获取投影分割点

通过对二值化图像进行水平/垂直投影分割有明显间隔的区域是图像处理中的常用技术，获取投影直方图很简单，但是要获得投影分割点并正确分割想要的区域或许并不是那么容易的事了，但多数情况下间隔明显，困难就减少了很多。

为了后面可以直接分割出一个个水瓶区域，并正确分割出部分瓶身，需要定义一个带属性的分割点类。

//带有属性的分割点
class BottleSegment
{
public:BottleSegment() = default;BottleSegment(int x_cord, bool state) : Segment_x(x_cord), SE_state(state) { }//垂直分割点int Segment_x;//true表示为瓶身开始点，false表示为瓶身结束点bool SE_state = false;
};

当SE_state = true表示该垂直分割点为瓶子的左侧，反之为瓶子的右侧。分割点取的是下降\上升趋势的中点，上升趋势对应瓶子左侧true，反之对应瓶子右侧false，原理很简单，具体可以参见代码。

在垂直投影处理过程中还有一个小技巧，当直方图的值大于某个值时，将其截断，那么非间隔区域的直方图值就不会影响到算法去寻找分割点，这在某些应用中还是挺好用的。

/** @brief 计算二值化原图的垂直投影，得到垂直投影直方图和带属性的瓶身垂直分割点@param binaryMat 输入原图的二值化图像
@param VerProj 输出垂直投影直方图
@param VerSegment_X 输出带属性的瓶身垂直分割点
*/
void calcVerSegment(const cv::Mat& binaryMat, std::vector<int>& VerProj, std::vector<BottleSegment>& VerSegment_X)
{CV_Assert(binaryMat.type() == CV_8UC1);const int MaxCount = 150;VerProj.resize(binaryMat.cols, 0);for (int x = 0; x < binaryMat.cols; ++x){for (int y = 0; y < binaryMat.rows; ++y){if (binaryMat.ptr<uchar>(y)[x] == 255 && VerProj[x] < MaxCount)VerProj[x]++;}}// 根据投影直方图得到分割点const int distThres = 10;for (int i = 1; i < VerProj.size() - 1; ++i){if (VerProj[i + 1] - VerProj[i - 1] > distThres && VerProj[i - 1] == 0)			//后一点大于前一点(瓶身起点){int StartIndex = i;while (VerProj[++i] != MaxCount);int Segment = (StartIndex + i) / 2;VerSegment_X.push_back(BottleSegment(Segment, true));}if (VerProj[i - 1] - VerProj[i + 1] > distThres && VerProj[i - 1] == MaxCount)	//前一点大于后一点(瓶身落点){int StartIndex = i;while (VerProj[++i] != 0);int Segment = (StartIndex + i) / 2;VerSegment_X.push_back(BottleSegment(Segment, false));}}
}

将分割点展示一下：效果很好，基本上是与肩部相切。

二、分离水瓶并求取肩部和颈部关键点

在根据垂直分割点分离出各个水瓶之前先建立一个瓶子的类，属性包括瓶身灰度图像Roi_image，瓶身水平投影直方图HorizontalProj_hist等。完整的类定义如下：

//瓶子类
class Bottle
{
friend void ResultVisualization(const Bottle& bottle, cv::Mat& inputOutput_ColorSrc, int shoulder, int neck, int water_line);public:Bottle() = default;Bottle(const cv::Mat& roi, bool ispart, cv::Point of);//将瓶身水平投影直方图进行可视化cv::Mat Horizon_HistMat() const{cv::Mat showMat = cv::Mat::zeros(Roi_image.size(), CV_8UC1);for (int i = 0; i < HorizontalProj_hist.size(); ++i){int x = HorizontalProj_hist[i] > Roi_image.cols ? Roi_image.cols : HorizontalProj_hist[i];cv::line(showMat, cv::Point(0, i), cv::Point(x, i), cv::Scalar::all(255), 1, 8);}return showMat;}//计算瓶身颈部关键点int neck_keyPoint() const;//计算肩部关键点int shoulder_keyPoint() const;//对瓶身做边缘检测以提取水位线图像(返回8位灰度图像)cv::Mat Water_level_line() const;//提取水位线区域掩模cv::Mat Water_line_region_mask() const;/成员属性定义
private://是否为部分瓶身bool is_partial = false;		//原图的灰度瓶身图像cv::Mat Roi_image;//瓶身水平投影直方图std::vector<int> HorizontalProj_hist;		//局部瓶身在原图中的偏移量(相对于左上角)cv::Point offset;
};

1. 分离水瓶形成一个个独立区域

根据提取的带属性的分割点来分离瓶身，处于第一个和最后一个的分割点要判断该水瓶是否是部分瓶身，分离结果保存在一个容器中。

/** @brief 根据投影分割点从原图中分离出瓶身区域@param GaussianSrc 输入原图的灰度图像
@param VerSegment_X 输入带属性的瓶身垂直分割点数组
@param bottles 输出瓶身区域数组
*/
void SplitBottles(const cv::Mat& GaussianSrc, std::vector<BottleSegment> VerProjSegment_X, std::vector<Bottle>& bottles)
{CV_Assert(GaussianSrc.type() == CV_8UC1);for (auto bottleSeg_it = VerProjSegment_X.begin(); bottleSeg_it != VerProjSegment_X.end(); bottleSeg_it++){if (bottleSeg_it == VerProjSegment_X.begin() && (*bottleSeg_it).SE_state == false)	//如果第一个分割点就是瓶身落点{cv::Mat temp = GaussianSrc(cv::Rect(0, 0, (*bottleSeg_it).Segment_x + 1, GaussianSrc.rows)).clone();bottles.push_back(Bottle(temp, true, cv::Point(0, 0)));}else if (bottleSeg_it == VerProjSegment_X.end() - 1 && (*bottleSeg_it).SE_state == true)	//如果最后一个分割点是瓶身起点{int x = (*bottleSeg_it).Segment_x;cv::Mat temp = GaussianSrc(cv::Rect(x, 0, (GaussianSrc.cols - x), GaussianSrc.rows)).clone();bottles.push_back(Bottle(temp, true, cv::Point(x, 0)));}else if ((*bottleSeg_it).SE_state == true && (*(bottleSeg_it + 1)).SE_state == false)	//前一个是起点后一个是落点{int Start_x = (*bottleSeg_it).Segment_x;int end_x = (*(bottleSeg_it + 1)).Segment_x;cv::Mat temp = GaussianSrc(cv::Rect(Start_x, 0, (end_x - Start_x), GaussianSrc.rows)).clone();bottles.push_back(Bottle(temp, false, cv::Point(Start_x, 0)));//已经完整分割一个瓶身，将迭代器更新到落点bottleSeg_it++;}}
}

分离的结果如下图所示，同时要对水瓶进行水平投影得到水平投影直方图，在本文解决方案中肩部和颈部关键点是通过水瓶的投影直方图来获取的，为了保证直方图的平滑性，还要将直方图进行一个平滑处理。

2. 肩部和颈部关键点提取

最开始尝试过几种方案，尝试过在水瓶图中拟合肩部斜直线和竖直直线的交点来得到关键点，但发现肩部斜直线并不是很直，存在误差，然后又试过角点检测提取，但是原图分辨率不高，角点提取的效果不是很好，鲁棒性不强。最后发现水平投影之后，肩部和颈部的特征就很明显了。

注意水平直方图中的凹坑区域就是颈部区域，只要定位到颈部区域上升趋势的开始点，和上升趋势的结束点，即为颈部关键点和肩部关键点，此方法简单易行，鲁棒性强，无论是否是部分瓶身均可有效提取(左侧部分瓶身出现问题，之后分析问题出现原因)。

肩部和颈部关键点提取代码如下：

/** @brief 通过瓶身水平投影直方图计算颈部关键点@note 先通过找最小值定位到颈部区域，然后从该区域遍历如果出现连续上升的三个点，则定位为颈部关键点
*/
int Bottle::neck_keyPoint() const
{int i = 0;while (HorizontalProj_hist[i++] < 10);int limit = Roi_image.rows / 2;				//遍历限制，不需要遍历整个直方图int neck_minValue = Roi_image.cols;int neck_minIndex = 0;//找到颈部最小值for (int k = i + 10; k < limit; ++k){if (HorizontalProj_hist[k] < neck_minValue){neck_minValue = HorizontalProj_hist[k];neck_minIndex = k;}}//上升次数计数，当满足3次则为颈部关键点const int INCREASE_TIME = 3;int increase = 0;//颈部关键点y坐标int keyPoint = 0;for (int j = neck_minIndex; j < limit; ++j){if (increase >= INCREASE_TIME)return keyPoint;else if (HorizontalProj_hist[j] < HorizontalProj_hist[j + 2]){if (increase == 0)keyPoint = j + 2;increase++;}else{increase = 0;keyPoint = 0;}}return keyPoint;
}/** @brief 通过瓶身水平投影直方图计算瓶身肩部关键点@note 直方图值等于瓶身宽度的最小索引值即为肩部关键点
*/
int Bottle::shoulder_keyPoint() const
{auto shoulder_it = std::find(HorizontalProj_hist.begin(), HorizontalProj_hist.end(), Roi_image.cols - 1);int temp = (int)(shoulder_it - HorizontalProj_hist.begin());int keyPoint = (shoulder_it != HorizontalProj_hist.end() ? temp : 0);return keyPoint;
}

三、水瓶液位曲线提取

1. 图像增强

为了更好的凸显出边缘细节以提取液位边缘，需要先对图像进行增强，图像增强之后必然会凸显噪点，所以还要再进行滤波，也有些视觉应用中是先滤波后增强，我在比较了两种效果之后选择了先增强后滤波。

• 伽马变换 + 高斯滤波(双边滤波)
• 高斯滤波 + 伽马变换

2. Y方向负边缘提取

边缘提取就是一个求导的过程，离散的求导就是差分，要提取水瓶液位1只需要进行Y方向的差分即可，差分的值有正有负，正数代表从黑到白的边缘，负数代表从白到黑的边缘。显然我们只需要从白到黑的液位边缘，所以需要将求得的边缘抹去正值，保留负数然后在归一化到0-255并转换到灰度图像。结合图像增强，完整的代码如下：

//伽马变换
void GammaTrans(const cv::Mat& inputimage, cv::Mat& outputimage, const float val)
{CV_Assert(inputimage.channels() == 1);cv::Mat normalImage = inputimage.clone();normalImage.convertTo(normalImage, CV_32FC1);cv::Mat TempImage = cv::Mat::zeros(inputimage.size(), CV_32FC1);cv::normalize(normalImage, normalImage, 1, 0, cv::NORM_MINMAX);cv::pow(normalImage, val, TempImage);cv::normalize(TempImage, TempImage, 255, 0, cv::NORM_MINMAX);TempImage.convertTo(TempImage, CV_8UC1);outputimage = TempImage.clone();
}//舍弃图像中的正值，只保留负值
void PositiveToZero(const cv::Mat& src, cv::Mat& Output)
{CV_Assert(src.type() == CV_16SC1);Output = src.clone();for (int i = 0; i < src.rows; ++i){for (int j = 0; j < src.cols; ++j){Output.ptr<short>(i)[j] = -(src.ptr<short>(i)[j] > 0 ? 0 : src.ptr<short>(i)[j]);}}
}/** @brief 瓶身边缘图求取以提取水位线区域@note 在边缘图求取，先增强再模糊效果好于先模糊再增强，获得的边缘图要先抹去正数然后再归一化
*/
cv::Mat Bottle::Water_level_line() const
{cv::Mat gamma;GammaTrans(Roi_image, gamma, 3);cv::Mat blur_roi_image;cv::GaussianBlur(gamma, blur_roi_image, cv::Size(7, 7), 0, 0);cv::Mat Sobel_Y;cv::Sobel(blur_roi_image, Sobel_Y, CV_16SC1, 0, 1, 3, 1, 0);//抹去正数，将图片归一化到0-255，并转化为CV_8UC1PositiveToZero(Sobel_Y, Sobel_Y);cv::normalize(Sobel_Y, Sobel_Y, 0, 255, cv::NORM_MINMAX);Sobel_Y.convertTo(Sobel_Y, CV_8UC1);return Sobel_Y;
}

提取到的边缘图如下所示：液位区域已经很明显了

3. 液位曲线提取

要提取液位边缘曲线根据精度要求有两种技术路线：

1. 边缘提取——阈值化——骨架化
2. 边缘提取——非极大抑制处理——阈值化——骨架化

此次练习简单处理，选择了第一条路线，首先对边缘图进行OTSU阈值处理，形成二值化图像，从上面边缘图中可以看出，阈值化后液位区域所占的面积最大，所以我们对二值化后的图进行轮廓提取，保留轮廓面积最大的轮廓作为掩模，也就是液位区域，对于液位区域，我们利用形态学细化提取区域骨架，形态学细化之前，也可执行形态学闭运算操作以平滑液位区域，形态学细化的细节和代码可以参考我之前的记录。

OpenCV实现形态学细化
实现的效果如下图所示：

提取液位区域部分的代码如下：

/** @brief 提取瓶身水位线区域掩模@note 先将边缘图阈值化，然后提取最大轮廓即为水位线区域
*/
cv::Mat Bottle::Water_line_region_mask() const
{cv::Mat Sobel_Y = Water_level_line();//OTSU阈值化cv::Mat Sobel_Y_thres;cv::threshold(Sobel_Y, Sobel_Y_thres, 0, 255, cv::THRESH_OTSU);//形态学闭运算cv::Mat element = cv::getStructuringElement(cv::MORPH_RECT, cv::Size(5, 5));cv::morphologyEx(Sobel_Y_thres, Sobel_Y_thres, cv::MORPH_CLOSE, element, cv::Point(-1, -1), 1);//提取水位线区域轮廓std::vector< std::vector<cv::Point> > contours;cv::findContours(Sobel_Y_thres, contours, cv::RETR_EXTERNAL, cv::CHAIN_APPROX_NONE);//提取最大轮廓面积所对应的区域int max_area = 0;int max_area_index = 0;for (int i = 0; i < contours.size(); ++i){int area = cv::contourArea(contours[i]);if (area > max_area){max_area = area;max_area_index = i;}}cv::Mat max_area_mask = cv::Mat::zeros(Roi_image.size(), CV_8UC1);cv::drawContours(max_area_mask, contours, max_area_index, cv::Scalar::all(255), cv::FILLED, 8);return max_area_mask;
}

四、水瓶液位合格判断并可视化

现在只差最后一步，判断液位线是否达到标准，液位标准就是肩部和颈部的中点，我们只需要将骨架化后形成的线段点对Y坐标求均值即可拟合液位线，在最终的可视化中，液位标准线用蓝色表示，液位合格时，液位线用绿色表示，反之用红色表示。

那么，最终的检测结果展示如下：

主程序代码如下：

int main(int argc, char** argv)
{std::string path = "F:\\NoteImage\\bottles.tif";cv::Mat ColorSrc = cv::imread(path, cv::IMREAD_COLOR);if (!ColorSrc.data) {std::cout << "Could not open or find the image" << std::endl;return -1;}// 灰度化并高斯平滑cv::Mat GaussianSrc, GraySrc;cv::cvtColor(ColorSrc, GraySrc, cv::COLOR_BGR2GRAY);cv::GaussianBlur(GraySrc, GaussianSrc, cv::Size(7, 7), 0, 0);cv::Mat binarySrc;cv::threshold(GaussianSrc, binarySrc, 50, 255, cv::THRESH_BINARY);//垂直投影直方图std::vector<int> VerProj(binarySrc.cols);//带有属性的垂直投影分割点std::vector<BottleSegment> VerProjSegment_X;				calcVerSegment(binarySrc, VerProj, VerProjSegment_X);//将每个瓶子单独分离开来std::vector<Bottle> bottles;SplitBottles(GraySrc, VerProjSegment_X, bottles);for (int i = 0; i < bottles.size(); ++i){int bottleNeckPoint = bottles[i].neck_keyPoint();int bottleShoulderPoint = bottles[i].shoulder_keyPoint();std::cout << "neckpoint " << bottleNeckPoint << "  shoulderpoint " << bottleShoulderPoint << std::endl;cv::Mat max_area_mask = bottles[i].Water_line_region_mask();//形态学细化提取区域骨架cv::Mat skeleton;Morph_Thinning(max_area_mask, skeleton, 100);//简单的利用平均值来估计水位线int water_line_row = 0;int total = 0;for (int y = 0; y < skeleton.rows; ++y){for (int x = 0; x < skeleton.cols; ++x){if (skeleton.ptr<uchar>(y)[x] == 255){water_line_row += y;total++;}}}water_line_row /= total;std::cout << "water line value = " << water_line_row << std::endl;//结果可视化ResultVisualization(bottles[i], ColorSrc, bottleShoulderPoint, bottleNeckPoint, water_line_row);}cv::imshow("src", ColorSrc);cv::waitKey(0);return 0;
}

结果分析

仔细观察会发现，左侧那个部分瓶身出现了很大的误差，其实问题就出现在，这个瓶子在肩部区域并不是对称的，左边肩部有一个缺口导致肩部下移，而该瓶身左侧缺失，所以检测到的肩部关键点为右侧那个偏上的关键点。
如果要处理，那就只能特殊处理，那么鲁棒的算法就会失去他所含的美感，但所幸这只是一个项目练习，那我们就让他将错就错吧！