Otsu阈值化

1979年日本大津展之提出了Otsu阈值法，他根据类间方差最大将灰度图分为目标和背景。

方法如下：

        若一张灰度图有n个像素， 取阈值为t, 小于t的个数为n1, 大于t的个数为n2, 小于t的像素平均值为m1, 大于t的像素平均值为m2, 权重w1 = n1 / n, w2 = n2 / n,

        类间方差 = w1 * w2 * (m1 - m2) ^ 2

        遍历每一个阈值，计算每一个类间方差，使类间方差最大的阈值，即为所求阈值。

        

#include "cv.h"
#include "highgui.h"
int Otsu(CvMat*gray);
int main()
{
IplImage* src = cvLoadImage("flower.jpg");

const int WIDTH = src->width;
const int HEIGHT = src->height;

CvMat* gray = cvCreateMat(HEIGHT, WIDTH, CV_8UC1);
cvCvtColor(src, gray, CV_BGR2GRAY);

int thres = Otsu(gray);
cvThreshold(gray, gray, thres, 255, CV_THRESH_BINARY);

cvShowImage("SRC", src);
cvShowImage("GRAY", gray);
cvWaitKey(0);

cvCvtColor(gray, src, CV_GRAY2BGR);
cvSaveImage("dst.bmp", src);
cvReleaseMat(&gray);
return 0;
}
int Otsu(CvMat*gray)
{
const int width = gray->width;
const int height = gray->height;

int histData[256] = {0};

for (int j = 0; j < height; j ++)
{
uchar*data = (uchar*)(gray->data.ptr + j * gray->step);
for(int i = 0; i < width; i ++)
{
histData[data[i]]++;
}
}

int total = width * height;

float sum = 0;

for (int i = 0 ; i < 256 ; i ++)
sum += i * histData[i];

float sumB = 0;
int wB = 0;
int wF = 0;

float varMax = 0;
int threshold = 0;

for (int i = 0 ; i < 256 ; i ++)
{
wB += histData[i];
if (wB == 0) continue;

wF = total - wB;
if (wF == 0) break;

sumB += (float) (i * histData[i]);

float mB = sumB / wB;
float mF = (sum - sumB) / wF;

float varBetween = (float)wB * (float)wF * (mB - mF) * (mB - mF);

if (varBetween > varMax)
{
varMax = varBetween;
threshold = i;
}
}

return threshold;

}

源图和效果图如下：