binarize.c

00001 /*Copyright (C) 2008-2009  Timothy B. Terriberry (tterribe@xiph.org)
00002   You can redistribute this library and/or modify it under the terms of the
00003    GNU Lesser General Public License as published by the Free Software
00004    Foundation; either version 2.1 of the License, or (at your option) any later
00005    version.*/
00006 #include <stdlib.h>
00007 #include <math.h>
00008 #include <string.h>
00009 #include "util.h"
00010 #include "image.h"
00011 #include "binarize.h"
00012 
00013 #if 0
00014 /*Binarization based on~\cite{GPP06}.
00015   @ARTICLE{GPP06,
00016     author="Basilios Gatos and Ioannis E. Pratikakis and Stavros J. Perantonis",
00017     title="Adaptive Degraded Document Image Binarization",
00018     journal="Pattern Recognition",
00019     volume=39,
00020     number=3,
00021     pages="317-327",
00022     month=Mar,
00023     year=2006
00024   }*/
00025 
00026 #if 0
00027 /*Applies a 5x5 Wiener filter to the image, in-place, emphasizing differences
00028    where the local variance is small, and de-emphasizing them where it is
00029    large.*/
00030 void qr_wiener_filter(unsigned char *_img,int _width,int _height){
00031   unsigned           *m_buf[8];
00032   unsigned           *sn2_buf[8];
00033   unsigned char       g;
00034   int                 x;
00035   int                 y;
00036   if(_width<=0||_height<=0)return;
00037   m_buf[0]=(unsigned *)malloc((_width+4<<3)*sizeof(*m_buf));
00038   sn2_buf[0]=(unsigned *)malloc((_width+4<<3)*sizeof(*sn2_buf));
00039   for(y=1;y<8;y++){
00040     m_buf[y]=m_buf[y-1]+_width+4;
00041     sn2_buf[y]=sn2_buf[y-1]+_width+4;
00042   }
00043   for(y=-4;y<_height;y++){
00044     unsigned *pm;
00045     unsigned *psn2;
00046     int       i;
00047     int       j;
00048     pm=m_buf[y+2&7];
00049     psn2=sn2_buf[y+2&7];
00050     for(x=-4;x<_width;x++){
00051       unsigned m;
00052       unsigned m2;
00053       m=m2=0;
00054       if(y>=0&&y<_height-4&&x>=0&&x<_width-4)for(i=0;i<5;i++)for(j=0;j<5;j++){
00055         g=_img[(y+i)*_width+x+j];
00056         m+=g;
00057         m2+=g*g;
00058       }
00059       else for(i=0;i<5;i++)for(j=0;j<5;j++){
00060         g=_img[QR_CLAMPI(0,y+i,_height-1)*_width+QR_CLAMPI(0,x+j,_width-1)];
00061         m+=g;
00062         m2+=g*g;
00063       }
00064       pm[x+4]=m;
00065       psn2[x+4]=(m2*25-m*m);
00066     }
00067     pm=m_buf[y&7];
00068     if(y>=0)for(x=0;x<_width;x++){
00069       int sn2;
00070       sn2=sn2_buf[y&7][x+2];
00071       if(sn2){
00072         int vn3;
00073         int m;
00074         /*Gatos et al. give the expression
00075             mu+(s2-v2)*(g-mu)/s2 ,
00076            which we reduce to
00077             mu+(s2-v2)*g/s2-(s2-v2)*mu/s2 ,
00078             g-(v2/s2)*g+(v2/s2)*mu ,
00079             g+(mu-g)*(v2/s2) .
00080            However, s2 is much noisier than v2, and dividing by it often gives
00081             extremely large adjustments, causing speckle near edges.
00082            Therefore we limit the ratio (v2/s2) to lie between 0 and 1.*/
00083         vn3=0;
00084         for(i=-2;i<3;i++){
00085           psn2=sn2_buf[y+i&7];
00086           for(j=0;j<5;j++)vn3+=psn2[x+j];
00087         }
00088         m=m_buf[y&7][x+2];
00089         vn3=vn3+1023>>10;
00090         sn2=25*sn2+1023>>10;
00091         if(vn3<sn2){
00092           int a;
00093           g=_img[y*_width+x];
00094           a=(m-25*g)*vn3;
00095           sn2*=25;
00096           _img[y*_width+x]=QR_CLAMP255(g+QR_DIVROUND(a,sn2));
00097         }
00098         else _img[y*_width+x]=(unsigned char)(((m<<1)+25)/50);
00099       }
00100     }
00101   }
00102   free(sn2_buf[0]);
00103   free(m_buf[0]);
00104 }
00105 
00106 #else
00107 /*Applies a 3x3 Wiener filter to the image, in-place, emphasizing differences
00108    where the local variance is small, and de-emphasizing them where it is
00109    large.*/
00110 void qr_wiener_filter(unsigned char *_img,int _width,int _height){
00111   unsigned           *m_buf[4];
00112   unsigned           *sn2_buf[4];
00113   unsigned char       g;
00114   int                 x;
00115   int                 y;
00116   if(_width<=0||_height<=0)return;
00117   m_buf[0]=(unsigned *)malloc((_width+2<<2)*sizeof(*m_buf));
00118   sn2_buf[0]=(unsigned *)malloc((_width+2<<2)*sizeof(*sn2_buf));
00119   for(y=1;y<4;y++){
00120     m_buf[y]=m_buf[y-1]+_width+2;
00121     sn2_buf[y]=sn2_buf[y-1]+_width+2;
00122   }
00123   for(y=-2;y<_height;y++){
00124     unsigned *pm;
00125     unsigned *psn2;
00126     int       i;
00127     int       j;
00128     pm=m_buf[y+1&3];
00129     psn2=sn2_buf[y+1&3];
00130     for(x=-2;x<_width;x++){
00131       unsigned m;
00132       unsigned m2;
00133       m=m2=0;
00134       if(y>=0&&y<_height-2&&x>=0&&x<_width-2)for(i=0;i<3;i++)for(j=0;j<3;j++){
00135         g=_img[(y+i)*_width+x+j];
00136         m+=g;
00137         m2+=g*g;
00138       }
00139       else for(i=0;i<3;i++)for(j=0;j<3;j++){
00140         g=_img[QR_CLAMPI(0,y+i,_height-1)*_width+QR_CLAMPI(0,x+j,_width-1)];
00141         m+=g;
00142         m2+=g*g;
00143       }
00144       pm[x+2]=m;
00145       psn2[x+2]=(m2*9-m*m);
00146     }
00147     pm=m_buf[y&3];
00148     if(y>=0)for(x=0;x<_width;x++){
00149       int sn2;
00150       sn2=sn2_buf[y&3][x+1];
00151       if(sn2){
00152         int m;
00153         int vn3;
00154         /*Gatos et al. give the expression
00155             mu+(s2-v2)*(g-mu)/s2 ,
00156            which we reduce to
00157             mu+(s2-v2)*g/s2-(s2-v2)*mu/s2 ,
00158             g-(v2/s2)*g+(v2/s2)*mu ,
00159             g+(mu-g)*(v2/s2) .
00160            However, s2 is much noisier than v2, and dividing by it often gives
00161             extremely large adjustments, causing speckle near edges.
00162            Therefore we limit the ratio (v2/s2) to lie between 0 and 1.*/
00163         vn3=0;
00164         for(i=-1;i<2;i++){
00165           psn2=sn2_buf[y+i&3];
00166           for(j=0;j<3;j++)vn3+=psn2[x+j];
00167         }
00168         m=m_buf[y&3][x+1];
00169         vn3=vn3+31>>5;
00170         sn2=9*sn2+31>>5;
00171         if(vn3<sn2){
00172           int a;
00173           g=_img[y*_width+x];
00174           a=m-9*g;
00175           sn2*=9;
00176           _img[y*_width+x]=QR_CLAMP255(g+QR_DIVROUND(a,sn2));
00177         }
00178         else _img[y*_width+x]=(unsigned char)(((m<<1)+9)/18);
00179       }
00180     }
00181   }
00182   free(sn2_buf[0]);
00183   free(m_buf[0]);
00184 }
00185 #endif
00186 
00187 /*Computes a (conservative) foreground mask using the adaptive binarization
00188    threshold given in~\cite{SP00}, but knocking the threshold parameter down to
00189    k=0.2.
00190   Note on dynamic range: we assume _width*_height<=0x1000000 (24 bits).
00191   Returns the average background value.
00192   @ARTICLE{SP00,
00193     author="Jaakko J. Sauvola and Matti Pietik\"{a}inen",
00194     title="Adaptive Document Image Binarization",
00195     volume=33,
00196     number=2,
00197     pages="225--236",
00198     month=Feb,
00199     year=2000
00200   }*/
00201 static void qr_sauvola_mask(unsigned char *_mask,unsigned *_b,int *_nb,
00202  const unsigned char *_img,int _width,int _height){
00203   unsigned  b;
00204   int       nb;
00205   b=0;
00206   nb=0;
00207   if(_width>0&&_height>0){
00208     unsigned *col_sums;
00209     unsigned *col2_sums;
00210     int       logwindw;
00211     int       logwindh;
00212     int       windw;
00213     int       windh;
00214     int       y0offs;
00215     int       y1offs;
00216     unsigned  g;
00217     unsigned  g2;
00218     int       x;
00219     int       y;
00220     /*We keep the window size fairly large to ensure it doesn't fit completely
00221        inside the center of a finder pattern of a version 1 QR code at full
00222        resolution.*/
00223     for(logwindw=4;logwindw<8&&(1<<logwindw)<(_width+7>>3);logwindw++);
00224     for(logwindh=4;logwindh<8&&(1<<logwindh)<(_height+7>>3);logwindh++);
00225     windw=1<<logwindw;
00226     windh=1<<logwindh;
00227     col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
00228     col2_sums=(unsigned *)malloc(_width*sizeof(*col2_sums));
00229     /*Initialize sums down each column.*/
00230     for(x=0;x<_width;x++){
00231       g=_img[x];
00232       g2=g*g;
00233       col_sums[x]=(g<<logwindh-1)+g;
00234       col2_sums[x]=(g2<<logwindh-1)+g2;
00235     }
00236     for(y=1;y<(windh>>1);y++){
00237       y1offs=QR_MINI(y,_height-1)*_width;
00238       for(x=0;x<_width;x++){
00239         g=_img[y1offs+x];
00240         col_sums[x]+=g;
00241         col2_sums[x]+=g*g;
00242       }
00243     }
00244     for(y=0;y<_height;y++){
00245       unsigned m;
00246       unsigned m2;
00247       int      x0;
00248       int      x1;
00249       /*Initialize the sums over the window.*/
00250       m=(col_sums[0]<<logwindw-1)+col_sums[0];
00251       m2=(col2_sums[0]<<logwindw-1)+col2_sums[0];
00252       for(x=1;x<(windw>>1);x++){
00253         x1=QR_MINI(x,_width-1);
00254         m+=col_sums[x1];
00255         m2+=col2_sums[x1];
00256       }
00257       for(x=0;x<_width;x++){
00258         int d;
00259         /*Perform the test against the threshold T = (m/n)*(1+k*(s/R-1)),
00260            where n=windw*windh, s=sqrt((m2-(m*m)/n)/n), and R=128.
00261           We don't actually compute the threshold directly, as that would
00262            require a square root.
00263           Instead we perform the equivalent test:
00264            (m/n)*(m/n)*(m2/n-(m/n)*(m/n))/16 > (((1/k)*g-((1-k)/k)*(m/n))*32)**2
00265           R is split up across each side of the inequality to maximize the
00266            dynamic range available for the right hand side, which requires
00267            31 bits in the worst case.*/
00268         /*(m/n)*(1+(1/5)*(sqrt((m2-m*m/n)/n)/128-1)) > g
00269           m*(1+(1/5)*(sqrt((m2-m*m/n)/n)/128-1)) > g*n
00270           m*sqrt((m2-m*m/n)/n) > 5*g*n-4*m<<7
00271           m*m*(m2*n-m*m) > (5*g*n-4*m<<7)**2*n*n || 5*g*n-4*m < 0 */
00272         g=_img[y*_width+x];
00273         d=(5*g<<logwindw+logwindh)-4*m;
00274         if(d>=0){
00275           unsigned mm;
00276           unsigned mms2;
00277           unsigned d2;
00278           mm=(m>>logwindw)*(m>>logwindh);
00279           mms2=(m2-mm>>logwindw+logwindh)*(mm>>logwindw+logwindh)+15>>4;
00280           d2=d>>logwindw+logwindh-5;
00281           d2*=d2;
00282           if(d2>=mms2){
00283             /*Update the background average.*/
00284             b+=g;
00285             nb++;
00286             _mask[y*_width+x]=0;
00287           }
00288           else _mask[y*_width+x]=0xFF;
00289         }
00290         else _mask[y*_width+x]=0xFF;
00291         /*Update the window sums.*/
00292         if(x+1<_width){
00293           x0=QR_MAXI(0,x-(windw>>1));
00294           x1=QR_MINI(x+(windw>>1),_width-1);
00295           m+=col_sums[x1]-col_sums[x0];
00296           m2+=col2_sums[x1]-col2_sums[x0];
00297         }
00298       }
00299       /*Update the column sums.*/
00300       if(y+1<_height){
00301         y0offs=QR_MAXI(0,y-(windh>>1))*_width;
00302         y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
00303         for(x=0;x<_width;x++){
00304           g=_img[y0offs+x];
00305           col_sums[x]-=g;
00306           col2_sums[x]-=g*g;
00307           g=_img[y1offs+x];
00308           col_sums[x]+=g;
00309           col2_sums[x]+=g*g;
00310         }
00311       }
00312     }
00313     free(col2_sums);
00314     free(col_sums);
00315   }
00316   *_b=b;
00317   *_nb=nb;
00318 }
00319 
00320 /*Interpolates a background image given the source and a conservative
00321    foreground mask.
00322   If the current window contains no foreground pixels, the average background
00323    value over the whole image is used.
00324   Note on dynamic range: we assume _width*_height<=0x8000000 (23 bits).
00325   Returns the average difference between the foreground and the interpolated
00326    background.*/
00327 static void qr_interpolate_background(unsigned char *_dst,
00328  int *_delta,int *_ndelta,const unsigned char *_img,const unsigned char *_mask,
00329  int _width,int _height,unsigned _b,int _nb){
00330   int delta;
00331   int ndelta;
00332   delta=ndelta=0;
00333   if(_width>0&&_height>0){
00334     unsigned *col_sums;
00335     unsigned *ncol_sums;
00336     int       logwindw;
00337     int       logwindh;
00338     int       windw;
00339     int       windh;
00340     int       y0offs;
00341     int       y1offs;
00342     unsigned  b;
00343     unsigned  g;
00344     int       x;
00345     int       y;
00346     b=_nb>0?((_b<<1)+_nb)/(_nb<<1):0xFF;
00347     for(logwindw=4;logwindw<8&&(1<<logwindw)<(_width+15>>4);logwindw++);
00348     for(logwindh=4;logwindh<8&&(1<<logwindh)<(_height+15>>4);logwindh++);
00349     windw=1<<logwindw;
00350     windh=1<<logwindh;
00351     col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
00352     ncol_sums=(unsigned *)malloc(_width*sizeof(*ncol_sums));
00353     /*Initialize sums down each column.*/
00354     for(x=0;x<_width;x++){
00355       if(!_mask[x]){
00356         g=_img[x];
00357         col_sums[x]=(g<<logwindh-1)+g;
00358         ncol_sums[x]=(1<<logwindh-1)+1;
00359       }
00360       else col_sums[x]=ncol_sums[x]=0;
00361     }
00362     for(y=1;y<(windh>>1);y++){
00363       y1offs=QR_MINI(y,_height-1)*_width;
00364       for(x=0;x<_width;x++)if(!_mask[y1offs+x]){
00365         col_sums[x]+=_img[y1offs+x];
00366         ncol_sums[x]++;
00367       }
00368     }
00369     for(y=0;y<_height;y++){
00370       unsigned n;
00371       unsigned m;
00372       int      x0;
00373       int      x1;
00374       /*Initialize the sums over the window.*/
00375       m=(col_sums[0]<<logwindw-1)+col_sums[0];
00376       n=(ncol_sums[0]<<logwindw-1)+ncol_sums[0];
00377       for(x=1;x<(windw>>1);x++){
00378         x1=QR_MINI(x,_width-1);
00379         m+=col_sums[x1];
00380         n+=ncol_sums[x1];
00381       }
00382       for(x=0;x<_width;x++){
00383         if(!_mask[y*_width+x])g=_img[y*_width+x];
00384         else{
00385           g=n>0?((m<<1)+n)/(n<<1):b;
00386           delta+=(int)g-_img[y*_width+x];
00387           ndelta++;
00388         }
00389         _dst[y*_width+x]=(unsigned char)g;
00390         /*Update the window sums.*/
00391         if(x+1<_width){
00392           x0=QR_MAXI(0,x-(windw>>1));
00393           x1=QR_MINI(x+(windw>>1),_width-1);
00394           m+=col_sums[x1]-col_sums[x0];
00395           n+=ncol_sums[x1]-ncol_sums[x0];
00396         }
00397       }
00398       /*Update the column sums.*/
00399       if(y+1<_height){
00400         y0offs=QR_MAXI(0,y-(windh>>1))*_width;
00401         y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
00402         for(x=0;x<_width;x++){
00403           if(!_mask[y0offs+x]){
00404             col_sums[x]-=_img[y0offs+x];
00405             ncol_sums[x]--;
00406           }
00407           if(!_mask[y1offs+x]){
00408             col_sums[x]+=_img[y1offs+x];
00409             ncol_sums[x]++;
00410           }
00411         }
00412       }
00413     }
00414     free(ncol_sums);
00415     free(col_sums);
00416   }
00417   *_delta=delta;
00418   *_ndelta=ndelta;
00419 }
00420 
00421 /*Parameters of the logistic sigmoid function that defines the threshold based
00422    on the background intensity.
00423   They should all be between 0 and 1.*/
00424 #define QR_GATOS_Q  (0.7)
00425 #define QR_GATOS_P1 (0.5)
00426 #define QR_GATOS_P2 (0.8)
00427 
00428 /*Compute the final binarization mask according to Gatos et al.'s
00429    method~\cite{GPP06}.*/
00430 static void qr_gatos_mask(unsigned char *_mask,const unsigned char *_img,
00431  const unsigned char *_background,int _width,int _height,
00432  unsigned _b,int _nb,int _delta,int _ndelta){
00433   unsigned thresh[256];
00434   unsigned g;
00435   double   delta;
00436   double   b;
00437   int      x;
00438   int      y;
00439   /*Construct a lookup table for the thresholds.
00440     This bit uses floating point, but doesn't need to do much calculation, so
00441      emulation should be fine.*/
00442   b=_nb>0?(_b+0.5)/_nb:0xFF;
00443   delta=_ndelta>0?(_delta+0.5)/_ndelta:0xFF;
00444   for(g=0;g<256;g++){
00445     double d;
00446     d=QR_GATOS_Q*delta*(QR_GATOS_P2+(1-QR_GATOS_P2)/
00447      (1+exp(2*(1+QR_GATOS_P1)/(1-QR_GATOS_P1)-4*g/(b*(1-QR_GATOS_P1)))));
00448     if(d<1)d=1;
00449     else if(d>0xFF)d=0xFF;
00450     thresh[g]=(unsigned)floor(d);
00451   }
00452   /*Apply the adaptive threshold.*/
00453   for(y=0;y<_height;y++)for(x=0;x<_width;x++){
00454     g=_background[y*_width+x];
00455     /*_background[y*_width+x]=thresh[g];*/
00456     _mask[y*_width+x]=(unsigned char)(-(g-_img[y*_width+x]>thresh[g])&0xFF);
00457   }
00458   /*{
00459     FILE *fout;
00460     fout=fopen("thresh.png","wb");
00461     image_write_png(_background,_width,_height,fout);
00462     fclose(fout);
00463   }*/
00464 }
00465 
00466 /*Binarizes a grayscale image.*/
00467 void qr_binarize(unsigned char *_img,int _width,int _height){
00468   unsigned char *mask;
00469   unsigned char *background;
00470   unsigned       b;
00471   int            nb;
00472   int            delta;
00473   int            ndelta;
00474   /*qr_wiener_filter(_img,_width,_height);
00475   {
00476     FILE *fout;
00477     fout=fopen("wiener.png","wb");
00478     image_write_png(_img,_width,_height,fout);
00479     fclose(fout);
00480   }*/
00481   mask=(unsigned char *)malloc(_width*_height*sizeof(*mask));
00482   qr_sauvola_mask(mask,&b,&nb,_img,_width,_height);
00483   /*{
00484     FILE *fout;
00485     fout=fopen("foreground.png","wb");
00486     image_write_png(mask,_width,_height,fout);
00487     fclose(fout);
00488   }*/
00489   background=(unsigned char *)malloc(_width*_height*sizeof(*mask));
00490   qr_interpolate_background(background,&delta,&ndelta,
00491    _img,mask,_width,_height,b,nb);
00492   /*{
00493     FILE *fout;
00494     fout=fopen("background.png","wb");
00495     image_write_png(background,_width,_height,fout);
00496     fclose(fout);
00497   }*/
00498   qr_gatos_mask(_img,_img,background,_width,_height,b,nb,delta,ndelta);
00499   free(background);
00500   free(mask);
00501 }
00502 
00503 #else
00504 /*The above algorithms are computationally expensive, and do not work as well
00505    as the simple algorithm below.
00506   Sauvola by itself does an excellent job of classifying regions outside the
00507    QR code as background, which greatly reduces the chance of false alarms.
00508   However, it also tends to over-shrink isolated black dots inside the code,
00509    making them easy to miss with even slight mis-alignment.
00510   Since the Gatos method uses Sauvola as input to its background interpolation
00511    method, it cannot possibly mark any pixels as foreground which Sauvola
00512    classified as background, and thus suffers from the same problem.
00513   The following simple adaptive threshold method does not have this problem,
00514    though it produces essentially random noise outside the QR code region.
00515   QR codes are structured well enough that this does not seem to lead to any
00516    actual false alarms in practice, and it allows many more codes to be
00517    detected and decoded successfully than the Sauvola or Gatos binarization
00518    methods.*/
00519 
00520 /*A simplified adaptive thresholder.
00521   This compares the current pixel value to the mean value of a (large) window
00522    surrounding it.*/
00523 unsigned char *qr_binarize(const unsigned char *_img,int _width,int _height){
00524   unsigned char *mask = NULL;
00525   if(_width>0&&_height>0){
00526     unsigned      *col_sums;
00527     int            logwindw;
00528     int            logwindh;
00529     int            windw;
00530     int            windh;
00531     int            y0offs;
00532     int            y1offs;
00533     unsigned       g;
00534     int            x;
00535     int            y;
00536     mask=(unsigned char *)malloc(_width*_height*sizeof(*mask));
00537     /*We keep the window size fairly large to ensure it doesn't fit completely
00538        inside the center of a finder pattern of a version 1 QR code at full
00539        resolution.*/
00540     for(logwindw=4;logwindw<8&&(1<<logwindw)<(_width+7>>3);logwindw++);
00541     for(logwindh=4;logwindh<8&&(1<<logwindh)<(_height+7>>3);logwindh++);
00542     windw=1<<logwindw;
00543     windh=1<<logwindh;
00544     col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
00545     /*Initialize sums down each column.*/
00546     for(x=0;x<_width;x++){
00547       g=_img[x];
00548       col_sums[x]=(g<<logwindh-1)+g;
00549     }
00550     for(y=1;y<(windh>>1);y++){
00551       y1offs=QR_MINI(y,_height-1)*_width;
00552       for(x=0;x<_width;x++){
00553         g=_img[y1offs+x];
00554         col_sums[x]+=g;
00555       }
00556     }
00557     for(y=0;y<_height;y++){
00558       unsigned m;
00559       int      x0;
00560       int      x1;
00561       /*Initialize the sum over the window.*/
00562       m=(col_sums[0]<<logwindw-1)+col_sums[0];
00563       for(x=1;x<(windw>>1);x++){
00564         x1=QR_MINI(x,_width-1);
00565         m+=col_sums[x1];
00566       }
00567       for(x=0;x<_width;x++){
00568         /*Perform the test against the threshold T = (m/n)-D,
00569            where n=windw*windh and D=3.*/
00570         g=_img[y*_width+x];
00571         mask[y*_width+x]=-(g+3<<logwindw+logwindh<m)&0xFF;
00572         /*Update the window sum.*/
00573         if(x+1<_width){
00574           x0=QR_MAXI(0,x-(windw>>1));
00575           x1=QR_MINI(x+(windw>>1),_width-1);
00576           m+=col_sums[x1]-col_sums[x0];
00577         }
00578       }
00579       /*Update the column sums.*/
00580       if(y+1<_height){
00581         y0offs=QR_MAXI(0,y-(windh>>1))*_width;
00582         y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
00583         for(x=0;x<_width;x++){
00584           col_sums[x]-=_img[y0offs+x];
00585           col_sums[x]+=_img[y1offs+x];
00586         }
00587       }
00588     }
00589     free(col_sums);
00590   }
00591 #if defined(QR_DEBUG)
00592   {
00593     FILE *fout;
00594     fout=fopen("binary.png","wb");
00595     image_write_png(_img,_width,_height,fout);
00596     fclose(fout);
00597   }
00598 #endif
00599   return(mask);
00600 }
00601 #endif
00602 
00603 #if defined(TEST_BINARIZE)
00604 #include <stdio.h>
00605 #include "image.c"
00606 
00607 int main(int _argc,char **_argv){
00608   unsigned char *img;
00609   int            width;
00610   int            height;
00611   int            x;
00612   int            y;
00613   if(_argc<2){
00614     fprintf(stderr,"usage: %s <image>.png\n",_argv[0]);
00615     return EXIT_FAILURE;
00616   }
00617   /*width=1182;
00618   height=1181;
00619   img=(unsigned char *)malloc(width*height*sizeof(*img));
00620   for(y=0;y<height;y++)for(x=0;x<width;x++){
00621     img[y*width+x]=(unsigned char)(-((x&1)^(y&1))&0xFF);
00622   }*/
00623   {
00624     FILE *fin;
00625     fin=fopen(_argv[1],"rb");
00626     image_read_png(&img,&width,&height,fin);
00627     fclose(fin);
00628   }
00629   qr_binarize(img,width,height);
00630   /*{
00631     FILE *fout;
00632     fout=fopen("binary.png","wb");
00633     image_write_png(img,width,height,fout);
00634     fclose(fout);
00635   }*/
00636   free(img);
00637   return EXIT_SUCCESS;
00638 }
00639 #endif
00640