c語言數字圖像處理(十):閾值處理


定義

全局閾值處理

假設某一副灰度圖有如下的直方圖,該圖像由暗色背景下的較亮物體組成,從背景中提取這一物體時,將閾值T作為分割點,分割后的圖像g(x, y)由下述公式給出,稱為全局閾值處理

 

多閾值處理

本文僅完成全局閾值處理的算法實現

基本全局閾值處理方法

1. 為全局閾值T選擇一個初始的估計值

2. 用T分割圖像,產生兩組像素:G1由大於T的像素組成,G2由小於T的像素組成

3. 對G1和G2的像素分別計算平均灰度值m1和m2

4. 計算新的閾值T = 1/2 * (m1 + m2)

5. 重復步驟2-4,直到連續迭代中的T值差小於一個預定義的參數ΔT

算法實現

 1 void threshold(short** in_array, short** out_array, long height, long width, int delt_t)
 2 {
 3     double T = 0xff / 2.0;
 4     double m1 = 0.0, m2 = 0.0;
 5     int m1_num = 0, m2_num = 0;
 6 
 7     while(dabs(T, 0.5*(m1 + m2)) > delt_t){
 8         T = 0.5 * (m1 + m2);
 9         m1 = 0.0;
10         m2 = 0.0;
11         m1_num = 0;
12         m2_num = 0;
13 
14         for (int i = 0; i < height; i++){
15             for (int j = 0; j < width; j++){
16                 if (in_array[i][j] <= T){
17                     m1 += in_array[i][j];
18                     m1_num++;
19                 }
20                 else{
21                     m2 += in_array[i][j];
22                     m2_num++;
23                 }            
24             }
25         }
26         if (m1_num != 0)
27             m1 /= m1_num;
28         if (m2_num != 0)
29             m2 /= m2_num;
30         printf("%lf\n", T);
31     }
32     for (int i = 0; i < height; i++){
33         for (int j = 0; j < width; j++){
34             if (in_array[i][j] <= T)
35                 out_array[i][j] = 0xff;
36             else
37                 out_array[i][j] = 0x00;
38         }
39     }
40 }

測試結果

 從實驗結果看出,第二組閾值處理的效果並不好,因此考慮更優的算法實現

Otsu方法進行最佳全局閾值處理

閾值處理可視為一種統計決策理論問題,其目的是在把像素分配給兩個或多個組的過程中引入的平均誤差最小。這一問題有個閉合形式的解,稱為貝葉斯決策規則。

Otsu方法在類間方差最大的情況下是最佳的

算法執行流程

代碼實現

  1 double dabs(double a, double b)
  2 {
  3     if (a < b)
  4         return b-a;
  5     else
  6         return a-b;
  7 }
  8 
  9 void calculate_histogram(long height, long width, short **image, unsigned long histogram[])
 10 {
 11     short k;
 12     for(int i=0; i < height; i++){
 13         for(int j=0; j < width; j++){
 14             k = image[i][j];
 15             histogram[k] = histogram[k] + 1;
 16         }
 17     }
 18 }
 19 
 20 void calculate_pi(long height, long width, unsigned long histogram[], double pi[])
 21 {
 22     for (int i = 0; i < GRAY_LEVELS; ++i){
 23         pi[i] = (double)histogram[i] / (double)(height * width);
 24     }
 25 }
 26 
 27 double p1(int k, double pi[])
 28 {
 29     double sum = 0.0;
 30 
 31     for (int i = 0; i <= k; i++){
 32         sum += pi[i];
 33     }
 34 
 35     return sum;
 36 }
 37 
 38 double m(int k, double pi[])
 39 {
 40     double sum = 0.0;
 41 
 42     for (int i = 0; i <= k; i++){
 43         sum += i * pi[i];
 44     }
 45 
 46     return sum;
 47 }
 48 
 49 double calculate_sigma(int k, double mg, double pi[])
 50 {
 51     double p1k = p1(k, pi);
 52     double mk = m(k, pi);
 53 
 54     if (p1k < 1e-10 || (1 - p1k) < 1e-10)
 55         return 0.0;
 56     else
 57         return pow(mg * p1k - mk, 2) / (p1k * (1 - p1k));
 58 }
 59 
 60 void otsu(short** in_array, short** out_array, long height, long width)
 61 {
 62     unsigned long histogram[GRAY_LEVELS] = {};
 63     double pi[GRAY_LEVELS] = {};
 64     double sigma[GRAY_LEVELS] = {};
 65     double mg;
 66     short max_count = 0;
 67     short k = 0;
 68     double max_value = 0.0;
 69     double k_star;
 70 
 71     calculate_histogram(height, width, in_array, histogram);
 72     calculate_pi(height, width, histogram, pi);
 73     mg = m(GRAY_LEVELS-1, pi);
 74 
 75     for (int i = 0; i < GRAY_LEVELS; i++)
 76         sigma[i] = calculate_sigma(i, mg, pi);
 77 
 78     max_value = sigma[0];
 79     max_count = 1;
 80     k = 0;
 81     for (int i = 1; i < GRAY_LEVELS; i++){
 82         if (dabs(sigma[i], max_value) < 1e-10){
 83             k += i;
 84             max_count++;
 85         }
 86         else if (sigma[i] > max_value){
 87             max_value = sigma[i];
 88             max_count = 1;
 89             k = i;
 90         }
 91     }
 92     k_star = k / max_count;
 93 
 94     printf("%lf\n", k_star);
 95     for (int i = 0; i < height; i++){
 96         for (int j = 0; j < width; j++){
 97             if (in_array[i][j] <= k_star)
 98                 out_array[i][j] = 0x00;
 99             else
100                 out_array[i][j] = 0xff;
101         }
102     }
103 }

結果

 


注意!

本站转载的文章为个人学习借鉴使用,本站对版权不负任何法律责任。如果侵犯了您的隐私权益,请联系我们删除。



 
粤ICP备14056181号  © 2014-2021 ITdaan.com