图像处理（六）：二值化

二值化：设置一个阈值，图像灰度值如果小于这个阈值，就将对应灰度值设置为0，反之，保留。

按照设置阈值的方式：

1、固定阈值：预先设置的固定阈值；

2、自适应阈值：通过某种算法对图像像素进行处理，自动计算出最佳的阈值；

固定阈值就不多说了，太简单了！这里说一下自适应阈值，最经典的要数大津法（OTSU），将图像直方图用某一灰度分成两组，将分割后两组间方差对应的灰度值设置为最佳阈值，具体原理可以参考博客：点击打开链接，该博客也有对应代码，代码写的更加简洁明了。

步骤如下：

1）直方图统计0~255对应像素的数目；

2）将0,1，...，244,255依次设置为阈值，计算分组后两组的均值及方差；

3）找寻最大方差对应的分割阈值即为最佳阈值；

代码如下：

#include <opencv2/opencv.hpp>
#include <iostream>

using namespace cv;
using namespace std;

int main()
{
	Mat src = imread("1.jpg");
	Mat gray_ori;
	cvtColor(src, gray_ori, CV_BGR2GRAY);

	Mat gray, gray_OTSU;
	gray_ori.copyTo(gray);
	gray_ori.copyTo(gray_OTSU);

	// 1、fixed threshold
	int thr = 30;
	int i, j;
	for (i = 0; i < gray.rows - 1; i++)
		for (j = 0; j < gray.cols - 1; j++)
		{
			if (gray.at<uchar>(i, j) < thr)
			{
				gray.at<uchar>(i, j) = 0;
			}
		}



	// 2、OTSU
	// 2.1、直方图统计
	int hist[256] = { 0 };
	for (i = 0; i < gray_ori.rows - 1; i++)
		for (j = 0; j < gray_ori.cols - 1; j++)
		{
			hist[gray_ori.at<uchar>(i, j)]++;

		}


	// 2.2、计算系数
	int s = 0;
	int hist_cum[256] = { 0 };
	float hist_coeff[256] = { 0 };
	for (i = 0; i < 256; i++)   // 累计统计
	{
		s += hist[i];
		hist_cum[i] = s;
		hist_coeff[i] = 1.0 * s / (gray_ori.rows * gray_ori.cols);

		//cout << hist_coeff[i] << endl;
	}

	// 2.3、计算加权值
	s = 0;
	int hist_cum_Add[256] = { 0 };	
	for (i = 0; i < 256; i++)   
	{
		s += i * hist[i];
		hist_cum_Add[i] = s;
		//cout << hist_cum_Add[i] << endl;
	}


	// 2.4、计算均值	
	float hist_mean[256] = { 0 };
	float hist_mean_l[256] = { 0 };
	float hist_mean_r[256] = { 0 };
	for (i = 0; i < 256; i++) 
	{
		if (hist_cum[i] == 0)
		{
			hist_mean_l[i] = 0;

		}
		else
		{
			hist_mean_l[i] = hist_cum_Add[i] / hist_cum[i];
		}
		
		hist_mean_r[i] = (s - hist_cum_Add[i]) / (gray.rows * gray.cols - hist_cum[i]);

		hist_mean[i] = hist_coeff[i] * hist_mean_l[i] + (1 - hist_coeff[i]) * hist_mean_r[i];
	
	}

	// 2.5、计算方差
	float sigma[256] = { 0 };
	for (i = 0; i < 256; i++)   
	{
		sigma[i] = hist_coeff[i] * (hist_mean_l[i] - hist_mean[i]) * (hist_mean_l[i] - hist_mean[i]) +
			(1 - hist_coeff[i]) * (hist_mean_r[i] - hist_mean[i]) * (hist_mean_r[i] - hist_mean[i]);
	}

	// 2.6、获取最佳阈值
	float temp = 0;
	int index = 0;
	for (i = 0; i < 256; i++)   
	{
		if (temp < sigma[i])
		{
			temp = sigma[i];
			index = i;
		}
	}

	cout << "best threshold is: " << index << endl;

	// 进行阈值滤波
	for (i = 0; i < gray_OTSU.rows - 1; i++)
		for (j = 0; j < gray_OTSU.cols - 1; j++)
		{
			if (gray_OTSU.at<uchar>(i, j) < index)
			{
				gray_OTSU.at<uchar>(i, j) = 0;
			}
		}

	imshow("gray", gray);
	imshow("OTSU", gray_OTSU);
	waitKey();

	return 0;





}

运行效果如下：

刚才又看了一种迭代阈值法，具体步骤如下：

1、以图像均值初始化 best_thr；

2、利用best_thr将图像中像素分为两组：Group1和Group2;

3、计算Group1和Group2对应的均值u1和u2；

4、重新选择阈值：best_thr = （u1 + u2）/ 2.0；

重复2~4，直到u1和u2不再发生变化。

代码如下：

#include <opencv2/opencv.hpp>
#include <vector>

using namespace cv;
using namespace std;

int main()
{
	Mat src = imread("1.jpg");
	Mat gray;
	cvtColor(src, gray, CV_BGR2GRAY);

	int i, j;
	double sum = 0;
	double best_thr = 0;
	for (i = 0; i < gray.rows; i++)
		for (j = 0; j < gray.cols; j++)
		{
			sum += gray.at<uchar>(i, j);
		}

	// 以图像像素均值来初始化阈值
	best_thr = sum / (gray.rows * gray.cols);

	vector<Point> Lows;
	vector<Point> Higs;
	bool flag = true;
	float u1 = 0, u2 = 0;
	float s1 = 0, s2 = 0;

	while (flag)
	{
		Lows.clear();
		Higs.clear();

		// 获取分割后的两个部分
		for (i = 0; i < gray.rows; i++)
			for (j = 0; j < gray.cols; j++)
			{
				if (gray.at<uchar>(i, j) < best_thr)
				{
					Lows.push_back(Point(i, j));
				}
				else
				{
					Higs.push_back(Point(i, j));
				}
			} // end for

		s1 = 0;
		u1 = 0;
		for (i = 0; i < Lows.size(); i++)
		{
			s1 += gray.at<uchar>(Lows[i].x, Lows[i].y);
		}
		u1 = s1 / Lows.size();

		s2 = 0;
		u2 = 0;
		for (i = 0; i < Higs.size(); i++)
		{
			s2 += gray.at<uchar>(Higs[i].x, Higs[i].y);
		}

		u2 = s2 / Higs.size();

		cout << "u1: " << u1 << " u2: " << u2 << endl;

		// 迭代的终止条件：u1和u2都不再变化
		if (abs(best_thr - (u1 + u2) / 2.0) < 0.1)
		{
			flag = false;
		}
		else
		{
			best_thr = (u1 + u2) / 2.0;
		}
		
		cout << best_thr << endl;
	}// end while

	cout << best_thr << endl;
	for (i = 0; i < gray.rows; i++)
		for (j = 0; j < gray.cols; j++)
		{
			if (gray.at<uchar>(i, j) < best_thr)
			{
				gray.at<uchar>(i, j) = 0;
			}
		}
	imshow("iteration", gray);
	waitKey();

	return 0;
}

效果如下：