详细及易读懂的 大津法(OTSU)原理 和 比opencv自带更快的算法实现

发布于:2021-10-20 19:55:54

OTSU算法原理简述:


最大类间方差是由日本学者大津(Nobuyuki Otsu)于1979年提出,是一种自适应的阈值确定方法。算法假设图像像素能够根据阈值,被分成背景[background]和目标[objects]两部分。然后,计算该最佳阈值来区分这两类像素,使得两类像素区分度最大。


公式:? 记?M = 256 单通道灰度分级 Sum = 像素总数


    背景像素占比???前景像素占比背景的*均灰度值前景的*均灰度值0~M灰度区间的灰度累计值类间方差:将公式3.4.5带入公式6 可得最终简化公式:?

?3、4的意思:每个灰度(0~255)比例 * 当前灰度的和。


代码:


#include
#include
#include
using namespace std;
using namespace cv;

int myOtsu(Mat & src, Mat &dst)
{

int th = 0;
const int GrayScale = 256; //单通道图像总灰度256级
int pixCount[GrayScale] = { 0 };//每个灰度值所占像素个数
int pixSum = src.cols * src.rows;//图像总像素点
float pixPro[GrayScale] = { 0 };//每个灰度值所占总像素比例
float SumpixPro[GrayScale] = { 0 }; // 比例的和
float WpixPro[GrayScale] = { 0 }; //比例 * 权重
float SumWpixPro[GrayScale] = { 0 };//比例 * 权重 的 和
float w0, w1, u0tmp, u1tmp, u0, u1, deltaTmp, deltaMax = 0;
double start = getTickCount(); //开始时间
for (int i = 0; i < src.rows; i++)
{
for (int j = 0; j < src.cols; j++)
{
int index = i * src.cols + j;
pixCount[src.data[index]]++;//统计每个灰度级中像素的个数
}
}
double c1 = getTickCount();
//cout << "c1 >> " << (c1 - start) / getTickFrequency() << endl;//输出时间
for (int i = 0; i < GrayScale; i++)
{
pixPro[i] = pixCount[i] * 1.0 / pixSum;//计算每个灰度级的像素数目占整幅图像的比例
WpixPro[i] = i * pixPro[i];
if (i == 0)
{
SumWpixPro[i] += WpixPro[i];
SumpixPro[i] += pixPro[i];
}
else
{
SumWpixPro[i] = WpixPro[i] + SumWpixPro[i - 1];
SumpixPro[i] = pixPro[i] + SumpixPro[i - 1];
}
}

for (int i = 0; i < GrayScale; i++)//遍历所有从0到255灰度级的阈值分割条件,测试哪一个的类间方差最大
{
w0 = w1 = u0tmp = u1tmp = u0 = u1 = deltaTmp = 0;
/*for (int j = 0; j < GrayScale; j++)
{
if (j <= i)//背景
{
w0 += pixPro[j];
u0tmp += WpixPro[j];
}
else//前景
{
w1 += pixPro[j];
u1tmp += WpixPro[j];
}
}*/
w0 = SumpixPro[i];
w1 = 1 - w0;

if (w0 == 0 || w1 == 0)
continue;
u0tmp = SumWpixPro[i];
u1tmp = SumWpixPro[255] - SumWpixPro[i];

u0 = u0tmp / w0;
u1 = u1tmp / w1;
deltaTmp = (float)(w0 *w1* pow((u0 - u1), 2)); //类间方差公式 g = w1 * w2 * (u1 - u2) ^ 2
if (deltaTmp > deltaMax)
{
deltaMax = deltaTmp;
th = i;
}
}
double c2 = getTickCount();
//cout << "c2 >> " << (c2 - c1) / getTickFrequency() << endl;//输出时间
for (int i = 0; i < src.rows; i++)
{
for (int j = 0; j < src.cols; j++)
{
int index = i * src.cols + j;
if (src.data[index] > th)
dst.data[index] = 255;
else
dst.data[index] = 0;
}
}
double c3 = getTickCount();
//cout << "c3 >> " << (c3 - c2) / getTickFrequency() << endl;//输出时间
return th;
}

