自适应阈值分割—大津法(OTSU算法)C++实现
2016-08-16 21:46
344 查看
大津法是一种图像灰度自适应的阈值分割算法,是1979年由日本学者大津提出,并由他的名字命名的。大津法按照图像上灰度值的分布,将图像分成背景和前景两部分看待,前景就是我们要按照阈值分割出来的部分。背景和前景的分界值就是我们要求出的阈值。遍历不同的阈值,计算不同阈值下对应的背景和前景之间的类内方差,当类内方差取得极大值时,此时对应的阈值就是大津法(OTSU算法)所求的阈值。
假设图像的背景较暗,并且图像的大小为M×N,图像中像素的灰度值小于阈值T的像素个数记作N0,像素灰度大于阈值T的像素个数记作N1,则有:
ω0=N0/ M×N (1)
ω1=N1/ M×N (2)
N0+N1=M×N (3)
ω0+ω1=1 (4)
μ=ω0*μ0+ω1*μ1 (5)
g=ω0(μ0-μ)^2+ω1(μ1-μ)^2 (6)
将式(5)代入式(6),得到等价公式:
g=ω0ω1(μ0-μ1)^2 (7) 这个就是类间方差的公式表述
采用遍历的方法得到使类间方差g最大的阈值T,即为所求。
2. 计算背景图像的平均灰度、背景图像像素数所占比例
3. 计算前景图像的平均灰度、前景图像像素数所占比例
4. 遍历0~255各灰阶,计算并寻找类间方差极大值
C++代码实现:
原图像:
该幅图像计算出来的大津阈值是104;
用这个阈值分割的图像:
跟Opencv threshold方法中使用CV_THRESH_OTSU参数计算出来的分割图像一致:
对上述代码稍加修改,增加画出直方图部分:
为显示清晰,本次使用一幅对比明显的灰度图:
OTSU分割效果:
对应阈值和直方图:
以上图像黑白对比度非常明显,从直方图上也可以看到只有两个波峰,求得的OTSU阈值为102。
上图中红色的竖线标识出了OTSu阈值分割线,显见,阈值大致位于两个波峰的低谷之间。
何为类间方差?
对于图像I(x,y),前景(即目标)和背景的分割阈值记作T,属于前景的像素点数占整幅图像的比例记为ω0,其平均灰度μ0;背景像素点数占整幅图像的比例为ω1,其平均灰度为μ1。图像的总平均灰度记为μ,类间方差记为g。假设图像的背景较暗,并且图像的大小为M×N,图像中像素的灰度值小于阈值T的像素个数记作N0,像素灰度大于阈值T的像素个数记作N1,则有:
ω0=N0/ M×N (1)
ω1=N1/ M×N (2)
N0+N1=M×N (3)
ω0+ω1=1 (4)
μ=ω0*μ0+ω1*μ1 (5)
g=ω0(μ0-μ)^2+ω1(μ1-μ)^2 (6)
将式(5)代入式(6),得到等价公式:
g=ω0ω1(μ0-μ1)^2 (7) 这个就是类间方差的公式表述
采用遍历的方法得到使类间方差g最大的阈值T,即为所求。
Otsu实现思路
1. 计算0~255各灰阶对应的像素个数,保存至一个数组中,该数组下标是灰度值,保存内容是当前灰度值对应像素数2. 计算背景图像的平均灰度、背景图像像素数所占比例
3. 计算前景图像的平均灰度、前景图像像素数所占比例
4. 遍历0~255各灰阶,计算并寻找类间方差极大值
C++代码实现:
#include <opencv2/highgui/highgui.hpp> #include <opencv2/imgproc/imgproc.hpp> #include <opencv2/core/core.hpp> #include <iostream> using namespace cv; using namespace std; //***************Otsu算法通过求类间方差极大值求自适应阈值****************** int OtsuAlgThreshold(const Mat image); int main(int argc,char *argv[]) { Mat image=imread(argv[1]); imshow("SoureImage",image); cvtColor(image,image,CV_RGB2GRAY); Mat imageOutput; Mat imageOtsu; int thresholdValue=OtsuAlgThreshold(image); cout<<"类间方差为: "<<thresholdValue<<endl; threshold(image,imageOutput,thresholdValue,255,CV_THRESH_BINARY); threshold(image,imageOtsu,0,255,CV_THRESH_OTSU); //Opencv Otsu算法 //imshow("SoureImage",image); imshow("Output Image",imageOutput); imshow("Opencv Otsu",imageOtsu); waitKey(); return 0; } int OtsuAlgThreshold(const Mat image) { if(image.channels()!=1) { cout<<"Please input Gray-image!"<<endl; return 0; } int T=0; //Otsu算法阈值 double varValue=0; //类间方差中间值保存 double w0=0; //前景像素点数所占比例 double w1=0; //背景像素点数所占比例 double u0=0; //前景平均灰度 double u1=0; //背景平均灰度 double Histogram[256]={0}; //灰度直方图,下标是灰度值,保存内容是灰度值对应的像素点总数 uchar *data=image.data; double totalNum=image.rows*image.cols; //像素总数 //计算灰度直方图分布,Histogram数组下标是灰度值,保存内容是灰度值对应像素点数 for(int i=0;i<image.rows;i++) //为表述清晰,并没有把rows和cols单独提出来 { for(int j=0;j<image.cols;j++) { Histogram[data[i*image.step+j]]++; } } for(int i=0;i<255;i++) { //每次遍历之前初始化各变量 w1=0; u1=0; w0=0; u0=0; //***********背景各分量值计算************************** for(int j=0;j<=i;j++) //背景部分各值计算 { w1+=Histogram[j]; //背景部分像素点总数 u1+=j*Histogram[j]; //背景部分像素总灰度和 } if(w1==0) //背景部分像素点数为0时退出 { break; } u1=u1/w1; //背景像素平均灰度 w1=w1/totalNum; // 背景部分像素点数所占比例 //***********背景各分量值计算************************** //***********前景各分量值计算************************** for(int k=i+1;k<255;k++) { w0+=Histogram[k]; //前景部分像素点总数 u0+=k*Histogram[k]; //前景部分像素总灰度和 } if(w0==0) //前景部分像素点数为0时退出 { break; } u0=u0/w0; //前景像素平均灰度 w0=w0/totalNum; // 前景部分像素点数所占比例 //***********前景各分量值计算************************** //***********类间方差计算****************************** double varValueI=w0*w1*(u1-u0)*(u1-u0); //当前类间方差计算 if(varValue<varValueI) { varValue=varValueI; T=i; } } return T; }
原图像:
该幅图像计算出来的大津阈值是104;
用这个阈值分割的图像:
跟Opencv threshold方法中使用CV_THRESH_OTSU参数计算出来的分割图像一致:
直方图直观理解
大津算法可以从图像直方图上有一个更为直观的理解:大津阈值大致上是直方图两个峰值之间低谷的值。对上述代码稍加修改,增加画出直方图部分:
#include <opencv2/highgui/highgui.hpp> #include <opencv2/imgproc/imgproc.hpp> #include <opencv2/core/core.hpp> #include <iostream> using namespace cv; using namespace std; //***************Otsu算法通过求类间方差极大值求自适应阈值****************** int OtsuAlgThreshold(const Mat image); int main(int argc,char *argv[]) { Mat image=imread(argv[1]); imshow("SoureImage",image); cvtColor(image,image,CV_RGB2GRAY); Mat imageOutput; Mat imageOtsu; int thresholdValue=OtsuAlgThreshold(image); cout<<"类间方差为: "<<thresholdValue<<endl; threshold(image,imageOutput,thresholdValue,255,CV_THRESH_BINARY); threshold(image,imageOtsu,0,255,CV_THRESH_OTSU); //Opencv Otsu算法 //imshow("SoureImage",image); imshow("Output Image",imageOutput); imshow("Opencv Otsu",imageOtsu); waitKey(); return 0; } int OtsuAlgThreshold(const Mat image) { if(image.channels()!=1) { cout<<"Please input Gray-image!"<<endl; return 0; } int T=0; //Otsu算法阈值 double varValue=0; //类间方差中间值保存 double w0=0; //前景像素点数所占比例 double w1=0; //背景像素点数所占比例 double u0=0; //前景平均灰度 double u1=0; //背景平均灰度 double Histogram[256]={0}; //灰度直方图,下标是灰度值,保存内容是灰度值对应的像素点总数 int Histogram1[256]={0}; uchar *data=image.data; double totalNum=image.rows*image.cols; //像素总数 //计算灰度直方图分布,Histogram数组下标是灰度值,保存内容是灰度值对应像素点数 for(int i=0;i<image.rows;i++) //为表述清晰,并没有把rows和cols单独提出来 { for(int j=0;j<image.cols;j++) { Histogram[data[i*image.step+j]]++; Histogram1[data[i*image.step+j]]++; } } //***********画出图像直方图******************************** Mat image1(255,255,CV_8UC3); for(int i=0;i<255;i++) { Histogram1[i]=Histogram1[i]%200; line(image1,Point(i,235),Point(i,235-Histogram1[i]),Scalar(255,0,0),1,8,0); if(i%50==0) { char ch[255]; sprintf(ch,"%d",i); string str=ch; putText(image1,str,Point(i,250),1,1,Scalar(0,0,255)); } } //***********画出图像直方图******************************** for(int i=0;i<255;i++) { //每次遍历之前初始化各变量 w1=0; u1=0; w0=0; u0=0; //***********背景各分量值计算************************** for(int j=0;j<=i;j++) //背景部分各值计算 { w1+=Histogram[j]; //背景部分像素点总数 u1+=j*Histogram[j]; //背景部分像素总灰度和 } if(w1==0) //背景部分像素点数为0时退出 { break; } u1=u1/w1; //背景像素平均灰度 w1=w1/totalNum; // 背景部分像素点数所占比例 //***********背景各分量值计算************************** //***********前景各分量值计算************************** for(int k=i+1;k<255;k++) { w0+=Histogram[k]; //前景部分像素点总数 u0+=k*Histogram[k]; //前景部分像素总灰度和 } if(w0==0) //前景部分像素点数为0时退出 { break; } u0=u0/w0; //前景像素平均灰度 w0=w0/totalNum; // 前景部分像素点数所占比例 //***********前景各分量值计算************************** //***********类间方差计算****************************** double varValueI=w0*w1*(u1-u0)*(u1-u0); //当前类间方差计算 if(varValue<varValueI) { varValue=varValueI; T=i; } } //画出以T为阈值的分割线 line(image1,Point(T,235),Point(T,0),Scalar(0,0,255),2,8); imshow("直方图",image1); return T; }
为显示清晰,本次使用一幅对比明显的灰度图:
OTSU分割效果:
对应阈值和直方图:
以上图像黑白对比度非常明显,从直方图上也可以看到只有两个波峰,求得的OTSU阈值为102。
上图中红色的竖线标识出了OTSu阈值分割线,显见,阈值大致位于两个波峰的低谷之间。
相关文章推荐
- 自适应阈值分割—大津法(OTSU算法)C++实现
- 自适应阈值分割—大津法(OTSU算法)C++实现
- opencv实现c++的otsu自适应阈值分割的算法描述
- 数字图像处理,自适应维纳滤波的C++实现
- 可变分区方式下的最优适应调度算法(C++环境下实现)
- C++实现otsu算法
- 大津法 c++实现
- 自适应Canny算法的C++实现[与Matlab的效果非常相近]
- Otsu算法在C/C++上的实现(从图片读取到像素值计算以及图片写回)
- Singleton模式的C++实现研究(示例代码)
- 在C++中实现“属性 (Property)”
- C++机理:虚拟机制的实现[兼谈对比于传统机制]
- 用C++ std::priority_queue 实现哈夫曼算法
- 用C++实现C#中的委托/事件(标准C++之升级版)
- 再探C++的单件实现
- 用 C++ 实现 C# 中的 委托/事件 (2-delegate event functor)
- 在C++中实现属性
- 用PHP实现通过Web执行C/C++程序
- 使用c++实现Format函数
- More Effective C++ Item 附2:一个auto_ptr的实现实例