设原始灰度级为M,灰度级为i的像素点个数为ni,对灰度直方图进行归一化:
opencv代码
// m_otsu.cpp : 定义控制台应用程序的入口点。
//
#include "stdafx.h"
#include "cv.h"
#include "highgui.h"
int _tmain(int argc, _TCHAR* argv[])
{
int i,j,nThresh;
int nHistogram[256] = {0};
double fStdHistogram[256] = {0.0};
double fGrayAccu[256] = {0.0};
double fGrayAve[256] = {0.0};
double fAverage = 0;
double fTemp = 0;
double fMax = 0;
IplImage *src,*dst;
src = cvLoadImage("test.jpg",CV_LOAD_IMAGE_GRAYSCALE);
dst = cvCreateImage(cvSize(src->width,src->height),IPL_DEPTH_8U,1);
//统计直方图
// 每行
for(i = 0; i < src->height; i++)
{
// 每列
for(j = 0; j < src->width; j++)
{
nHistogram[(unsigned char)src->imageData[i*src->width+j]] ++;
}
}
//归一化直方图
for(i = 0; i <= 255;i++)
{
fStdHistogram[i] = nHistogram[i]/(double)(src->width * src->height); //Pi
//printf("%f\n",fStdHistogram[i]);
}
for(i=0;i<=255;i++)
{
for(j=0;j<=i;j++)
{
fGrayAccu[i] += fStdHistogram[j]; //所有灰度级,关于w0的数组
fGrayAve[i] += j*fStdHistogram[j]; //所有灰度级,关于u(t)的数组
}
fAverage += i*fStdHistogram[i]; //uT
//printf("%f\n",fAverage);
}
//计算OSTU
for(i=0;i<=255;i++)
{
fTemp=(fAverage*fGrayAccu[i]-fGrayAve[i])*(fAverage*fGrayAccu[i]-fGrayAve[i])/(fGrayAccu[i]*(1-fGrayAccu[i]));
if(fTemp>fMax)
{
fMax=fTemp;
nThresh=i;
}
}
//计算二值图像
for (i=0;i<src->height;i++)
{
for (j=0;j<src->width;j++)
{
if ((unsigned char)src->imageData[i*src->width+j]<nThresh)
{
dst->imageData[i*src->width+j] = 0;
}else{
dst->imageData[i*src->width+j] = 255;
}
}
}
printf("%d",nThresh);
cvNamedWindow("otsu",0);
cvShowImage("otsu",dst);
cvSaveImage("otsu_result.jpg",dst);
cvWaitKey(0);
return 0;
}