C++实现Otsu(大津法）

C++实现：

#include <iostream>
#include <opencv2/core.hpp>
#include <opencv2/highgui.hpp>
#include <opencv2/imgproc.hpp>
using namespace std;
using namespace cv;
int Otsu(Mat& src, Mat& dst, int thresh) {
	const int Grayscale = 256;
	int graynum[Grayscale] = { 0 };   //数组下标0～255的数字全部写为0;
	int r = src.rows;
	int c = src.cols;
	for (int i = 0; i < r; ++i) {
		const uchar* ptr = src.ptr<uchar>(i); //返回图像第i行的头指针
		for (int j = 0; j < c; ++j) {        //直方图统计
			graynum[ptr[j]]++;  //ptr[j]返回该点像素的像素值，graynum【】存储着每个像素级的像素数
		}
	}

	double P[Grayscale] = { 0 };
	double PK[Grayscale] = { 0 };
	double MK[Grayscale] = { 0 };
	double srcpixnum = (double)r * c, sumtmpPK = 0, sumtmpMK = 0; //scripixnum是像素数，
	for (int i = 0; i < Grayscale; ++i) {
		P[i] = graynum[i] / srcpixnum;   //每个灰度级出现的概率

		PK[i] = sumtmpPK + P[i];         //概率累计和 
		sumtmpPK = PK[i];              //它没有什么特别的意义，只是累计和

		MK[i] = sumtmpMK + i * P[i];       //灰度级的累加均值，平均灰度级                                                                                                                                                                                                                                                                                                                                                                                                        
		sumtmpMK = MK[i];             //平均灰度级,就是放眼望去，哪个灰度级处在中间的位置，哪个颜色最占据人眼
	}

	//计算类间方差
	double Var = 0;
	for (int k = 0; k < Grayscale; ++k) {
		if ((sumtmpMK * PK[k] - MK[k]) * (sumtmpMK * PK[k] - MK[k]) / (PK[k] * (1 - PK[k])) > Var) {
			Var = (sumtmpMK * PK[k] - MK[k]) * (sumtmpMK * PK[k] - MK[k]) / (PK[k] * (1 - PK[k]));
			thresh = k;
		}
	}

	//阈值处理
	src.copyTo(dst);
	for (int i = 0; i < r; ++i) {
		uchar* ptr = dst.ptr<uchar>(i);   //该点的像素值大于阈值
		for (int j = 0; j < c; ++j) {
			if (ptr[j] > thresh)
				ptr[j] = 255;
			else
				ptr[j] = 0;
		}
	}
	return thresh;
}


int main() {
	Mat src = imread("D:/PT/许昕.jpg");
	if (src.empty()) {
		return -1;
	}
	double scale = 0.2;
	Size dsize = Size(src.cols * scale, src.rows * scale);
	//Mat src2 = Mat(dsize, CV_32S);
	resize(src, src, dsize);

	if (src.channels() > 1)
		cvtColor(src, src, CV_RGB2GRAY);

	Mat dst, dst2;
	int thresh = 0;
	double t2 = (double)getTickCount();
	thresh = Otsu(src, dst, thresh); //Otsu阈值
	cout << "Mythresh=" << thresh << endl;
	t2 = (double)getTickCount() - t2;
	double time2 = (t2 * 1000.) / ((double)getTickFrequency());
	cout << "my_process=" << time2 << " ms. " << endl << endl;
	double  Otsu = 0;

	Otsu = threshold(src, dst2, Otsu, 255, CV_THRESH_OTSU + CV_THRESH_BINARY);
	cout << "OpenCVthresh=" << Otsu << endl;

	namedWindow("src", CV_WINDOW_AUTOSIZE);
	imshow("src", src);
	namedWindow("dst", CV_WINDOW_AUTOSIZE);
	imshow("dst", dst);
	namedWindow("dst2", CV_WINDOW_AUTOSIZE);
	imshow("dst2", dst2);

	imwrite("D:/PT/xuxin2.jpg", dst);
	imwrite("D:/PT/xuxin.jpg", src);
	waitKey(0);
}