アルゴリズム

オトス法(最大クラス間分散法、大津アルゴリズムとも呼ばれる場合があります)は、クラスターの思想を用いて、画像の階調数をグレーレベルで2つの部分に分けて、2つの部分の間の階調値の差が最大となり、各部分の間の階調差が最小となり、分散の計算により適切な階調レベルを求めて分割します. したがって、二値化の際には、otsアルゴリズムを用いて、閾値を自動的に選択して二値化することができます.otsアルゴリズムは画像分割における閾値選択の最適アルゴリズムとして考えられ,計算は単純であり,画像の輝度とコントラストに影響されない.したがって、クラス間の分散を最大にする分割は、エラー確率が最小となることを意味する.
tを設定の閾値とする.
私:別れてから 前景ピクセル数がイメージに占める比率
uo:  別れてから 前景ピクセルポイントの平均階調
w 1:別れてから 被景ピクセル数が画像に占める割合
u 1:  別れてから 被景画素点の平均階調数
u=w 0*u 0+w 1*u 1:画像総平均階調数
L個のグレーレベルからtを遍歴して、tがある値である場合、見通しと背景の分散が一番大きい場合、このtは 値は私たちが要求している閾値です.
ここで、分散の計算式は以下の通りである.
g=wo*(uo-u)*(uo-u)+w 1*(u 1-u)*(u 1-u)
［            この数式の計算量は大きいです.     g=wo***1*(uo-u 1)*(uo-u 1)               ]
otsアルゴリズムは画像の階調をクラスター化するために、so otsアルゴリズムを実行する前に、この画像の階調ヒストグラムを計算する必要がある.
上記の説明に従って参照コードは以下の通りです.

#include "stdafx.h"
#include "stdio.h"
#include "cv.h"
#include "highgui.h"
#include "Math.h"

int Otsu(IplImage* src);

int _tmain(int argc, _TCHAR* argv[])
{
	IplImage* img = cvLoadImage("c:\\aSa.jpg",0);
	IplImage* dst = cvCreateImage(cvGetSize(img), 8, 1);
	int threshold = Otsu(img);

	cvThreshold(img, dst, threshold, 255, CV_THRESH_BINARY);


	cvNamedWindow( "img", 1 );
	cvShowImage("img", dst);


	cvWaitKey(-1);

	cvReleaseImage(&img);
	cvReleaseImage(&dst);

	cvDestroyWindow( "dst" );
	return 0;
}

int Otsu(IplImage* src)  
{  
	int height=src->height;  
	int width=src->width;      
	long size = height * width; 

	//histogram  
	float histogram[256] = {0};  
	for(int m=0; m < height; m++)
	{  
		unsigned char* p=(unsigned char*)src->imageData + src->widthStep * m;  
		for(int n = 0; n < width; n++) 
		{  
			histogram[int(*p++)]++;  
		}  
	}  

	int threshold;    
	long sum0 = 0, sum1 = 0; //                
	long cnt0 = 0, cnt1 = 0; //             
	double w0 = 0, w1 = 0; //              
	double u0 = 0, u1 = 0;  //          
	double variance = 0; //      
	int i, j;
	double u = 0;
	double maxVariance = 0;
	for(i = 1; i < 256; i++) //        
	{  
		sum0 = 0;
		sum1 = 0; 
		cnt0 = 0;
		cnt1 = 0;
		w0 = 0;
	    w1 = 0;
		for(j = 0; j < i; j++)
		{
			cnt0 += histogram[j];
			sum0 += j * histogram[j];
		}

		u0 = (double)sum0 /  cnt0; 
		w0 = (double)cnt0 / size;

		for(j = i ; j <= 255; j++)
		{
			cnt1 += histogram[j];
			sum1 += j * histogram[j];
		}

		u1 = (double)sum1 / cnt1;
		w1 = 1 - w0; // (double)cnt1 / size;

		u = u0 * w0 + u1 * w1; //       
		printf("u = %f
", u);
		//variance =  w0 * pow((u0 - u), 2) + w1 * pow((u1 - u), 2);
		variance =  w0 * w1 *  (u0 - u1) * (u0 - u1);
		if(variance > maxVariance) 
		{  
			maxVariance = variance;  
			threshold = i;  
		} 
	}  

	printf("threshold = %d
", threshold);
	return threshold;  
}  

 w 1をw 0・・私のデバッグにして長い間～いつも不真面目で、頭がふらふらしている・・・これも見えません...
=====================================================================
そうだ、前に集めたotsのアルゴリズムのコードは次の通りです.

#include "stdafx.h"
#include "stdio.h"
#include "cv.h"
#include "highgui.h"
#include "Math.h"

int Otsu(IplImage* src);

int _tmain(int argc, _TCHAR* argv[])
{
	IplImage* img = cvLoadImage("c:\\aSa.jpg",0);
	IplImage* dst = cvCreateImage(cvGetSize(img), 8, 1);
	int threshold = Otsu(img);
	printf("threshold = %d
", threshold);
	cvThreshold(img, dst, threshold, 255, CV_THRESH_BINARY);

	cvNamedWindow( "img", 1 );
	cvShowImage("img", dst);


	cvWaitKey(-1);

	cvReleaseImage(&img);
	cvReleaseImage(&dst);
	
	cvDestroyWindow( "dst" );
	return 0;
}

int Otsu(IplImage* src)  
{  
	int height=src->height;  
	int width=src->width;      

	//histogram  
	float histogram[256] = {0};  
	for(int i=0; i < height; i++)
	{  
		unsigned char* p=(unsigned char*)src->imageData + src->widthStep * i;  
		for(int j = 0; j < width; j++) 
		{  
			histogram[*p++]++;  
		}  
	}  
	//normalize histogram  
	int size = height * width;  
	for(int i = 0; i < 256; i++)
	{  
		histogram[i] = histogram[i] / size;  
	}  

	//average pixel value  
	float avgValue=0;  
	for(int i=0; i < 256; i++)
	{  
		avgValue += i * histogram[i];  //         
	}   

	int threshold;    
	float maxVariance=0;  
	float w = 0, u = 0;  
	for(int i = 0; i < 256; i++) 
	{  
		w += histogram[i];  //      i   , 0~i      (                  )          
		u += i * histogram[i];  //   i      (0~i)      ：           

		float t = avgValue * w - u;  
		float variance = t * t / (w * (1 - w) );  
		if(variance > maxVariance) 
		{  
			maxVariance = variance;  
			threshold = i;  
		}  
	}  

	return threshold;  
}