int main()
{
Mat src = imread("scene.jpg", IMREAD_GRAYSCALE);//单通道读取图像
/*my_dst: 自己实现的大津法 得到的处理图像
otsu_dst:opencv自带的大津法 得到的处理图像
sub:两个处理图像相差图
*/
Mat my_dst, otsu_dst, sub;
/*my_th: 自己实现的大津法 得到的最大类件方差 即阈值
th:opencv自带的大津法 得到的最大类件方差 即阈值
*/
int my_th, th;

/*计算开销时间,对比两个算法效率*/
int times = 100;

double start = getTickCount(); //开始时间
my_dst = Mat::zeros(src.size(), CV_8UC1);
while(times--)
{

my_th = myOtsu(src, my_dst);
}
double end = getTickCount(); //结束时间
cout << "myOtsu end algorithm >> " << (end - start) / getTickFrequency() << endl; //输出时间
cout << "myOtsu threshold >> " << my_th << endl;

start = getTickCount(); //开始时间
times = 100;
while(times--)
{
th = threshold(src, otsu_dst, 0, 255, THRESH_OTSU);
}
end = getTickCount(); //结束时间
cout << "Otsu end algorithm >> " << (end - start) / getTickFrequency() << endl;//输出时间
cout << "Otsu threshold >> " << th << endl;


waitKey();
system("pause");
return 0;
}

原图:用的是1920*1080的图像,循环100次得到测试结果和图像。



结论:可以看到自己实现的otsu和opencv自带的自适应阈值算法效率上相差大约3倍左右。


分析:遍历图片没有加速,opencv内部有加速。


分析时间慢的原因:第一次遍历整张图片,统计每个灰度级像素个数,第二次遍历整张图片进行阈值操作。


似曾相识第一次遍历不就是灰度直方图吗,然后第二次遍历不就是LUT图像压缩。


改进后源码:


#if 1
#include
#include
#include
using namespace std;
using namespace cv;

void getHistogram(Mat &src, int *dst)
{
Mat hist;
int channels[1] = { 0 };
int histSize[1] = { 256 };
float hranges[2] = { 0, 256.0 };
const float *ranges[1] = { hranges };
calcHist(&src, 1, channels, Mat(), hist, 1, histSize, ranges);
//cout << hist ;
for (int i = 0; i < 256; i++)
{
float binVal = hist.at(i);
dst[i] = int(binVal);
}
}


