(patched) gocr-0.44
[swftools.git] / lib / gocr / otsu.c
diff --git a/lib/gocr/otsu.c b/lib/gocr/otsu.c
new file mode 100644 (file)
index 0000000..50c754a
--- /dev/null
@@ -0,0 +1,284 @@
+/*
+This is a Optical-Character-Recognition program
+Copyright (C) 2000-2007  Joerg Schulenburg
+
+This program is free software; you can redistribute it and/or
+modify it under the terms of the GNU General Public License
+as published by the Free Software Foundation; either version 2
+of the License, or (at your option) any later version.
+
+This program is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU General Public License for more details.
+
+You should have received a copy of the GNU General Public License
+along with this program; if not, write to the Free Software
+Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
+
+ see README for EMAIL-address
+
+ the following code was send by Ryan Dibble <dibbler@umich.edu>
+  The algorithm is very simple but works good hopefully.
+  Compare the grayscale histogram with a mass density diagram:
+  I think the algorithm is a kind of
+  divide a body into two parts in a way that the mass 
+  centers have the largest distance from each other,
+  the function is weighted in a way that same masses have a advantage
+  
+  - otsu algorithm is failing on diskrete multi color images
+  
+  TODO:
+    RGB: do the same with all colors (CMYG?) seperately 
+
+    test: hardest case = two colors
+       bbg: test done, using a two color gray file. Output:
+       # threshold: Value = 43 gmin=43 gmax=188
+
+ my changes:
+   - float -> double
+   - debug option added (vvv & 1..2)
+   - **image => *image,  &image[i][1] => &image[i*cols+1]
+   - do only count pixels near contrast regions
+     this makes otsu much better for shadowed fonts or multi colored text 
+     on white background
+
+ (m) Joerg Schulenburg (see README for email address)
+
+ ToDo:
+   - measure contrast
+   - detect low-contrast regions
+
+ */
+
+#include <stdio.h>
+#include <string.h>
+
+#define Abs(x) ((x<0)?-(x):x)
+
+/*======================================================================*/
+/* global thresholding routine                                     */
+/*   takes a 2D unsigned char array pointer, number of rows, and        */
+/*   number of cols in the array. returns the value of the threshold    */
+/*======================================================================*/
+int
+otsu (unsigned char *image, int rows, int cols,
+      int x0, int y0, int dx, int dy, int vvv) {
+
+  unsigned char *np;    // pointer to position in the image we are working with
+  unsigned char op1, op2;   // predecessor of pixel *np (start value)
+  int maxc=0;           // maximum contrast (start value) 
+  int thresholdValue=1; // value we will threshold at
+  int ihist[256];       // image histogram
+  int chist[256];       // contrast histogram
+
+  int i, j, k;          // various counters
+  int is, i1, i2, ns, n1, n2, gmin, gmax;
+  double m1, m2, sum, csum, fmax, sb;
+
+  // zero out histogram ...
+  memset(ihist, 0, sizeof(ihist));
+  memset(chist, 0, sizeof(chist));
+  op1=op2=0;
+
+  gmin=255; gmax=0; k=dy/512+1;
+  // v0.43 first get max contrast, dont do it together with next step
+  //  because it failes if we have pattern as background (on top)
+  for (i =  0; i <  dy ; i+=k) {
+    np = &image[(y0+i)*cols+x0];
+    for (j = 0; j < dx ; j++) {
+      ihist[*np]++;
+      if(*np > gmax) gmax=*np;
+      if(*np < gmin) gmin=*np;
+      if (Abs(*np-op1)>maxc) maxc=Abs(*np-op1); /* new maximum contrast */
+      if (Abs(*np-op2)>maxc) maxc=Abs(*np-op2); /* new maximum contrast */
+      /* we hope that maxc will be find its maximum very fast */
+      op2=op1;    /* shift old pixel to next older */
+      op1=*np;     /* store old pixel for contrast check */
+      np++;       /* next pixel */
+    }
+  }
+
+  // generate the histogram
+  // Aug06 images with large white or black homogeneous
+  //   areas give bad results, so we only add pixels on contrast edges
+  for (i =  0; i <  dy ; i+=k) {
+    np = &image[(y0+i)*cols+x0];
+    for (j = 0; j < dx ; j++) {
+      if (Abs(*np-op1)>maxc/4
+       || Abs(*np-op2)>maxc/4)
+         chist[*np]++; // count only relevant pixels
+      op2=op1;    /* shift old pixel to next older */
+      op1=*np;     /* store old pixel for contrast check */
+      np++;       /* next pixel */
+    }
+  }
+
+  // set up everything
+  sum = csum = 0.0;
+  ns = 0;
+  is = 0;
+  
+  for (k = 0; k <= 255; k++) {
+    sum += (double) k * (double) chist[k];  /* x*f(x) cmass moment */
+    ns  += chist[k];                        /*  f(x)    cmass      */
+    is  += ihist[k];                        /*  f(x)    imass      */
+    // Debug: output to out_hist.dat?
+    // fprintf(stderr,"\chistogram %3d %6d (brightness weight)", k, ihist[k]);
+  }
+
+  if (!ns) {
+    // if n has no value we have problems...
+    fprintf (stderr, "NOT NORMAL, thresholdValue = 160\n");
+    return (160);
+  }
+
+  // ToDo: only care about extremas in a 3 pixel environment
+  //       check if there are more than 2 mass centers (more colors)
+  //       return object colors and color radius instead of threshold value
+  //        also the reagion, where colored objects are found
+  //       what if more than one background color? no otsu at all?
+  //       whats background? box with lot of other boxes in it
+  //       threshold each box (examples/invers.png,colors.png) 
+  //  get maximum white and minimum black pixel color (possible range)
+  //    check range between them for low..high contrast ???
+  // typical scenes (which must be covered): 
+  //    - white page with text of different colors (gray values)
+  //    - binear page: background (gray=1) + black text (gray=0)
+  //    - text mixed with big (dark) images
+  //  ToDo: recursive clustering for maximum multipol moments?
+  //  idea: normalize ihist to max=1024 before otsu?
+  
+  // do the otsu global thresholding method
+
+  if ((vvv&1)) // Debug
+     fprintf(stderr,"# threshold: value ihist chist mass_dipol_moment\n");
+  fmax = -1.0;
+  n1 = 0;
+  for (k = 0; k < 255; k++) {
+    n1 += chist[k];          // left  mass (integration)
+    if (!n1) continue;       // we need at least one foreground pixel
+    n2 = ns - n1;            // right mass (num pixels - left mass)
+    if (n2 == 0) break;      // we need at least one background pixel
+    csum += (double) k *chist[k];  // left mass moment
+    m1 =        csum  / n1;        // left  mass center (black chars)
+    m2 = (sum - csum) / n2;        // right mass center (white background)
+    // max. dipol moment?
+    // orig: sb = (double) n1 *(double) n2 * (m1 - m2) * (m1 - m2);
+    sb = (double) n1 *(double) n2 * (m2 - m1); // seems to be better Aug06
+    /* bbg: note: can be optimized. */
+    if (sb > fmax) {
+      fmax = sb;
+      thresholdValue = k + 1;
+      // thresholdValue = (m1 + 3 * m2) / 4;
+    }
+    if ((vvv&1) && ihist[k]) // Debug
+     fprintf(stderr,"# threshold: %3d %6d %6d %8.2f\n",
+       k, ihist[k], chist[k],
+       sb/(dx*dy));  /* normalized dipol moment */
+  }
+  // ToDo: error = left/right point where sb is 90% of maximum?
+  // now we count all pixels for background detection
+  i1 = 0;
+  for (k = 0; k < thresholdValue; k++) {
+    i1 += ihist[k];          // left  mass (integration)
+  }
+  i2 = is - i1;            // right mass (num pixels - left mass)
+
+  // at this point we have our thresholding value
+  // black_char: value<cs,  white_background: value>=cs
+  
+  // can it happen? check for sureness
+  if (thresholdValue >  gmax) {
+    fprintf(stderr,"# threshold: Value >gmax\n");
+    thresholdValue = gmax;
+  }
+  if (thresholdValue <= gmin) {
+    fprintf(stderr,"# threshold: Value<=gmin\n");
+    thresholdValue = gmin+1;
+  }
+
+  // debug code to display thresholding values
+  if ( vvv & 1 )
+  fprintf(stderr,"# threshold: Value = %d gmin=%d gmax=%d cmax=%d"
+                 " i= %d %d\n",
+     thresholdValue, gmin, gmax, maxc, i1, i2);
+
+  if (i1>=4*i2) { // black>=4*white, obviously black is background
+    if ( vvv & 1 )
+      fprintf(stderr,"# threshold: invert the image\n");
+    // we do inversion here (no data lost)
+    for (i =  0; i <  dy ; i++) {
+      np = &image[(y0+i)*cols+x0];
+      for (j = 0; j < dx ; j++) {
+        *np=255-*np;
+        np++;       /* next pixel */
+      }
+    }
+    thresholdValue=255-thresholdValue+1;
+  }
+
+  return(thresholdValue); 
+  /* range: 0 < thresholdValue <= 255, example: 1 on b/w images */
+  /* 0..threshold-1 is foreground */
+  /* threshold..255 is background */
+  /* ToDo:  min=blackmasscenter/2,thresh,max=(whitemasscenter+255)/2 */
+}
+
+/*======================================================================*/
+/* thresholding the image  (set threshold to 128+32=160=0xA0)           */
+/*   now we have a fixed thresholdValue good to recognize on gray image */
+/*   - so lower bits can used for other things (bad design?)            */
+/* ToDo: different foreground colors, gray on black/white background    */
+/*======================================================================*/
+int
+thresholding (unsigned char *image, int rows, int cols,
+  int x0, int y0, int dx, int dy, int thresholdValue) {
+
+  unsigned char *np; // pointer to position in the image we are working with
+
+  int i, j;          // various counters
+  int gmin=255,gmax=0;
+  int nmin=255,nmax=0;
+
+  // calculate min/max (twice?)
+  for (i = y0 + 1; i < y0 + dy - 1; i++) {
+    np = &image[i*cols+x0+1];
+    for (j = x0 + 1; j < x0 + dx - 1; j++) {
+      if(*np > gmax) gmax=*np;
+      if(*np < gmin) gmin=*np;
+      np++; /* next pixel */
+    }
+  }
+  
+  /* allowed_threshold=gmin+1..gmax v0.43 */
+  if (thresholdValue<=gmin || thresholdValue>gmax){
+    thresholdValue=(gmin+gmax+1)/2; /* range=0..1 -> threshold=1 */
+    fprintf(stderr,"# thresholdValue out of range %d..%d, reset to %d\n",
+     gmin, gmax, thresholdValue);
+  }
+
+  /* b/w: min=0,tresh=1,max=1 v0.43 */
+  // actually performs the thresholding of the image...
+  //  later: grayvalues should also be used, only rescaling threshold=160=0xA0
+  for (i = y0; i < y0+dy; i++) {
+    np = &image[i*cols+x0];
+    for (j = x0; j < x0+dx; j++) {
+      *np = (unsigned char) (*np >= thresholdValue ?
+         (255-(gmax - *np)* 80/(gmax - thresholdValue + 1)) :
+         (  0+(*np - gmin)*150/(thresholdValue - gmin    )) );
+      if(*np > nmax) nmax=*np;
+      if(*np < nmin) nmin=*np;
+      np++;
+    }
+  }
+
+  // fprintf(stderr,"# thresholding:  nmin=%d nmax=%d\n", nmin, nmax);
+
+  return(128+32); // return the new normalized threshold value
+  /*   0..159 is foreground */
+  /* 160..255 is background */
+}
+