utilities/CxImage/ximadsp.cpp

   1 // xImaDsp.cpp : DSP functions
   2 /* 07/08/2001 v1.00 - Davide Pizzolato - www.xdp.it
   3  * CxImage version 6.0.0 02/Feb/2008
   4  */
   5
   6 #include "ximage.h"
   7
   8 #include "ximaiter.h"
   9
  10 #if CXIMAGE_SUPPORT_DSP
  11
  12 ////////////////////////////////////////////////////////////////////////////////
  13 /**
  14  * Converts the image to B&W.
  15  * The OptimalThreshold() function can be used for calculating the optimal threshold.
  16  * \param level: the lightness threshold.
  17  * \return true if everything is ok
  18  */
  19 bool CxImage::Threshold(BYTE level)
  20 {
  21         if (!pDib) return false;
  22         if (head.biBitCount == 1) return true;
  23
  24         GrayScale();
  25
  26         CxImage tmp(head.biWidth,head.biHeight,1);
  27         if (!tmp.IsValid()){
  28                 strcpy(info.szLastError,tmp.GetLastError());
  29                 return false;
  30         }
  31
  32         for (long y=0;y<head.biHeight;y++){
  33                 info.nProgress = (long)(100*y/head.biHeight);
  34                 if (info.nEscape) break;
  35                 for (long x=0;x<head.biWidth;x++){
  36                         if (BlindGetPixelIndex(x,y)>level)
  37                                 tmp.BlindSetPixelIndex(x,y,1);
  38                         else
  39                                 tmp.BlindSetPixelIndex(x,y,0);
  40                 }
  41         }
  42         tmp.SetPaletteColor(0,0,0,0);
  43         tmp.SetPaletteColor(1,255,255,255);
  44         Transfer(tmp);
  45         return true;
  46 }
  47 ////////////////////////////////////////////////////////////////////////////////
  48 /**
  49  * Converts the image to B&W, using a threshold mask
  50  * \param pThresholdMask: the lightness threshold mask.
  51  * the pThresholdMask image must be grayscale with same with and height of the current image
  52  * \return true if everything is ok
  53  */
  54 bool CxImage::Threshold(CxImage* pThresholdMask)
  55 {
  56         if (!pDib) return false;
  57         if (head.biBitCount == 1) return true;
  58
  59         if (!pThresholdMask) return false;
  60
  61         if (!pThresholdMask->IsValid() ||
  62                 !pThresholdMask->IsGrayScale() ||
  63                 pThresholdMask->GetWidth() != GetWidth() ||
  64                 pThresholdMask->GetHeight() != GetHeight()){
  65                 strcpy(info.szLastError,"invalid ThresholdMask");
  66                 return false;
  67         }
  68
  69         GrayScale();
  70
  71         CxImage tmp(head.biWidth,head.biHeight,1);
  72         if (!tmp.IsValid()){
  73                 strcpy(info.szLastError,tmp.GetLastError());
  74                 return false;
  75         }
  76
  77         for (long y=0;y<head.biHeight;y++){
  78                 info.nProgress = (long)(100*y/head.biHeight);
  79                 if (info.nEscape) break;
  80                 for (long x=0;x<head.biWidth;x++){
  81                         if (BlindGetPixelIndex(x,y)>pThresholdMask->BlindGetPixelIndex(x,y))
  82                                 tmp.BlindSetPixelIndex(x,y,1);
  83                         else
  84                                 tmp.BlindSetPixelIndex(x,y,0);
  85                 }
  86         }
  87         tmp.SetPaletteColor(0,0,0,0);
  88         tmp.SetPaletteColor(1,255,255,255);
  89         Transfer(tmp);
  90         return true;
  91 }
  92 ////////////////////////////////////////////////////////////////////////////////
  93 /**
  94  * Filters only the pixels with a lightness less (or more) than the threshold level,
  95  * and preserves the colors for the unfiltered pixels.
  96  * \param level = the lightness threshold.
  97  * \param bDirection = false: filter dark pixels, true: filter light pixels
  98  * \param nBkgndColor =  filtered pixels are set to nBkgndColor color
  99  * \param bSetAlpha = if true, sets also the alpha component for the filtered pixels, with nBkgndColor.rgbReserved
 100  * \return true if everything is ok
 101  * \author [DP], [wangsongtao]
 102  */
 103 ////////////////////////////////////////////////////////////////////////////////
 104 bool CxImage::Threshold2(BYTE level, bool bDirection, RGBQUAD nBkgndColor, bool bSetAlpha)
 105 {
 106         if (!pDib) return false;
 107         if (head.biBitCount == 1) return true;
 108
 109         CxImage tmp(*this, true, false, false);
 110         if (!tmp.IsValid()){
 111                 strcpy(info.szLastError,tmp.GetLastError());
 112                 return false;
 113         }
 114
 115         tmp.GrayScale();
 116
 117         long xmin,xmax,ymin,ymax;
 118         if (pSelection){
 119                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
 120                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
 121         } else {
 122                 xmin = ymin = 0;
 123                 xmax = head.biWidth; ymax=head.biHeight;
 124         }
 125
 126         for(long y=ymin; y<ymax; y++){
 127                 info.nProgress = (long)(100*y/head.biHeight);
 128                 if (info.nEscape) break;
 129                 for(long x=xmin; x<xmax; x++){
 130 #if CXIMAGE_SUPPORT_SELECTION
 131                         if (BlindSelectionIsInside(x,y))
 132 #endif //CXIMAGE_SUPPORT_SELECTION
 133                         {
 134                                 BYTE i = tmp.BlindGetPixelIndex(x,y);
 135                                 if (!bDirection && i<level) BlindSetPixelColor(x,y,nBkgndColor,bSetAlpha);
 136                                 if (bDirection && i>=level) BlindSetPixelColor(x,y,nBkgndColor,bSetAlpha);
 137                         }
 138                 }
 139         }
 140
 141         return true;
 142 }
 143 ////////////////////////////////////////////////////////////////////////////////
 144 /**
 145  * Extract RGB channels from the image. Each channel is an 8 bit grayscale image.
 146  * \param r,g,b: pointers to CxImage objects, to store the splited channels
 147  * \return true if everything is ok
 148  */
 149 bool CxImage::SplitRGB(CxImage* r,CxImage* g,CxImage* b)
 150 {
 151         if (!pDib) return false;
 152         if (r==NULL && g==NULL && b==NULL) return false;
 153
 154         CxImage tmpr(head.biWidth,head.biHeight,8);
 155         CxImage tmpg(head.biWidth,head.biHeight,8);
 156         CxImage tmpb(head.biWidth,head.biHeight,8);
 157
 158         RGBQUAD color;
 159         for(long y=0; y<head.biHeight; y++){
 160                 for(long x=0; x<head.biWidth; x++){
 161                         color = BlindGetPixelColor(x,y);
 162                         if (r) tmpr.BlindSetPixelIndex(x,y,color.rgbRed);
 163                         if (g) tmpg.BlindSetPixelIndex(x,y,color.rgbGreen);
 164                         if (b) tmpb.BlindSetPixelIndex(x,y,color.rgbBlue);
 165                 }
 166         }
 167
 168         if (r) tmpr.SetGrayPalette();
 169         if (g) tmpg.SetGrayPalette();
 170         if (b) tmpb.SetGrayPalette();
 171
 172         /*for(long j=0; j<256; j++){
 173                 BYTE i=(BYTE)j;
 174                 if (r) tmpr.SetPaletteColor(i,i,0,0);
 175                 if (g) tmpg.SetPaletteColor(i,0,i,0);
 176                 if (b) tmpb.SetPaletteColor(i,0,0,i);
 177         }*/
 178
 179         if (r) r->Transfer(tmpr);
 180         if (g) g->Transfer(tmpg);
 181         if (b) b->Transfer(tmpb);
 182
 183         return true;
 184 }
 185 ////////////////////////////////////////////////////////////////////////////////
 186 /**
 187  * Extract CMYK channels from the image. Each channel is an 8 bit grayscale image.
 188  * \param c,m,y,k: pointers to CxImage objects, to store the splited channels
 189  * \return true if everything is ok
 190  */
 191 bool CxImage::SplitCMYK(CxImage* c,CxImage* m,CxImage* y,CxImage* k)
 192 {
 193         if (!pDib) return false;
 194         if (c==NULL && m==NULL && y==NULL && k==NULL) return false;
 195
 196         CxImage tmpc(head.biWidth,head.biHeight,8);
 197         CxImage tmpm(head.biWidth,head.biHeight,8);
 198         CxImage tmpy(head.biWidth,head.biHeight,8);
 199         CxImage tmpk(head.biWidth,head.biHeight,8);
 200
 201         RGBQUAD color;
 202         for(long yy=0; yy<head.biHeight; yy++){
 203                 for(long xx=0; xx<head.biWidth; xx++){
 204                         color = BlindGetPixelColor(xx,yy);
 205                         if (c) tmpc.BlindSetPixelIndex(xx,yy,(BYTE)(255-color.rgbRed));
 206                         if (m) tmpm.BlindSetPixelIndex(xx,yy,(BYTE)(255-color.rgbGreen));
 207                         if (y) tmpy.BlindSetPixelIndex(xx,yy,(BYTE)(255-color.rgbBlue));
 208                         if (k) tmpk.BlindSetPixelIndex(xx,yy,(BYTE)RGB2GRAY(color.rgbRed,color.rgbGreen,color.rgbBlue));
 209                 }
 210         }
 211
 212         if (c) tmpc.SetGrayPalette();
 213         if (m) tmpm.SetGrayPalette();
 214         if (y) tmpy.SetGrayPalette();
 215         if (k) tmpk.SetGrayPalette();
 216
 217         if (c) c->Transfer(tmpc);
 218         if (m) m->Transfer(tmpm);
 219         if (y) y->Transfer(tmpy);
 220         if (k) k->Transfer(tmpk);
 221
 222         return true;
 223 }
 224 ////////////////////////////////////////////////////////////////////////////////
 225 /**
 226  * Extract YUV channels from the image. Each channel is an 8 bit grayscale image.
 227  * \param y,u,v: pointers to CxImage objects, to store the splited channels
 228  * \return true if everything is ok
 229  */
 230 bool CxImage::SplitYUV(CxImage* y,CxImage* u,CxImage* v)
 231 {
 232         if (!pDib) return false;
 233         if (y==NULL && u==NULL && v==NULL) return false;
 234
 235         CxImage tmpy(head.biWidth,head.biHeight,8);
 236         CxImage tmpu(head.biWidth,head.biHeight,8);
 237         CxImage tmpv(head.biWidth,head.biHeight,8);
 238
 239         RGBQUAD color;
 240         for(long yy=0; yy<head.biHeight; yy++){
 241                 for(long x=0; x<head.biWidth; x++){
 242                         color = RGBtoYUV(BlindGetPixelColor(x,yy));
 243                         if (y) tmpy.BlindSetPixelIndex(x,yy,color.rgbRed);
 244                         if (u) tmpu.BlindSetPixelIndex(x,yy,color.rgbGreen);
 245                         if (v) tmpv.BlindSetPixelIndex(x,yy,color.rgbBlue);
 246                 }
 247         }
 248
 249         if (y) tmpy.SetGrayPalette();
 250         if (u) tmpu.SetGrayPalette();
 251         if (v) tmpv.SetGrayPalette();
 252
 253         if (y) y->Transfer(tmpy);
 254         if (u) u->Transfer(tmpu);
 255         if (v) v->Transfer(tmpv);
 256
 257         return true;
 258 }
 259 ////////////////////////////////////////////////////////////////////////////////
 260 /**
 261  * Extract YIQ channels from the image. Each channel is an 8 bit grayscale image.
 262  * \param y,i,q: pointers to CxImage objects, to store the splited channels
 263  * \return true if everything is ok
 264  */
 265 bool CxImage::SplitYIQ(CxImage* y,CxImage* i,CxImage* q)
 266 {
 267         if (!pDib) return false;
 268         if (y==NULL && i==NULL && q==NULL) return false;
 269
 270         CxImage tmpy(head.biWidth,head.biHeight,8);
 271         CxImage tmpi(head.biWidth,head.biHeight,8);
 272         CxImage tmpq(head.biWidth,head.biHeight,8);
 273
 274         RGBQUAD color;
 275         for(long yy=0; yy<head.biHeight; yy++){
 276                 for(long x=0; x<head.biWidth; x++){
 277                         color = RGBtoYIQ(BlindGetPixelColor(x,yy));
 278                         if (y) tmpy.BlindSetPixelIndex(x,yy,color.rgbRed);
 279                         if (i) tmpi.BlindSetPixelIndex(x,yy,color.rgbGreen);
 280                         if (q) tmpq.BlindSetPixelIndex(x,yy,color.rgbBlue);
 281                 }
 282         }
 283
 284         if (y) tmpy.SetGrayPalette();
 285         if (i) tmpi.SetGrayPalette();
 286         if (q) tmpq.SetGrayPalette();
 287
 288         if (y) y->Transfer(tmpy);
 289         if (i) i->Transfer(tmpi);
 290         if (q) q->Transfer(tmpq);
 291
 292         return true;
 293 }
 294 ////////////////////////////////////////////////////////////////////////////////
 295 /**
 296  * Extract XYZ channels from the image. Each channel is an 8 bit grayscale image.
 297  * \param x,y,z: pointers to CxImage objects, to store the splited channels
 298  * \return true if everything is ok
 299  */
 300 bool CxImage::SplitXYZ(CxImage* x,CxImage* y,CxImage* z)
 301 {
 302         if (!pDib) return false;
 303         if (x==NULL && y==NULL && z==NULL) return false;
 304
 305         CxImage tmpx(head.biWidth,head.biHeight,8);
 306         CxImage tmpy(head.biWidth,head.biHeight,8);
 307         CxImage tmpz(head.biWidth,head.biHeight,8);
 308
 309         RGBQUAD color;
 310         for(long yy=0; yy<head.biHeight; yy++){
 311                 for(long xx=0; xx<head.biWidth; xx++){
 312                         color = RGBtoXYZ(BlindGetPixelColor(xx,yy));
 313                         if (x) tmpx.BlindSetPixelIndex(xx,yy,color.rgbRed);
 314                         if (y) tmpy.BlindSetPixelIndex(xx,yy,color.rgbGreen);
 315                         if (z) tmpz.BlindSetPixelIndex(xx,yy,color.rgbBlue);
 316                 }
 317         }
 318
 319         if (x) tmpx.SetGrayPalette();
 320         if (y) tmpy.SetGrayPalette();
 321         if (z) tmpz.SetGrayPalette();
 322
 323         if (x) x->Transfer(tmpx);
 324         if (y) y->Transfer(tmpy);
 325         if (z) z->Transfer(tmpz);
 326
 327         return true;
 328 }
 329 ////////////////////////////////////////////////////////////////////////////////
 330 /**
 331  * Extract HSL channels from the image. Each channel is an 8 bit grayscale image.
 332  * \param h,s,l: pointers to CxImage objects, to store the splited channels
 333  * \return true if everything is ok
 334  */
 335 bool CxImage::SplitHSL(CxImage* h,CxImage* s,CxImage* l)
 336 {
 337         if (!pDib) return false;
 338         if (h==NULL && s==NULL && l==NULL) return false;
 339
 340         CxImage tmph(head.biWidth,head.biHeight,8);
 341         CxImage tmps(head.biWidth,head.biHeight,8);
 342         CxImage tmpl(head.biWidth,head.biHeight,8);
 343
 344         RGBQUAD color;
 345         for(long y=0; y<head.biHeight; y++){
 346                 for(long x=0; x<head.biWidth; x++){
 347                         color = RGBtoHSL(BlindGetPixelColor(x,y));
 348                         if (h) tmph.BlindSetPixelIndex(x,y,color.rgbRed);
 349                         if (s) tmps.BlindSetPixelIndex(x,y,color.rgbGreen);
 350                         if (l) tmpl.BlindSetPixelIndex(x,y,color.rgbBlue);
 351                 }
 352         }
 353
 354         if (h) tmph.SetGrayPalette();
 355         if (s) tmps.SetGrayPalette();
 356         if (l) tmpl.SetGrayPalette();
 357
 358         /* pseudo-color generator for hue channel (visual debug)
 359         if (h) for(long j=0; j<256; j++){
 360                 BYTE i=(BYTE)j;
 361                 RGBQUAD hsl={120,240,i,0};
 362                 tmph.SetPaletteColor(i,HSLtoRGB(hsl));
 363         }*/
 364
 365         if (h) h->Transfer(tmph);
 366         if (s) s->Transfer(tmps);
 367         if (l) l->Transfer(tmpl);
 368
 369         return true;
 370 }
 371 ////////////////////////////////////////////////////////////////////////////////
 372 #define  HSLMAX   255   /* H,L, and S vary over 0-HSLMAX */
 373 #define  RGBMAX   255   /* R,G, and B vary over 0-RGBMAX */
 374                         /* HSLMAX BEST IF DIVISIBLE BY 6 */
 375                         /* RGBMAX, HSLMAX must each fit in a BYTE. */
 376 /* Hue is undefined if Saturation is 0 (grey-scale) */
 377 /* This value determines where the Hue scrollbar is */
 378 /* initially set for achromatic colors */
 379 #define HSLUNDEFINED (HSLMAX*2/3)
 380 ////////////////////////////////////////////////////////////////////////////////
 381 RGBQUAD CxImage::RGBtoHSL(RGBQUAD lRGBColor)
 382 {
 383         BYTE R,G,B;                                     /* input RGB values */
 384         BYTE H,L,S;                                     /* output HSL values */
 385         BYTE cMax,cMin;                         /* max and min RGB values */
 386         WORD Rdelta,Gdelta,Bdelta;      /* intermediate value: % of spread from max*/
 387
 388         R = lRGBColor.rgbRed;   /* get R, G, and B out of DWORD */
 389         G = lRGBColor.rgbGreen;
 390         B = lRGBColor.rgbBlue;
 391
 392         cMax = max( max(R,G), B);       /* calculate lightness */
 393         cMin = min( min(R,G), B);
 394         L = (BYTE)((((cMax+cMin)*HSLMAX)+RGBMAX)/(2*RGBMAX));
 395
 396         if (cMax==cMin){                        /* r=g=b --> achromatic case */
 397                 S = 0;                                  /* saturation */
 398                 H = HSLUNDEFINED;               /* hue */
 399         } else {                                        /* chromatic case */
 400                 if (L <= (HSLMAX/2))    /* saturation */
 401                         S = (BYTE)((((cMax-cMin)*HSLMAX)+((cMax+cMin)/2))/(cMax+cMin));
 402                 else
 403                         S = (BYTE)((((cMax-cMin)*HSLMAX)+((2*RGBMAX-cMax-cMin)/2))/(2*RGBMAX-cMax-cMin));
 404                 /* hue */
 405                 Rdelta = (WORD)((((cMax-R)*(HSLMAX/6)) + ((cMax-cMin)/2) ) / (cMax-cMin));
 406                 Gdelta = (WORD)((((cMax-G)*(HSLMAX/6)) + ((cMax-cMin)/2) ) / (cMax-cMin));
 407                 Bdelta = (WORD)((((cMax-B)*(HSLMAX/6)) + ((cMax-cMin)/2) ) / (cMax-cMin));
 408
 409                 if (R == cMax)
 410                         H = (BYTE)(Bdelta - Gdelta);
 411                 else if (G == cMax)
 412                         H = (BYTE)((HSLMAX/3) + Rdelta - Bdelta);
 413                 else /* B == cMax */
 414                         H = (BYTE)(((2*HSLMAX)/3) + Gdelta - Rdelta);
 415
 416 //              if (H < 0) H += HSLMAX;     //always false
 417                 if (H > HSLMAX) H -= HSLMAX;
 418         }
 419         RGBQUAD hsl={L,S,H,0};
 420         return hsl;
 421 }
 422 ////////////////////////////////////////////////////////////////////////////////
 423 float CxImage::HueToRGB(float n1,float n2, float hue)
 424 {
 425         //<F. Livraghi> fixed implementation for HSL2RGB routine
 426         float rValue;
 427
 428         if (hue > 360)
 429                 hue = hue - 360;
 430         else if (hue < 0)
 431                 hue = hue + 360;
 432
 433         if (hue < 60)
 434                 rValue = n1 + (n2-n1)*hue/60.0f;
 435         else if (hue < 180)
 436                 rValue = n2;
 437         else if (hue < 240)
 438                 rValue = n1+(n2-n1)*(240-hue)/60;
 439         else
 440                 rValue = n1;
 441
 442         return rValue;
 443 }
 444 ////////////////////////////////////////////////////////////////////////////////
 445 RGBQUAD CxImage::HSLtoRGB(COLORREF cHSLColor)
 446 {
 447         return HSLtoRGB(RGBtoRGBQUAD(cHSLColor));
 448 }
 449 ////////////////////////////////////////////////////////////////////////////////
 450 RGBQUAD CxImage::HSLtoRGB(RGBQUAD lHSLColor)
 451 {
 452         //<F. Livraghi> fixed implementation for HSL2RGB routine
 453         float h,s,l;
 454         float m1,m2;
 455         BYTE r,g,b;
 456
 457         h = (float)lHSLColor.rgbRed * 360.0f/255.0f;
 458         s = (float)lHSLColor.rgbGreen/255.0f;
 459         l = (float)lHSLColor.rgbBlue/255.0f;
 460
 461         if (l <= 0.5)   m2 = l * (1+s);
 462         else                    m2 = l + s - l*s;
 463
 464         m1 = 2 * l - m2;
 465
 466         if (s == 0) {
 467                 r=g=b=(BYTE)(l*255.0f);
 468         } else {
 469                 r = (BYTE)(HueToRGB(m1,m2,h+120) * 255.0f);
 470                 g = (BYTE)(HueToRGB(m1,m2,h) * 255.0f);
 471                 b = (BYTE)(HueToRGB(m1,m2,h-120) * 255.0f);
 472         }
 473
 474         RGBQUAD rgb = {b,g,r,0};
 475         return rgb;
 476 }
 477 ////////////////////////////////////////////////////////////////////////////////
 478 RGBQUAD CxImage::YUVtoRGB(RGBQUAD lYUVColor)
 479 {
 480         int U,V,R,G,B;
 481         float Y = lYUVColor.rgbRed;
 482         U = lYUVColor.rgbGreen - 128;
 483         V = lYUVColor.rgbBlue - 128;
 484
 485 //      R = (int)(1.164 * Y + 2.018 * U);
 486 //      G = (int)(1.164 * Y - 0.813 * V - 0.391 * U);
 487 //      B = (int)(1.164 * Y + 1.596 * V);
 488         R = (int)( Y + 1.403f * V);
 489         G = (int)( Y - 0.344f * U - 0.714f * V);
 490         B = (int)( Y + 1.770f * U);
 491
 492         R= min(255,max(0,R));
 493         G= min(255,max(0,G));
 494         B= min(255,max(0,B));
 495         RGBQUAD rgb={(BYTE)B,(BYTE)G,(BYTE)R,0};
 496         return rgb;
 497 }
 498 ////////////////////////////////////////////////////////////////////////////////
 499 RGBQUAD CxImage::RGBtoYUV(RGBQUAD lRGBColor)
 500 {
 501         int Y,U,V,R,G,B;
 502         R = lRGBColor.rgbRed;
 503         G = lRGBColor.rgbGreen;
 504         B = lRGBColor.rgbBlue;
 505
 506 //      Y = (int)( 0.257 * R + 0.504 * G + 0.098 * B);
 507 //      U = (int)( 0.439 * R - 0.368 * G - 0.071 * B + 128);
 508 //      V = (int)(-0.148 * R - 0.291 * G + 0.439 * B + 128);
 509         Y = (int)(0.299f * R + 0.587f * G + 0.114f * B);
 510         U = (int)((B-Y) * 0.565f + 128);
 511         V = (int)((R-Y) * 0.713f + 128);
 512
 513         Y= min(255,max(0,Y));
 514         U= min(255,max(0,U));
 515         V= min(255,max(0,V));
 516         RGBQUAD yuv={(BYTE)V,(BYTE)U,(BYTE)Y,0};
 517         return yuv;
 518 }
 519 ////////////////////////////////////////////////////////////////////////////////
 520 RGBQUAD CxImage::YIQtoRGB(RGBQUAD lYIQColor)
 521 {
 522         int I,Q,R,G,B;
 523         float Y = lYIQColor.rgbRed;
 524         I = lYIQColor.rgbGreen - 128;
 525         Q = lYIQColor.rgbBlue - 128;
 526
 527         R = (int)( Y + 0.956f * I + 0.621f * Q);
 528         G = (int)( Y - 0.273f * I - 0.647f * Q);
 529         B = (int)( Y - 1.104f * I + 1.701f * Q);
 530
 531         R= min(255,max(0,R));
 532         G= min(255,max(0,G));
 533         B= min(255,max(0,B));
 534         RGBQUAD rgb={(BYTE)B,(BYTE)G,(BYTE)R,0};
 535         return rgb;
 536 }
 537 ////////////////////////////////////////////////////////////////////////////////
 538 RGBQUAD CxImage::RGBtoYIQ(RGBQUAD lRGBColor)
 539 {
 540         int Y,I,Q,R,G,B;
 541         R = lRGBColor.rgbRed;
 542         G = lRGBColor.rgbGreen;
 543         B = lRGBColor.rgbBlue;
 544
 545         Y = (int)( 0.2992f * R + 0.5868f * G + 0.1140f * B);
 546         I = (int)( 0.5960f * R - 0.2742f * G - 0.3219f * B + 128);
 547         Q = (int)( 0.2109f * R - 0.5229f * G + 0.3120f * B + 128);
 548
 549         Y= min(255,max(0,Y));
 550         I= min(255,max(0,I));
 551         Q= min(255,max(0,Q));
 552         RGBQUAD yiq={(BYTE)Q,(BYTE)I,(BYTE)Y,0};
 553         return yiq;
 554 }
 555 ////////////////////////////////////////////////////////////////////////////////
 556 RGBQUAD CxImage::XYZtoRGB(RGBQUAD lXYZColor)
 557 {
 558         int X,Y,Z,R,G,B;
 559         X = lXYZColor.rgbRed;
 560         Y = lXYZColor.rgbGreen;
 561         Z = lXYZColor.rgbBlue;
 562         double k=1.088751;
 563
 564         R = (int)(  3.240479f * X - 1.537150f * Y - 0.498535f * Z * k);
 565         G = (int)( -0.969256f * X + 1.875992f * Y + 0.041556f * Z * k);
 566         B = (int)(  0.055648f * X - 0.204043f * Y + 1.057311f * Z * k);
 567
 568         R= min(255,max(0,R));
 569         G= min(255,max(0,G));
 570         B= min(255,max(0,B));
 571         RGBQUAD rgb={(BYTE)B,(BYTE)G,(BYTE)R,0};
 572         return rgb;
 573 }
 574 ////////////////////////////////////////////////////////////////////////////////
 575 RGBQUAD CxImage::RGBtoXYZ(RGBQUAD lRGBColor)
 576 {
 577         int X,Y,Z,R,G,B;
 578         R = lRGBColor.rgbRed;
 579         G = lRGBColor.rgbGreen;
 580         B = lRGBColor.rgbBlue;
 581
 582         X = (int)( 0.412453f * R + 0.357580f * G + 0.180423f * B);
 583         Y = (int)( 0.212671f * R + 0.715160f * G + 0.072169f * B);
 584         Z = (int)((0.019334f * R + 0.119193f * G + 0.950227f * B)*0.918483657f);
 585
 586         //X= min(255,max(0,X));
 587         //Y= min(255,max(0,Y));
 588         //Z= min(255,max(0,Z));
 589         RGBQUAD xyz={(BYTE)Z,(BYTE)Y,(BYTE)X,0};
 590         return xyz;
 591 }
 592 ////////////////////////////////////////////////////////////////////////////////
 593 /**
 594  * Generates a "rainbow" palette with saturated colors
 595  * \param correction: 1 generates a single hue spectrum. 0.75 is nice for scientific applications.
 596  */
 597 void CxImage::HuePalette(float correction)
 598 {
 599         if (head.biClrUsed==0) return;
 600
 601         for(DWORD j=0; j<head.biClrUsed; j++){
 602                 BYTE i=(BYTE)(j*correction*(255/(head.biClrUsed-1)));
 603                 RGBQUAD hsl={120,240,i,0};
 604                 SetPaletteColor((BYTE)j,HSLtoRGB(hsl));
 605         }
 606 }
 607 ////////////////////////////////////////////////////////////////////////////////
 608 /**
 609  * Replaces the original hue and saturation values.
 610  * \param hue: hue
 611  * \param sat: saturation
 612  * \param blend: can be from 0 (no effect) to 1 (full effect)
 613  * \return true if everything is ok
 614  */
 615 bool CxImage::Colorize(BYTE hue, BYTE sat, float blend)
 616 {
 617         if (!pDib) return false;
 618
 619         if (blend < 0.0f) blend = 0.0f;
 620         if (blend > 1.0f) blend = 1.0f;
 621         int a0 = (int)(256*blend);
 622         int a1 = 256 - a0;
 623
 624         bool bFullBlend = false;
 625         if (blend > 0.999f)     bFullBlend = true;
 626
 627         RGBQUAD color,hsl;
 628         if (head.biClrUsed==0){
 629
 630                 long xmin,xmax,ymin,ymax;
 631                 if (pSelection){
 632                         xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
 633                         ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
 634                 } else {
 635                         xmin = ymin = 0;
 636                         xmax = head.biWidth; ymax=head.biHeight;
 637                 }
 638
 639                 for(long y=ymin; y<ymax; y++){
 640                         info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
 641                         if (info.nEscape) break;
 642                         for(long x=xmin; x<xmax; x++){
 643 #if CXIMAGE_SUPPORT_SELECTION
 644                                 if (BlindSelectionIsInside(x,y))
 645 #endif //CXIMAGE_SUPPORT_SELECTION
 646                                 {
 647                                         if (bFullBlend){
 648                                                 color = RGBtoHSL(BlindGetPixelColor(x,y));
 649                                                 color.rgbRed=hue;
 650                                                 color.rgbGreen=sat;
 651                                                 BlindSetPixelColor(x,y,HSLtoRGB(color));
 652                                         } else {
 653                                                 color = BlindGetPixelColor(x,y);
 654                                                 hsl.rgbRed=hue;
 655                                                 hsl.rgbGreen=sat;
 656                                                 hsl.rgbBlue = (BYTE)RGB2GRAY(color.rgbRed,color.rgbGreen,color.rgbBlue);
 657                                                 hsl = HSLtoRGB(hsl);
 658                                                 //BlendPixelColor(x,y,hsl,blend);
 659                                                 //color.rgbRed = (BYTE)(hsl.rgbRed * blend + color.rgbRed * (1.0f - blend));
 660                                                 //color.rgbBlue = (BYTE)(hsl.rgbBlue * blend + color.rgbBlue * (1.0f - blend));
 661                                                 //color.rgbGreen = (BYTE)(hsl.rgbGreen * blend + color.rgbGreen * (1.0f - blend));
 662                                                 color.rgbRed = (BYTE)((hsl.rgbRed * a0 + color.rgbRed * a1)>>8);
 663                                                 color.rgbBlue = (BYTE)((hsl.rgbBlue * a0 + color.rgbBlue * a1)>>8);
 664                                                 color.rgbGreen = (BYTE)((hsl.rgbGreen * a0 + color.rgbGreen * a1)>>8);
 665                                                 BlindSetPixelColor(x,y,color);
 666                                         }
 667                                 }
 668                         }
 669                 }
 670         } else {
 671                 for(DWORD j=0; j<head.biClrUsed; j++){
 672                         if (bFullBlend){
 673                                 color = RGBtoHSL(GetPaletteColor((BYTE)j));
 674                                 color.rgbRed=hue;
 675                                 color.rgbGreen=sat;
 676                                 SetPaletteColor((BYTE)j,HSLtoRGB(color));
 677                         } else {
 678                                 color = GetPaletteColor((BYTE)j);
 679                                 hsl.rgbRed=hue;
 680                                 hsl.rgbGreen=sat;
 681                                 hsl.rgbBlue = (BYTE)RGB2GRAY(color.rgbRed,color.rgbGreen,color.rgbBlue);
 682                                 hsl = HSLtoRGB(hsl);
 683                                 color.rgbRed = (BYTE)(hsl.rgbRed * blend + color.rgbRed * (1.0f - blend));
 684                                 color.rgbBlue = (BYTE)(hsl.rgbBlue * blend + color.rgbBlue * (1.0f - blend));
 685                                 color.rgbGreen = (BYTE)(hsl.rgbGreen * blend + color.rgbGreen * (1.0f - blend));
 686                                 SetPaletteColor((BYTE)j,color);
 687                         }
 688                 }
 689         }
 690
 691         return true;
 692 }
 693 ////////////////////////////////////////////////////////////////////////////////
 694 /**
 695  * Changes the brightness and the contrast of the image.
 696  * \param brightness: can be from -255 to 255, if brightness is negative, the image becomes dark.
 697  * \param contrast: can be from -100 to 100, the neutral value is 0.
 698  * \return true if everything is ok
 699  */
 700 bool CxImage::Light(long brightness, long contrast)
 701 {
 702         if (!pDib) return false;
 703         float c=(100 + contrast)/100.0f;
 704         brightness+=128;
 705
 706         BYTE cTable[256]; //<nipper>
 707         for (int i=0;i<256;i++) {
 708                 cTable[i] = (BYTE)max(0,min(255,(int)((i-128)*c + brightness + 0.5f)));
 709         }
 710
 711         return Lut(cTable);
 712 }
 713 ////////////////////////////////////////////////////////////////////////////////
 714 /**
 715  * \return mean lightness of the image. Useful with Threshold() and Light()
 716  */
 717 float CxImage::Mean()
 718 {
 719         if (!pDib) return 0;
 720
 721         CxImage tmp(*this,true);
 722         if (!tmp.IsValid()){
 723                 strcpy(info.szLastError,tmp.GetLastError());
 724                 return false;
 725         }
 726
 727         tmp.GrayScale();
 728         float sum=0;
 729
 730         long xmin,xmax,ymin,ymax;
 731         if (pSelection){
 732                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
 733                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
 734         } else {
 735                 xmin = ymin = 0;
 736                 xmax = head.biWidth; ymax=head.biHeight;
 737         }
 738         if (xmin==xmax || ymin==ymax) return (float)0.0;
 739
 740         BYTE *iSrc=tmp.info.pImage;
 741         iSrc += tmp.info.dwEffWidth*ymin; // necessary for selections <Admir Hodzic>
 742
 743         for(long y=ymin; y<ymax; y++){
 744                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin)); //<zhanghk><Anatoly Ivasyuk>
 745                 for(long x=xmin; x<xmax; x++){
 746                         sum+=iSrc[x];
 747                 }
 748                 iSrc+=tmp.info.dwEffWidth;
 749         }
 750         return sum/(xmax-xmin)/(ymax-ymin);
 751 }
 752 ////////////////////////////////////////////////////////////////////////////////
 753 /**
 754  * 2D linear filter
 755  * \param kernel: convolving matrix, in row format.
 756  * \param Ksize: size of the kernel.
 757  * \param Kfactor: normalization constant.
 758  * \param Koffset: bias.
 759  * \verbatim Example: the "soften" filter uses this kernel:
 760         1 1 1
 761         1 8 1
 762         1 1 1
 763  the function needs: kernel={1,1,1,1,8,1,1,1,1}; Ksize=3; Kfactor=16; Koffset=0; \endverbatim
 764  * \return true if everything is ok
 765  */
 766 bool CxImage::Filter(long* kernel, long Ksize, long Kfactor, long Koffset)
 767 {
 768         if (!pDib) return false;
 769
 770         long k2 = Ksize/2;
 771         long kmax= Ksize-k2;
 772         long r,g,b,i;
 773         long ksumcur,ksumtot;
 774         RGBQUAD c;
 775
 776         CxImage tmp(*this);
 777         if (!tmp.IsValid()){
 778                 strcpy(info.szLastError,tmp.GetLastError());
 779                 return false;
 780         }
 781
 782         long xmin,xmax,ymin,ymax;
 783         if (pSelection){
 784                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
 785                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
 786         } else {
 787                 xmin = ymin = 0;
 788                 xmax = head.biWidth; ymax=head.biHeight;
 789         }
 790
 791         ksumtot = 0;
 792         for(long j=-k2;j<kmax;j++){
 793                 for(long k=-k2;k<kmax;k++){
 794                         ksumtot += kernel[(j+k2)+Ksize*(k+k2)];
 795                 }
 796         }
 797
 798         if ((head.biBitCount==8) && IsGrayScale())
 799         {
 800                 unsigned char* cPtr;
 801                 unsigned char* cPtr2;
 802                 int iCount;
 803                 int iY, iY2, iY1;
 804                 cPtr = info.pImage;
 805                 cPtr2 = (unsigned char *)tmp.info.pImage;
 806                 for(long y=ymin; y<ymax; y++){
 807                         info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
 808                         if (info.nEscape) break;
 809                         iY1 = y*info.dwEffWidth+xmin;
 810                         for(long x=xmin; x<xmax; x++, iY1++){
 811 #if CXIMAGE_SUPPORT_SELECTION
 812                                 if (BlindSelectionIsInside(x,y))
 813 #endif //CXIMAGE_SUPPORT_SELECTION
 814                                 {
 815                                         b=ksumcur=0;
 816                                         iCount = 0;
 817                                         iY2 = ((y-k2)*info.dwEffWidth);
 818                                         for(long j=-k2;j<kmax;j++, iY2+=info.dwEffWidth)
 819                                         {
 820                                                 if (0>(y+j) || (y+j)>=head.biHeight) continue;
 821                                                 iY = iY2+x;
 822                                                 for(long k=-k2;k<kmax;k++, iCount++)
 823                                                 {
 824                                                         if (0>(x+k) || (x+k)>=head.biWidth) continue;
 825                                                         i=kernel[iCount];
 826                                                         b += cPtr[iY+k] * i;
 827                                                         ksumcur += i;
 828                                                 }
 829                                         }
 830                                         if (Kfactor==0 || ksumcur==0){
 831                                                 cPtr2[iY1] = (BYTE)min(255, max(0,(int)(b + Koffset)));
 832                                         } else if (ksumtot == ksumcur) {
 833                                                 cPtr2[iY1] = (BYTE)min(255, max(0,(int)(b/Kfactor + Koffset)));
 834                                         } else {
 835                                                 cPtr2[iY1] = (BYTE)min(255, max(0,(int)((b*ksumtot)/(ksumcur*Kfactor) + Koffset)));
 836                                         }
 837                                 }
 838                         }
 839                 }
 840         }
 841         else
 842         {
 843                 for(long y=ymin; y<ymax; y++){
 844                         info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
 845                         if (info.nEscape) break;
 846                         for(long x=xmin; x<xmax; x++){
 847         #if CXIMAGE_SUPPORT_SELECTION
 848                                 if (BlindSelectionIsInside(x,y))
 849         #endif //CXIMAGE_SUPPORT_SELECTION
 850                                         {
 851                                         r=b=g=ksumcur=0;
 852                                         for(long j=-k2;j<kmax;j++){
 853                                                 for(long k=-k2;k<kmax;k++){
 854                                                         if (!IsInside(x+j,y+k)) continue;
 855                                                         c = BlindGetPixelColor(x+j,y+k);
 856                                                         i = kernel[(j+k2)+Ksize*(k+k2)];
 857                                                         r += c.rgbRed * i;
 858                                                         g += c.rgbGreen * i;
 859                                                         b += c.rgbBlue * i;
 860                                                         ksumcur += i;
 861                                                 }
 862                                         }
 863                                         if (Kfactor==0 || ksumcur==0){
 864                                                 c.rgbRed   = (BYTE)min(255, max(0,(int)(r + Koffset)));
 865                                                 c.rgbGreen = (BYTE)min(255, max(0,(int)(g + Koffset)));
 866                                                 c.rgbBlue  = (BYTE)min(255, max(0,(int)(b + Koffset)));
 867                                         } else if (ksumtot == ksumcur) {
 868                                                 c.rgbRed   = (BYTE)min(255, max(0,(int)(r/Kfactor + Koffset)));
 869                                                 c.rgbGreen = (BYTE)min(255, max(0,(int)(g/Kfactor + Koffset)));
 870                                                 c.rgbBlue  = (BYTE)min(255, max(0,(int)(b/Kfactor + Koffset)));
 871                                         } else {
 872                                                 c.rgbRed   = (BYTE)min(255, max(0,(int)((r*ksumtot)/(ksumcur*Kfactor) + Koffset)));
 873                                                 c.rgbGreen = (BYTE)min(255, max(0,(int)((g*ksumtot)/(ksumcur*Kfactor) + Koffset)));
 874                                                 c.rgbBlue  = (BYTE)min(255, max(0,(int)((b*ksumtot)/(ksumcur*Kfactor) + Koffset)));
 875                                         }
 876                                         tmp.BlindSetPixelColor(x,y,c);
 877                                 }
 878                         }
 879                 }
 880         }
 881         Transfer(tmp);
 882         return true;
 883 }
 884 ////////////////////////////////////////////////////////////////////////////////
 885 /**
 886  * Enhance the dark areas of the image
 887  * \param Ksize: size of the kernel.
 888  * \return true if everything is ok
 889  */
 890 bool CxImage::Erode(long Ksize)
 891 {
 892         if (!pDib) return false;
 893
 894         long k2 = Ksize/2;
 895         long kmax= Ksize-k2;
 896         BYTE r,g,b;
 897         RGBQUAD c;
 898
 899         CxImage tmp(*this);
 900         if (!tmp.IsValid()){
 901                 strcpy(info.szLastError,tmp.GetLastError());
 902                 return false;
 903         }
 904
 905         long xmin,xmax,ymin,ymax;
 906         if (pSelection){
 907                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
 908                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
 909         } else {
 910                 xmin = ymin = 0;
 911                 xmax = head.biWidth; ymax=head.biHeight;
 912         }
 913
 914         for(long y=ymin; y<ymax; y++){
 915                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
 916                 if (info.nEscape) break;
 917                 for(long x=xmin; x<xmax; x++){
 918 #if CXIMAGE_SUPPORT_SELECTION
 919                         if (BlindSelectionIsInside(x,y))
 920 #endif //CXIMAGE_SUPPORT_SELECTION
 921                         {
 922                                 r=b=g=255;
 923                                 for(long j=-k2;j<kmax;j++){
 924                                         for(long k=-k2;k<kmax;k++){
 925                                                 if (!IsInside(x+j,y+k)) continue;
 926                                                 c = BlindGetPixelColor(x+j,y+k);
 927                                                 if (c.rgbRed < r) r=c.rgbRed;
 928                                                 if (c.rgbGreen < g) g=c.rgbGreen;
 929                                                 if (c.rgbBlue < b) b=c.rgbBlue;
 930                                         }
 931                                 }
 932                                 c.rgbRed   = r;
 933                                 c.rgbGreen = g;
 934                                 c.rgbBlue  = b;
 935                                 tmp.BlindSetPixelColor(x,y,c);
 936                         }
 937                 }
 938         }
 939         Transfer(tmp);
 940         return true;
 941 }
 942 ////////////////////////////////////////////////////////////////////////////////
 943 /**
 944  * Enhance the light areas of the image
 945  * \param Ksize: size of the kernel.
 946  * \return true if everything is ok
 947  */
 948 bool CxImage::Dilate(long Ksize)
 949 {
 950         if (!pDib) return false;
 951
 952         long k2 = Ksize/2;
 953         long kmax= Ksize-k2;
 954         BYTE r,g,b;
 955         RGBQUAD c;
 956
 957         CxImage tmp(*this);
 958         if (!tmp.IsValid()){
 959                 strcpy(info.szLastError,tmp.GetLastError());
 960                 return false;
 961         }
 962
 963         long xmin,xmax,ymin,ymax;
 964         if (pSelection){
 965                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
 966                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
 967         } else {
 968                 xmin = ymin = 0;
 969                 xmax = head.biWidth; ymax=head.biHeight;
 970         }
 971
 972         for(long y=ymin; y<ymax; y++){
 973                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
 974                 if (info.nEscape) break;
 975                 for(long x=xmin; x<xmax; x++){
 976 #if CXIMAGE_SUPPORT_SELECTION
 977                         if (BlindSelectionIsInside(x,y))
 978 #endif //CXIMAGE_SUPPORT_SELECTION
 979                         {
 980                                 r=b=g=0;
 981                                 for(long j=-k2;j<kmax;j++){
 982                                         for(long k=-k2;k<kmax;k++){
 983                                                 if (!IsInside(x+j,y+k)) continue;
 984                                                 c = BlindGetPixelColor(x+j,y+k);
 985                                                 if (c.rgbRed > r) r=c.rgbRed;
 986                                                 if (c.rgbGreen > g) g=c.rgbGreen;
 987                                                 if (c.rgbBlue > b) b=c.rgbBlue;
 988                                         }
 989                                 }
 990                                 c.rgbRed   = r;
 991                                 c.rgbGreen = g;
 992                                 c.rgbBlue  = b;
 993                                 tmp.BlindSetPixelColor(x,y,c);
 994                         }
 995                 }
 996         }
 997         Transfer(tmp);
 998         return true;
 999 }
