判別分析二値化プログラム

誰かの役に立つかもと思うので、ソースをさらす。
判別分析二値化法というのは、画像の濃度ヒストグラムをとってやって、その分離度が最大となるように、二値化のためのしきい値を決めてやるというもの。

詳しくは、ここら辺を読んでみて。
http://ipr20.cs.ehime-u.ac.jp/column/gazo_syori/chapter4.html#3rd
http://imagingsolution.blog107.fc2.com/blog-entry-113.html

// 判別分析法による閾値の決定(現実装では１次元ヒストグラムのみ)
float getThreshold(CvHistogram* histStruct)
{
	assert(!CV_IS_SPARSE_HIST(histStruct));
	CvMatND *mat;
	mat = &(histStruct->mat);

	assert(mat->dims == 1);
	int hist_size = (int)(mat->dim[0].size);
	int type = CV_MAT_TYPE(mat->type);
	assert(type == CV_32FC1);

	float* histgram = mat->data.fl;
	int* integral = new int[hist_size];	// ヒストグラムの累積
	int* average = new int[hist_size];
	int threshold = hist_size / 2;

	integral[0] = histgram[0];
	average[0] = 0;
	int i;
	for(i = 1; i < hist_size; i++){
		integral[i] = integral[i-1] + histgram[i];
		average[i] = average[i-1] + i * histgram[i];
	}
	
	int allnum = integral[hist_size - 1];
	int allavg = average[hist_size - 1];
	int w1, w2;
	float avg1, avg2;
	float tmp;
	double val;
	double max = 0;
	for(i=0; i 0 && w2 > 0){
			avg1 = (float)average[i] / w1;
			avg2 = (float)(allavg - average[i]) / w2;
			tmp = avg1 - avg2;
			val = w1 * w2 * tmp * tmp;
			if(val > max){
				max = val;
				threshold = i;
			}
		}
	}
	delete integral;
	delete average;

	// 最終的な閾値を求める
	float minf = histStruct->thresh[0][0];
	float maxf = histStruct->thresh[0][1];

	float bin_size = (maxf - minf) / hist_size;
	float ret = bin_size * (0.5 + threshold) + minf;

	return ret;
}

int _tmain(int argc, char* argv){
	char filename[512];
	printf("Img File Name: ");
	scanf("%s",filename);

	//　画像をグレースケールで読み込み
	IplImage* src_img = cvLoadImage(filename, CV_LOAD_IMAGE_GRAYSCALE);

	//　ヒストグラムの作成
	int histsize = 256;
	float range_0 = { 0, 256 };
	float *ranges[] = { range_0 };
	CvHistogram* hist = cvCreateHist(1, &histsize,CV_HIST_ARRAY, ranges, 1);
	cvCalcHist(&src_img, hist);

	//　判別分析法による閾値算出
	float threshold = getThreshold(hist);
	printf("threshold: %f\n",threshold);
	cvReleaseHist(&hist);

	//　結果の表示
	IplImage *dest_img = cvCreateImage(cvGetSize(src_img), IPL_DEPTH_8U, 1);
	cvThreshold(src_img, dest_img, threshold, 255, CV_THRESH_BINARY);
	cvNamedWindow ("result", CV_WINDOW_AUTOSIZE);
	cvShowImage ("result", src_img);
	cvWaitKey (0);
	cvShowImage ("result", dest_img);
	cvWaitKey (0);
	cvDestroyWindow("result");
	cvReleaseImage(&dest_img);

	cvReleaseImage(&src_img);
}

このプログラムを使うと、こんな感じ。

基画像：

二値化画像：