OTSU算法（也称最大类间差法，有时也称之为大津算法）

在图像处理应用中二值化操作是一个很常用的处理方式，例如零器件图片的处理、文本图片和验证码图片中字符的提取、车牌识别中的字符分割，以及视频图像中的运动目标检测中的前景分割，等等。

   较为常用的图像二值化方法有：1）全局固定阈值；2）局部自适应阈值；3）OTSU等。此篇讲述的是OTSU算法。

OTSU算法也称最大类间差法，有时也称之为大津算法，由大津于1979年提出，被认为是图像分割中阈值选取的最佳算法，计算简单，不受图像亮度和对比度的影响，因此在数字图像处理上得到了广泛的应用。它是按图像的灰度特性,将图像分成背景和前景两部分。因方差是灰度分布均匀性的一种度量,背景和前景之间的类间方差越大,说明构成图像的两部分的差别越大,当部分前景错分为背景或部分背景错分为前景都会导致两部分差别变小。因此,使类间方差最大的分割意味着错分概率最小。

  otsu使用的是聚类的思想，把图像的灰度数按灰度级分成2个部分，使得两个部分之间的灰度值差异最大，每个部分之间的灰度差异最小，通过方差的计算来寻找一个合适的灰度级别 来划分。  所以 可以在二值化的时候 采用otsu算法来自动选取阈值进行二值化。otsu算法被认为是图像分割中阈值选取的最佳算法，计算简单，不受图像亮度和对比度的影响。因此,使类间方差最大的分割意味着错分概率最小。

使用otsu大津法实现直方图的单阈值分割，可用于图像的分割。大津是用最大类间方差来判读阈值的，也可以用最小类内方差判读，最大类间方差+最小类内方差=总方差

OPENCV的二值化操作中，有一种“大津阈值处理”的方法，使用函数cvThreshold(image,image2,0,255,CV_THRESH_OTSU) 实现，该函数就会使用大律法OTSU得到的全局自适应阈值来进行二值化图片，而参数中的threshold不再起作用。

设t为设定的阈值。

wo： 分开后  前景像素点数占图像的比例

uo：  分开后  前景像素点的平均灰度

w1：分开后  被景像素点数占图像的比例

u1：  分开后  被景像素点的平均灰度

u=w0*u0 + w1*u1 :图像总平均灰度

从L个灰度级遍历t，使得t为某个值的时候，前景和背景的方差最大， 则 这个 t  值便是我们要求得的阈值。

其中，方差的计算公式如下：

g=wo * (uo - u) * (uo - u) + w1 * (u1 - u) * (u1 - u)

[             此公式计算量较大，可以采用：      g = wo * w1 * (uo - u1) * (uo - u1)                ]

由于otsu算法是对图像的灰度级进行聚类，so  在执行otsu算法之前，需要计算该图像的灰度直方图。

【第一套代码】