1000 ////////////////////////////////////////////////////////////////////////////////
1001 /**
1002  * Enhance the variations between adjacent pixels.
1003  * Similar results can be achieved using Filter(),
1004  * but the algorithms are different both in Edge() and in Contour().
1005  * \param Ksize: size of the kernel.
1006  * \return true if everything is ok
1007  */
1008 bool CxImage::Edge(long Ksize)
1009 {
1010         if (!pDib) return false;
1011
1012         long k2 = Ksize/2;
1013         long kmax= Ksize-k2;
1014         BYTE r,g,b,rr,gg,bb;
1015         RGBQUAD c;
1016
1017         CxImage tmp(*this);
1018         if (!tmp.IsValid()){
1019                 strcpy(info.szLastError,tmp.GetLastError());
1020                 return false;
1021         }
1022
1023         long xmin,xmax,ymin,ymax;
1024         if (pSelection){
1025                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
1026                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
1027         } else {
1028                 xmin = ymin = 0;
1029                 xmax = head.biWidth; ymax=head.biHeight;
1030         }
1031
1032         for(long y=ymin; y<ymax; y++){
1033                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
1034                 if (info.nEscape) break;
1035                 for(long x=xmin; x<xmax; x++){
1036 #if CXIMAGE_SUPPORT_SELECTION
1037                         if (BlindSelectionIsInside(x,y))
1038 #endif //CXIMAGE_SUPPORT_SELECTION
1039                         {
1040                                 r=b=g=0;
1041                                 rr=bb=gg=255;
1042                                 for(long j=-k2;j<kmax;j++){
1043                                         for(long k=-k2;k<kmax;k++){
1044                                                 if (!IsInside(x+j,y+k)) continue;
1045                                                 c = BlindGetPixelColor(x+j,y+k);
1046                                                 if (c.rgbRed > r) r=c.rgbRed;
1047                                                 if (c.rgbGreen > g) g=c.rgbGreen;
1048                                                 if (c.rgbBlue > b) b=c.rgbBlue;
1049
1050                                                 if (c.rgbRed < rr) rr=c.rgbRed;
1051                                                 if (c.rgbGreen < gg) gg=c.rgbGreen;
1052                                                 if (c.rgbBlue < bb) bb=c.rgbBlue;
1053                                         }
1054                                 }
1055                                 c.rgbRed   = (BYTE)(255-abs(r-rr));
1056                                 c.rgbGreen = (BYTE)(255-abs(g-gg));
1057                                 c.rgbBlue  = (BYTE)(255-abs(b-bb));
1058                                 tmp.BlindSetPixelColor(x,y,c);
1059                         }
1060                 }
1061         }
1062         Transfer(tmp);
1063         return true;
1064 }
1065 ////////////////////////////////////////////////////////////////////////////////
1066 /**
1067  * Blends two images
1068  * \param imgsrc2: image to be mixed with this
1069  * \param op: blending method; see ImageOpType
1070  * \param lXOffset, lYOffset: image displacement
1071  * \param bMixAlpha: if true and imgsrc2 has a valid alpha layer, it will be mixed in the destination image.
1072  * \return true if everything is ok
1073  *
1074  * thanks to Mwolski
1075  */
1076 //
1077 void CxImage::Mix(CxImage & imgsrc2, ImageOpType op, long lXOffset, long lYOffset, bool bMixAlpha)
1078 {
1079     long lWide = min(GetWidth(),imgsrc2.GetWidth()-lXOffset);
1080     long lHeight = min(GetHeight(),imgsrc2.GetHeight()-lYOffset);
1081
1082         bool bEditAlpha = imgsrc2.AlphaIsValid() & bMixAlpha;
1083
1084         if (bEditAlpha && AlphaIsValid()==false){
1085                 AlphaCreate();
1086         }
1087
1088     RGBQUAD rgbBackgrnd1 = GetTransColor();
1089     RGBQUAD rgb1, rgb2, rgbDest;
1090
1091     for(long lY=0;lY<lHeight;lY++)
1092     {
1093                 info.nProgress = (long)(100*lY/head.biHeight);
1094                 if (info.nEscape) break;
1095
1096         for(long lX=0;lX<lWide;lX++)
1097         {
1098 #if CXIMAGE_SUPPORT_SELECTION
1099                         if (SelectionIsInside(lX,lY) && imgsrc2.SelectionIsInside(lX+lXOffset,lY+lYOffset))
1100 #endif //CXIMAGE_SUPPORT_SELECTION
1101                         {
1102                                 rgb1 = GetPixelColor(lX,lY);
1103                                 rgb2 = imgsrc2.GetPixelColor(lX+lXOffset,lY+lYOffset);
1104                                 switch(op)
1105                                 {
1106                                         case OpAvg:
1107                                                 rgbDest.rgbBlue =  (BYTE)((rgb1.rgbBlue+rgb2.rgbBlue)/2);
1108                                                 rgbDest.rgbGreen = (BYTE)((rgb1.rgbGreen+rgb2.rgbGreen)/2);
1109                                                 rgbDest.rgbRed =   (BYTE)((rgb1.rgbRed+rgb2.rgbRed)/2);
1110                                                 if (bEditAlpha) rgbDest.rgbReserved = (BYTE)((rgb1.rgbReserved+rgb2.rgbReserved)/2);
1111                                         break;
1112                                         case OpAdd:
1113                                                 rgbDest.rgbBlue = (BYTE)max(0,min(255,rgb1.rgbBlue+rgb2.rgbBlue));
1114                                                 rgbDest.rgbGreen = (BYTE)max(0,min(255,rgb1.rgbGreen+rgb2.rgbGreen));
1115                                                 rgbDest.rgbRed = (BYTE)max(0,min(255,rgb1.rgbRed+rgb2.rgbRed));
1116                                                 if (bEditAlpha) rgbDest.rgbReserved = (BYTE)max(0,min(255,rgb1.rgbReserved+rgb2.rgbReserved));
1117                                         break;
1118                                         case OpSub:
1119                                                 rgbDest.rgbBlue = (BYTE)max(0,min(255,rgb1.rgbBlue-rgb2.rgbBlue));
1120                                                 rgbDest.rgbGreen = (BYTE)max(0,min(255,rgb1.rgbGreen-rgb2.rgbGreen));
1121                                                 rgbDest.rgbRed = (BYTE)max(0,min(255,rgb1.rgbRed-rgb2.rgbRed));
1122                                                 if (bEditAlpha) rgbDest.rgbReserved = (BYTE)max(0,min(255,rgb1.rgbReserved-rgb2.rgbReserved));
1123                                         break;
1124                                         case OpAnd:
1125                                                 rgbDest.rgbBlue = (BYTE)(rgb1.rgbBlue&rgb2.rgbBlue);
1126                                                 rgbDest.rgbGreen = (BYTE)(rgb1.rgbGreen&rgb2.rgbGreen);
1127                                                 rgbDest.rgbRed = (BYTE)(rgb1.rgbRed&rgb2.rgbRed);
1128                                                 if (bEditAlpha) rgbDest.rgbReserved = (BYTE)(rgb1.rgbReserved&rgb2.rgbReserved);
1129                                         break;
1130                                         case OpXor:
1131                                                 rgbDest.rgbBlue = (BYTE)(rgb1.rgbBlue^rgb2.rgbBlue);
1132                                                 rgbDest.rgbGreen = (BYTE)(rgb1.rgbGreen^rgb2.rgbGreen);
1133                                                 rgbDest.rgbRed = (BYTE)(rgb1.rgbRed^rgb2.rgbRed);
1134                                                 if (bEditAlpha) rgbDest.rgbReserved = (BYTE)(rgb1.rgbReserved^rgb2.rgbReserved);
1135                                         break;
1136                                         case OpOr:
1137                                                 rgbDest.rgbBlue = (BYTE)(rgb1.rgbBlue|rgb2.rgbBlue);
1138                                                 rgbDest.rgbGreen = (BYTE)(rgb1.rgbGreen|rgb2.rgbGreen);
1139                                                 rgbDest.rgbRed = (BYTE)(rgb1.rgbRed|rgb2.rgbRed);
1140                                                 if (bEditAlpha) rgbDest.rgbReserved = (BYTE)(rgb1.rgbReserved|rgb2.rgbReserved);
1141                                         break;
1142                                         case OpMask:
1143                                                 if(rgb2.rgbBlue==0 && rgb2.rgbGreen==0 && rgb2.rgbRed==0)
1144                                                         rgbDest = rgbBackgrnd1;
1145                                                 else
1146                                                         rgbDest = rgb1;
1147                                                 break;
1148                                         case OpSrcCopy:
1149                                                 if(IsTransparent(lX,lY))
1150                                                         rgbDest = rgb2;
1151                                                 else // copy straight over
1152                                                         rgbDest = rgb1;
1153                                                 break;
1154                                         case OpDstCopy:
1155                                                 if(imgsrc2.IsTransparent(lX+lXOffset,lY+lYOffset))
1156                                                         rgbDest = rgb1;
1157                                                 else // copy straight over
1158                                                         rgbDest = rgb2;
1159                                                 break;
1160                                         case OpScreen:
1161                                                 {
1162                                                         BYTE a,a1;
1163
1164                                                         if (imgsrc2.IsTransparent(lX+lXOffset,lY+lYOffset)){
1165                                                                 a=0;
1166                                                         } else if (imgsrc2.AlphaIsValid()){
1167                                                                 a=imgsrc2.AlphaGet(lX+lXOffset,lY+lYOffset);
1168                                                                 a =(BYTE)((a*imgsrc2.info.nAlphaMax)/255);
1169                                                         } else {
1170                                                                 a=255;
1171                                                         }
1172
1173                                                         if (a==0){ //transparent
1174                                                                 rgbDest = rgb1;
1175                                                         } else if (a==255){ //opaque
1176                                                                 rgbDest = rgb2;
1177                                                         } else { //blend
1178                                                                 a1 = (BYTE)~a;
1179                                                                 rgbDest.rgbBlue = (BYTE)((rgb1.rgbBlue*a1+rgb2.rgbBlue*a)/255);
1180                                                                 rgbDest.rgbGreen = (BYTE)((rgb1.rgbGreen*a1+rgb2.rgbGreen*a)/255);
1181                                                                 rgbDest.rgbRed = (BYTE)((rgb1.rgbRed*a1+rgb2.rgbRed*a)/255);
1182                                                         }
1183
1184                                                         if (bEditAlpha) rgbDest.rgbReserved = (BYTE)((rgb1.rgbReserved*a)/255);
1185                                                 }
1186                                                 break;
1187                                         case OpSrcBlend:
1188                                                 if(IsTransparent(lX,lY))
1189                                                         rgbDest = rgb2;
1190                                                 else
1191                                                 {
1192                                                         long lBDiff = abs(rgb1.rgbBlue - rgbBackgrnd1.rgbBlue);
1193                                                         long lGDiff = abs(rgb1.rgbGreen - rgbBackgrnd1.rgbGreen);
1194                                                         long lRDiff = abs(rgb1.rgbRed - rgbBackgrnd1.rgbRed);
1195
1196                                                         double lAverage = (lBDiff+lGDiff+lRDiff)/3;
1197                                                         double lThresh = 16;
1198                                                         double dLarge = lAverage/lThresh;
1199                                                         double dSmall = (lThresh-lAverage)/lThresh;
1200                                                         double dSmallAmt = dSmall*((double)rgb2.rgbBlue);
1201
1202                                                         if( lAverage < lThresh+1){
1203                                                                 rgbDest.rgbBlue = (BYTE)max(0,min(255,(int)(dLarge*((double)rgb1.rgbBlue) +
1204                                                                                                 dSmallAmt)));
1205                                                                 rgbDest.rgbGreen = (BYTE)max(0,min(255,(int)(dLarge*((double)rgb1.rgbGreen) +
1206                                                                                                 dSmallAmt)));
1207                                                                 rgbDest.rgbRed = (BYTE)max(0,min(255,(int)(dLarge*((double)rgb1.rgbRed) +
1208                                                                                                 dSmallAmt)));
1209                                                         }
1210                                                         else
1211                                                                 rgbDest = rgb1;
1212                                                 }
1213                                                 break;
1214                                                 default:
1215                                                 return;
1216                                 }
1217                                 SetPixelColor(lX,lY,rgbDest,bEditAlpha);
1218                         }
1219                 }
1220         }
1221 }
1222 ////////////////////////////////////////////////////////////////////////////////
1223 // thanks to Kenneth Ballard
1224 void CxImage::MixFrom(CxImage & imagesrc2, long lXOffset, long lYOffset)
1225 {
1226     long width = imagesrc2.GetWidth();
1227     long height = imagesrc2.GetHeight();
1228
1229     int x, y;
1230
1231         if (imagesrc2.IsTransparent()) {
1232                 for(x = 0; x < width; x++) {
1233                         for(y = 0; y < height; y++) {
1234                                 if(!imagesrc2.IsTransparent(x,y)){
1235                                         SetPixelColor(x + lXOffset, y + lYOffset, imagesrc2.BlindGetPixelColor(x, y));
1236                                 }
1237                         }
1238                 }
1239         } else { //no transparency so just set it <Matt>
1240                 for(x = 0; x < width; x++) {
1241                         for(y = 0; y < height; y++) {
1242                                 SetPixelColor(x + lXOffset, y + lYOffset, imagesrc2.BlindGetPixelColor(x, y));
1243                         }
1244                 }
1245         }
1246 }
1247 ////////////////////////////////////////////////////////////////////////////////
1248 /**
1249  * Adjusts separately the red, green, and blue values in the image.
1250  * \param r, g, b: can be from -255 to +255.
1251  * \return true if everything is ok
1252  */
1253 bool CxImage::ShiftRGB(long r, long g, long b)
1254 {
1255         if (!pDib) return false;
1256         RGBQUAD color;
1257         if (head.biClrUsed==0){
1258
1259                 long xmin,xmax,ymin,ymax;
1260                 if (pSelection){
1261                         xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
1262                         ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
1263                 } else {
1264                         xmin = ymin = 0;
1265                         xmax = head.biWidth; ymax=head.biHeight;
1266                 }
1267
1268                 for(long y=ymin; y<ymax; y++){
1269                         for(long x=xmin; x<xmax; x++){
1270 #if CXIMAGE_SUPPORT_SELECTION
1271                                 if (BlindSelectionIsInside(x,y))
1272 #endif //CXIMAGE_SUPPORT_SELECTION
1273                                 {
1274                                         color = BlindGetPixelColor(x,y);
1275                                         color.rgbRed = (BYTE)max(0,min(255,(int)(color.rgbRed + r)));
1276                                         color.rgbGreen = (BYTE)max(0,min(255,(int)(color.rgbGreen + g)));
1277                                         color.rgbBlue = (BYTE)max(0,min(255,(int)(color.rgbBlue + b)));
1278                                         BlindSetPixelColor(x,y,color);
1279                                 }
1280                         }
1281                 }
1282         } else {
1283                 for(DWORD j=0; j<head.biClrUsed; j++){
1284                         color = GetPaletteColor((BYTE)j);
1285                         color.rgbRed = (BYTE)max(0,min(255,(int)(color.rgbRed + r)));
1286                         color.rgbGreen = (BYTE)max(0,min(255,(int)(color.rgbGreen + g)));
1287                         color.rgbBlue = (BYTE)max(0,min(255,(int)(color.rgbBlue + b)));
1288                         SetPaletteColor((BYTE)j,color);
1289                 }
1290         }
1291         return true;
1292 }
1293 ////////////////////////////////////////////////////////////////////////////////
1294 /**
1295  * Adjusts the color balance of the image
1296  * \param gamma can be from 0.1 to 5.
1297  * \return true if everything is ok
1298  * \sa GammaRGB
1299  */
1300 bool CxImage::Gamma(float gamma)
1301 {
1302         if (!pDib) return false;
1303
1304         if (gamma <= 0.0f) return false;
1305
1306         double dinvgamma = 1/gamma;
1307         double dMax = pow(255.0, dinvgamma) / 255.0;
1308
1309         BYTE cTable[256]; //<nipper>
1310         for (int i=0;i<256;i++) {
1311                 cTable[i] = (BYTE)max(0,min(255,(int)( pow((double)i, dinvgamma) / dMax)));
1312         }
1313
1314         return Lut(cTable);
1315 }
1316 ////////////////////////////////////////////////////////////////////////////////
1317 /**
1318  * Adjusts the color balance indipendent for each color channel
1319  * \param gammaR, gammaG, gammaB  can be from 0.1 to 5.
1320  * \return true if everything is ok
1321  * \sa Gamma
1322  */
1323 bool CxImage::GammaRGB(float gammaR, float gammaG, float gammaB)
1324 {
1325         if (!pDib) return false;
1326
1327         if (gammaR <= 0.0f) return false;
1328         if (gammaG <= 0.0f) return false;
1329         if (gammaB <= 0.0f) return false;
1330
1331         double dinvgamma, dMax;
1332         int i;
1333
1334         dinvgamma = 1/gammaR;
1335         dMax = pow(255.0, dinvgamma) / 255.0;
1336         BYTE cTableR[256];
1337         for (i=0;i<256;i++)     {
1338                 cTableR[i] = (BYTE)max(0,min(255,(int)( pow((double)i, dinvgamma) / dMax)));
1339         }
1340
1341         dinvgamma = 1/gammaG;
1342         dMax = pow(255.0, dinvgamma) / 255.0;
1343         BYTE cTableG[256];
1344         for (i=0;i<256;i++)     {
1345                 cTableG[i] = (BYTE)max(0,min(255,(int)( pow((double)i, dinvgamma) / dMax)));
1346         }
1347
1348         dinvgamma = 1/gammaB;
1349         dMax = pow(255.0, dinvgamma) / 255.0;
1350         BYTE cTableB[256];
1351         for (i=0;i<256;i++)     {
1352                 cTableB[i] = (BYTE)max(0,min(255,(int)( pow((double)i, dinvgamma) / dMax)));
1353         }
1354
1355         return Lut(cTableR, cTableG, cTableB);
1356 }
1357 ////////////////////////////////////////////////////////////////////////////////
1358
1359 //#if !defined (_WIN32_WCE)
1360 /**
1361  * Adjusts the intensity of each pixel to the median intensity of its surrounding pixels.
1362  * \param Ksize: size of the kernel.
1363  * \return true if everything is ok
1364  */
1365 bool CxImage::Median(long Ksize)
1366 {
1367         if (!pDib) return false;
1368
1369         long k2 = Ksize/2;
1370         long kmax= Ksize-k2;
1371         long i,j,k;
1372
1373         RGBQUAD* kernel = (RGBQUAD*)malloc(Ksize*Ksize*sizeof(RGBQUAD));
1374
1375         CxImage tmp(*this);
1376         if (!tmp.IsValid()){
1377                 strcpy(info.szLastError,tmp.GetLastError());
1378                 return false;
1379         }
1380
1381         long xmin,xmax,ymin,ymax;
1382         if (pSelection){
1383                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
1384                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
1385         } else {
1386                 xmin = ymin = 0;
1387                 xmax = head.biWidth; ymax=head.biHeight;
1388         }
1389
1390         for(long y=ymin; y<ymax; y++){
1391                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
1392                 if (info.nEscape) break;
1393                 for(long x=xmin; x<xmax; x++){
1394 #if CXIMAGE_SUPPORT_SELECTION
1395                         if (BlindSelectionIsInside(x,y))
1396 #endif //CXIMAGE_SUPPORT_SELECTION
1397                                 {
1398                                 for(j=-k2, i=0;j<kmax;j++)
1399                                         for(k=-k2;k<kmax;k++)
1400                                                 if (IsInside(x+j,y+k))
1401                                                         kernel[i++]=BlindGetPixelColor(x+j,y+k);
1402
1403                                 qsort(kernel, i, sizeof(RGBQUAD), CompareColors);
1404                                 tmp.SetPixelColor(x,y,kernel[i/2]);
1405                         }
1406                 }
1407         }
1408         free(kernel);
1409         Transfer(tmp);
1410         return true;
1411 }
1412 //#endif //_WIN32_WCE
1413 ////////////////////////////////////////////////////////////////////////////////
1414 /**
1415  * Adds an uniform noise to the image
1416  * \param level: can be from 0 (no noise) to 255 (lot of noise).
1417  * \return true if everything is ok
1418  */
1419 bool CxImage::Noise(long level)
1420 {
1421         if (!pDib) return false;
1422         RGBQUAD color;
1423
1424         long xmin,xmax,ymin,ymax,n;
1425         if (pSelection){
1426                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
1427                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
1428         } else {
1429                 xmin = ymin = 0;
1430                 xmax = head.biWidth; ymax=head.biHeight;
1431         }
1432
1433         for(long y=ymin; y<ymax; y++){
1434                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin)); //<zhanghk><Anatoly Ivasyuk>
1435                 for(long x=xmin; x<xmax; x++){
1436 #if CXIMAGE_SUPPORT_SELECTION
1437                         if (BlindSelectionIsInside(x,y))
1438 #endif //CXIMAGE_SUPPORT_SELECTION
1439                         {
1440                                 color = BlindGetPixelColor(x,y);
1441                                 n=(long)((rand()/(float)RAND_MAX - 0.5)*level);
1442                                 color.rgbRed = (BYTE)max(0,min(255,(int)(color.rgbRed + n)));
1443                                 n=(long)((rand()/(float)RAND_MAX - 0.5)*level);
1444                                 color.rgbGreen = (BYTE)max(0,min(255,(int)(color.rgbGreen + n)));
1445                                 n=(long)((rand()/(float)RAND_MAX - 0.5)*level);
1446                                 color.rgbBlue = (BYTE)max(0,min(255,(int)(color.rgbBlue + n)));
1447                                 BlindSetPixelColor(x,y,color);
1448                         }
1449                 }
1450         }
1451         return true;
1452 }
1453 ////////////////////////////////////////////////////////////////////////////////
1454 /**
1455  * Computes the bidimensional FFT or DFT of the image.
1456  * - The images are processed as grayscale
1457  * - If the dimensions of the image are a power of, 2 the FFT is performed automatically.
1458  * - If dstReal and/or dstImag are NULL, the resulting images replaces the original(s).
1459  * - Note: with 8 bits there is a HUGE loss in the dynamics. The function tries
1460  *   to keep an acceptable SNR, but 8bit = 48dB...
1461  *
1462  * \param srcReal, srcImag: source images: One can be NULL, but not both
1463  * \param dstReal, dstImag: destination images. Can be NULL.
1464  * \param direction: 1 = forward, -1 = inverse.
1465  * \param bForceFFT: if true, the images are resampled to make the dimensions a power of 2.
1466  * \param bMagnitude: if true, the real part returns the magnitude, the imaginary part returns the phase
1467  * \return true if everything is ok
1468  */
1469 bool CxImage::FFT2(CxImage* srcReal, CxImage* srcImag, CxImage* dstReal, CxImage* dstImag,
1470                                    long direction, bool bForceFFT, bool bMagnitude)
1471 {
1472         //check if there is something to convert
1473         if (srcReal==NULL && srcImag==NULL) return false;
1474
1475         long w,h;
1476         //get width and height
1477         if (srcReal) {
1478                 w=srcReal->GetWidth();
1479                 h=srcReal->GetHeight();
1480         } else {
1481                 w=srcImag->GetWidth();
1482                 h=srcImag->GetHeight();
1483         }
1484
1485         bool bXpow2 = IsPowerof2(w);
1486         bool bYpow2 = IsPowerof2(h);
1487         //if bForceFFT, width AND height must be powers of 2
1488         if (bForceFFT && !(bXpow2 && bYpow2)) {
1489                 long i;
1490
1491                 i=0;
1492                 while((1<<i)<w) i++;
1493                 w=1<<i;
1494                 bXpow2=true;
1495
1496                 i=0;
1497                 while((1<<i)<h) i++;
1498                 h=1<<i;
1499                 bYpow2=true;
1500         }
1501
1502         // I/O images for FFT
1503         CxImage *tmpReal,*tmpImag;
1504
1505         // select output
1506         tmpReal = (dstReal) ? dstReal : srcReal;
1507         tmpImag = (dstImag) ? dstImag : srcImag;
1508
1509         // src!=dst -> copy the image
1510         if (srcReal && dstReal) tmpReal->Copy(*srcReal,true,false,false);
1511         if (srcImag && dstImag) tmpImag->Copy(*srcImag,true,false,false);
1512
1513         // dst&&src are empty -> create new one, else turn to GrayScale
1514         if (srcReal==0 && dstReal==0){
1515                 tmpReal = new CxImage(w,h,8);
1516                 tmpReal->Clear(0);
1517                 tmpReal->SetGrayPalette();
1518         } else {
1519                 if (!tmpReal->IsGrayScale()) tmpReal->GrayScale();
1520         }
1521         if (srcImag==0 && dstImag==0){
1522                 tmpImag = new CxImage(w,h,8);
1523                 tmpImag->Clear(0);
1524                 tmpImag->SetGrayPalette();
1525         } else {
1526                 if (!tmpImag->IsGrayScale()) tmpImag->GrayScale();
1527         }
1528
1529         if (!(tmpReal->IsValid() && tmpImag->IsValid())){
1530                 if (srcReal==0 && dstReal==0) delete tmpReal;
1531                 if (srcImag==0 && dstImag==0) delete tmpImag;
1532                 return false;
1533         }
1534
1535         //resample for FFT, if necessary
1536         tmpReal->Resample(w,h,0);
1537         tmpImag->Resample(w,h,0);
1538
1539         //ok, here we have 2 (w x h), grayscale images ready for a FFT
1540
1541         double* real;
1542         double* imag;
1543         long j,k,m;
1544
1545         _complex **grid;
1546         //double mean = tmpReal->Mean();
1547         /* Allocate memory for the grid */
1548         grid = (_complex **)malloc(w * sizeof(_complex));
1549         for (k=0;k<w;k++) {
1550                 grid[k] = (_complex *)malloc(h * sizeof(_complex));
1551         }
1552         for (j=0;j<h;j++) {
1553                 for (k=0;k<w;k++) {
1554                         grid[k][j].x = tmpReal->GetPixelIndex(k,j)-128;
1555                         grid[k][j].y = tmpImag->GetPixelIndex(k,j)-128;
1556                 }
1557         }
1558
1559         //DFT buffers
1560         double *real2,*imag2;
1561         real2 = (double*)malloc(max(w,h) * sizeof(double));
1562         imag2 = (double*)malloc(max(w,h) * sizeof(double));
1563
1564         /* Transform the rows */
1565         real = (double *)malloc(w * sizeof(double));
1566         imag = (double *)malloc(w * sizeof(double));
1567
1568         m=0;
1569         while((1<<m)<w) m++;
1570
1571         for (j=0;j<h;j++) {
1572                 for (k=0;k<w;k++) {
1573                         real[k] = grid[k][j].x;
1574                         imag[k] = grid[k][j].y;
1575                 }
1576
1577                 if (bXpow2) FFT(direction,m,real,imag);
1578                 else            DFT(direction,w,real,imag,real2,imag2);
1579
1580                 for (k=0;k<w;k++) {
1581                         grid[k][j].x = real[k];
1582                         grid[k][j].y = imag[k];
1583                 }
1584         }
1585         free(real);
1586         free(imag);
1587
1588         /* Transform the columns */
1589         real = (double *)malloc(h * sizeof(double));
1590         imag = (double *)malloc(h * sizeof(double));
1591
1592         m=0;
1593         while((1<<m)<h) m++;
1594
1595         for (k=0;k<w;k++) {
1596                 for (j=0;j<h;j++) {
1597                         real[j] = grid[k][j].x;
1598                         imag[j] = grid[k][j].y;
1599                 }
1600
1601                 if (bYpow2) FFT(direction,m,real,imag);
1602                 else            DFT(direction,h,real,imag,real2,imag2);
1603
1604                 for (j=0;j<h;j++) {
1605                         grid[k][j].x = real[j];
1606                         grid[k][j].y = imag[j];
1607                 }
1608         }
1609         free(real);
1610         free(imag);
1611
1612         free(real2);
1613         free(imag2);
1614
1615         /* converting from double to byte, there is a HUGE loss in the dynamics
1616           "nn" tries to keep an acceptable SNR, but 8bit=48dB: don't ask more */
1617         double nn=pow((double)2,(double)log((double)max(w,h))/(double)log((double)2)-4);
1618         //reversed gain for reversed transform
1619         if (direction==-1) nn=1/nn;
1620         //bMagnitude : just to see it on the screen
1621         if (bMagnitude) nn*=4;
1622
1623         for (j=0;j<h;j++) {
1624                 for (k=0;k<w;k++) {
1625                         if (bMagnitude){
1626                                 tmpReal->SetPixelIndex(k,j,(BYTE)max(0,min(255,(nn*(3+log(_cabs(grid[k][j])))))));
1627                                 if (grid[k][j].x==0){
1628                                         tmpImag->SetPixelIndex(k,j,(BYTE)max(0,min(255,(128+(atan(grid[k][j].y/0.0000000001)*nn)))));
1629                                 } else {
1630                                         tmpImag->SetPixelIndex(k,j,(BYTE)max(0,min(255,(128+(atan(grid[k][j].y/grid[k][j].x)*nn)))));
1631                                 }
1632                         } else {
1633                                 tmpReal->SetPixelIndex(k,j,(BYTE)max(0,min(255,(128 + grid[k][j].x*nn))));
1634                                 tmpImag->SetPixelIndex(k,j,(BYTE)max(0,min(255,(128 + grid[k][j].y*nn))));
1635                         }
1636                 }
1637         }
1638
1639         for (k=0;k<w;k++) free (grid[k]);
1640         free (grid);
1641
1642         if (srcReal==0 && dstReal==0) delete tmpReal;
1643         if (srcImag==0 && dstImag==0) delete tmpImag;
1644
1645         return true;
1646 }
1647 ////////////////////////////////////////////////////////////////////////////////
1648 bool CxImage::IsPowerof2(long x)
1649 {
1650         long i=0;
1651         while ((1<<i)<x) i++;
1652         if (x==(1<<i)) return true;
1653         return false;
1654 }
1655 ////////////////////////////////////////////////////////////////////////////////
1656 /**
1657    This computes an in-place complex-to-complex FFT
1658    x and y are the real and imaginary arrays of n=2^m points.
1659    o(n)=n*log2(n)
1660    dir =  1 gives forward transform
1661    dir = -1 gives reverse transform
1662    Written by Paul Bourke, July 1998
1663    FFT algorithm by Cooley and Tukey, 1965
1664 */
1665 bool CxImage::FFT(int dir,int m,double *x,double *y)
1666 {
1667         long nn,i,i1,j,k,i2,l,l1,l2;
1668         double c1,c2,tx,ty,t1,t2,u1,u2,z;
1669
1670         /* Calculate the number of points */
1671         nn = 1<<m;
1672
1673         /* Do the bit reversal */
1674         i2 = nn >> 1;
1675         j = 0;
1676         for (i=0;i<nn-1;i++) {
1677                 if (i < j) {
1678                         tx = x[i];
1679                         ty = y[i];
1680                         x[i] = x[j];
1681                         y[i] = y[j];
1682                         x[j] = tx;
1683                         y[j] = ty;
1684                 }
1685                 k = i2;
1686                 while (k <= j) {
1687                         j -= k;
1688                         k >>= 1;
1689                 }
1690                 j += k;
1691         }
1692
1693         /* Compute the FFT */
1694         c1 = -1.0;
1695         c2 = 0.0;
1696         l2 = 1;
1697         for (l=0;l<m;l++) {
1698                 l1 = l2;
1699                 l2 <<= 1;
1700                 u1 = 1.0;
1701                 u2 = 0.0;
1702                 for (j=0;j<l1;j++) {
1703                         for (i=j;i<nn;i+=l2) {
1704                                 i1 = i + l1;
1705                                 t1 = u1 * x[i1] - u2 * y[i1];
1706                                 t2 = u1 * y[i1] + u2 * x[i1];
1707                                 x[i1] = x[i] - t1;
1708                                 y[i1] = y[i] - t2;
1709                                 x[i] += t1;
1710                                 y[i] += t2;
1711                         }
1712                         z =  u1 * c1 - u2 * c2;
1713                         u2 = u1 * c2 + u2 * c1;
1714                         u1 = z;
1715                 }
1716                 c2 = sqrt((1.0 - c1) / 2.0);
1717                 if (dir == 1)
1718                         c2 = -c2;
1719                 c1 = sqrt((1.0 + c1) / 2.0);
1720         }
1721
1722         /* Scaling for forward transform */
1723         if (dir == 1) {
1724                 for (i=0;i<nn;i++) {
1725                         x[i] /= (double)nn;
1726                         y[i] /= (double)nn;
1727                 }
1728         }
1729
1730    return true;
1731 }
1732 ////////////////////////////////////////////////////////////////////////////////
1733 /**
1734    Direct fourier transform o(n)=n^2
1735    Written by Paul Bourke, July 1998
1736 */
1737 bool CxImage::DFT(int dir,long m,double *x1,double *y1,double *x2,double *y2)
1738 {
1739    long i,k;
1740    double arg;
1741    double cosarg,sinarg;
1742
1743    for (i=0;i<m;i++) {
1744       x2[i] = 0;
1745       y2[i] = 0;
1746       arg = - dir * 2.0 * PI * i / (double)m;
1747       for (k=0;k<m;k++) {
1748          cosarg = cos(k * arg);
1749          sinarg = sin(k * arg);
1750          x2[i] += (x1[k] * cosarg - y1[k] * sinarg);
1751          y2[i] += (x1[k] * sinarg + y1[k] * cosarg);
1752       }
1753    }
1754
1755    /* Copy the data back */
1756    if (dir == 1) {
1757       for (i=0;i<m;i++) {
1758          x1[i] = x2[i] / m;
1759          y1[i] = y2[i] / m;
1760       }
1761    } else {
1762       for (i=0;i<m;i++) {
1763          x1[i] = x2[i];
1764          y1[i] = y2[i];
1765       }
1766    }
1767
1768    return true;
1769 }
1770 ////////////////////////////////////////////////////////////////////////////////
1771 /**
1772  * Combines different color components into a single image
1773  * \param r,g,b: color channels
1774  * \param a: alpha layer, can be NULL
1775  * \param colorspace: 0 = RGB, 1 = HSL, 2 = YUV, 3 = YIQ, 4 = XYZ
1776  * \return true if everything is ok
1777  */
1778 bool CxImage::Combine(CxImage* r,CxImage* g,CxImage* b,CxImage* a, long colorspace)
1779 {
1780         if (r==0 || g==0 || b==0) return false;
1781
1782         long w = r->GetWidth();
1783         long h = r->GetHeight();
1784
1785         Create(w,h,24);
1786
1787         g->Resample(w,h);
1788         b->Resample(w,h);
1789
1790         if (a) {
1791                 a->Resample(w,h);
1792 #if CXIMAGE_SUPPORT_ALPHA
1793                 AlphaCreate();
1794 #endif //CXIMAGE_SUPPORT_ALPHA
1795         }
1796
1797         RGBQUAD c;
1798         for (long y=0;y<h;y++){
1799                 info.nProgress = (long)(100*y/h); //<Anatoly Ivasyuk>
1800                 for (long x=0;x<w;x++){
1801                         c.rgbRed=r->GetPixelIndex(x,y);
1802                         c.rgbGreen=g->GetPixelIndex(x,y);
1803                         c.rgbBlue=b->GetPixelIndex(x,y);
1804                         switch (colorspace){
1805                         case 1:
1806                                 BlindSetPixelColor(x,y,HSLtoRGB(c));
1807                                 break;
1808                         case 2:
1809                                 BlindSetPixelColor(x,y,YUVtoRGB(c));
1810                                 break;
1811                         case 3:
1812                                 BlindSetPixelColor(x,y,YIQtoRGB(c));
1813                                 break;
1814                         case 4:
1815                                 BlindSetPixelColor(x,y,XYZtoRGB(c));
1816                                 break;
1817                         default:
1818                                 BlindSetPixelColor(x,y,c);
1819                         }
1820 #if CXIMAGE_SUPPORT_ALPHA
1821                         if (a) AlphaSet(x,y,a->GetPixelIndex(x,y));
1822 #endif //CXIMAGE_SUPPORT_ALPHA
1823                 }
1824         }
1825
1826         return true;
1827 }
1828 ////////////////////////////////////////////////////////////////////////////////
1829 /**
1830  * Smart blurring to remove small defects, dithering or artifacts.
1831  * \param radius: normally between 0.01 and 0.5
1832  * \param niterations: should be trimmed with radius, to avoid blurring should be (radius*niterations)<1
1833  * \param colorspace: 0 = RGB, 1 = HSL, 2 = YUV, 3 = YIQ, 4 = XYZ
1834  * \return true if everything is ok
1835  */
1836 bool CxImage::Repair(float radius, long niterations, long colorspace)
1837 {
1838         if (!IsValid()) return false;
1839
1840         long w = GetWidth();
1841         long h = GetHeight();
1842
1843         CxImage r,g,b;
1844
1845         r.Create(w,h,8);
1846         g.Create(w,h,8);
1847         b.Create(w,h,8);
1848
1849         switch (colorspace){
1850         case 1:
1851                 SplitHSL(&r,&g,&b);
1852                 break;
1853         case 2:
1854                 SplitYUV(&r,&g,&b);
1855                 break;
1856         case 3:
1857                 SplitYIQ(&r,&g,&b);
1858                 break;
1859         case 4:
1860                 SplitXYZ(&r,&g,&b);
1861                 break;
1862         default:
1863                 SplitRGB(&r,&g,&b);
1864         }
1865
1866         for (int i=0; i<niterations; i++){
1867                 RepairChannel(&r,radius);
1868                 RepairChannel(&g,radius);
1869                 RepairChannel(&b,radius);
1870         }
1871
1872         CxImage* a=NULL;
1873 #if CXIMAGE_SUPPORT_ALPHA
1874         if (AlphaIsValid()){
1875                 a = new CxImage();
1876                 AlphaSplit(a);
1877         }
1878 #endif
1879
1880         Combine(&r,&g,&b,a,colorspace);
1881
1882         delete a;
1883
1884         return true;
1885 }
1886 ////////////////////////////////////////////////////////////////////////////////
1887 bool CxImage::RepairChannel(CxImage *ch, float radius)
1888 {
1889         if (ch==NULL) return false;
1890
1891         CxImage tmp(*ch);
1892         if (!tmp.IsValid()){
1893                 strcpy(info.szLastError,tmp.GetLastError());
1894                 return false;
1895         }
1896
1897         long w = ch->GetWidth()-1;
1898         long h = ch->GetHeight()-1;
1899
1900         double correction,ix,iy,ixx,ixy,iyy;
1901         int x,y,xy0,xp1,xm1,yp1,ym1;
1902
1903         for(x=1; x<w; x++){
1904                 for(y=1; y<h; y++){
1905
1906                         xy0 = ch->BlindGetPixelIndex(x,y);
1907                         xm1 = ch->BlindGetPixelIndex(x-1,y);
1908                         xp1 = ch->BlindGetPixelIndex(x+1,y);
1909                         ym1 = ch->BlindGetPixelIndex(x,y-1);
1910                         yp1 = ch->BlindGetPixelIndex(x,y+1);
1911
1912                         ix= (xp1-xm1)/2.0;
1913                         iy= (yp1-ym1)/2.0;
1914                         ixx= xp1 - 2.0 * xy0 + xm1;
1915                         iyy= yp1 - 2.0 * xy0 + ym1;
1916                         ixy=(ch->BlindGetPixelIndex(x+1,y+1) + ch->BlindGetPixelIndex(x-1,y-1) -
1917                                  ch->BlindGetPixelIndex(x-1,y+1) - ch->BlindGetPixelIndex(x+1,y-1))/4.0;
1918
1919                         correction = ((1.0+iy*iy)*ixx - ix*iy*ixy + (1.0+ix*ix)*iyy)/(1.0+ix*ix+iy*iy);
1920
1921                         tmp.BlindSetPixelIndex(x,y,(BYTE)min(255,max(0,(xy0 + radius * correction + 0.5))));
1922                 }
1923         }
1924
1925         for (x=0;x<=w;x++){
1926                 for(y=0; y<=h; y+=h){
1927                         xy0 = ch->BlindGetPixelIndex(x,y);
1928                         xm1 = ch->GetPixelIndex(x-1,y);
1929                         xp1 = ch->GetPixelIndex(x+1,y);
1930                         ym1 = ch->GetPixelIndex(x,y-1);
1931                         yp1 = ch->GetPixelIndex(x,y+1);
1932
1933                         ix= (xp1-xm1)/2.0;
1934                         iy= (yp1-ym1)/2.0;
1935                         ixx= xp1 - 2.0 * xy0 + xm1;
1936                         iyy= yp1 - 2.0 * xy0 + ym1;
1937                         ixy=(ch->GetPixelIndex(x+1,y+1) + ch->GetPixelIndex(x-1,y-1) -
1938                                  ch->GetPixelIndex(x-1,y+1) - ch->GetPixelIndex(x+1,y-1))/4.0;
1939
1940                         correction = ((1.0+iy*iy)*ixx - ix*iy*ixy + (1.0+ix*ix)*iyy)/(1.0+ix*ix+iy*iy);
1941
1942                         tmp.BlindSetPixelIndex(x,y,(BYTE)min(255,max(0,(xy0 + radius * correction + 0.5))));
1943                 }
1944         }
1945         for (x=0;x<=w;x+=w){
1946                 for (y=0;y<=h;y++){
1947                         xy0 = ch->BlindGetPixelIndex(x,y);
1948                         xm1 = ch->GetPixelIndex(x-1,y);
1949                         xp1 = ch->GetPixelIndex(x+1,y);
1950                         ym1 = ch->GetPixelIndex(x,y-1);
1951                         yp1 = ch->GetPixelIndex(x,y+1);
1952
1953                         ix= (xp1-xm1)/2.0;
1954                         iy= (yp1-ym1)/2.0;
1955                         ixx= xp1 - 2.0 * xy0 + xm1;
1956                         iyy= yp1 - 2.0 * xy0 + ym1;
1957                         ixy=(ch->GetPixelIndex(x+1,y+1) + ch->GetPixelIndex(x-1,y-1) -
1958                                  ch->GetPixelIndex(x-1,y+1) - ch->GetPixelIndex(x+1,y-1))/4.0;
1959
1960                         correction = ((1.0+iy*iy)*ixx - ix*iy*ixy + (1.0+ix*ix)*iyy)/(1.0+ix*ix+iy*iy);
1961
1962                         tmp.BlindSetPixelIndex(x,y,(BYTE)min(255,max(0,(xy0 + radius * correction + 0.5))));
1963                 }
1964         }
1965
1966         ch->Transfer(tmp);
1967         return true;
1968 }
1969 ////////////////////////////////////////////////////////////////////////////////
1970 /**
1971  * Enhance the variations between adjacent pixels.
1972  * Similar results can be achieved using Filter(),
1973  * but the algorithms are different both in Edge() and in Contour().
1974  * \return true if everything is ok
1975  */
1976 bool CxImage::Contour()
1977 {
1978         if (!pDib) return false;
1979
1980         long Ksize = 3;
1981         long k2 = Ksize/2;
1982         long kmax= Ksize-k2;
1983         long i,j,k;
1984         BYTE maxr,maxg,maxb;
1985         RGBQUAD pix1,pix2;
1986
1987         CxImage tmp(*this);
1988         if (!tmp.IsValid()){
1989                 strcpy(info.szLastError,tmp.GetLastError());
1990                 return false;
1991         }
1992
1993         long xmin,xmax,ymin,ymax;
1994         if (pSelection){
1995                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
1996                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
1997         } else {
1998                 xmin = ymin = 0;
1999                 xmax = head.biWidth; ymax=head.biHeight;
2000         }
2001
2002         for(long y=ymin; y<ymax; y++){
2003                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
2004                 if (info.nEscape) break;
2005                 for(long x=xmin; x<xmax; x++){
2006 #if CXIMAGE_SUPPORT_SELECTION
2007                         if (BlindSelectionIsInside(x,y))
2008 #endif //CXIMAGE_SUPPORT_SELECTION
2009                                 {
2010                                 pix1 = BlindGetPixelColor(x,y);
2011                                 maxr=maxg=maxb=0;
2012                                 for(j=-k2, i=0;j<kmax;j++){
2013                                         for(k=-k2;k<kmax;k++, i++){
2014                                                 if (!IsInside(x+j,y+k)) continue;
2015                                                 pix2 = BlindGetPixelColor(x+j,y+k);
2016                                                 if ((pix2.rgbBlue-pix1.rgbBlue)>maxb) maxb = pix2.rgbBlue;
2017                                                 if ((pix2.rgbGreen-pix1.rgbGreen)>maxg) maxg = pix2.rgbGreen;
2018                                                 if ((pix2.rgbRed-pix1.rgbRed)>maxr) maxr = pix2.rgbRed;
2019                                         }
2020                                 }
2021                                 pix1.rgbBlue=(BYTE)(255-maxb);
2022                                 pix1.rgbGreen=(BYTE)(255-maxg);
2023                                 pix1.rgbRed=(BYTE)(255-maxr);
2024                                 tmp.BlindSetPixelColor(x,y,pix1);
2025                         }
2026                 }
2027         }
2028         Transfer(tmp);
2029         return true;
2030 }
2031 ////////////////////////////////////////////////////////////////////////////////
2032 /**
2033  * Adds a random offset to each pixel in the image
2034  * \param radius: maximum pixel displacement
2035  * \return true if everything is ok
2036  */
2037 bool CxImage::Jitter(long radius)
2038 {
2039         if (!pDib) return false;
2040
2041         long nx,ny;
2042
2043         CxImage tmp(*this);
2044         if (!tmp.IsValid()){
2045                 strcpy(info.szLastError,tmp.GetLastError());
2046                 return false;
2047         }
2048
2049         long xmin,xmax,ymin,ymax;
2050         if (pSelection){
2051                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2052                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2053         } else {
2054                 xmin = ymin = 0;
2055                 xmax = head.biWidth; ymax=head.biHeight;
2056         }
2057
2058         for(long y=ymin; y<ymax; y++){
2059                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
2060                 if (info.nEscape) break;
2061                 for(long x=xmin; x<xmax; x++){
2062 #if CXIMAGE_SUPPORT_SELECTION
2063                         if (BlindSelectionIsInside(x,y))
2064 #endif //CXIMAGE_SUPPORT_SELECTION
2065                         {
2066                                 nx=x+(long)((rand()/(float)RAND_MAX - 0.5)*(radius*2));
2067                                 ny=y+(long)((rand()/(float)RAND_MAX - 0.5)*(radius*2));
2068                                 if (!IsInside(nx,ny)) {
2069                                         nx=x;
2070                                         ny=y;
2071                                 }
2072                                 if (head.biClrUsed==0){
2073                                         tmp.BlindSetPixelColor(x,y,BlindGetPixelColor(nx,ny));
2074                                 } else {
2075                                         tmp.BlindSetPixelIndex(x,y,BlindGetPixelIndex(nx,ny));
2076                                 }
2077 #if CXIMAGE_SUPPORT_ALPHA
2078                                 tmp.AlphaSet(x,y,AlphaGet(nx,ny));
2079 #endif //CXIMAGE_SUPPORT_ALPHA
2080                         }
2081                 }
2082         }
2083         Transfer(tmp);
2084         return true;
2085 }
2086 ////////////////////////////////////////////////////////////////////////////////
2087 /**
2088  * generates a 1-D convolution matrix to be used for each pass of
2089  * a two-pass gaussian blur.  Returns the length of the matrix.
2090  * \author [nipper]
2091  */
2092 int CxImage::gen_convolve_matrix (float radius, float **cmatrix_p)
2093 {
2094         int matrix_length;
2095         int matrix_midpoint;
2096         float* cmatrix;
2097         int i,j;
2098         float std_dev;
2099         float sum;
2100
2101         /* we want to generate a matrix that goes out a certain radius
2102         * from the center, so we have to go out ceil(rad-0.5) pixels,
2103         * inlcuding the center pixel.  Of course, that's only in one direction,
2104         * so we have to go the same amount in the other direction, but not count
2105         * the center pixel again.  So we double the previous result and subtract
2106         * one.
2107         * The radius parameter that is passed to this function is used as
2108         * the standard deviation, and the radius of effect is the
2109         * standard deviation * 2.  It's a little confusing.
2110         * <DP> modified scaling, so that matrix_lenght = 1+2*radius parameter
2111         */
2112         radius = (float)fabs(0.5*radius) + 0.25f;
2113
2114         std_dev = radius;
2115         radius = std_dev * 2;
2116
2117         /* go out 'radius' in each direction */
2118         matrix_length = int (2 * ceil(radius-0.5) + 1);
2119         if (matrix_length <= 0) matrix_length = 1;
2120         matrix_midpoint = matrix_length/2 + 1;
2121         *cmatrix_p = new float[matrix_length];
2122         cmatrix = *cmatrix_p;
2123
2124         /*  Now we fill the matrix by doing a numeric integration approximation
2125         * from -2*std_dev to 2*std_dev, sampling 50 points per pixel.
2126         * We do the bottom half, mirror it to the top half, then compute the
2127         * center point.  Otherwise asymmetric quantization errors will occur.
2128         *  The formula to integrate is e^-(x^2/2s^2).
2129         */
2130
2131         /* first we do the top (right) half of matrix */
2132         for (i = matrix_length/2 + 1; i < matrix_length; i++)
2133     {
2134                 float base_x = i - (float)floor((float)(matrix_length/2)) - 0.5f;
2135                 sum = 0;
2136                 for (j = 1; j <= 50; j++)
2137                 {
2138                         if ( base_x+0.02*j <= radius )
2139                                 sum += (float)exp (-(base_x+0.02*j)*(base_x+0.02*j) /
2140                                 (2*std_dev*std_dev));
2141                 }
2142                 cmatrix[i] = sum/50;
2143     }
2144
2145         /* mirror the thing to the bottom half */
2146         for (i=0; i<=matrix_length/2; i++) {
2147                 cmatrix[i] = cmatrix[matrix_length-1-i];
2148         }
2149
2150         /* find center val -- calculate an odd number of quanta to make it symmetric,
2151         * even if the center point is weighted slightly higher than others. */
2152         sum = 0;
2153         for (j=0; j<=50; j++)
2154     {
2155                 sum += (float)exp (-(0.5+0.02*j)*(0.5+0.02*j) /
2156                         (2*std_dev*std_dev));
2157     }
2158         cmatrix[matrix_length/2] = sum/51;
2159
2160         /* normalize the distribution by scaling the total sum to one */
2161         sum=0;
2162         for (i=0; i<matrix_length; i++) sum += cmatrix[i];
2163         for (i=0; i<matrix_length; i++) cmatrix[i] = cmatrix[i] / sum;
2164
2165         return matrix_length;
2166 }
2167 ////////////////////////////////////////////////////////////////////////////////
2168 /**
2169  * generates a lookup table for every possible product of 0-255 and
2170  * each value in the convolution matrix.  The returned array is
2171  * indexed first by matrix position, then by input multiplicand (?)
2172  * value.
2173  * \author [nipper]
2174  */
2175 float* CxImage::gen_lookup_table (float *cmatrix, int cmatrix_length)
2176 {
2177         float* lookup_table = new float[cmatrix_length * 256];
2178         float* lookup_table_p = lookup_table;
2179         float* cmatrix_p      = cmatrix;
2180
2181         for (int i=0; i<cmatrix_length; i++)
2182     {
2183                 for (int j=0; j<256; j++)
2184                 {
2185                         *(lookup_table_p++) = *cmatrix_p * (float)j;
2186                 }
2187                 cmatrix_p++;
2188     }
2189
2190         return lookup_table;
2191 }
2192 ////////////////////////////////////////////////////////////////////////////////
2193 /**
2194  * this function is written as if it is blurring a column at a time,
2195  * even though it can operate on rows, too.  There is no difference
2196  * in the processing of the lines, at least to the blur_line function.
2197  * \author [nipper]
2198  */
2199 void CxImage::blur_line (float *ctable, float *cmatrix, int cmatrix_length, BYTE* cur_col, BYTE* dest_col, int y, long bytes)
2200 {
2201         float scale;
2202         float sum;
2203         int i=0, j=0;
2204         int row;
2205         int cmatrix_middle = cmatrix_length/2;
2206
2207         float *cmatrix_p;
2208         BYTE  *cur_col_p;
2209         BYTE  *cur_col_p1;
2210         BYTE  *dest_col_p;
2211         float *ctable_p;
2212
2213         /* this first block is the same as the non-optimized version --
2214         * it is only used for very small pictures, so speed isn't a
2215         * big concern.
2216         */
2217         if (cmatrix_length > y)
2218     {
2219                 for (row = 0; row < y ; row++)
2220                 {
2221                         scale=0;
2222                         /* find the scale factor */
2223                         for (j = 0; j < y ; j++)
2224                         {
2225                                 /* if the index is in bounds, add it to the scale counter */
2226                                 if ((j + cmatrix_middle - row >= 0) &&
2227                                         (j + cmatrix_middle - row < cmatrix_length))
2228                                         scale += cmatrix[j + cmatrix_middle - row];
2229                         }
2230                         for (i = 0; i<bytes; i++)
2231                         {
2232                                 sum = 0;
2233                                 for (j = 0; j < y; j++)
2234                                 {
2235                                         if ((j >= row - cmatrix_middle) &&
2236                                                 (j <= row + cmatrix_middle))
2237                                                 sum += cur_col[j*bytes + i] * cmatrix[j];
2238                                 }
2239                                 dest_col[row*bytes + i] = (BYTE)(0.5f + sum / scale);
2240                         }
2241                 }
2242     }
2243         else
2244     {
2245                 /* for the edge condition, we only use available info and scale to one */
2246                 for (row = 0; row < cmatrix_middle; row++)
2247                 {
2248                         /* find scale factor */
2249                         scale=0;
2250                         for (j = cmatrix_middle - row; j<cmatrix_length; j++)
2251                                 scale += cmatrix[j];
2252                         for (i = 0; i<bytes; i++)
2253                         {
2254                                 sum = 0;
2255                                 for (j = cmatrix_middle - row; j<cmatrix_length; j++)
2256                                 {
2257                                         sum += cur_col[(row + j-cmatrix_middle)*bytes + i] * cmatrix[j];
2258                                 }
2259                                 dest_col[row*bytes + i] = (BYTE)(0.5f + sum / scale);
2260                         }
2261                 }
2262                 /* go through each pixel in each col */
2263                 dest_col_p = dest_col + row*bytes;
2264                 for (; row < y-cmatrix_middle; row++)
2265                 {
2266                         cur_col_p = (row - cmatrix_middle) * bytes + cur_col;
2267                         for (i = 0; i<bytes; i++)
2268                         {
2269                                 sum = 0;
2270                                 cmatrix_p = cmatrix;
2271                                 cur_col_p1 = cur_col_p;
2272                                 ctable_p = ctable;
2273                                 for (j = cmatrix_length; j>0; j--)
2274                                 {
2275                                         sum += *(ctable_p + *cur_col_p1);
2276                                         cur_col_p1 += bytes;
2277                                         ctable_p += 256;
2278                                 }
2279                                 cur_col_p++;
2280                                 *(dest_col_p++) = (BYTE)(0.5f + sum);
2281                         }
2282                 }
2283
2284                 /* for the edge condition , we only use available info, and scale to one */
2285                 for (; row < y; row++)
2286                 {
2287                         /* find scale factor */
2288                         scale=0;
2289                         for (j = 0; j< y-row + cmatrix_middle; j++)
2290                                 scale += cmatrix[j];
2291                         for (i = 0; i<bytes; i++)
2292                         {
2293                                 sum = 0;
2294                                 for (j = 0; j<y-row + cmatrix_middle; j++)
2295                                 {
2296                                         sum += cur_col[(row + j-cmatrix_middle)*bytes + i] * cmatrix[j];
2297                                 }
2298                                 dest_col[row*bytes + i] = (BYTE) (0.5f + sum / scale);
2299                         }
2300                 }
2301     }
2302 }
2303 ////////////////////////////////////////////////////////////////////////////////
2304 /**
2305  * \author [DP]
2306  */
2307 void CxImage::blur_text (BYTE threshold, BYTE decay, BYTE max_depth, CxImage* iSrc, CxImage* iDst, BYTE bytes)
2308 {
2309         long x,y,z,m;
2310         BYTE *pSrc, *pSrc2, *pSrc3, *pDst;
2311         BYTE step,n;
2312         int pivot;
2313
2314         if (max_depth<1) max_depth = 1;
2315
2316         long nmin,nmax,xmin,xmax,ymin,ymax;
2317         xmin = ymin = 0;
2318         xmax = iSrc->head.biWidth;
2319         ymax = iSrc->head.biHeight;
2320
2321         if (xmin==xmax || ymin==ymax) return;
2322
2323         nmin = xmin * bytes;
2324         nmax = xmax * bytes;
2325
2326         CImageIterator itSrc(iSrc);
2327         CImageIterator itTmp(iDst);
2328
2329         double dbScaler = 100.0f/(ymax-ymin)/bytes;
2330
2331         for (n=0; n<bytes; n++){
2332                 for (y=ymin+1;y<(ymax-1);y++)
2333                 {
2334                         if (info.nEscape) break;
2335                         info.nProgress = (long)((y-ymin)*dbScaler*(1+n));
2336
2337                         pSrc  = itSrc.GetRow(y);
2338                         pSrc2 = itSrc.GetRow(y+1);
2339                         pSrc3 = itSrc.GetRow(y-1);
2340                         pDst  = itTmp.GetRow(y);
2341
2342                         //scan left to right
2343                         for (x=n+nmin /*,i=xmin*/; x<(nmax-1); x+=bytes /*,i++*/)
2344                         {
2345                                 z=x+bytes;
2346                                 pivot = pSrc[z]-threshold;
2347                                 //find upper corner
2348                                 if (pSrc[x]<pivot && pSrc2[z]<pivot && pSrc3[x]>=pivot){
2349                                         while (z<nmax && pSrc2[z]<pSrc[x+bytes] && pSrc[x+bytes]<=pSrc[z]){
2350                                                 z+=bytes;
2351                                         }
2352                                         m = z-x;
2353                                         m = (decay>1) ? ((m/bytes)/decay+1) : m/bytes;
2354                                         if (m>max_depth) m = max_depth;
2355                                         step = (BYTE)((pSrc[x+bytes]-pSrc[x])/(m+1));
2356                                         while (m-->1){
2357                                                 pDst[x+m*bytes] = (BYTE)(pDst[x]+(step*(m+1)));
2358                                         }
2359                                 }
2360                                 //find lower corner
2361                                 z=x+bytes;
2362                                 if (pSrc[x]<pivot && pSrc3[z]<pivot && pSrc2[x]>=pivot){
2363                                         while (z<nmax && pSrc3[z]<pSrc[x+bytes] && pSrc[x+bytes]<=pSrc[z]){
2364                                                 z+=bytes;
2365                                         }
2366                                         m = z-x;
2367                                         m = (decay>1) ? ((m/bytes)/decay+1) : m/bytes;
2368                                         if (m>max_depth) m = max_depth;
2369                                         step = (BYTE)((pSrc[x+bytes]-pSrc[x])/(m+1));
2370                                         while (m-->1){
2371                                                 pDst[x+m*bytes] = (BYTE)(pDst[x]+(step*(m+1)));
2372                                         }
2373                                 }
2374                         }
2375                         //scan right to left
2376                         for (x=nmax-1-n /*,i=(xmax-1)*/; x>0; x-=bytes /*,i--*/)
2377                         {
2378                                 z=x-bytes;
2379                                 pivot = pSrc[z]-threshold;
2380                                 //find upper corner
2381                                 if (pSrc[x]<pivot && pSrc2[z]<pivot && pSrc3[x]>=pivot){
2382                                         while (z>n && pSrc2[z]<pSrc[x-bytes] && pSrc[x-bytes]<=pSrc[z]){
2383                                                 z-=bytes;
2384                                         }
2385                                         m = x-z;
2386                                         m = (decay>1) ? ((m/bytes)/decay+1) : m/bytes;
2387                                         if (m>max_depth) m = max_depth;
2388                                         step = (BYTE)((pSrc[x-bytes]-pSrc[x])/(m+1));
2389                                         while (m-->1){
2390                                                 pDst[x-m*bytes] = (BYTE)(pDst[x]+(step*(m+1)));
2391                                         }
2392                                 }
2393                                 //find lower corner
2394                                 z=x-bytes;
2395                                 if (pSrc[x]<pivot && pSrc3[z]<pivot && pSrc2[x]>=pivot){
2396                                         while (z>n && pSrc3[z]<pSrc[x-bytes] && pSrc[x-bytes]<=pSrc[z]){
2397                                                 z-=bytes;
2398                                         }
2399                                         m = x-z;
2400                                         m = (decay>1) ? ((m/bytes)/decay+1) : m/bytes;
2401                                         if (m>max_depth) m = max_depth;
2402                                         step = (BYTE)((pSrc[x-bytes]-pSrc[x])/(m+1));
2403                                         while (m-->1){
2404                                                 pDst[x-m*bytes] = (BYTE)(pDst[x]+(step*(m+1)));
2405                                         }
2406                                 }
2407                         }
2408                 }
2409         }
2410 }
2411 ////////////////////////////////////////////////////////////////////////////////
2412 /**
2413  * \author [DP]
2414  */
2415 bool CxImage::TextBlur(BYTE threshold, BYTE decay, BYTE max_depth, bool bBlurHorizontal, bool bBlurVertical, CxImage* iDst)
2416 {
2417         if (!pDib) return false;
2418
2419         RGBQUAD* pPalette=NULL;
2420         WORD bpp = GetBpp();
2421
2422         //the routine is optimized for RGB or GrayScale images
2423         if (!(head.biBitCount == 24 || IsGrayScale())){
2424                 pPalette = new RGBQUAD[head.biClrUsed];
2425                 memcpy(pPalette, GetPalette(),GetPaletteSize());
2426                 if (!IncreaseBpp(24))
2427                         return false;
2428         }
2429
2430         CxImage tmp(*this);
2431         if (!tmp.IsValid()){
2432                 strcpy(info.szLastError,tmp.GetLastError());
2433                 return false;
2434         }
2435
2436         if (bBlurHorizontal)
2437                 blur_text(threshold, decay, max_depth, this, &tmp, head.biBitCount>>3);
2438
2439         if (bBlurVertical){
2440                 CxImage src2(*this);
2441                 src2.RotateLeft();
2442                 tmp.RotateLeft();
2443                 blur_text(threshold, decay, max_depth, &src2, &tmp, head.biBitCount>>3);
2444                 tmp.RotateRight();
2445         }
2446
2447 #if CXIMAGE_SUPPORT_SELECTION
2448         //restore the non selected region
2449         if (pSelection){
2450                 for(long y=0; y<head.biHeight; y++){
2451                         for(long x=0; x<head.biWidth; x++){
2452                                 if (!BlindSelectionIsInside(x,y)){
2453                                         tmp.BlindSetPixelColor(x,y,BlindGetPixelColor(x,y));
2454                                 }
2455                         }
2456                 }
2457         }
2458 #endif //CXIMAGE_SUPPORT_SELECTION
2459
2460         //if necessary, restore the original BPP and palette
2461         if (pPalette){
2462                 tmp.DecreaseBpp(bpp, true, pPalette);
2463                 delete [] pPalette;
2464         }
2465
2466         if (iDst) iDst->Transfer(tmp);
2467         else Transfer(tmp);
2468
2469         return true;
2470 }
2471 ////////////////////////////////////////////////////////////////////////////////
2472 /**
2473  * \author [nipper]; changes [DP]
2474  */
2475 bool CxImage::GaussianBlur(float radius /*= 1.0f*/, CxImage* iDst /*= 0*/)
2476 {
2477         if (!pDib) return false;
2478
2479         RGBQUAD* pPalette=NULL;
2480         WORD bpp = GetBpp();
2481
2482         //the routine is optimized for RGB or GrayScale images
2483         if (!(head.biBitCount == 24 || IsGrayScale())){
2484                 pPalette = new RGBQUAD[head.biClrUsed];
2485                 memcpy(pPalette, GetPalette(),GetPaletteSize());
2486                 if (!IncreaseBpp(24))
2487                         return false;
2488         }
2489
2490         CxImage tmp_x(*this, false, true, true);
2491         if (!tmp_x.IsValid()){
2492                 strcpy(info.szLastError,tmp_x.GetLastError());
2493                 return false;
2494         }
2495
2496         // generate convolution matrix and make sure it's smaller than each dimension
2497         float *cmatrix = NULL;
2498         int cmatrix_length = gen_convolve_matrix(radius, &cmatrix);
2499         // generate lookup table
2500         float *ctable = gen_lookup_table(cmatrix, cmatrix_length);
2501
2502         long x,y;
2503         int bypp = head.biBitCount>>3;
2504
2505         CImageIterator itSrc(this);
2506         CImageIterator itTmp(&tmp_x);
2507
2508         double dbScaler = 50.0f/head.biHeight;
2509
2510         // blur the rows
2511     for (y=0;y<head.biHeight;y++)
2512         {
2513                 if (info.nEscape) break;
2514                 info.nProgress = (long)(y*dbScaler);
2515
2516                 blur_line(ctable, cmatrix, cmatrix_length, itSrc.GetRow(y), itTmp.GetRow(y), head.biWidth, bypp);
2517         }
2518
2519         CxImage tmp_y(tmp_x, false, true, true);
2520         if (!tmp_y.IsValid()){
2521                 strcpy(info.szLastError,tmp_y.GetLastError());
2522                 return false;
2523         }
2524
2525         CImageIterator itDst(&tmp_y);
2526
2527         // blur the cols
2528     BYTE* cur_col = (BYTE*)malloc(bypp*head.biHeight);
2529     BYTE* dest_col = (BYTE*)malloc(bypp*head.biHeight);
2530
2531         dbScaler = 50.0f/head.biWidth;
2532
2533         for (x=0;x<head.biWidth;x++)
2534         {
2535                 if (info.nEscape) break;
2536                 info.nProgress = (long)(50.0f+x*dbScaler);
2537
2538                 itTmp.GetCol(cur_col, x);
2539                 itDst.GetCol(dest_col, x);
2540                 blur_line(ctable, cmatrix, cmatrix_length, cur_col, dest_col, head.biHeight, bypp);
2541                 itDst.SetCol(dest_col, x);
2542         }
2543
2544         free(cur_col);
2545         free(dest_col);
2546
2547         delete [] cmatrix;
2548         delete [] ctable;
2549
2550 #if CXIMAGE_SUPPORT_SELECTION
2551         //restore the non selected region
2552         if (pSelection){
2553                 for(y=0; y<head.biHeight; y++){
2554                         for(x=0; x<head.biWidth; x++){
2555                                 if (!BlindSelectionIsInside(x,y)){
2556                                         tmp_y.BlindSetPixelColor(x,y,BlindGetPixelColor(x,y));
2557                                 }
2558                         }
2559                 }
2560         }
2561 #endif //CXIMAGE_SUPPORT_SELECTION
2562
2563         //if necessary, restore the original BPP and palette
2564         if (pPalette){
2565                 tmp_y.DecreaseBpp(bpp, false, pPalette);
2566                 if (iDst) DecreaseBpp(bpp, false, pPalette);
2567                 delete [] pPalette;
2568         }
2569
2570         if (iDst) iDst->Transfer(tmp_y);
2571         else Transfer(tmp_y);
2572
2573         return true;
2574 }
2575 ////////////////////////////////////////////////////////////////////////////////
2576 /**
2577  * \author [DP],[nipper]
2578  */
2579 bool CxImage::SelectiveBlur(float radius, BYTE threshold, CxImage* iDst)
2580 {
2581         if (!pDib) return false;
2582
2583         RGBQUAD* pPalette=NULL;
2584         WORD bpp = GetBpp();
2585
2586         CxImage Tmp(*this, true, true, true);
2587         if (!Tmp.IsValid()){
2588                 strcpy(info.szLastError,Tmp.GetLastError());
2589                 return false;
2590         }
2591
2592         //the routine is optimized for RGB or GrayScale images
2593         if (!(head.biBitCount == 24 || IsGrayScale())){
2594                 pPalette = new RGBQUAD[head.biClrUsed];
2595                 memcpy(pPalette, GetPalette(),GetPaletteSize());
2596                 if (!Tmp.IncreaseBpp(24))
2597                         return false;
2598         }
2599
2600         CxImage Dst(Tmp, true, true, true);
2601         if (!Dst.IsValid()){
2602                 strcpy(info.szLastError,Dst.GetLastError());
2603                 return false;
2604         }
2605
2606         //build the difference mask
2607         BYTE thresh_dw = (BYTE)max( 0 ,(int)(128 - threshold));
2608         BYTE thresh_up = (BYTE)min(255,(int)(128 + threshold));
2609         long kernel[]={-100,-100,-100,-100,801,-100,-100,-100,-100};
2610         if (!Tmp.Filter(kernel,3,800,128)){
2611                 strcpy(info.szLastError,Tmp.GetLastError());
2612                 return false;
2613         }
2614
2615         //if the image has no selection, build a selection for the whole image
2616         if (!Tmp.SelectionIsValid()){
2617                 Tmp.SelectionCreate();
2618                 Tmp.SelectionClear(255);
2619         }
2620
2621         long xmin,xmax,ymin,ymax;
2622         xmin = Tmp.info.rSelectionBox.left;
2623         xmax = Tmp.info.rSelectionBox.right;
2624         ymin = Tmp.info.rSelectionBox.bottom;
2625         ymax = Tmp.info.rSelectionBox.top;
2626
2627         //modify the selection where the difference mask is over the threshold
2628         for(long y=ymin; y<ymax; y++){
2629                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
2630                 if (info.nEscape) break;
2631                 for(long x=xmin; x<xmax; x++){
2632                         if(Tmp.BlindSelectionIsInside(x,y)){
2633                                 RGBQUAD c = Tmp.BlindGetPixelColor(x,y);
2634                                 if ((c.rgbRed   < thresh_dw || c.rgbRed   > thresh_up) ||
2635                                         (c.rgbGreen < thresh_dw || c.rgbGreen > thresh_up) ||
2636                                         (c.rgbBlue  < thresh_dw || c.rgbBlue  > thresh_up))
2637                                 {
2638                                         Tmp.SelectionSet(x,y,0);
2639                                 }
2640                         }
2641                 }
2642         }
2643
2644         //blur the image (only in the selected pixels)
2645         Dst.SelectionCopy(Tmp);
2646         if (!Dst.GaussianBlur(radius)){
2647                 strcpy(info.szLastError,Dst.GetLastError());
2648                 return false;
2649         }
2650
2651         //restore the original selection
2652         Dst.SelectionCopy(*this);
2653
2654         //if necessary, restore the original BPP and palette
2655         if (pPalette){
2656                 Dst.DecreaseBpp(bpp, false, pPalette);
2657                 delete [] pPalette;
2658         }
2659
2660         if (iDst) iDst->Transfer(Dst);
2661         else Transfer(Dst);
2662
2663         return true;
2664 }
2665 ////////////////////////////////////////////////////////////////////////////////
2666 /**
2667  * sharpen the image by subtracting a blurred copy from the original image.
2668  * \param radius: width in pixels of the blurring effect. Range: >0; default = 5.
2669  * \param amount: strength of the filter. Range: 0.0 (none) to 1.0 (max); default = 0.5
2670  * \param threshold: difference, between blurred and original pixel, to trigger the filter
2671  *                   Range: 0 (always triggered) to 255 (never triggered); default = 0.
2672  * \return true if everything is ok
2673  * \author [nipper]; changes [DP]
2674  */
2675 bool CxImage::UnsharpMask(float radius /*= 5.0*/, float amount /*= 0.5*/, int threshold /*= 0*/)
2676 {
2677         if (!pDib) return false;
2678
2679         RGBQUAD* pPalette=NULL;
2680         WORD bpp = GetBpp();
2681
2682         //the routine is optimized for RGB or GrayScale images
2683         if (!(head.biBitCount == 24 || IsGrayScale())){
2684                 pPalette = new RGBQUAD[head.biClrUsed];
2685                 memcpy(pPalette, GetPalette(),GetPaletteSize());
2686                 if (!IncreaseBpp(24))
2687                         return false;
2688         }
2689
2690         CxImage iDst;
2691         if (!GaussianBlur(radius,&iDst))
2692                 return false;
2693
2694         CImageIterator itSrc(this);
2695         CImageIterator itDst(&iDst);
2696
2697         long xmin,xmax,ymin,ymax;
2698         if (pSelection){
2699                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2700                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2701         } else {
2702                 xmin = ymin = 0;
2703                 xmax = head.biWidth; ymax=head.biHeight;
2704         }
2705
2706         if (xmin==xmax || ymin==ymax)
2707                 return false;
2708
2709         double dbScaler = 100.0/(ymax-ymin);
2710         int bypp = head.biBitCount>>3;
2711
2712         // merge the source and destination (which currently contains
2713         // the blurred version) images
2714     for (long y=ymin; y<ymax; y++)
2715         {
2716                 if (info.nEscape) break;
2717                 info.nProgress = (long)((y-ymin)*dbScaler);
2718
2719                 // get source row
2720                 BYTE* cur_row = itSrc.GetRow(y);
2721                 // get dest row
2722                 BYTE* dest_row = itDst.GetRow(y);
2723                 // combine the two
2724                 for (long x=xmin; x<xmax; x++) {
2725 #if CXIMAGE_SUPPORT_SELECTION
2726                         if (BlindSelectionIsInside(x,y))
2727 #endif //CXIMAGE_SUPPORT_SELECTION
2728                         {
2729                                 for (long b=0, z=x*bypp; b<bypp; b++, z++){
2730                                         int diff = cur_row[z] - dest_row[z];
2731
2732                                         // do tresholding
2733                                         if (abs(diff) < threshold){
2734                                                 dest_row[z] = cur_row[z];
2735                                         } else {
2736                                                 dest_row[z] = (BYTE)min(255, max(0,(int)(cur_row[z] + amount * diff)));
2737                                         }
2738                                 }
2739                         }
2740                 }
2741         }
2742
2743         //if necessary, restore the original BPP and palette
2744         if (pPalette){
2745                 iDst.DecreaseBpp(bpp, false, pPalette);
2746                 delete [] pPalette;
2747         }
2748
2749         Transfer(iDst);
2750
2751         return true;
2752 }
2753 ////////////////////////////////////////////////////////////////////////////////
2754 /**
2755  * Apply a look up table to the image.
2756  * \param pLut: BYTE[256] look up table
2757  * \return true if everything is ok
2758  */
2759 bool CxImage::Lut(BYTE* pLut)
2760 {
2761         if (!pDib || !pLut) return false;
2762         RGBQUAD color;
2763
2764         double dbScaler;
2765         if (head.biClrUsed==0){
2766
2767                 long xmin,xmax,ymin,ymax;
2768                 if (pSelection){
2769                         xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2770                         ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2771                 } else {
2772                         // faster loop for full image
2773                         BYTE *iSrc=info.pImage;
2774                         for(unsigned long i=0; i < head.biSizeImage ; i++){
2775                                 *iSrc++ = pLut[*iSrc];
2776                         }
2777                         return true;
2778                 }
2779
2780                 if (xmin==xmax || ymin==ymax)
2781                         return false;
2782
2783                 dbScaler = 100.0/(ymax-ymin);
2784
2785                 for(long y=ymin; y<ymax; y++){
2786                         info.nProgress = (long)((y-ymin)*dbScaler); //<Anatoly Ivasyuk>
2787                         for(long x=xmin; x<xmax; x++){
2788 #if CXIMAGE_SUPPORT_SELECTION
2789                                 if (BlindSelectionIsInside(x,y))
2790 #endif //CXIMAGE_SUPPORT_SELECTION
2791                                 {
2792                                         color = BlindGetPixelColor(x,y);
2793                                         color.rgbRed = pLut[color.rgbRed];
2794                                         color.rgbGreen = pLut[color.rgbGreen];
2795                                         color.rgbBlue = pLut[color.rgbBlue];
2796                                         BlindSetPixelColor(x,y,color);
2797                                 }
2798                         }
2799                 }
2800 #if CXIMAGE_SUPPORT_SELECTION
2801         } else if (pSelection && (head.biBitCount==8) && IsGrayScale()){
2802                 long xmin,xmax,ymin,ymax;
2803                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2804                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2805
2806                 if (xmin==xmax || ymin==ymax)
2807                         return false;
2808
2809                 dbScaler = 100.0/(ymax-ymin);
2810                 for(long y=ymin; y<ymax; y++){
2811                         info.nProgress = (long)((y-ymin)*dbScaler);
2812                         for(long x=xmin; x<xmax; x++){
2813                                 if (BlindSelectionIsInside(x,y))
2814                                 {
2815                                         BlindSetPixelIndex(x,y,pLut[BlindGetPixelIndex(x,y)]);
2816                                 }
2817                         }
2818                 }
2819 #endif //CXIMAGE_SUPPORT_SELECTION
2820         } else {
2821                 bool bIsGrayScale = IsGrayScale();
2822                 for(DWORD j=0; j<head.biClrUsed; j++){
2823                         color = GetPaletteColor((BYTE)j);
2824                         color.rgbRed = pLut[color.rgbRed];
2825                         color.rgbGreen = pLut[color.rgbGreen];
2826                         color.rgbBlue = pLut[color.rgbBlue];
2827                         SetPaletteColor((BYTE)j,color);
2828                 }
2829                 if (bIsGrayScale) GrayScale();
2830         }
2831         return true;
2832
2833 }
2834 ////////////////////////////////////////////////////////////////////////////////
2835 /**
2836  * Apply an indipendent look up table for each channel
2837  * \param pLutR, pLutG, pLutB, pLutA: BYTE[256] look up tables
2838  * \return true if everything is ok
2839  */
2840 bool CxImage::Lut(BYTE* pLutR, BYTE* pLutG, BYTE* pLutB, BYTE* pLutA)
2841 {
2842         if (!pDib || !pLutR || !pLutG || !pLutB) return false;
2843         RGBQUAD color;
2844
2845         double dbScaler;
2846         if (head.biClrUsed==0){
2847
2848                 long xmin,xmax,ymin,ymax;
2849                 if (pSelection){
2850                         xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2851                         ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2852                 } else {
2853                         xmin = ymin = 0;
2854                         xmax = head.biWidth; ymax=head.biHeight;
2855                 }
2856
2857                 if (xmin==xmax || ymin==ymax)
2858                         return false;
2859
2860                 dbScaler = 100.0/(ymax-ymin);
2861
2862                 for(long y=ymin; y<ymax; y++){
2863                         info.nProgress = (long)((y-ymin)*dbScaler);
2864                         for(long x=xmin; x<xmax; x++){
2865 #if CXIMAGE_SUPPORT_SELECTION
2866                                 if (BlindSelectionIsInside(x,y))
2867 #endif //CXIMAGE_SUPPORT_SELECTION
2868                                 {
2869                                         color = BlindGetPixelColor(x,y);
2870                                         color.rgbRed =   pLutR[color.rgbRed];
2871                                         color.rgbGreen = pLutG[color.rgbGreen];
2872                                         color.rgbBlue =  pLutB[color.rgbBlue];
2873                                         if (pLutA) color.rgbReserved=pLutA[color.rgbReserved];
2874                                         BlindSetPixelColor(x,y,color,true);
2875                                 }
2876                         }
2877                 }
2878         } else {
2879                 bool bIsGrayScale = IsGrayScale();
2880                 for(DWORD j=0; j<head.biClrUsed; j++){
2881                         color = GetPaletteColor((BYTE)j);
2882                         color.rgbRed =   pLutR[color.rgbRed];
2883                         color.rgbGreen = pLutG[color.rgbGreen];
2884                         color.rgbBlue =  pLutB[color.rgbBlue];
2885                         SetPaletteColor((BYTE)j,color);
2886                 }
2887                 if (bIsGrayScale) GrayScale();
2888         }
2889
2890         return true;
2891
2892 }
2893 ////////////////////////////////////////////////////////////////////////////////
2894 /**
2895  * Use the RedEyeRemove function to remove the red-eye effect that frequently
2896  * occurs in photographs of humans and animals. You must select the region
2897  * where the function will filter the red channel.
2898  * \param strength: range from 0.0f (no effect) to 1.0f (full effect). Default = 0.8
2899  * \return true if everything is ok
2900  */
2901 bool CxImage::RedEyeRemove(float strength)
2902 {
2903         if (!pDib) return false;
2904         RGBQUAD color;
2905
2906         long xmin,xmax,ymin,ymax;
2907         if (pSelection){
2908                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2909                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2910         } else {
2911                 xmin = ymin = 0;
2912                 xmax = head.biWidth; ymax=head.biHeight;
2913         }
2914
2915         if (xmin==xmax || ymin==ymax)
2916                 return false;
2917
2918         if (strength<0.0f) strength = 0.0f;
2919         if (strength>1.0f) strength = 1.0f;
2920
2921         for(long y=ymin; y<ymax; y++){
2922                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
2923                 if (info.nEscape) break;
2924                 for(long x=xmin; x<xmax; x++){
2925 #if CXIMAGE_SUPPORT_SELECTION
2926                         if (BlindSelectionIsInside(x,y))
2927 #endif //CXIMAGE_SUPPORT_SELECTION
2928                         {
2929                                 float a = 1.0f-5.0f*((float)((x-0.5f*(xmax+xmin))*(x-0.5f*(xmax+xmin))+(y-0.5f*(ymax+ymin))*(y-0.5f*(ymax+ymin))))/((float)((xmax-xmin)*(ymax-ymin)));
2930                                 if (a<0) a=0;
2931                                 color = BlindGetPixelColor(x,y);
2932                                 color.rgbRed = (BYTE)(a*min(color.rgbGreen,color.rgbBlue)+(1.0f-a)*color.rgbRed);
2933                                 BlindSetPixelColor(x,y,color);
2934                         }
2935                 }
2936         }
2937         return true;
2938 }
2939 ////////////////////////////////////////////////////////////////////////////////
2940 /**
2941  * Changes the saturation of the image.
2942  * \param saturation: can be from -100 to 100, positive values increase the saturation.
2943  * \param colorspace: can be 1 (HSL) or 2 (YUV).
2944  * \return true if everything is ok
2945  */
2946 bool CxImage::Saturate(const long saturation, const long colorspace)
2947 {
2948         if (!pDib)
2949                 return false;
2950
2951         long xmin,xmax,ymin,ymax;
2952         if (pSelection){
2953                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
2954                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
2955         } else {
2956                 xmin = ymin = 0;
2957                 xmax = head.biWidth; ymax=head.biHeight;
2958         }
2959
2960         if (xmin==xmax || ymin==ymax)
2961                 return false;
2962
2963         BYTE cTable[256];
2964
2965         switch(colorspace)
2966         {
2967         case 1:
2968                 {
2969                         for (int i=0;i<256;i++) {
2970                                 cTable[i] = (BYTE)max(0,min(255,(int)(i + saturation)));
2971                         }
2972                         for(long y=ymin; y<ymax; y++){
2973                                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
2974                                 if (info.nEscape) break;
2975                                 for(long x=xmin; x<xmax; x++){
2976 #if CXIMAGE_SUPPORT_SELECTION
2977                                         if (BlindSelectionIsInside(x,y))
2978 #endif //CXIMAGE_SUPPORT_SELECTION
2979                                         {
2980                                                 RGBQUAD c = RGBtoHSL(BlindGetPixelColor(x,y));
2981                                                 c.rgbGreen  = cTable[c.rgbGreen];
2982                                                 c = HSLtoRGB(c);
2983                                                 BlindSetPixelColor(x,y,c);
2984                                         }
2985                                 }
2986                         }
2987                 }
2988                 break;
2989         case 2:
2990                 {
2991                         for (int i=0;i<256;i++) {
2992                                 cTable[i] = (BYTE)max(0,min(255,(int)((i-128)*(100 + saturation)/100.0f + 128.5f)));
2993                         }
2994                         for(long y=ymin; y<ymax; y++){
2995                                 info.nProgress = (long)(100*(y-ymin)/(ymax-ymin));
2996                                 if (info.nEscape) break;
2997                                 for(long x=xmin; x<xmax; x++){
2998 #if CXIMAGE_SUPPORT_SELECTION
2999                                         if (BlindSelectionIsInside(x,y))
3000 #endif //CXIMAGE_SUPPORT_SELECTION
3001                                         {
3002                                                 RGBQUAD c = RGBtoYUV(BlindGetPixelColor(x,y));
3003                                                 c.rgbGreen  = cTable[c.rgbGreen];
3004                                                 c.rgbBlue = cTable[c.rgbBlue];
3005                                                 c = YUVtoRGB(c);
3006                                                 BlindSetPixelColor(x,y,c);
3007                                         }
3008                                 }
3009                         }
3010                 }
3011                 break;
3012         default:
3013                 strcpy(info.szLastError,"Saturate: wrong colorspace");
3014                 return false;
3015         }
3016         return true;
3017 }
3018
3019 ////////////////////////////////////////////////////////////////////////////////
3020 /**
3021  * Solarize: convert all colors above a given lightness level into their negative
3022  * \param  level : lightness threshold. Range = 0 to 255; default = 128.
3023  * \param  bLinkedChannels: true = compare with luminance, preserve colors (default)
3024  *                         false = compare with independent R,G,B levels
3025  * \return true if everything is ok
3026  * \author [Priyank Bolia] (priyank_bolia(at)yahoo(dot)com); changes [DP]
3027  */
3028 bool CxImage::Solarize(BYTE level, bool bLinkedChannels)
3029 {
3030         if (!pDib) return false;
3031
3032         long xmin,xmax,ymin,ymax;
3033         if (pSelection){
3034                 xmin = info.rSelectionBox.left; xmax = info.rSelectionBox.right;
3035                 ymin = info.rSelectionBox.bottom; ymax = info.rSelectionBox.top;
3036         } else {
3037                 xmin = ymin = 0;
3038                 xmax = head.biWidth; ymax=head.biHeight;
3039         }
3040
3041         if (head.biBitCount<=8){
3042                 if (IsGrayScale()){ //GRAYSCALE, selection
3043                         for(long y=ymin; y<ymax; y++){
3044                                 for(long x=xmin; x<xmax; x++){
3045 #if CXIMAGE_SUPPORT_SELECTION
3046                                         if (BlindSelectionIsInside(x,y))
3047 #endif //CXIMAGE_SUPPORT_SELECTION
3048                                         {
3049                                                 BYTE index = BlindGetPixelIndex(x,y);
3050                                                 RGBQUAD color = GetPaletteColor(index);
3051                                                 if ((BYTE)RGB2GRAY(color.rgbRed,color.rgbGreen,color.rgbBlue)>level){
3052                                                         BlindSetPixelIndex(x,y,255-index);
3053                                                 }
3054                                         }
3055                                 }
3056                         }
3057                 } else { //PALETTE, full image
3058                         RGBQUAD* ppal=GetPalette();
3059                         for(DWORD i=0;i<head.biClrUsed;i++){
3060                                 RGBQUAD color = GetPaletteColor((BYTE)i);
3061                                 if (bLinkedChannels){
3062                                         if ((BYTE)RGB2GRAY(color.rgbRed,color.rgbGreen,color.rgbBlue)>level){
3063                                                 ppal[i].rgbBlue =(BYTE)(255-ppal[i].rgbBlue);
3064                                                 ppal[i].rgbGreen =(BYTE)(255-ppal[i].rgbGreen);
3065                                                 ppal[i].rgbRed =(BYTE)(255-ppal[i].rgbRed);
3066                                         }
3067                                 } else {
3068                                         if (color.rgbBlue>level)        ppal[i].rgbBlue =(BYTE)(255-ppal[i].rgbBlue);
3069                                         if (color.rgbGreen>level)       ppal[i].rgbGreen =(BYTE)(255-ppal[i].rgbGreen);
3070                                         if (color.rgbRed>level)         ppal[i].rgbRed =(BYTE)(255-ppal[i].rgbRed);
3071                                 }
3072                         }
3073                 }
3074         } else { //RGB, selection
3075                 for(long y=ymin; y<ymax; y++){
3076                         for(long x=xmin; x<xmax; x++){
3077 #if CXIMAGE_SUPPORT_SELECTION
3078                                 if (BlindSelectionIsInside(x,y))
3079 #endif //CXIMAGE_SUPPORT_SELECTION
3080                                 {
3081                                         RGBQUAD color = BlindGetPixelColor(x,y);
3082                                         if (bLinkedChannels){
3083                                                 if ((BYTE)RGB2GRAY(color.rgbRed,color.rgbGreen,color.rgbBlue)>level){
3084                                                         color.rgbRed = (BYTE)(255-color.rgbRed);
3085                                                         color.rgbGreen = (BYTE)(255-color.rgbGreen);
3086                                                         color.rgbBlue = (BYTE)(255-color.rgbBlue);
3087                                                 }
3088                                         } else {
3089                                                 if (color.rgbBlue>level)        color.rgbBlue =(BYTE)(255-color.rgbBlue);
3090                                                 if (color.rgbGreen>level)       color.rgbGreen =(BYTE)(255-color.rgbGreen);
3091                                                 if (color.rgbRed>level)         color.rgbRed =(BYTE)(255-color.rgbRed);
3092                                         }
3093                                         BlindSetPixelColor(x,y,color);
3094                                 }
3095                         }
3096                 }
3097         }
3098
3099         //invert transparent color only in case of full image processing
3100         if (pSelection==0 || (!IsGrayScale() && IsIndexed())){
3101                 if (bLinkedChannels){
3102                         if ((BYTE)RGB2GRAY(info.nBkgndColor.rgbRed,info.nBkgndColor.rgbGreen,info.nBkgndColor.rgbBlue)>level){
3103                                 info.nBkgndColor.rgbBlue = (BYTE)(255-info.nBkgndColor.rgbBlue);
3104                                 info.nBkgndColor.rgbGreen = (BYTE)(255-info.nBkgndColor.rgbGreen);
3105                                 info.nBkgndColor.rgbRed = (BYTE)(255-info.nBkgndColor.rgbRed);
3106                         }
3107                 } else {
3108                         if (info.nBkgndColor.rgbBlue>level)      info.nBkgndColor.rgbBlue = (BYTE)(255-info.nBkgndColor.rgbBlue);
3109                         if (info.nBkgndColor.rgbGreen>level) info.nBkgndColor.rgbGreen = (BYTE)(255-info.nBkgndColor.rgbGreen);
3110                         if (info.nBkgndColor.rgbRed>level)       info.nBkgndColor.rgbRed = (BYTE)(255-info.nBkgndColor.rgbRed);
3111                 }
3112         }
3113
3114         return true;
3115 }
3116
3117 ////////////////////////////////////////////////////////////////////////////////
3118 /**
3119  * Converts the RGB triplets to and from different colorspace
3120  * \param dstColorSpace: destination colorspace; 0 = RGB, 1 = HSL, 2 = YUV, 3 = YIQ, 4 = XYZ
3121  * \param srcColorSpace: source colorspace; 0 = RGB, 1 = HSL, 2 = YUV, 3 = YIQ, 4 = XYZ
3122  * \return true if everything is ok
3123  */
3124 bool CxImage::ConvertColorSpace(const long dstColorSpace, const long srcColorSpace)
3125 {
3126         if (!pDib)
3127                 return false;
3128
3129         if (dstColorSpace == srcColorSpace)
3130                 return true;
3131
3132         long w = GetWidth();
3133         long h = GetHeight();
3134
3135         for (long y=0;y<h;y++){
3136                 info.nProgress = (long)(100*y/h);
3137                 if (info.nEscape) break;
3138                 for (long x=0;x<w;x++){
3139                         RGBQUAD c = BlindGetPixelColor(x,y);
3140                         switch (srcColorSpace){
3141                         case 0:
3142                                 break;
3143                         case 1:
3144                                 c = HSLtoRGB(c);
3145                                 break;
3146                         case 2:
3147                                 c = YUVtoRGB(c);
3148                                 break;
3149                         case 3:
3150                                 c = YIQtoRGB(c);
3151                                 break;
3152                         case 4:
3153                                 c = XYZtoRGB(c);
3154                                 break;
3155                         default:
3156                                 strcpy(info.szLastError,"ConvertColorSpace: unknown source colorspace");
3157                                 return false;
3158                         }
3159                         switch (dstColorSpace){
3160                         case 0:
3161                                 break;
3162                         case 1:
3163                                 c = RGBtoHSL(c);
3164                                 break;
3165                         case 2:
3166                                 c = RGBtoYUV(c);
3167                                 break;
3168                         case 3:
3169                                 c = RGBtoYIQ(c);
3170                                 break;
3171                         case 4:
3172                                 c = RGBtoXYZ(c);
3173                                 break;
3174                         default:
3175                                 strcpy(info.szLastError,"ConvertColorSpace: unknown destination colorspace");
3176                                 return false;
3177                         }
3178                         BlindSetPixelColor(x,y,c);
3179                 }
3180         }
3181         return true;
3182 }
3183 ////////////////////////////////////////////////////////////////////////////////
3184 /**
3185  * Finds the optimal (global or local) treshold for image binarization
3186  * \param method: 0 = average all methods (default); 1 = Otsu; 2 = Kittler & Illingworth; 3 = max entropy; 4 = potential difference;
3187  * \param pBox: region from where the threshold is computed; 0 = full image (default).
3188  * \param pContrastMask: limit the computation only in regions with contrasted (!=0) pixels; default = 0.
3189  * the pContrastMask image must be grayscale with same with and height of the current image,
3190  * can be obtained from the current image with a filter:
3191  *   CxImage iContrastMask(*image,true,false,false);
3192  *   iContrastMask.GrayScale();
3193  *   long edge[]={-1,-1,-1,-1,8,-1,-1,-1,-1};
3194  *   iContrastMask.Filter(edge,3,1,0);
3195  *   long blur[]={1,1,1,1,1,1,1,1,1};
3196  *   iContrastMask.Filter(blur,3,9,0);
3197  * \return optimal threshold; -1 = error.
3198  * \sa AdaptiveThreshold
3199  */
3200 int  CxImage::OptimalThreshold(long method, RECT * pBox, CxImage* pContrastMask)
3201 {
3202         if (!pDib)
3203                 return false;
3204
3205         if (head.biBitCount!=8){
3206                 strcpy(info.szLastError,"OptimalThreshold works only on 8 bit images");
3207                 return -1;
3208         }
3209
3210         if (pContrastMask){
3211                 if (!pContrastMask->IsValid() ||
3212                         !pContrastMask->IsGrayScale() ||
3213                         pContrastMask->GetWidth() != GetWidth() ||
3214                         pContrastMask->GetHeight() != GetHeight()){
3215                         strcpy(info.szLastError,"OptimalThreshold invalid ContrastMask");
3216                         return -1;
3217                 }
3218         }
3219
3220         long xmin,xmax,ymin,ymax;
3221         if (pBox){
3222                 xmin = max(pBox->left,0);
3223                 xmax = min(pBox->right,head.biWidth);
3224                 ymin = max(pBox->bottom,0);
3225                 ymax = min(pBox->top,head.biHeight);
3226         } else {
3227                 xmin = ymin = 0;
3228                 xmax = head.biWidth; ymax=head.biHeight;
3229         }
3230
3231         if (xmin>=xmax || ymin>=ymax)
3232                 return -1;
3233
3234         double p[256];
3235         memset(p,  0, 256*sizeof(double));
3236         //build histogram
3237         for (long y = ymin; y<ymax; y++){
3238                 BYTE* pGray = GetBits(y) + xmin;
3239                 BYTE* pContr = 0;
3240                 if (pContrastMask) pContr = pContrastMask->GetBits(y) + xmin;
3241                 for (long x = xmin; x<xmax; x++){
3242                         BYTE n = *pGray++;
3243                         if (pContr){
3244                                 if (*pContr) p[n]++;
3245                                 pContr++;
3246                         } else {
3247                                 p[n]++;
3248                         }
3249                 }
3250         }
3251
3252         //find histogram limits
3253         int gray_min = 0;
3254         while (gray_min<255 && p[gray_min]==0) gray_min++;
3255         int gray_max = 255;
3256         while (gray_max>0 && p[gray_max]==0) gray_max--;
3257         if (gray_min > gray_max)
3258                 return -1;
3259         if (gray_min == gray_max){
3260                 if (gray_min == 0)
3261                         return 0;
3262                 else
3263                         return gray_max-1;
3264         }
3265
3266         //compute total moments 0th,1st,2nd order
3267         int i,k;
3268         double w_tot = 0;
3269         double m_tot = 0;
3270         double q_tot = 0;
3271         for (i = gray_min; i <= gray_max; i++){
3272                 w_tot += p[i];
3273                 m_tot += i*p[i];
3274                 q_tot += i*i*p[i];
3275         }
3276
3277         double L, L1max, L2max, L3max, L4max; //objective functions
3278         int th1,th2,th3,th4; //optimal thresholds
3279         L1max = L2max = L3max = L4max = 0;
3280         th1 = th2 = th3 = th4 = -1;
3281
3282         double w1, w2, m1, m2, q1, q2, s1, s2;
3283         w1 = m1 = q1 = 0;
3284         for (i = gray_min; i < gray_max; i++){
3285                 w1 += p[i];
3286                 w2 = w_tot - w1;
3287                 m1 += i*p[i];
3288                 m2 = m_tot - m1;
3289                 q1 += i*i*p[i];
3290                 q2 = q_tot - q1;
3291                 s1 = q1/w1-m1*m1/w1/w1; //s1 = q1/w1-pow(m1/w1,2);
3292                 s2 = q2/w2-m2*m2/w2/w2; //s2 = q2/w2-pow(m2/w2,2);
3293
3294                 //Otsu
3295                 L = -(s1*w1 + s2*w2); //implemented as definition
3296                 //L = w1 * w2 * (m2/w2 - m1/w1)*(m2/w2 - m1/w1); //implementation that doesn't need s1 & s2
3297                 if (L1max < L || th1<0){
3298                         L1max = L;
3299                         th1 = i;
3300                 }
3301
3302                 //Kittler and Illingworth
3303                 if (s1>0 && s2>0){
3304                         L = w1*log(w1/sqrt(s1))+w2*log(w2/sqrt(s2));
3305                         //L = w1*log(w1*w1/s1)+w2*log(w2*w2/s2);
3306                         if (L2max < L || th2<0){
3307                                 L2max = L;
3308                                 th2 = i;
3309                         }
3310                 }
3311
3312                 //max entropy
3313                 L = 0;
3314                 for (k=gray_min;k<=i;k++) if (p[k] > 0) L -= p[k]*log(p[k]/w1)/w1;
3315                 for (k;k<=gray_max;k++) if (p[k] > 0)   L -= p[k]*log(p[k]/w2)/w2;
3316                 if (L3max < L || th3<0){
3317                         L3max = L;
3318                         th3 = i;
3319                 }
3320
3321                 //potential difference (based on Electrostatic Binarization method by J. Acharya & G. Sreechakra)
3322                 // L=-fabs(vdiff/vsum); è molto selettivo, sembra che L=-fabs(vdiff) o L=-(vsum)
3323                 // abbiano lo stesso valore di soglia... il che semplificherebbe molto la routine
3324                 double vdiff = 0;
3325                 for (k=gray_min;k<=i;k++)
3326                         vdiff += p[k]*(i-k)*(i-k);
3327                 double vsum = vdiff;
3328                 for (k;k<=gray_max;k++){
3329                         double dv = p[k]*(k-i)*(k-i);
3330                         vdiff -= dv;
3331                         vsum += dv;
3332                 }
3333                 if (vsum>0) L = -fabs(vdiff/vsum); else L = 0;
3334                 if (L4max < L || th4<0){
3335                         L4max = L;
3336                         th4 = i;
3337                 }
3338         }
3339
3340         int threshold;
3341         switch (method){
3342         case 1: //Otsu
3343                 threshold = th1;
3344                 break;
3345         case 2: //Kittler and Illingworth
3346                 threshold = th2;
3347                 break;
3348         case 3: //max entropy
3349                 threshold = th3;
3350                 break;
3351         case 4: //potential difference
3352                 threshold = th4;
3353                 break;
3354         default: //auto
3355                 {
3356                         int nt = 0;
3357                         threshold = 0;
3358                         if (th1>=0) { threshold += th1; nt++;}
3359                         if (th2>=0) { threshold += th2; nt++;}
3360                         if (th3>=0) { threshold += th3; nt++;}
3361                         if (th4>=0) { threshold += th4; nt++;}
3362                         if (nt)
3363                                 threshold /= nt;
3364                         else
3365                                 threshold = (gray_min+gray_max)/2;
3366
3367                         /*better(?) but really expensive alternative:
3368                         n = 0:255;
3369                         pth1 = c1(th1)/sqrt(2*pi*s1(th1))*exp(-((n - m1(th1)).^2)/2/s1(th1)) + c2(th1)/sqrt(2*pi*s2(th1))*exp(-((n - m2(th1)).^2)/2/s2(th1));
3370                         pth2 = c1(th2)/sqrt(2*pi*s1(th2))*exp(-((n - m1(th2)).^2)/2/s1(th2)) + c2(th2)/sqrt(2*pi*s2(th2))*exp(-((n - m2(th2)).^2)/2/s2(th2));
3371                         ...
3372                         mse_th1 = sum((p-pth1).^2);
3373                         mse_th2 = sum((p-pth2).^2);
3374                         ...
3375                         select th# that gives minimum mse_th#
3376                         */
3377
3378                 }
3379         }
3380
3381         if (threshold <= gray_min || threshold >= gray_max)
3382                 threshold = (gray_min+gray_max)/2;
3383
3384         return threshold;
3385 }
3386 ///////////////////////////////////////////////////////////////////////////////
3387 /**
3388  * Converts the image to B&W, using an optimal threshold mask
3389  * \param method: 0 = average all methods (default); 1 = Otsu; 2 = Kittler & Illingworth; 3 = max entropy; 4 = potential difference;
3390  * \param nBoxSize: the image is divided into "nBoxSize x nBoxSize" blocks, from where the threshold is computed; min = 8; default = 64.
3391  * \param pContrastMask: limit the computation only in regions with contrasted (!=0) pixels; default = 0.
3392  * \param nBias: global offset added to the threshold mask; default = 0.
3393  * \param fGlobalLocalBalance: balance between local and global threshold. default = 0.5
3394  * fGlobalLocalBalance can be from 0.0 (use only local threshold) to 1.0 (use only global threshold)
3395  * the pContrastMask image must be grayscale with same with and height of the current image,
3396  * \return true if everything is ok.
3397  * \sa OptimalThreshold
3398  */
3399 bool CxImage::AdaptiveThreshold(long method, long nBoxSize, CxImage* pContrastMask, long nBias, float fGlobalLocalBalance)
3400 {
3401         if (!pDib)
3402                 return false;
3403
3404         if (pContrastMask){
3405                 if (!pContrastMask->IsValid() ||
3406                         !pContrastMask->IsGrayScale() ||
3407                         pContrastMask->GetWidth() != GetWidth() ||
3408                         pContrastMask->GetHeight() != GetHeight()){
3409                         strcpy(info.szLastError,"AdaptiveThreshold invalid ContrastMask");
3410                         return false;
3411                 }
3412         }
3413
3414         if (nBoxSize<8) nBoxSize = 8;
3415         if (fGlobalLocalBalance<0.0f) fGlobalLocalBalance = 0.0f;
3416         if (fGlobalLocalBalance>1.0f) fGlobalLocalBalance = 1.0f;
3417
3418         long mw = (head.biWidth + nBoxSize - 1)/nBoxSize;
3419         long mh = (head.biHeight + nBoxSize - 1)/nBoxSize;
3420
3421         CxImage mask(mw,mh,8);
3422         if(!mask.GrayScale())
3423                 return false;
3424
3425         if(!GrayScale())
3426                 return false;
3427
3428         int globalthreshold = OptimalThreshold(method, 0, pContrastMask);
3429         if (globalthreshold <0)
3430                 return false;
3431
3432         for (long y=0; y<mh; y++){
3433                 for (long x=0; x<mw; x++){
3434                         info.nProgress = (long)(100*(x+y*mw)/(mw*mh));
3435                         if (info.nEscape) break;
3436                         RECT r;
3437                         r.left = x*nBoxSize;
3438                         r.right = r.left + nBoxSize;
3439                         r.bottom = y*nBoxSize;
3440                         r.top = r.bottom + nBoxSize;
3441                         int threshold = OptimalThreshold(method, &r, pContrastMask);
3442                         if (threshold <0) return false;
3443                         mask.SetPixelIndex(x,y,(BYTE)max(0,min(255,nBias+((1.0f-fGlobalLocalBalance)*threshold + fGlobalLocalBalance*globalthreshold))));
3444                 }
3445         }
3446
3447         mask.Resample(mw*nBoxSize,mh*nBoxSize,0);
3448         mask.Crop(0,head.biHeight,head.biWidth,0);
3449
3450         if(!Threshold(&mask))
3451                 return false;
3452
3453         return true;
3454 }
3455
3456 ////////////////////////////////////////////////////////////////////////////////
3457 #include <queue>
3458 ////////////////////////////////////////////////////////////////////////////////
3459 /**
3460  * Flood Fill
3461  * \param xStart, yStart: starting point
3462  * \param cFillColor: filling color
3463  * \param nTolerance: deviation from the starting point color
3464  * \param nOpacity: can be from 0 (transparent) to 255 (opaque, default)
3465  * \param bSelectFilledArea: if true, the pixels in the region are also set in the selection layer; default = false
3466  * \param nSelectionLevel: if bSelectFilledArea is true, the selected pixels are set to nSelectionLevel; default = 255
3467  * Note: nOpacity=0 && bSelectFilledArea=true act as a "magic wand"
3468  * \return true if everything is ok
3469  */
3470 bool CxImage::FloodFill(const long xStart, const long yStart, const RGBQUAD cFillColor, const BYTE nTolerance,
3471                                                 BYTE nOpacity, const bool bSelectFilledArea, const BYTE nSelectionLevel)
3472 {
3473         if (!pDib)
3474                 return false;
3475
3476         if (!IsInside(xStart,yStart))
3477                 return true;
3478
3479 #if CXIMAGE_SUPPORT_SELECTION
3480         if (!SelectionIsInside(xStart,yStart))
3481                 return true;
3482 #endif //CXIMAGE_SUPPORT_SELECTION
3483
3484         RGBQUAD* pPalette=NULL;
3485         WORD bpp = GetBpp();
3486         //nTolerance or nOpacity implemented only for grayscale or 24bpp images
3487         if ((nTolerance || nOpacity != 255) &&  !(head.biBitCount == 24 || IsGrayScale())){
3488                 pPalette = new RGBQUAD[head.biClrUsed];
3489                 memcpy(pPalette, GetPalette(),GetPaletteSize());
3490                 if (!IncreaseBpp(24))
3491                         return false;
3492         }
3493
3494         BYTE* pFillMask = (BYTE*)calloc(head.biWidth * head.biHeight,1);
3495         if (!pFillMask)
3496                 return false;
3497
3498 //------------------------------------- Begin of Flood Fill
3499         POINT offset[4] = {{-1,0},{0,-1},{1,0},{0,1}};
3500         std::queue<POINT> q;
3501         POINT point = {xStart,yStart};
3502         q.push(point);
3503
3504         if (IsIndexed()){ //--- Generic indexed image, no tolerance OR Grayscale image with tolerance
3505                 BYTE idxRef = GetPixelIndex(xStart,yStart);
3506                 BYTE idxFill = GetNearestIndex(cFillColor);
3507                 BYTE idxMin = (BYTE)min(255, max(0,(int)(idxRef - nTolerance)));
3508                 BYTE idxMax = (BYTE)min(255, max(0,(int)(idxRef + nTolerance)));
3509
3510                 while(!q.empty())
3511                 {
3512                         point = q.front();
3513                         q.pop();
3514
3515                         for (int z=0; z<4; z++){
3516                                 int x = point.x + offset[z].x;
3517                                 int y = point.y + offset[z].y;
3518                                 if(IsInside(x,y)){
3519 #if CXIMAGE_SUPPORT_SELECTION
3520                                   if (BlindSelectionIsInside(x,y))
3521 #endif //CXIMAGE_SUPPORT_SELECTION
3522                                   {
3523                                         BYTE idx = BlindGetPixelIndex(x, y);
3524                                         BYTE* pFill = pFillMask + x + y * head.biWidth;
3525                                         if (*pFill==0 && idxMin <= idx && idx <= idxMax )
3526                                         {
3527                                                 if (nOpacity>0){
3528                                                         if (nOpacity == 255)
3529                                                                 BlindSetPixelIndex(x, y, idxFill);
3530                                                         else
3531                                                                 BlindSetPixelIndex(x, y, (BYTE)((idxFill * nOpacity + idx * (255-nOpacity))>>8));
3532                                                 }
3533                                                 POINT pt = {x,y};
3534                                                 q.push(pt);
3535                                                 *pFill = 1;
3536                                         }
3537                                   }
3538                                 }
3539                         }
3540                 }
3541         } else { //--- RGB image
3542                 RGBQUAD cRef = GetPixelColor(xStart,yStart);
3543                 RGBQUAD cRefMin, cRefMax;
3544                 cRefMin.rgbRed   = (BYTE)min(255, max(0,(int)(cRef.rgbRed   - nTolerance)));
3545                 cRefMin.rgbGreen = (BYTE)min(255, max(0,(int)(cRef.rgbGreen - nTolerance)));
3546                 cRefMin.rgbBlue  = (BYTE)min(255, max(0,(int)(cRef.rgbBlue  - nTolerance)));
3547                 cRefMax.rgbRed   = (BYTE)min(255, max(0,(int)(cRef.rgbRed   + nTolerance)));
3548                 cRefMax.rgbGreen = (BYTE)min(255, max(0,(int)(cRef.rgbGreen + nTolerance)));
3549                 cRefMax.rgbBlue  = (BYTE)min(255, max(0,(int)(cRef.rgbBlue  + nTolerance)));
3550
3551                 while(!q.empty())
3552                 {
3553                         point = q.front();
3554                         q.pop();
3555
3556                         for (int z=0; z<4; z++){
3557                                 int x = point.x + offset[z].x;
3558                                 int y = point.y + offset[z].y;
3559                                 if(IsInside(x,y)){
3560 #if CXIMAGE_SUPPORT_SELECTION
3561                                   if (BlindSelectionIsInside(x,y))
3562 #endif //CXIMAGE_SUPPORT_SELECTION
3563                                   {
3564                                         RGBQUAD cc = BlindGetPixelColor(x, y);
3565                                         BYTE* pFill = pFillMask + x + y * head.biWidth;
3566                                         if (*pFill==0 &&
3567                                                 cRefMin.rgbRed   <= cc.rgbRed   && cc.rgbRed   <= cRefMax.rgbRed   &&
3568                                                 cRefMin.rgbGreen <= cc.rgbGreen && cc.rgbGreen <= cRefMax.rgbGreen &&
3569                                                 cRefMin.rgbBlue  <= cc.rgbBlue  && cc.rgbBlue  <= cRefMax.rgbBlue )
3570                                         {
3571                                                 if (nOpacity>0){
3572                                                         if (nOpacity == 255)
3573                                                                 BlindSetPixelColor(x, y, cFillColor);
3574                                                         else
3575                                                         {
3576                                                                 cc.rgbRed   = (BYTE)((cFillColor.rgbRed   * nOpacity + cc.rgbRed   * (255-nOpacity))>>8);
3577                                                                 cc.rgbGreen = (BYTE)((cFillColor.rgbGreen * nOpacity + cc.rgbGreen * (255-nOpacity))>>8);
3578                                                                 cc.rgbBlue  = (BYTE)((cFillColor.rgbBlue  * nOpacity + cc.rgbBlue  * (255-nOpacity))>>8);
3579                                                                 BlindSetPixelColor(x, y, cc);
3580                                                         }
3581                                                 }
3582                                                 POINT pt = {x,y};
3583                                                 q.push(pt);
3584                                                 *pFill = 1;
3585                                         }
3586                                   }
3587                                 }
3588                         }
3589                 }
3590         }
3591         if (pFillMask[xStart+yStart*head.biWidth] == 0 && nOpacity>0){
3592                 if (nOpacity == 255)
3593                         BlindSetPixelColor(xStart, yStart, cFillColor);
3594                 else
3595                 {
3596                         RGBQUAD cc = BlindGetPixelColor(xStart, yStart);
3597                         cc.rgbRed   = (BYTE)((cFillColor.rgbRed   * nOpacity + cc.rgbRed   * (255-nOpacity))>>8);
3598                         cc.rgbGreen = (BYTE)((cFillColor.rgbGreen * nOpacity + cc.rgbGreen * (255-nOpacity))>>8);
3599                         cc.rgbBlue  = (BYTE)((cFillColor.rgbBlue  * nOpacity + cc.rgbBlue  * (255-nOpacity))>>8);
3600                         BlindSetPixelColor(xStart, yStart, cc);
3601                 }
3602         }
3603         pFillMask[xStart+yStart*head.biWidth] = 1;
3604 //------------------------------------- End of Flood Fill
3605
3606         //if necessary, restore the original BPP and palette
3607         if (pPalette){
3608                 DecreaseBpp(bpp, false, pPalette);
3609                 delete [] pPalette;
3610         }
3611
3612 #if CXIMAGE_SUPPORT_SELECTION
3613         if (bSelectFilledArea){
3614                 if (!SelectionIsValid()){
3615                         if (!SelectionCreate()){
3616                                 return false;
3617                         }
3618                         SelectionClear();
3619                         info.rSelectionBox.right = head.biWidth;
3620                         info.rSelectionBox.top = head.biHeight;
3621                         info.rSelectionBox.left = info.rSelectionBox.bottom = 0;
3622                 }
3623                 RECT r;
3624                 SelectionGetBox(r);
3625                 for (long y = r.bottom; y < r.top; y++){
3626                         BYTE* pFill = pFillMask + r.left + y * head.biWidth;
3627                         for (long x = r.left; x<r.right; x++){
3628                                 if (*pFill)     SelectionSet(x,y,nSelectionLevel);
3629                                 pFill++;
3630                         }
3631                 }
3632                 SelectionRebuildBox();
3633         }
3634 #endif //CXIMAGE_SUPPORT_SELECTION
3635
3636         free(pFillMask);
3637
3638         return true;
3639 }
3640
3641 ////////////////////////////////////////////////////////////////////////////////
3642 #endif //CXIMAGE_SUPPORT_DSP