sauvola二值化算法研究
sauvola是一种考虑局部均值亮度的图像二值化方法, 以局部均值为基准在根据标准差做些微调.算法实现上一般用积分图方法
来实现.这个方法能很好的解决全局阈值方法的短板—关照不均图像二值化不好的问题.先贴代码
| //************************************ // 函数名称: sauvola // 函数说明: 局部均值二值化 // 参 数: // const unsigned char * grayImage [in] 输入图像数据 // const unsigned char * biImage [out] 输出图像数据 // const int w [in] 输入输出图像数据宽 // const int h [in] 输入输出图像数据高 // const int k [in] threshold = mean*(1 + k*((std / 128) - 1)) // const int windowSize [in] 处理区域宽高 // 返 回 值: void //************************************ void sauvola(const unsigned char * grayImage, const unsigned char * biImage, const int w, const int h, const int k, const int windowSize){ int whalf = windowSize >> 1; int i, j; int IMAGE_WIDTH = w; int IMAGE_HEIGHT = h; // create the integral image unsigned long * integralImg = (unsigned long*)malloc(IMAGE_WIDTH*IMAGE_HEIGHT*sizeof(unsigned long*)); unsigned long * integralImgSqrt = (unsigned long*)malloc(IMAGE_WIDTH*IMAGE_HEIGHT*sizeof(unsigned long*)); int sum = 0; int sqrtsum = 0; int index; //收集数据 integralImg像素和积分图 integralImgSqrt像素平方和积分图 for (i = 0; i < IMAGE_HEIGHT; i++){ // reset this column sum sum = 0; sqrtsum = 0; for (j = 0; j < IMAGE_WIDTH; j++) { index = i*IMAGE_WIDTH + j; sum += grayImage[index]; sqrtsum += grayImage[index] * grayImage[index]; if (i == 0){ integralImg[index] = sum; integralImgSqrt[index] = sqrtsum; } else{ integralImgSqrt[index] = integralImgSqrt[(i - 1)*IMAGE_WIDTH + j] + sqrtsum; integralImg[index] = integralImg[(i - 1)*IMAGE_WIDTH + j] + sum; } } } //Calculate the mean and standard deviation using the integral image int xmin, ymin, xmax, ymax; double mean, std, threshold; double diagsum, idiagsum, diff, sqdiagsum, sqidiagsum, sqdiff, area; for (i = 0; i < IMAGE_WIDTH; i++){ for (j = 0; j < IMAGE_HEIGHT; j++){ xmin = max(0, i - whalf); ymin = max(0, j - whalf); xmax = min(IMAGE_WIDTH - 1, i + whalf); ymax = min(IMAGE_HEIGHT - 1, j + whalf); area = (xmax - xmin + 1) * (ymax - ymin + 1); if (area <= 0){ biImage[i * IMAGE_WIDTH + j] = 255; continue; } if (xmin == 0 && ymin == 0){ diff = integralImg[ymax * IMAGE_WIDTH + xmax]; sqdiff = integralImgSqrt[ymax * IMAGE_WIDTH + xmax]; } else if (xmin > 0 && ymin == 0){ diff = integralImg[ymax * IMAGE_WIDTH + xmax] - integralImg[ymax * IMAGE_WIDTH + xmin - 1]; sqdiff = integralImgSqrt[ymax * IMAGE_WIDTH + xmax] - integralImgSqrt[ymax * IMAGE_WIDTH + xmin - 1]; } else if (xmin == 0 && ymin > 0){ diff = integralImg[ymax * IMAGE_WIDTH + xmax] - integralImg[(ymin - 1) * IMAGE_WIDTH + xmax]; sqdiff = integralImgSqrt[ymax * IMAGE_WIDTH + xmax] - integralImgSqrt[(ymin - 1) * IMAGE_WIDTH + xmax];; } else{ diagsum = integralImg[ymax * IMAGE_WIDTH + xmax] + integralImg[(ymin - 1) * IMAGE_WIDTH + xmin - 1]; idiagsum = integralImg[(ymin - 1) * IMAGE_WIDTH + xmax] + integralImg[ymax * IMAGE_WIDTH + xmin - 1]; diff = diagsum - idiagsum; sqdiagsum = integralImgSqrt[ymax * IMAGE_WIDTH + xmax] + integralImgSqrt[(ymin - 1) * IMAGE_WIDTH + xmin - 1]; sqidiagsum = integralImgSqrt[(ymin - 1) * IMAGE_WIDTH + xmax] + integralImgSqrt[ymax * IMAGE_WIDTH + xmin - 1]; sqdiff = sqdiagsum - sqidiagsum; } mean = diff / area; std = sqrt((sqdiff - diff*diff / area) / (area - 1)); threshold = mean*(1 + k*((std / 128) - 1)); if (grayImage[j*IMAGE_WIDTH + i] < threshold) biImage[j*IMAGE_WIDTH + i] = 0; else biImage[j*IMAGE_WIDTH + i] = 255; } } free(integralImg); free(integralImgSqrt); }
|
代码要注意下面几点:
1 计算区域像素和,几乎使用积分图技术是必然的选择.
2 标准差的表示方法: std = sqrt((sqdiff - diff*diff / area) / (area - 1)) 终于感到高等代数没有白学,
可以看百度百科关于方差的说明
3 判定方程 threshold = mean*(1 + k*((std / 128) - 1)). 首先均值是基础, 如果标准差大写,阈值就会大些,标准差小些,阈值就会小些.
这个方法对一些不是光照不均的图片有时候效果不好,现在还在找较好的方法,初步打算先用全局均值做二值化,如何效果不好再用局部均值的方法.
----为什么我的博客没人看啊……………………………………