开源社区 
二值图像是指每个像素亮度只有两个可能值的数字图像，其亮度函数可以描述为：
P(x,y) \in \left \{ 0, 1\right \}
其中：P 表示像素亮度， x,y 表示像素坐标，\left \{ 0, 1\right \} 表示像素亮度值的可能取值。
图像二值化方法在图像处理中时最常见的操作，在工业领域是许多应用的前置操作，例如在缺陷处理，轮廓识别，特征处理等应用中起到至关重要的作用。
二值化图像方法一般可以分为两种：1. 全局阈值二值化；2. 局部阈值二值化。在全局阈值二值化中，整幅图像使用一个固定的阈值T对图像像素进行二值化操作。最常见的一种情况是当像素亮度高于阈值，新亮度取1，否则取0。目前最常见的全局阈值方法有：OTSU，Triangle，其中OTSU适用于直方图为双峰的图像，Triangle适用于直方图为单峰的图像。
\begin{aligned}P'(x, y) = \begin{cases} 1, & \text{if } \text{p(x,y)} \ge T \\ 0, & \text{otherwise} \end{cases} \end{aligned}
全局阈值二值化方法，简单易行，算法成本低，然而在图像光照条件不佳时难以应用。局部阈值二值化则是根据图像不同部分的亮度计算出不同的阈值进行处理，其阈值实际上是一个随像素改变的二维函数曲面，可以有效的解决光照问题。局部阈值效果较好，然后算法成本高，速度慢，综合来看，目前使用最为广泛的依然是使用OTSU的全局二值化方法。本文将介绍OTSU的主要流程以及cpp实现的核心代码。
1. 流程
对于一副RGB图像，欲将其二值化，通常可分为2步：1)灰度化：将RGB图像转变为灰度图像；2）二值化：选取二值化方法将灰度图像转化为二值图像
对于OTSU方法可分为以下6个步骤：
1. 计算每个强度级的直方图和概率
2. 设置类概率 \omega_i(0) 和类均值 \mu_i(0) 的初始值
3. 遍历所有可能的阈值 t = 1 \ldots 最大强度
更新 \omega_i 和 \mu_i
计算类的方差 \sigma^2_b(t)
4. 所需的阈值对应于最大的 \sigma^2_b(t)
5. 计算两个最大值（和两个对应的）。 \sigma^2_{b1}(t) 是最大值而 \sigma^2_{b2}(t) 是更大的或相等的最大值
6. 所需的阈值 = \frac{T_1 + T_2}{2}
c伪代码如下：
| unsigned char otsu_threshold( int* histogram, int pixel_total ){ unsigned int sumB = 0; unsigned int sum1 = 0; float wB = 0.0f; float wF = 0.0f; float mF = 0.0f; float max_var = 0.0f; float inter_var = 0.0f; unsigned char threshold = 0; unsigned short index_histo = 0; for ( index_histo = 1; index_histo < 256; ++index_histo ) { sum1 += index_histo * histogram[ index_histo ]; } for (index_histo = 1; index_histo < 256; ++index_histo) { wB = wB + histogram[ index_histo ]; wF = pixel_total - wB; if ( wB == 0 || wF == 0 ) { continue; } sumB = sumB + index_histo * histogram[ index_histo ]; mF = ( sum1 - sumB ) / wF; inter_var = wB * wF * ( ( sumB / wB ) - mF ) * ( ( sumB / wB ) - mF ); if ( inter_var >= max_var ) { threshold = index_histo; max_var = inter_var; } } return threshold;} |
对于Opencv而言使用cvtColor函数可以将RGB/BGR图像转化为灰度图像，使用threshold函数可以将灰度图像转化为二值化图像，并返回计算的阈值。
OpenCV使用两个宏来控制阈值计算，分别表示OTSU，Triangle方法：
· THRESH_OTSU
· THRESH_TRIANGLE
函数原型如下：
| double cv::threshold (InputArray _src_, OutputArray _dst_, double _thresh_, double _maxval_, int _type_ ) |
OpenCV函数说明参考: OpenCV: Miscellaneous Image Transformations
2. 代码时间
为了比对OpenCV中实现的threshold性能，根据前文所述的OTSU流程，可以实现一个最基础的OTSU函数，看一看到底有多大差距。除此之外，为了便于测算时间，使用OpenCV中的cv::getCPUTickCount()与cv::getTickFrequency()函数构建一个timer类，顾名思义，cv::getCPUTickCount()与cv::getTickFrequency()类似与Linux中的Time与Clock可以返回系统的CPU时钟数与频率，从而计算出运行时间。
| // timer.hppclass Timer{public: Timer(const std::string &label) : label_(label) { this->cpu_count = cv::getCPUTickCount(); } ~Timer() { for (auto &ele : time_point) { printf("%s: %s elapsed: %.3f ms.\n", label_.c_str(), ele.first.c_str(), ele.second); } } double operator()(const std::string &point) { int64_t now = cv::getCPUTickCount(); double time = 1000 * static_cast<double>((now - this->cpu_count)) / cv::getTickFrequency(); time_point.emplace_back(std::make_pair(point, time)); update_timestamp(); return time; } private: const std::string label_; std::vector<std::pair<std::string, double>> time_point; int64_t cpu_count; void update_timestamp() { this->cpu_count = cv::getCPUTickCount(); }}; |
| // otsu.hppint myOtsu(Mat & src){ int th; const int GrayScale = 256; //单通道图像总灰度256级 int pixCount[GrayScale] = {0};//每个灰度值所占像素个数 int pixSum = src.cols * src.rows;//图像总像素点 float pixPro[GrayScale] = {0};//每个灰度值所占总像素比例 float w0, w1, u0tmp, u1tmp, u0, u1, deltaTmp, deltaMax = 0; for(int i = 0; i < src.cols; i++) { for(int j = 0; j < src.rows; j++) { pixCount[src.at<uchar>(j,i)]++;//统计每个灰度级中像素的个数 } } for(int i = 0; i < GrayScale; i++) { pixPro[i] = pixCount[i] * 1.0 / pixSum;//计算每个灰度级的像素数目占整幅图像的比例 } 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 += j * pixPro[j]; } else//前景 { w1 += pixPro[j]; u1tmp += j * pixPro[j]; } } 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; } } return th;} |
万事俱备，现在可以运行了，测试代码如下（注意 OpenCV默认读取的图像是BGR格式，故转换灰度时，需要使用cv::COLOR_BGR2GRAY宏）：
| //main.cpp#include <iostream>#include <opencv2/imgcodecs.hpp>#include <opencv2/opencv.hpp> using cv::Mat;using std::string; #include "timer.hpp"#include "otsu.hpp" int main(const int argc, const char *argv[]){ const string file_name = std::string(argv[1]); Timer tracker("imbinary"); Mat img_ = cv::imread(file_name, cv::IMREAD_COLOR); Mat img; cv::cvtColor(img_, img, cv::COLOR_BGR2GRAY); Mat dst_c; Mat dst_; tracker("imread"); if (img.empty()) { printf("Failed read image from %s\n", file_name.c_str()); } else { printf("Read Image [%d x %d] from %s\n", img.cols, img.rows, file_name.c_str()); uint32_t th = myOtsu(img); cv::threshold(img, dst_c, th, 255, cv::THRESH_BINARY); tracker("otsu"); uint32_t th_cv = cv::threshold(img, dst_, 0, 255, cv::THRESH_OTSU); tracker("otsu with opencv"); printf("th: %u, th_cv: %u\n", th, th_cv); cv::imwrite("color.bmp", img_); cv::imwrite("gray.bmp", img); cv::imwrite("result.bmp", dst_c); cv::imwrite("result_cv.bmp", dst_); } return 0;} |
使用前篇文章的CMakeLists.txt文件，编译运行：
| $ cmake -Bbuild -S .$ cmake --build build |
| user@PhytiumPi-Pro:~/Documents/opencv_binary_cross/build$ ./image_demo lena512rgb.tifRead Image [513 x 513] from lena512rgb.tifth: 92, th_cv: 92imbinary: imread elapsed: 22.401 ms.imbinary: otsu elapsed: 5.752 ms.imbinary: otsu with opencv elapsed: 0.914 ms. |
可以看到myOTSU函数计算了和OpenCV相同的结果，然而OpenCV的速度几乎是myOTSU的16倍，这可能是OpenCV使用了一些类似于并行，查表之类的黑魔法用于算法的加速。
3.小结
本文针对飞腾派PRO中OTSU算法的原理与使用做了初步的介绍，同时测试了OpenCV中二值化函数再飞腾派PRO上运行的性能。在接下来的文章中，将继续介绍经典机器视觉算法在飞腾派PRO的使用。
4. 论文（深入了解）
1. Nobuyuki Otsu. A threshold selection method from gray-level histograms. IEEE Trans. Sys., Man., Cyber. 1979, 9 (1): 62–66. doi:10.1109/TSMC.1979.4310076
