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)). 首先均值是基礎, 如果標准差大寫,閾值就會大些,標准差小些,閾值就會小些.
這個方法對一些不是光照不均的圖片有時候效果不好,現在還在找較好的方法,初步打算先用全局均值做二值化,如何效果不好再用局部均值的方法.
----為什么我的博客沒人看啊……………………………………
