阈值分割-Otsu法

 算法实现（不调用函数）

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

using namespace std;
using namespace cv;

//实现灰度直方图的绘制
void drawpicture(Mat &inpicture, Mat &outpicture)
{
MatND hist;

const int bins = 256;
int hist_size[] = { bins };

float range[] = { 0, 256 };
const float* ranges[] = { range };

int channels[] = { 0 };

calcHist(&inpicture, 1, channels, Mat(), hist, 1, hist_size, ranges, true, false);

double maxnum;
minMaxLoc(hist, 0, &maxnum, 0, 0);
int scale = 1;
int histheight = 256;

for (int i = 0; i < bins; i++)
{
float binvalue = hist.at<float>(i);
int height = cvRound(binvalue*histheight / maxnum);
// 绘制长方形
rectangle(outpicture, Point(i*scale, histheight), Point((i + 1)*scale, histheight - height), Scalar(255));
namedWindow("windown", CV_WINDOW_NORMAL);
imshow("windown", outpicture);
}
}

int main()
{
//以灰度方式读入图片
Mat img = imread("F:\work_three_grade\DigitalImage\d.jpg", IMREAD_GRAYSCALE);

double max = 0.0; //方差
int otsu = 0; //阈值

//得到像素点总个数
int h = img.rows;
int w = img.cols;
int num = h * w;

//声明一个一维数组，用来记录各像素级的个数
int m[256];
for (int i = 0; i < 256; i++)
{
m[i] = 0;
}
//Mat M(256, CV_8U, m);
for (int i = 0; i < h; i++)
{
for (int j = 0; j < w; j++)
{
int temp = img.at<uchar>(i, j);
m[temp]++;
}
}

//小于等于t的和大于t的像素出现的概率和平均灰度
int t;//将阈值从0-255依次遍历，寻找结果最好的
double outcome[256];//存放阈值为i的类间方差

for (t = 0; t < 256; t++)
{
double w0=0.0;//小于等于t的像素出现的概率
double w1= 0.0;//大于t的像素出现的概率
double u0= 0.0;// 平均灰度
double u1= 0.0;// 平均灰度
double sum1 = 0.0;
double sum2 = 0.0;
for (int i = 0; i <= t; i++)
{
sum1 += m[i];
}
w0 = sum1 / num; //小于等于t的像素出现的概率,应该是double类型，后面再说
w1 = 1 - w0; //大于t的像素出现的概率,应该是double类型，后面再说
sum1 = 0;
for (int i = 0; i <= t; i++)
{
sum1 += i * m[i];
}
u0 = sum1 / num / w0;
for (int i = t + 1; i < 256; i++)
{
sum2 += i * m[i];
}
u1 = sum2 / num / w1;

//暂存结果
outcome[t] = w0 * w1*(u0 - u1)*(u0 - u1);
}

//求最大类间方差

for (int i = 0; i < 256; i++)
{
if (max < outcome[i])
{
max = outcome[i];
otsu = i;
}
}

//输出原图像
namedWindow("img",CV_WINDOW_NORMAL);
imshow("img", img);
//输出灰度直方图
Mat grayimg = Mat::zeros(256, 256, CV_8UC1);
std::string grayimgwin;
drawpicture(img, grayimg);
cout << "阈值： " << otsu << endl;
cout << "最大类间方差： " << max << endl;

waitKey();
return 0;

}