1. #include "stdafx.h"  
2. #include "stdio.h"  
3. #include "cv.h"  
4. #include "highgui.h"  
5. #include "Math.h"  
6.   
7. int Otsu(IplImage* src);  
8.   
9. int _tmain(int argc, _TCHAR* argv[])  
10. {  
11.     IplImage* img = cvLoadImage("c:\\aSa.jpg",0);  
12.     IplImage* dst = cvCreateImage(cvGetSize(img), 8, 1);  
13.     int threshold = Otsu(img);  
14.   
15.     cvThreshold(img, dst, threshold, 255, CV_THRESH_BINARY);  
16.   
17.   
18.     cvNamedWindow( "img", 1 );  
19.     cvShowImage("img", dst);  
20.   
21.   
22.     cvWaitKey(-1);  
23.   
24.     cvReleaseImage(&img);  
25.     cvReleaseImage(&dst);  
26.   
27.     cvDestroyWindow( "dst" );  
28.     return 0;  
29. }  
30.   
31. int Otsu(IplImage* src)    
32. {    
33.     int height=src->height;    
34.     int width=src->width;        
35.     long size = height * width;   
36.   
37.     //histogram    
38.     float histogram[256] = {0};    
39.     for(int m=0; m < height; m++)  
40.     {    
41.         unsigned char* p=(unsigned char*)src->imageData + src->widthStep * m;    
42.         for(int n = 0; n < width; n++)   
43.         {    
44.             histogram[int(*p++)]++;    
45.         }    
46.     }    
47.   
48.     int threshold;      
49.     long sum0 = 0, sum1 = 0; //存储前景的灰度总和和背景灰度总和  
50.     long cnt0 = 0, cnt1 = 0; //前景的总个数和背景的总个数  
51.     double w0 = 0, w1 = 0; //前景和背景所占整幅图像的比例  
52.     double u0 = 0, u1 = 0;  //前景和背景的平均灰度  
53.     double variance = 0; //最大类间方差  
54.     int i, j;  
55.     double u = 0;  
56.     double maxVariance = 0;  
57.     for(i = 1; i < 256; i++) //一次遍历每个像素  
58.     {    
59.         sum0 = 0;  
60.         sum1 = 0;   
61.         cnt0 = 0;  
62.         cnt1 = 0;  
63.         w0 = 0;  
64.         w1 = 0;  
65.         for(j = 0; j < i; j++)  
66.         {  
67.             cnt0 += histogram[j];  
68.             sum0 += j * histogram[j];  
69.         }  
70.   
71.         u0 = (double)sum0 /  cnt0;   
72.         w0 = (double)cnt0 / size;  
73.   
74.         for(j = i ; j <= 255; j++)  
75.         {  
76.             cnt1 += histogram[j];  
77.             sum1 += j * histogram[j];  
78.         }  
79.   
80.         u1 = (double)sum1 / cnt1;  
81.         w1 = 1 - w0; // (double)cnt1 / size;  
82.   
83.         u = u0 * w0 + u1 * w1; //图像的平均灰度  
84.         printf("u = %f\n", u);  
85.         //variance =  w0 * pow((u0 - u), 2) + w1 * pow((u1 - u), 2);  
86.         variance =  w0 * w1 *  (u0 - u1) * (u0 - u1);  
87.         if(variance > maxVariance)   
88.         {    
89.             maxVariance = variance;    
90.             threshold = i;    
91.         }   
92.     }    
93.   
94.     printf("threshold = %d\n", threshold);  
95.     return threshold;    
96. }    

【第二套代码】

1. #include "stdafx.h"  
2. #include "stdio.h"  
3. #include "cv.h"  
4. #include "highgui.h"  
5. #include "Math.h"  
6.   
7. int Otsu(IplImage* src);  
8.   
9. int _tmain(int argc, _TCHAR* argv[])  
10. {  
11.     IplImage* img = cvLoadImage("c:\\aSa.jpg",0);  
12.     IplImage* dst = cvCreateImage(cvGetSize(img), 8, 1);  
13.     int threshold = Otsu(img);  
14.     printf("threshold = %d\n", threshold);  
15.     cvThreshold(img, dst, threshold, 255, CV_THRESH_BINARY);  
16.   
17.     cvNamedWindow( "img", 1 );  
18.     cvShowImage("img", dst);  
19.   
20.   
21.     cvWaitKey(-1);  
22.   
23.     cvReleaseImage(&img);  
24.     cvReleaseImage(&dst);  
25.       
26.     cvDestroyWindow( "dst" );  
27.     return 0;  
28. }  
29.   
30. int Otsu(IplImage* src)    
31. {    
32.     int height=src->height;    
33.     int width=src->width;        
34.   
35.     //histogram    
36.     float histogram[256] = {0};    
37.     for(int i=0; i < height; i++)  
38.     {    
39.         unsigned char* p=(unsigned char*)src->imageData + src->widthStep * i;    
40.         for(int j = 0; j < width; j++)   
41.         {    
42.             histogram[*p++]++;    
43.         }    
44.     }    
45.     //normalize histogram    
46.     int size = height * width;    
47.     for(int i = 0; i < 256; i++)  
48.     {    
49.         histogram[i] = histogram[i] / size;    
50.     }    
51.   
52.     //average pixel value    
53.     float avgValue=0;    
54.     for(int i=0; i < 256; i++)  
55.     {    
56.         avgValue += i * histogram[i];  //整幅图像的平均灰度  
57.     }     
58.   
59.     int threshold;      
60.     float maxVariance=0;    
61.     float w = 0, u = 0;    
62.     for(int i = 0; i < 256; i++)   
63.     {    
64.         w += histogram[i];  //假设当前灰度i为阈值, 0~i 灰度的像素(假设像素值在此范围的像素叫做前景像素) 所占整幅图像的比例  
65.         u += i * histogram[i];  // 灰度i 之前的像素(0~i)的平均灰度值： 前景像素的平均灰度值  
66.   
67.         float t = avgValue * w - u;    
68.         float variance = t * t / (w * (1 - w) );    
69.         if(variance > maxVariance)   
70.         {    
71.             maxVariance = variance;    
72.             threshold = i;    
73.         }    
74.     }    
75.   
76.     return threshold;    
77. }