int myOtsu(Mat & src, Mat &dst)
{

int th = 0;
const int GrayScale = 256; //单通道图像总灰度256级
int pixCount[GrayScale] = { 0 };//每个灰度值所占像素个数
int pixSum = src.cols * src.rows;//图像总像素点
float pixPro[GrayScale] = { 0 };//每个灰度值所占总像素比例
float SumpixPro[GrayScale] = { 0 }; // 比例的和
float WpixPro[GrayScale] = { 0 }; //比例 * 权重
float SumWpixPro[GrayScale] = { 0 };//比例 * 权重 的 和
float w0, w1, u0tmp, u1tmp, u0, u1, deltaTmp, deltaMax = 0;
double start = getTickCount(); //开始时间

/*for (int i = 0; i < src.rows; i++)
{
for (int j = 0; j < src.cols; j++)
{
int index = i * src.cols + j;
pixCount[src.data[index]]++;//统计每个灰度级中像素的个数
}
}*/
getHistogram(src, pixCount);

double c1 = getTickCount();
//cout << "c1 >> " << (c1 - start) / getTickFrequency() << endl;//输出时间
for (int i = 0; i < GrayScale; i++)
{
pixPro[i] = pixCount[i] * 1.0 / pixSum;//计算每个灰度级的像素数目占整幅图像的比例
WpixPro[i] = i * pixPro[i];
if (i == 0)
{
SumWpixPro[i] += WpixPro[i];
SumpixPro[i] += pixPro[i];
}
else
{
SumWpixPro[i] = WpixPro[i] + SumWpixPro[i - 1];
SumpixPro[i] = pixPro[i] + SumpixPro[i - 1];
}
}

for (int i = 0; i < GrayScale; i++)//遍历所有从0到255灰度级的阈值分割条件,测试哪一个的类间方差最大
{
w0 = w1 = u0tmp = u1tmp = u0 = u1 = deltaTmp = 0;
/*for (int j = 0; j < GrayScale; j++)
{
if (j <= i)//背景
{
w0 += pixPro[j];
u0tmp += WpixPro[j];
}
else//前景
{
w1 += pixPro[j];
u1tmp += WpixPro[j];
}
}*/
w0 = SumpixPro[i];
w1 = 1 - w0;

if (w0 == 0 || w1 == 0)
continue;
u0tmp = SumWpixPro[i];
u1tmp = SumWpixPro[255] - SumWpixPro[i];

u0 = u0tmp / w0;
u1 = u1tmp / w1;
deltaTmp = (float)(w0 *w1* pow((u0 - u1), 2)); //类间方差公式 g = w1 * w2 * (u1 - u2) ^ 2
if (deltaTmp > deltaMax)
{
deltaMax = deltaTmp;
th = i;
}
}
double c2 = getTickCount();
//cout << "c2 >> " << (c2 - c1) / getTickFrequency() << endl;//输出时间
uchar lutData[256];
for (int i = 0; i < 256; i++)
{
if (i > th)
lutData[i] = 255;
else
lutData[i] = 0;
}
Mat lut(1, 256, CV_8UC1, lutData);
LUT(src, lut, dst);
/*for (int i = 0; i < src.rows; i++)
{
for (int j = 0; j < src.cols; j++)
{
int index = i * src.cols + j;
if (src.data[index] > th)
dst.data[index] = 255;
else
dst.data[index] = 0;
}
}*/
double c3 = getTickCount();
//cout << "c3 >> " << (c3 - c2) / getTickFrequency() << endl;//输出时间
return th;
}

int main()
{
Mat src = imread("scene.jpg", IMREAD_GRAYSCALE);//单通道读取图像
/*my_dst: 自己实现的大津法 得到的处理图像
otsu_dst:opencv自带的大津法 得到的处理图像
sub:两个处理图像相差图
*/
//Mat src = Mat(Size(3000, 3000), CV_8UC1, Scalar(127));
Mat my_dst, otsu_dst, sub;
/*my_th: 自己实现的大津法 得到的最大类件方差 即阈值
th:opencv自带的大津法 得到的最大类件方差 即阈值
*/
int my_th, th;

/*计算开销时间,对比两个算法效率*/
int times = 100;

double start = getTickCount(); //开始时间
my_dst = Mat::zeros(src.size(), CV_8UC1);
while(times--)
{

my_th = myOtsu(src, my_dst);
}
double end = getTickCount(); //结束时间
cout << "myOtsu end algorithm >> " << (end - start) / getTickFrequency() << endl; //输出时间
cout << "myOtsu threshold >> " << my_th << endl;

start = getTickCount(); //开始时间
times = 100;
while(times--)
{
th = threshold(src, otsu_dst, 0, 255, THRESH_OTSU);
}
end = getTickCount(); //结束时间
cout << "Otsu end algorithm >> " << (end - start) / getTickFrequency() << endl;//输出时间
cout << "Otsu threshold >> " << th << endl;


waitKey();
system("pause");
return 0;
}
#endif


发现优化后的效率比opencv自带的还快。


最后:欢迎大家的批评,很高心与大家分享,谢谢大家。


?


?

相关推荐

最新更新

猜你喜欢