granicus.if.org Git - imagemagick/blob - magick/gem.c

   1 /*
   2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   3 %                                                                             %
   4 %                                                                             %
   5 %                                                                             %
   6 %                              GGGG  EEEEE  M   M                             %
   7 %                             G      E      MM MM                             %
   8 %                             G GG   EEE    M M M                             %
   9 %                             G   G  E      M   M                             %
  10 %                              GGGG  EEEEE  M   M                             %
  11 %                                                                             %
  12 %                                                                             %
  13 %                    Graphic Gems - Graphic Support Methods                   %
  14 %                                                                             %
  15 %                               Software Design                               %
  16 %                                 John Cristy                                 %
  17 %                                 August 1996                                 %
  18 %                                                                             %
  19 %                                                                             %
  20 %  Copyright 1999-2009 ImageMagick Studio LLC, a non-profit organization      %
  21 %  dedicated to making software imaging solutions freely available.           %
  22 %                                                                             %
  23 %  You may not use this file except in compliance with the License.  You may  %
  24 %  obtain a copy of the License at                                            %
  25 %                                                                             %
  26 %    http://www.imagemagick.org/script/license.php                            %
  27 %                                                                             %
  28 %  Unless required by applicable law or agreed to in writing, software        %
  29 %  distributed under the License is distributed on an "AS IS" BASIS,          %
  30 %  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.   %
  31 %  See the License for the specific language governing permissions and        %
  32 %  limitations under the License.                                             %
  33 %                                                                             %
  34 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  35 %
  36 %
  37 %
  38 */
  39 \f
  40 /*
  41   Include declarations.
  42 */
  43 #include "magick/studio.h"
  44 #include "magick/color-private.h"
  45 #include "magick/draw.h"
  46 #include "magick/gem.h"
  47 #include "magick/image.h"
  48 #include "magick/image-private.h"
  49 #include "magick/log.h"
  50 #include "magick/memory_.h"
  51 #include "magick/pixel-private.h"
  52 #include "magick/quantum.h"
  53 #include "magick/random_.h"
  54 #include "magick/resize.h"
  55 #include "magick/transform.h"
  56 #include "magick/signature-private.h"
  57 \f
  58 /*
  59 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  60 %                                                                             %
  61 %                                                                             %
  62 %                                                                             %
  63 %   C o n v e r t H S B T o R G B                                             %
  64 %                                                                             %
  65 %                                                                             %
  66 %                                                                             %
  67 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  68 %
  69 %  ConvertHSBToRGB() transforms a (hue, saturation, brightness) to a (red,
  70 %  green, blue) triple.
  71 %
  72 %  The format of the ConvertHSBToRGBImage method is:
  73 %
  74 %      void ConvertHSBToRGB(const double hue,const double saturation,
  75 %        const double brightness,Quantum *red,Quantum *green,Quantum *blue)
  76 %
  77 %  A description of each parameter follows:
  78 %
  79 %    o hue, saturation, brightness: A double value representing a
  80 %      component of the HSB color space.
  81 %
  82 %    o red, green, blue: A pointer to a pixel component of type Quantum.
  83 %
  84 */
  85 MagickExport void ConvertHSBToRGB(const double hue,const double saturation,
  86   const double brightness,Quantum *red,Quantum *green,Quantum *blue)
  87 {
  88   MagickRealType
  89     f,
  90     h,
  91     p,
  92     q,
  93     t;
  94
  95   /*
  96     Convert HSB to RGB colorspace.
  97   */
  98   assert(red != (Quantum *) NULL);
  99   assert(green != (Quantum *) NULL);
 100   assert(blue != (Quantum *) NULL);
 101   if (saturation == 0.0)
 102     {
 103       *red=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 104       *green=(*red);
 105       *blue=(*red);
 106       return;
 107     }
 108   h=6.0*(hue-floor(hue));
 109   f=h-floor((double) h);
 110   p=brightness*(1.0-saturation);
 111   q=brightness*(1.0-saturation*f);
 112   t=brightness*(1.0-(saturation*(1.0-f)));
 113   switch ((int) h)
 114   {
 115     case 0:
 116     default:
 117     {
 118       *red=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 119       *green=RoundToQuantum((MagickRealType) QuantumRange*t);
 120       *blue=RoundToQuantum((MagickRealType) QuantumRange*p);
 121       break;
 122     }
 123     case 1:
 124     {
 125       *red=RoundToQuantum((MagickRealType) QuantumRange*q);
 126       *green=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 127       *blue=RoundToQuantum((MagickRealType) QuantumRange*p);
 128       break;
 129     }
 130     case 2:
 131     {
 132       *red=RoundToQuantum((MagickRealType) QuantumRange*p);
 133       *green=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 134       *blue=RoundToQuantum((MagickRealType) QuantumRange*t);
 135       break;
 136     }
 137     case 3:
 138     {
 139       *red=RoundToQuantum((MagickRealType) QuantumRange*p);
 140       *green=RoundToQuantum((MagickRealType) QuantumRange*q);
 141       *blue=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 142       break;
 143     }
 144     case 4:
 145     {
 146       *red=RoundToQuantum((MagickRealType) QuantumRange*t);
 147       *green=RoundToQuantum((MagickRealType) QuantumRange*p);
 148       *blue=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 149       break;
 150     }
 151     case 5:
 152     {
 153       *red=RoundToQuantum((MagickRealType) QuantumRange*brightness);
 154       *green=RoundToQuantum((MagickRealType) QuantumRange*p);
 155       *blue=RoundToQuantum((MagickRealType) QuantumRange*q);
 156       break;
 157     }
 158   }
 159 }
 160 \f
 161 /*
 162 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 163 %                                                                             %
 164 %                                                                             %
 165 %                                                                             %
 166 %   C o n v e r t H S L T o R G B                                             %
 167 %                                                                             %
 168 %                                                                             %
 169 %                                                                             %
 170 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 171 %
 172 %  ConvertHSLToRGB() transforms a (hue, saturation, lightness) to a (red,
 173 %  green, blue) triple.
 174 %
 175 %  The format of the ConvertHSLToRGBImage method is:
 176 %
 177 %      void ConvertHSLToRGB(const double hue,const double saturation,
 178 %        const double lightness,Quantum *red,Quantum *green,Quantum *blue)
 179 %
 180 %  A description of each parameter follows:
 181 %
 182 %    o hue, saturation, lightness: A double value representing a
 183 %      component of the HSL color space.
 184 %
 185 %    o red, green, blue: A pointer to a pixel component of type Quantum.
 186 %
 187 */
 188
 189 static inline MagickRealType ConvertHueToRGB(MagickRealType m1,
 190   MagickRealType m2,MagickRealType hue)
 191 {
 192   if (hue < 0.0)
 193     hue+=1.0;
 194   if (hue > 1.0)
 195     hue-=1.0;
 196   if ((6.0*hue) < 1.0)
 197     return(m1+6.0*(m2-m1)*hue);
 198   if ((2.0*hue) < 1.0)
 199     return(m2);
 200   if ((3.0*hue) < 2.0)
 201     return(m1+6.0*(m2-m1)*(2.0/3.0-hue));
 202   return(m1);
 203 }
 204
 205 MagickExport void ConvertHSLToRGB(const double hue,const double saturation,
 206   const double lightness,Quantum *red,Quantum *green,Quantum *blue)
 207 {
 208   MagickRealType
 209     b,
 210     g,
 211     r,
 212     m1,
 213     m2;
 214
 215   /*
 216     Convert HSL to RGB colorspace.
 217   */
 218   assert(red != (Quantum *) NULL);
 219   assert(green != (Quantum *) NULL);
 220   assert(blue != (Quantum *) NULL);
 221   if (saturation == 0)
 222     {
 223       *red=RoundToQuantum((MagickRealType) QuantumRange*lightness);
 224       *green=(*red);
 225       *blue=(*red);
 226       return;
 227     }
 228   if (lightness <= 0.5)
 229     m2=lightness*(saturation+1.0);
 230   else
 231     m2=(lightness+saturation)-(lightness*saturation);
 232   m1=2.0*lightness-m2;
 233   r=ConvertHueToRGB(m1,m2,hue+1.0/3.0);
 234   g=ConvertHueToRGB(m1,m2,hue);
 235   b=ConvertHueToRGB(m1,m2,hue-1.0/3.0);
 236   *red=RoundToQuantum((MagickRealType) QuantumRange*r);
 237   *green=RoundToQuantum((MagickRealType) QuantumRange*g);
 238   *blue=RoundToQuantum((MagickRealType) QuantumRange*b);
 239 }
 240 \f
 241 /*
 242 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 243 %                                                                             %
 244 %                                                                             %
 245 %                                                                             %
 246 %   C o n v e r t H W B T o R G B                                             %
 247 %                                                                             %
 248 %                                                                             %
 249 %                                                                             %
 250 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 251 %
 252 %  ConvertHWBToRGB() transforms a (hue, whiteness, blackness) to a (red, green,
 253 %  blue) triple.
 254 %
 255 %  The format of the ConvertHWBToRGBImage method is:
 256 %
 257 %      void ConvertHWBToRGB(const double hue,const double whiteness,
 258 %        const double blackness,Quantum *red,Quantum *green,Quantum *blue)
 259 %
 260 %  A description of each parameter follows:
 261 %
 262 %    o hue, whiteness, blackness: A double value representing a
 263 %      component of the HWB color space.
 264 %
 265 %    o red, green, blue: A pointer to a pixel component of type Quantum.
 266 %
 267 */
 268 MagickExport void ConvertHWBToRGB(const double hue,const double whiteness,
 269   const double blackness,Quantum *red,Quantum *green,Quantum *blue)
 270 {
 271   MagickRealType
 272     b,
 273     f,
 274     g,
 275     n,
 276     r,
 277     v;
 278
 279   register long
 280     i;
 281
 282   /*
 283     Convert HWB to RGB colorspace.
 284   */
 285   assert(red != (Quantum *) NULL);
 286   assert(green != (Quantum *) NULL);
 287   assert(blue != (Quantum *) NULL);
 288   v=1.0-blackness;
 289   if (hue == 0.0)
 290     {
 291       *red=RoundToQuantum((MagickRealType) QuantumRange*v);
 292       *green=RoundToQuantum((MagickRealType) QuantumRange*v);
 293       *blue=RoundToQuantum((MagickRealType) QuantumRange*v);
 294       return;
 295     }
 296   i=(long) floor(6.0*hue);
 297   f=6.0*hue-i;
 298   if ((i & 0x01) != 0)
 299     f=1.0-f;
 300   n=whiteness+f*(v-whiteness);  /* linear interpolation */
 301   switch (i)
 302   {
 303     default:
 304     case 6:
 305     case 0: r=v; g=n; b=whiteness; break;
 306     case 1: r=n; g=v; b=whiteness; break;
 307     case 2: r=whiteness; g=v; b=n; break;
 308     case 3: r=whiteness; g=n; b=v; break;
 309     case 4: r=n; g=whiteness; b=v; break;
 310     case 5: r=v; g=whiteness; b=n; break;
 311   }
 312   *red=RoundToQuantum((MagickRealType) QuantumRange*r);
 313   *green=RoundToQuantum((MagickRealType) QuantumRange*g);
 314   *blue=RoundToQuantum((MagickRealType) QuantumRange*b);
 315 }
 316 \f
 317 /*
 318 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 319 %                                                                             %
 320 %                                                                             %
 321 %                                                                             %
 322 %   C o n v e r t R G B T o H S B                                             %
 323 %                                                                             %
 324 %                                                                             %
 325 %                                                                             %
 326 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 327 %
 328 %  ConvertRGBToHSB() transforms a (red, green, blue) to a (hue, saturation,
 329 %  brightness) triple.
 330 %
 331 %  The format of the ConvertRGBToHSB method is:
 332 %
 333 %      void ConvertRGBToHSB(const Quantum red,const Quantum green,
 334 %        const Quantum blue,double *hue,double *saturation,double *brightness)
 335 %
 336 %  A description of each parameter follows:
 337 %
 338 %    o red, green, blue: A Quantum value representing the red, green, and
 339 %      blue component of a pixel..
 340 %
 341 %    o hue, saturation, brightness: A pointer to a double value representing a
 342 %      component of the HSB color space.
 343 %
 344 */
 345 MagickExport void ConvertRGBToHSB(const Quantum red,const Quantum green,
 346   const Quantum blue,double *hue,double *saturation,double *brightness)
 347 {
 348   MagickRealType
 349     delta,
 350     max,
 351     min;
 352
 353   /*
 354     Convert RGB to HSB colorspace.
 355   */
 356   assert(hue != (double *) NULL);
 357   assert(saturation != (double *) NULL);
 358   assert(brightness != (double *) NULL);
 359   *hue=0.0;
 360   *saturation=0.0;
 361   *brightness=0.0;
 362   min=(MagickRealType) (red < green ? red : green);
 363   if ((MagickRealType) blue < min)
 364     min=(MagickRealType) blue;
 365   max=(MagickRealType) (red > green ? red : green);
 366   if ((MagickRealType) blue > max)
 367     max=(MagickRealType) blue;
 368   if (max == 0.0)
 369     return;
 370   delta=max-min;
 371   *saturation=(double) (delta/max);
 372   *brightness=(double) (QuantumScale*max);
 373   if (delta == 0.0)
 374     return;
 375   if ((MagickRealType) red == max)
 376     *hue=(double) ((green-(MagickRealType) blue)/delta);
 377   else
 378     if ((MagickRealType) green == max)
 379       *hue=(double) (2.0+(blue-(MagickRealType) red)/delta);
 380     else
 381       *hue=(double) (4.0+(red-(MagickRealType) green)/delta);
 382   *hue/=6.0;
 383   if (*hue < 0.0)
 384     *hue+=1.0;
 385 }
 386 \f
 387 /*
 388 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 389 %                                                                             %
 390 %                                                                             %
 391 %                                                                             %
 392 %   C o n v e r t R G B T o H S L                                             %
 393 %                                                                             %
 394 %                                                                             %
 395 %                                                                             %
 396 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 397 %
 398 %  ConvertRGBToHSL() transforms a (red, green, blue) to a (hue, saturation,
 399 %  lightness) triple.
 400 %
 401 %  The format of the ConvertRGBToHSL method is:
 402 %
 403 %      void ConvertRGBToHSL(const Quantum red,const Quantum green,
 404 %        const Quantum blue,double *hue,double *saturation,double *lightness)
 405 %
 406 %  A description of each parameter follows:
 407 %
 408 %    o red, green, blue: A Quantum value representing the red, green, and
 409 %      blue component of a pixel..
 410 %
 411 %    o hue, saturation, lightness: A pointer to a double value representing a
 412 %      component of the HSL color space.
 413 %
 414 */
 415
 416 static inline double MagickMax(const double x,const double y)
 417 {
 418   if (x > y)
 419     return(x);
 420   return(y);
 421 }
 422
 423 static inline double MagickMin(const double x,const double y)
 424 {
 425   if (x < y)
 426     return(x);
 427   return(y);
 428 }
 429
 430 MagickExport void ConvertRGBToHSL(const Quantum red,const Quantum green,
 431   const Quantum blue,double *hue,double *saturation,double *lightness)
 432 {
 433   MagickRealType
 434     b,
 435     delta,
 436     g,
 437     max,
 438     min,
 439     r;
 440
 441   /*
 442     Convert RGB to HSL colorspace.
 443   */
 444   assert(hue != (double *) NULL);
 445   assert(saturation != (double *) NULL);
 446   assert(lightness != (double *) NULL);
 447   r=QuantumScale*red;
 448   g=QuantumScale*green;
 449   b=QuantumScale*blue;
 450   max=MagickMax(r,MagickMax(g,b));
 451   min=MagickMin(r,MagickMin(g,b));
 452   *lightness=(double) ((min+max)/2.0);
 453   delta=max-min;
 454   if (delta == 0.0)
 455     {
 456       *hue=0.0;
 457       *saturation=0.0;
 458       return;
 459     }
 460   if (*lightness < 0.5)
 461     *saturation=(double) (delta/(min+max));
 462   else
 463     *saturation=(double) (delta/(2.0-max-min));
 464   if (r == max)
 465     *hue=((((max-b)/6.0)+(delta/2.0))-(((max-g)/6.0)+(delta/2.0)))/delta;
 466   else
 467     if (g == max)
 468       *hue=(1.0/3.0)+((((max-r)/6.0)+(delta/2.0))-(((max-b)/6.0)+(delta/2.0)))/
 469         delta;
 470     else
 471       if (b == max)
 472         *hue=(2.0/3.0)+((((max-g)/6.0)+(delta/2.0))-(((max-r)/6.0)+
 473           (delta/2.0)))/delta;
 474   if (*hue < 0.0)
 475     *hue+=1.0;
 476   if (*hue > 1.0)
 477     *hue-=1.0;
 478 }
 479 \f
 480 /*
 481 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 482 %                                                                             %
 483 %                                                                             %
 484 %                                                                             %
 485 %   C o n v e r t R G B T o H W B                                             %
 486 %                                                                             %
 487 %                                                                             %
 488 %                                                                             %
 489 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 490 %
 491 %  ConvertRGBToHWB() transforms a (red, green, blue) to a (hue, whiteness,
 492 %  blackness) triple.
 493 %
 494 %  The format of the ConvertRGBToHWB method is:
 495 %
 496 %      void ConvertRGBToHWB(const Quantum red,const Quantum green,
 497 %        const Quantum blue,double *hue,double *whiteness,double *blackness)
 498 %
 499 %  A description of each parameter follows:
 500 %
 501 %    o red, green, blue: A Quantum value representing the red, green, and
 502 %      blue component of a pixel.
 503 %
 504 %    o hue, whiteness, blackness: A pointer to a double value representing a
 505 %      component of the HWB color space.
 506 %
 507 */
 508 MagickExport void ConvertRGBToHWB(const Quantum red,const Quantum green,
 509   const Quantum blue,double *hue,double *whiteness,double *blackness)
 510 {
 511   MagickRealType
 512     f,
 513     v,
 514     w;
 515
 516   register long
 517     i;
 518
 519   /*
 520     Convert RGB to HWB colorspace.
 521   */
 522   assert(hue != (double *) NULL);
 523   assert(whiteness != (double *) NULL);
 524   assert(blackness != (double *) NULL);
 525   w=(MagickRealType) MagickMin((double) red,MagickMin((double) green,(double)
 526     blue));
 527   v=(MagickRealType) MagickMax((double) red,MagickMax((double) green,(double)
 528     blue));
 529   *blackness=1.0-QuantumScale*v;
 530   *whiteness=QuantumScale*w;
 531   if (v == w)
 532     {
 533       *hue=0.0;
 534       return;
 535     }
 536   f=((MagickRealType) red == w) ? green-(MagickRealType) blue :
 537     (((MagickRealType) green == w) ? blue-(MagickRealType) red : red-
 538     (MagickRealType) green);
 539   i=((MagickRealType) red == w) ? 3 : (((MagickRealType) green == w) ? 5 : 1);
 540   *hue=((double) i-f/(v-1.0*w))/6.0;
 541 }
 542 \f
 543 /*
 544 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 545 %                                                                             %
 546 %                                                                             %
 547 %                                                                             %
 548 %   E x p a n d A f f i n e                                                   %
 549 %                                                                             %
 550 %                                                                             %
 551 %                                                                             %
 552 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 553 %
 554 %  ExpandAffine() computes the affine's expansion factor, i.e. the square root
 555 %  of the factor by which the affine transform affects area. In an affine
 556 %  transform composed of scaling, rotation, shearing, and translation, returns
 557 %  the amount of scaling.
 558 %
 559 %  The format of the ExpandAffine method is:
 560 %
 561 %      double ExpandAffine(const AffineMatrix *affine)
 562 %
 563 %  A description of each parameter follows:
 564 %
 565 %    o expansion: Method ExpandAffine returns the affine's expansion factor.
 566 %
 567 %    o affine: A pointer the affine transform of type AffineMatrix.
 568 %
 569 */
 570 MagickExport double ExpandAffine(const AffineMatrix *affine)
 571 {
 572   assert(affine != (const AffineMatrix *) NULL);
 573   return(sqrt(fabs(affine->sx*affine->sy-affine->rx*affine->ry)));
 574 }
 575 \f
 576 /*
 577 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 578 %                                                                             %
 579 %                                                                             %
 580 %                                                                             %
 581 %   G e n e r a t e D i f f e r e n t i a l N o i s e                         %
 582 %                                                                             %
 583 %                                                                             %
 584 %                                                                             %
 585 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 586 %
 587 %  GenerateDifferentialNoise()
 588 %
 589 %  The format of the GenerateDifferentialNoise method is:
 590 %
 591 %      double GenerateDifferentialNoise(RandomInfo *random_info,
 592 %        const Quantum pixel,const NoiseType noise_type,
 593 %        const MagickRealType attenuate)
 594 %
 595 %  A description of each parameter follows:
 596 %
 597 %    o random_info: the random info.
 598 %
 599 %    o pixel: noise is relative to this pixel value.
 600 %
 601 %    o noise_type: the type of noise.
 602 %
 603 %    o attenuate:  attenuate the noise.
 604 %
 605 */
 606 MagickExport double GenerateDifferentialNoise(RandomInfo *random_info,
 607   const Quantum pixel,const NoiseType noise_type,const MagickRealType attenuate)
 608 {
 609 #define NoiseEpsilon  (attenuate*1.0e-5)
 610 #define SigmaUniform  ScaleCharToQuantum((unsigned char) (attenuate*4.0+0.5))
 611 #define SigmaGaussian  ScaleCharToQuantum((unsigned char) (attenuate*4.0+0.5))
 612 #define SigmaImpulse  (attenuate*0.10)
 613 #define SigmaLaplacian ScaleCharToQuantum((unsigned char) (attenuate*10.0+0.5))
 614 #define SigmaMultiplicativeGaussian \
 615   ScaleCharToQuantum((unsigned char) (attenuate*1.0+0.5))
 616 #define SigmaPoisson  (attenuate*0.05)
 617 #define TauGaussian  ScaleCharToQuantum((unsigned char) (attenuate*20.0+0.5))
 618
 619   MagickRealType
 620     alpha,
 621     beta,
 622     noise,
 623     sigma;
 624
 625   alpha=GetPseudoRandomValue(random_info);
 626   if (alpha == 0.0)
 627     alpha=1.0;
 628   switch (noise_type)
 629   {
 630     case UniformNoise:
 631     default:
 632     {
 633       noise=(MagickRealType) pixel+SigmaUniform*(alpha-0.5);
 634       break;
 635     }
 636     case GaussianNoise:
 637     {
 638       MagickRealType
 639         tau;
 640
 641       beta=GetPseudoRandomValue(random_info);
 642       sigma=sqrt(-2.0*log((double) alpha))*cos((double) (2.0*MagickPI*beta));
 643       tau=sqrt(-2.0*log((double) alpha))*sin((double) (2.0*MagickPI*beta));
 644       noise=(MagickRealType) pixel+sqrt((double) pixel)*SigmaGaussian*sigma+
 645         TauGaussian*tau;
 646       break;
 647     }
 648     case MultiplicativeGaussianNoise:
 649     {
 650       if (alpha <= NoiseEpsilon)
 651         sigma=(MagickRealType) QuantumRange;
 652       else
 653         sigma=sqrt(-2.0*log((double) alpha));
 654       beta=GetPseudoRandomValue(random_info);
 655       noise=(MagickRealType) pixel+pixel*SigmaMultiplicativeGaussian*sigma/2.0*
 656         cos((double) (2.0*MagickPI*beta));
 657       break;
 658     }
 659     case ImpulseNoise:
 660     {
 661       if (alpha < (SigmaImpulse/2.0))
 662         noise=0.0;
 663        else
 664          if (alpha >= (1.0-(SigmaImpulse/2.0)))
 665            noise=(MagickRealType) QuantumRange;
 666          else
 667            noise=(MagickRealType) pixel;
 668       break;
 669     }
 670     case LaplacianNoise:
 671     {
 672       if (alpha <= 0.5)
 673         {
 674           if (alpha <= NoiseEpsilon)
 675             noise=(MagickRealType) pixel-(MagickRealType) QuantumRange;
 676           else
 677             noise=(MagickRealType) pixel+ScaleCharToQuantum((unsigned char)
 678               (SigmaLaplacian*log((double) (2.0*alpha))+0.5));
 679           break;
 680         }
 681       beta=1.0-alpha;
 682       if (beta <= (0.5*NoiseEpsilon))
 683         noise=(MagickRealType) (pixel+QuantumRange);
 684       else
 685         noise=(MagickRealType) pixel-ScaleCharToQuantum((unsigned char)
 686           (SigmaLaplacian*log((double) (2.0*beta))+0.5));
 687       break;
 688     }
 689     case PoissonNoise:
 690     {
 691       MagickRealType
 692         poisson;
 693
 694       register long
 695         i;
 696
 697       poisson=exp(-SigmaPoisson*(double) ScaleQuantumToChar(pixel));
 698       for (i=0; alpha > poisson; i++)
 699       {
 700         beta=GetPseudoRandomValue(random_info);
 701         alpha*=beta;
 702       }
 703       noise=(MagickRealType) ScaleCharToQuantum((unsigned char)
 704         (i/SigmaPoisson));
 705       break;
 706     }
 707     case RandomNoise:
 708     {
 709       noise=(MagickRealType) QuantumRange*GetPseudoRandomValue(random_info);
 710       break;
 711     }
 712   }
 713   return(noise);
 714 }
 715 \f
 716 /*
 717 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 718 %                                                                             %
 719 %                                                                             %
 720 %                                                                             %
 721 %   G e t O p t i m a l K e r n e l W i d t h                                 %
 722 %                                                                             %
 723 %                                                                             %
 724 %                                                                             %
 725 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 726 %
 727 %  GetOptimalKernelWidth() computes the optimal kernel radius for a convolution
 728 %  filter.  Start with the minimum value of 3 pixels and walk out until we drop
 729 %  below the threshold of one pixel numerical accuracy.
 730 %
 731 %  The format of the GetOptimalKernelWidth method is:
 732 %
 733 %      unsigned long GetOptimalKernelWidth(const double radius,
 734 %        const double sigma)
 735 %
 736 %  A description of each parameter follows:
 737 %
 738 %    o width: Method GetOptimalKernelWidth returns the optimal width of
 739 %      a convolution kernel.
 740 %
 741 %    o radius: the radius of the Gaussian, in pixels, not counting the center
 742 %      pixel.
 743 %
 744 %    o sigma: the standard deviation of the Gaussian, in pixels.
 745 %
 746 */
 747 MagickExport unsigned long GetOptimalKernelWidth1D(const double radius,
 748   const double sigma)
 749 {
 750   long
 751     width;
 752
 753   MagickRealType
 754     normalize,
 755     value;
 756
 757   register long
 758     u;
 759
 760   (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
 761   if (radius > 0.0)
 762     return((unsigned long) (2.0*ceil(radius)+1.0));
 763   if (fabs(sigma) <= MagickEpsilon)
 764     return(1);
 765   for (width=5; ; )
 766   {
 767     normalize=0.0;
 768     for (u=(-width/2); u <= (width/2); u++)
 769       normalize+=exp(-((double) u*u)/(2.0*sigma*sigma))/(MagickSQ2PI*sigma);
 770     u=width/2;
 771     value=exp(-((double) u*u)/(2.0*sigma*sigma))/(MagickSQ2PI*sigma)/normalize;
 772     if ((long) (QuantumRange*value) <= 0L)
 773       break;
 774     width+=2;
 775   }
 776   return((unsigned long) (width-2));
 777 }
 778
 779 MagickExport unsigned long GetOptimalKernelWidth2D(const double radius,
 780   const double sigma)
 781 {
 782
 783   long
 784     width;
 785
 786   MagickRealType
 787     alpha,
 788     normalize,
 789     value;
 790
 791   register long
 792     u,
 793     v;
 794
 795   (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
 796   if (radius > 0.0)
 797     return((unsigned long) (2.0*ceil(radius)+1.0));
 798   if (fabs(sigma) <= MagickEpsilon)
 799     return(1);
 800   for (width=5; ; )
 801   {
 802     normalize=0.0;
 803     for (v=(-width/2); v <= (width/2); v++)
 804     {
 805       for (u=(-width/2); u <= (width/2); u++)
 806       {
 807         alpha=exp(-((double) u*u+v*v)/(2.0*sigma*sigma));
 808         normalize+=alpha/(2.0*MagickPI*sigma*sigma);
 809       }
 810     }
 811     v=width/2;
 812     value=exp(-((double) v*v)/(2.0*sigma*sigma))/normalize;
 813     if ((long) (QuantumRange*value) <= 0L)
 814       break;
 815     width+=2;
 816   }
 817   return((unsigned long) (width-2));
 818 }
 819
 820 MagickExport unsigned long  GetOptimalKernelWidth(const double radius,
 821   const double sigma)
 822 {
 823   return(GetOptimalKernelWidth1D(radius,sigma));
 824 }