Tesseract  3.02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
otsuthr.cpp
Go to the documentation of this file.
1 /**********************************************************************
2  * File: otsuthr.cpp
3  * Description: Simple Otsu thresholding for binarizing images.
4  * Author: Ray Smith
5  * Created: Fri Mar 07 12:31:01 PST 2008
6  *
7  * (C) Copyright 2008, Google Inc.
8  ** Licensed under the Apache License, Version 2.0 (the "License");
9  ** you may not use this file except in compliance with the License.
10  ** You may obtain a copy of the License at
11  ** http://www.apache.org/licenses/LICENSE-2.0
12  ** Unless required by applicable law or agreed to in writing, software
13  ** distributed under the License is distributed on an "AS IS" BASIS,
14  ** WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
15  ** See the License for the specific language governing permissions and
16  ** limitations under the License.
17  *
18  **********************************************************************/
19 
20 #include <string.h>
21 #include "otsuthr.h"
22 
23 namespace tesseract {
24 
25 // Compute the Otsu threshold(s) for the given image rectangle, making one
26 // for each channel. Each channel is always one byte per pixel.
27 // Returns an array of threshold values and an array of hi_values, such
28 // that a pixel value >threshold[channel] is considered foreground if
29 // hi_values[channel] is 0 or background if 1. A hi_value of -1 indicates
30 // that there is no apparent foreground. At least one hi_value will not be -1.
31 // Delete thresholds and hi_values with delete [] after use.
32 void OtsuThreshold(const unsigned char* imagedata,
33  int bytes_per_pixel, int bytes_per_line,
34  int left, int top, int width, int height,
35  int** thresholds, int** hi_values) {
36  // Of all channels with no good hi_value, keep the best so we can always
37  // produce at least one answer.
38  int best_hi_value = 1;
39  int best_hi_index = 0;
40  bool any_good_hivalue = false;
41  double best_hi_dist = 0.0;
42  *thresholds = new int[bytes_per_pixel];
43  *hi_values = new int[bytes_per_pixel];
44 
45  for (int ch = 0; ch < bytes_per_pixel; ++ch) {
46  (*thresholds)[ch] = -1;
47  (*hi_values)[ch] = -1;
48  // Compute the histogram of the image rectangle.
49  int histogram[kHistogramSize];
50  HistogramRect(imagedata + ch, bytes_per_pixel, bytes_per_line,
51  left, top, width, height, histogram);
52  int H;
53  int best_omega_0;
54  int best_t = OtsuStats(histogram, &H, &best_omega_0);
55  if (best_omega_0 == 0 || best_omega_0 == H) {
56  // This channel is empty.
57  continue;
58  }
59  // To be a convincing foreground we must have a small fraction of H
60  // or to be a convincing background we must have a large fraction of H.
61  // In between we assume this channel contains no thresholding information.
62  int hi_value = best_omega_0 < H * 0.5;
63  (*thresholds)[ch] = best_t;
64  if (best_omega_0 > H * 0.75) {
65  any_good_hivalue = true;
66  (*hi_values)[ch] = 0;
67  } else if (best_omega_0 < H * 0.25) {
68  any_good_hivalue = true;
69  (*hi_values)[ch] = 1;
70  } else {
71  // In case all channels are like this, keep the best of the bad lot.
72  double hi_dist = hi_value ? (H - best_omega_0) : best_omega_0;
73  if (hi_dist > best_hi_dist) {
74  best_hi_dist = hi_dist;
75  best_hi_value = hi_value;
76  best_hi_index = ch;
77  }
78  }
79  }
80  if (!any_good_hivalue) {
81  // Use the best of the ones that were not good enough.
82  (*hi_values)[best_hi_index] = best_hi_value;
83  }
84 }
85 
86 // Compute the histogram for the given image rectangle, and the given
87 // channel. (Channel pointed to by imagedata.) Each channel is always
88 // one byte per pixel.
89 // Bytes per pixel is used to skip channels not being
90 // counted with this call in a multi-channel (pixel-major) image.
91 // Histogram is always a kHistogramSize(256) element array to count
92 // occurrences of each pixel value.
93 void HistogramRect(const unsigned char* imagedata,
94  int bytes_per_pixel, int bytes_per_line,
95  int left, int top, int width, int height,
96  int* histogram) {
97  int bottom = top + height;
98  memset(histogram, 0, sizeof(*histogram) * kHistogramSize);
99  const unsigned char* pixels = imagedata +
100  top * bytes_per_line +
101  left * bytes_per_pixel;
102  for (int y = top; y < bottom; ++y) {
103  for (int x = 0; x < width; ++x) {
104  ++histogram[pixels[x * bytes_per_pixel]];
105  }
106  pixels += bytes_per_line;
107  }
108 }
109 
110 // Compute the Otsu threshold(s) for the given histogram.
111 // Also returns H = total count in histogram, and
112 // omega0 = count of histogram below threshold.
113 int OtsuStats(const int* histogram, int* H_out, int* omega0_out) {
114  int H = 0;
115  double mu_T = 0.0;
116  for (int i = 0; i < kHistogramSize; ++i) {
117  H += histogram[i];
118  mu_T += static_cast<double>(i) * histogram[i];
119  }
120 
121  // Now maximize sig_sq_B over t.
122  // http://www.ctie.monash.edu.au/hargreave/Cornall_Terry_328.pdf
123  int best_t = -1;
124  int omega_0, omega_1;
125  int best_omega_0 = 0;
126  double best_sig_sq_B = 0.0;
127  double mu_0, mu_1, mu_t;
128  omega_0 = 0;
129  mu_t = 0.0;
130  for (int t = 0; t < kHistogramSize - 1; ++t) {
131  omega_0 += histogram[t];
132  mu_t += t * static_cast<double>(histogram[t]);
133  if (omega_0 == 0)
134  continue;
135  omega_1 = H - omega_0;
136  if (omega_1 == 0)
137  break;
138  mu_0 = mu_t / omega_0;
139  mu_1 = (mu_T - mu_t) / omega_1;
140  double sig_sq_B = mu_1 - mu_0;
141  sig_sq_B *= sig_sq_B * omega_0 * omega_1;
142  if (best_t < 0 || sig_sq_B > best_sig_sq_B) {
143  best_sig_sq_B = sig_sq_B;
144  best_t = t;
145  best_omega_0 = omega_0;
146  }
147  }
148  if (H_out != NULL) *H_out = H;
149  if (omega0_out != NULL) *omega0_out = best_omega_0;
150  return best_t;
151 }
152 
153 } // namespace tesseract.