自动得到二值化阀值 ostu大津法

一、前言

Ostu方法又名最大类间差方法，通过统计整个图像的直方图特性来实现全局阈值T的自动选取，其算法步骤为：

1)  先计算图像的直方图，即将图像所有的像素点按照0~255共256个bin，统计落在每个bin的像素点数量

2)  归一化直方图，也即将每个bin中像素点数量除以总的像素点

3)  i表示分类的阈值，也即一个灰度级，从0开始迭代

4)  通过归一化的直方图，统计0~i 灰度级的像素(假设像素值在此范围的像素叫做前景像素) 所占整幅图像的比例w0，并统计前景像素的平均灰度u0；统计i~255灰度级的像素(假设像素值在此范围的像素叫做背景像素) 所占整幅图像的比例w1，并统计背景像素的平均灰度u1；

5)  计算前景像素和背景像素的方差 g = w0*w1*(u0-u1) (u0-u1)

6)  i++；转到4)，直到i为256时结束迭代

7）将最大g相应的i值作为图像的全局阈值

二、实现代码

Ostu算法实现函数为下面中的int getOstu(const Mat& in)，其实在opencv中已经在threshold中实现了，下面代码也比较了两者的结果

[cpp] view plain copy

 print?

1. #include "opencv2/core/core.hpp"  
2. #include "opencv2/highgui/highgui.hpp"  
3. #include "opencv2/imgproc/imgproc.hpp"  
4.   
5. using namespace cv;  
6.   
7. #include <iostream>  
8.   
9. int getOstu(const Mat & in);  
10.   
11. int main()  
12. {  
13.     Mat img = imread("id_card.jpg" ,0);  
14.     Mat img_high_Light = imread("id_card2.jpg" ,0);  
15.     Mat  dst , dst_HL;  
16.   
17.     if(img.empty()  | img_high_Light.empty())  
18.     {  
19.         std::cout<<"Error!!";  
20.         return -1;  
21.     }  
22.   
23.     std::cout<<"The return value of getOstu is: "<<getOstu(img);  
24.     std::cout<<"\n"<<"The return value of opencv threshold is: "<<threshold(img , dst ,0,255,CV_THRESH_OTSU);//opencv已实现的大津法  
25.   
26.     imshow("光照均匀的原图像" ,img);  
27.     imshow("图像光照均匀的处理结果" , dst);  
28.   
29.     threshold(img_high_Light , dst_HL ,0,255,CV_THRESH_OTSU);  
30.     imshow("光照不均匀的原图像" ,img_high_Light );  
31.     imshow("图像光照不均匀的处理结果",  dst_HL);  
32.   
33.   
34.     waitKey(0);  
35.   
36.     return 0;  
37. }  
38.   
39. int getOstu(const Mat & in)  
40. {  
41.     int rows = in.rows;  
42.     int cols = in.cols;  
43.     long size = rows * cols;  
44.   
45.     float histogram[256] = {0};    
46.     for( int i = 0; i < rows; ++i)   
47.     {   
48.         //获取第 i行首像素指针   
49.         const uchar * p = in.ptr<uchar>(i);   
50.         //对第i 行的每个像素(byte)操作   
51.         for( int j = 0; j < cols; ++j )   
52.         {  
53.             histogram[int(*p++)]++;   
54.         }  
55.     }  
56.     int threshold;        
57.     long sum0 = 0, sum1 = 0; //存储前景的灰度总和及背景灰度总和    
58.     long cnt0 = 0, cnt1 = 0; //前景的总个数及背景的总个数    
59.     double w0 = 0, w1 = 0; //前景及背景所占整幅图像的比例    
60.     double u0 = 0, u1 = 0;  //前景及背景的平均灰度    
61.     double variance = 0; //最大类间方差    
62.   
63.   
64.     double maxVariance = 0;    
65.     for(int i = 1; i < 256; i++) //一次遍历每个像素    
66.     {      
67.         sum0 = 0;    
68.         sum1 = 0;     
69.         cnt0 = 0;    
70.         cnt1 = 0;    
71.         w0 = 0;    
72.         w1 = 0;    
73.         for(int j = 0; j < i; j++)    
74.         {    
75.             cnt0 += histogram[j];    
76.             sum0 += j * histogram[j];    
77.         }    
78.   
79.         u0 = (double)sum0 /  cnt0;     
80.         w0 = (double)cnt0 / size;    
81.   
82.         for(int j = i ; j <= 255; j++)    
83.         {    
84.             cnt1 += histogram[j];    
85.             sum1 += j * histogram[j];    
86.         }    
87.   
88.         u1 = (double)sum1 / cnt1;    
89.         w1 = 1 - w0; // (double)cnt1 / size;    
90.   
91.         variance =  w0 * w1 *  (u0 - u1) * (u0 - u1);    
92.         if(variance > maxVariance)     
93.         {      
94.             maxVariance = variance;      
95.             threshold = i;      
96.         }     
97.     }      
98.   
99.     return threshold;      
100. }