granicus.if.org Git - imagemagick/blob - MagickCore/effect.c

   1 /*
   2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   3 %                                                                             %
   4 %                                                                             %
   5 %                                                                             %
   6 %                   EEEEE  FFFFF  FFFFF  EEEEE  CCCC  TTTTT                   %
   7 %                   E      F      F      E     C        T                     %
   8 %                   EEE    FFF    FFF    EEE   C        T                     %
   9 %                   E      F      F      E     C        T                     %
  10 %                   EEEEE  F      F      EEEEE  CCCC    T                     %
  11 %                                                                             %
  12 %                                                                             %
  13 %                       MagickCore Image Effects Methods                      %
  14 %                                                                             %
  15 %                               Software Design                               %
  16 %                                 John Cristy                                 %
  17 %                                 October 1996                                %
  18 %                                                                             %
  19 %                                                                             %
  20 %  Copyright 1999-2012 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 "MagickCore/studio.h"
  44 #include "MagickCore/accelerate.h"
  45 #include "MagickCore/blob.h"
  46 #include "MagickCore/cache-view.h"
  47 #include "MagickCore/color.h"
  48 #include "MagickCore/color-private.h"
  49 #include "MagickCore/colorspace.h"
  50 #include "MagickCore/constitute.h"
  51 #include "MagickCore/decorate.h"
  52 #include "MagickCore/distort.h"
  53 #include "MagickCore/draw.h"
  54 #include "MagickCore/enhance.h"
  55 #include "MagickCore/exception.h"
  56 #include "MagickCore/exception-private.h"
  57 #include "MagickCore/effect.h"
  58 #include "MagickCore/fx.h"
  59 #include "MagickCore/gem.h"
  60 #include "MagickCore/gem-private.h"
  61 #include "MagickCore/geometry.h"
  62 #include "MagickCore/image-private.h"
  63 #include "MagickCore/list.h"
  64 #include "MagickCore/log.h"
  65 #include "MagickCore/memory_.h"
  66 #include "MagickCore/monitor.h"
  67 #include "MagickCore/monitor-private.h"
  68 #include "MagickCore/montage.h"
  69 #include "MagickCore/morphology.h"
  70 #include "MagickCore/paint.h"
  71 #include "MagickCore/pixel-accessor.h"
  72 #include "MagickCore/property.h"
  73 #include "MagickCore/quantize.h"
  74 #include "MagickCore/quantum.h"
  75 #include "MagickCore/quantum-private.h"
  76 #include "MagickCore/random_.h"
  77 #include "MagickCore/random-private.h"
  78 #include "MagickCore/resample.h"
  79 #include "MagickCore/resample-private.h"
  80 #include "MagickCore/resize.h"
  81 #include "MagickCore/resource_.h"
  82 #include "MagickCore/segment.h"
  83 #include "MagickCore/shear.h"
  84 #include "MagickCore/signature-private.h"
  85 #include "MagickCore/statistic.h"
  86 #include "MagickCore/string_.h"
  87 #include "MagickCore/thread-private.h"
  88 #include "MagickCore/transform.h"
  89 #include "MagickCore/threshold.h"
  90 \f
  91 /*
  92 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  93 %                                                                             %
  94 %                                                                             %
  95 %                                                                             %
  96 %     A d a p t i v e B l u r I m a g e                                       %
  97 %                                                                             %
  98 %                                                                             %
  99 %                                                                             %
 100 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 101 %
 102 %  AdaptiveBlurImage() adaptively blurs the image by blurring less
 103 %  intensely near image edges and more intensely far from edges.  We blur the
 104 %  image with a Gaussian operator of the given radius and standard deviation
 105 %  (sigma).  For reasonable results, radius should be larger than sigma.  Use a
 106 %  radius of 0 and AdaptiveBlurImage() selects a suitable radius for you.
 107 %
 108 %  The format of the AdaptiveBlurImage method is:
 109 %
 110 %      Image *AdaptiveBlurImage(const Image *image,const double radius,
 111 %        const double sigma,const double bias,ExceptionInfo *exception)
 112 %
 113 %  A description of each parameter follows:
 114 %
 115 %    o image: the image.
 116 %
 117 %    o radius: the radius of the Gaussian, in pixels, not counting the center
 118 %      pixel.
 119 %
 120 %    o sigma: the standard deviation of the Laplacian, in pixels.
 121 %
 122 %    o bias: the bias.
 123 %
 124 %    o exception: return any errors or warnings in this structure.
 125 %
 126 */
 127
 128 MagickExport MagickBooleanType AdaptiveLevelImage(Image *image,
 129   const char *levels,ExceptionInfo *exception)
 130 {
 131   double
 132     black_point,
 133     gamma,
 134     white_point;
 135
 136   GeometryInfo
 137     geometry_info;
 138
 139   MagickBooleanType
 140     status;
 141
 142   MagickStatusType
 143     flags;
 144
 145   /*
 146     Parse levels.
 147   */
 148   if (levels == (char *) NULL)
 149     return(MagickFalse);
 150   flags=ParseGeometry(levels,&geometry_info);
 151   black_point=geometry_info.rho;
 152   white_point=(double) QuantumRange;
 153   if ((flags & SigmaValue) != 0)
 154     white_point=geometry_info.sigma;
 155   gamma=1.0;
 156   if ((flags & XiValue) != 0)
 157     gamma=geometry_info.xi;
 158   if ((flags & PercentValue) != 0)
 159     {
 160       black_point*=(double) image->columns*image->rows/100.0;
 161       white_point*=(double) image->columns*image->rows/100.0;
 162     }
 163   if ((flags & SigmaValue) == 0)
 164     white_point=(double) QuantumRange-black_point;
 165   if ((flags & AspectValue ) == 0)
 166     status=LevelImage(image,black_point,white_point,gamma,exception);
 167   else
 168     status=LevelizeImage(image,black_point,white_point,gamma,exception);
 169   return(status);
 170 }
 171
 172 MagickExport Image *AdaptiveBlurImage(const Image *image,const double radius,
 173   const double sigma,const double bias,ExceptionInfo *exception)
 174 {
 175 #define AdaptiveBlurImageTag  "Convolve/Image"
 176 #define MagickSigma  (fabs(sigma) <= MagickEpsilon ? 1.0 : sigma)
 177
 178   CacheView
 179     *blur_view,
 180     *edge_view,
 181     *image_view;
 182
 183   double
 184     **kernel,
 185     normalize;
 186
 187   Image
 188     *blur_image,
 189     *edge_image,
 190     *gaussian_image;
 191
 192   MagickBooleanType
 193     status;
 194
 195   MagickOffsetType
 196     progress;
 197
 198   register ssize_t
 199     i;
 200
 201   size_t
 202     width;
 203
 204   ssize_t
 205     j,
 206     k,
 207     u,
 208     v,
 209     y;
 210
 211   assert(image != (const Image *) NULL);
 212   assert(image->signature == MagickSignature);
 213   if (image->debug != MagickFalse)
 214     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 215   assert(exception != (ExceptionInfo *) NULL);
 216   assert(exception->signature == MagickSignature);
 217   blur_image=CloneImage(image,image->columns,image->rows,MagickTrue,exception);
 218   if (blur_image == (Image *) NULL)
 219     return((Image *) NULL);
 220   if (fabs(sigma) <= MagickEpsilon)
 221     return(blur_image);
 222   if (SetImageStorageClass(blur_image,DirectClass,exception) == MagickFalse)
 223     {
 224       blur_image=DestroyImage(blur_image);
 225       return((Image *) NULL);
 226     }
 227   /*
 228     Edge detect the image brighness channel, level, blur, and level again.
 229   */
 230   edge_image=EdgeImage(image,radius,sigma,exception);
 231   if (edge_image == (Image *) NULL)
 232     {
 233       blur_image=DestroyImage(blur_image);
 234       return((Image *) NULL);
 235     }
 236   (void) AdaptiveLevelImage(edge_image,"20%,95%",exception);
 237   gaussian_image=GaussianBlurImage(edge_image,radius,sigma,bias,exception);
 238   if (gaussian_image != (Image *) NULL)
 239     {
 240       edge_image=DestroyImage(edge_image);
 241       edge_image=gaussian_image;
 242     }
 243   (void) AdaptiveLevelImage(edge_image,"10%,95%",exception);
 244   /*
 245     Create a set of kernels from maximum (radius,sigma) to minimum.
 246   */
 247   width=GetOptimalKernelWidth2D(radius,sigma);
 248   kernel=(double **) AcquireQuantumMemory((size_t) width,sizeof(*kernel));
 249   if (kernel == (double **) NULL)
 250     {
 251       edge_image=DestroyImage(edge_image);
 252       blur_image=DestroyImage(blur_image);
 253       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
 254     }
 255   (void) ResetMagickMemory(kernel,0,(size_t) width*sizeof(*kernel));
 256   for (i=0; i < (ssize_t) width; i+=2)
 257   {
 258     kernel[i]=(double *) AcquireQuantumMemory((size_t) (width-i),(width-i)*
 259       sizeof(**kernel));
 260     if (kernel[i] == (double *) NULL)
 261       break;
 262     normalize=0.0;
 263     j=(ssize_t) (width-i)/2;
 264     k=0;
 265     for (v=(-j); v <= j; v++)
 266     {
 267       for (u=(-j); u <= j; u++)
 268       {
 269         kernel[i][k]=(double) (exp(-((double) u*u+v*v)/(2.0*MagickSigma*
 270           MagickSigma))/(2.0*MagickPI*MagickSigma*MagickSigma));
 271         normalize+=kernel[i][k];
 272         k++;
 273       }
 274     }
 275     if (fabs(normalize) <= MagickEpsilon)
 276       normalize=1.0;
 277     normalize=1.0/normalize;
 278     for (k=0; k < (j*j); k++)
 279       kernel[i][k]=normalize*kernel[i][k];
 280   }
 281   if (i < (ssize_t) width)
 282     {
 283       for (i-=2; i >= 0; i-=2)
 284         kernel[i]=(double *) RelinquishMagickMemory(kernel[i]);
 285       kernel=(double **) RelinquishMagickMemory(kernel);
 286       edge_image=DestroyImage(edge_image);
 287       blur_image=DestroyImage(blur_image);
 288       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
 289     }
 290   /*
 291     Adaptively blur image.
 292   */
 293   status=MagickTrue;
 294   progress=0;
 295   image_view=AcquireCacheView(image);
 296   edge_view=AcquireCacheView(edge_image);
 297   blur_view=AcquireCacheView(blur_image);
 298 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 299   #pragma omp parallel for schedule(static,4) shared(progress,status)
 300 #endif
 301   for (y=0; y < (ssize_t) blur_image->rows; y++)
 302   {
 303     register const Quantum
 304       *restrict r;
 305
 306     register Quantum
 307       *restrict q;
 308
 309     register ssize_t
 310       x;
 311
 312     if (status == MagickFalse)
 313       continue;
 314     r=GetCacheViewVirtualPixels(edge_view,0,y,edge_image->columns,1,exception);
 315     q=QueueCacheViewAuthenticPixels(blur_view,0,y,blur_image->columns,1,
 316       exception);
 317     if ((r == (const Quantum *) NULL) || (q == (Quantum *) NULL))
 318       {
 319         status=MagickFalse;
 320         continue;
 321       }
 322     for (x=0; x < (ssize_t) blur_image->columns; x++)
 323     {
 324       register const Quantum
 325         *restrict p;
 326
 327       register ssize_t
 328         i;
 329
 330       ssize_t
 331         center,
 332         j;
 333
 334       j=(ssize_t) ceil((double) width*QuantumScale*
 335         GetPixelIntensity(edge_image,r)-0.5);
 336       if (j < 0)
 337         j=0;
 338       else
 339         if (j > (ssize_t) width)
 340           j=(ssize_t) width;
 341       if ((j & 0x01) != 0)
 342         j--;
 343       p=GetCacheViewVirtualPixels(image_view,x-((ssize_t) (width-j)/2L),y-
 344         (ssize_t) ((width-j)/2L),width-j,width-j,exception);
 345       if (p == (const Quantum *) NULL)
 346         break;
 347       center=(ssize_t) GetPixelChannels(image)*(width-j)*
 348         ((width-j)/2L)+GetPixelChannels(image)*((width-j)/2L);
 349       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
 350       {
 351         MagickRealType
 352           alpha,
 353           gamma,
 354           pixel;
 355
 356         PixelChannel
 357           channel;
 358
 359         PixelTrait
 360           blur_traits,
 361           traits;
 362
 363         register const double
 364           *restrict k;
 365
 366         register const Quantum
 367           *restrict pixels;
 368
 369         register ssize_t
 370           u;
 371
 372         ssize_t
 373           v;
 374
 375         channel=GetPixelChannelMapChannel(image,i);
 376         traits=GetPixelChannelMapTraits(image,channel);
 377         blur_traits=GetPixelChannelMapTraits(blur_image,channel);
 378         if ((traits == UndefinedPixelTrait) ||
 379             (blur_traits == UndefinedPixelTrait))
 380           continue;
 381         if ((blur_traits & CopyPixelTrait) != 0)
 382           {
 383             SetPixelChannel(blur_image,channel,p[center+i],q);
 384             continue;
 385           }
 386         k=kernel[j];
 387         pixels=p;
 388         pixel=bias;
 389         gamma=0.0;
 390         if ((blur_traits & BlendPixelTrait) == 0)
 391           {
 392             /*
 393               No alpha blending.
 394             */
 395             for (v=0; v < (ssize_t) (width-j); v++)
 396             {
 397               for (u=0; u < (ssize_t) (width-j); u++)
 398               {
 399                 pixel+=(*k)*pixels[i];
 400                 gamma+=(*k);
 401                 k++;
 402                 pixels+=GetPixelChannels(image);
 403               }
 404             }
 405             gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
 406             SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
 407             continue;
 408           }
 409         /*
 410           Alpha blending.
 411         */
 412         for (v=0; v < (ssize_t) (width-j); v++)
 413         {
 414           for (u=0; u < (ssize_t) (width-j); u++)
 415           {
 416             alpha=(MagickRealType) (QuantumScale*GetPixelAlpha(image,pixels));
 417             pixel+=(*k)*alpha*pixels[i];
 418             gamma+=(*k)*alpha;
 419             k++;
 420             pixels+=GetPixelChannels(image);
 421           }
 422         }
 423         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
 424         SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
 425       }
 426       q+=GetPixelChannels(blur_image);
 427       r+=GetPixelChannels(edge_image);
 428     }
 429     if (SyncCacheViewAuthenticPixels(blur_view,exception) == MagickFalse)
 430       status=MagickFalse;
 431     if (image->progress_monitor != (MagickProgressMonitor) NULL)
 432       {
 433         MagickBooleanType
 434           proceed;
 435
 436 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 437   #pragma omp critical (MagickCore_AdaptiveBlurImage)
 438 #endif
 439         proceed=SetImageProgress(image,AdaptiveBlurImageTag,progress++,
 440           image->rows);
 441         if (proceed == MagickFalse)
 442           status=MagickFalse;
 443       }
 444   }
 445   blur_image->type=image->type;
 446   blur_view=DestroyCacheView(blur_view);
 447   edge_view=DestroyCacheView(edge_view);
 448   image_view=DestroyCacheView(image_view);
 449   edge_image=DestroyImage(edge_image);
 450   for (i=0; i < (ssize_t) width;  i+=2)
 451     kernel[i]=(double *) RelinquishMagickMemory(kernel[i]);
 452   kernel=(double **) RelinquishMagickMemory(kernel);
 453   if (status == MagickFalse)
 454     blur_image=DestroyImage(blur_image);
 455   return(blur_image);
 456 }
 457 \f
 458 /*
 459 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 460 %                                                                             %
 461 %                                                                             %
 462 %                                                                             %
 463 %     A d a p t i v e S h a r p e n I m a g e                                 %
 464 %                                                                             %
 465 %                                                                             %
 466 %                                                                             %
 467 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 468 %
 469 %  AdaptiveSharpenImage() adaptively sharpens the image by sharpening more
 470 %  intensely near image edges and less intensely far from edges. We sharpen the
 471 %  image with a Gaussian operator of the given radius and standard deviation
 472 %  (sigma).  For reasonable results, radius should be larger than sigma.  Use a
 473 %  radius of 0 and AdaptiveSharpenImage() selects a suitable radius for you.
 474 %
 475 %  The format of the AdaptiveSharpenImage method is:
 476 %
 477 %      Image *AdaptiveSharpenImage(const Image *image,const double radius,
 478 %        const double sigma,const double bias,ExceptionInfo *exception)
 479 %
 480 %  A description of each parameter follows:
 481 %
 482 %    o image: the image.
 483 %
 484 %    o radius: the radius of the Gaussian, in pixels, not counting the center
 485 %      pixel.
 486 %
 487 %    o sigma: the standard deviation of the Laplacian, in pixels.
 488 %
 489 %    o bias: the bias.
 490 %
 491 %    o exception: return any errors or warnings in this structure.
 492 %
 493 */
 494 MagickExport Image *AdaptiveSharpenImage(const Image *image,const double radius,
 495   const double sigma,const double bias,ExceptionInfo *exception)
 496 {
 497 #define AdaptiveSharpenImageTag  "Convolve/Image"
 498 #define MagickSigma  (fabs(sigma) <= MagickEpsilon ? 1.0 : sigma)
 499
 500   CacheView
 501     *sharp_view,
 502     *edge_view,
 503     *image_view;
 504
 505   double
 506     **kernel,
 507     normalize;
 508
 509   Image
 510     *sharp_image,
 511     *edge_image,
 512     *gaussian_image;
 513
 514   MagickBooleanType
 515     status;
 516
 517   MagickOffsetType
 518     progress;
 519
 520   register ssize_t
 521     i;
 522
 523   size_t
 524     width;
 525
 526   ssize_t
 527     j,
 528     k,
 529     u,
 530     v,
 531     y;
 532
 533   assert(image != (const Image *) NULL);
 534   assert(image->signature == MagickSignature);
 535   if (image->debug != MagickFalse)
 536     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 537   assert(exception != (ExceptionInfo *) NULL);
 538   assert(exception->signature == MagickSignature);
 539   sharp_image=CloneImage(image,0,0,MagickTrue,exception);
 540   if (sharp_image == (Image *) NULL)
 541     return((Image *) NULL);
 542   if (fabs(sigma) <= MagickEpsilon)
 543     return(sharp_image);
 544   if (SetImageStorageClass(sharp_image,DirectClass,exception) == MagickFalse)
 545     {
 546       sharp_image=DestroyImage(sharp_image);
 547       return((Image *) NULL);
 548     }
 549   /*
 550     Edge detect the image brighness channel, level, sharp, and level again.
 551   */
 552   edge_image=EdgeImage(image,radius,sigma,exception);
 553   if (edge_image == (Image *) NULL)
 554     {
 555       sharp_image=DestroyImage(sharp_image);
 556       return((Image *) NULL);
 557     }
 558   (void) AdaptiveLevelImage(edge_image,"20%,95%",exception);
 559   gaussian_image=GaussianBlurImage(edge_image,radius,sigma,bias,exception);
 560   if (gaussian_image != (Image *) NULL)
 561     {
 562       edge_image=DestroyImage(edge_image);
 563       edge_image=gaussian_image;
 564     }
 565   (void) AdaptiveLevelImage(edge_image,"10%,95%",exception);
 566   /*
 567     Create a set of kernels from maximum (radius,sigma) to minimum.
 568   */
 569   width=GetOptimalKernelWidth2D(radius,sigma);
 570   kernel=(double **) AcquireQuantumMemory((size_t) width,sizeof(*kernel));
 571   if (kernel == (double **) NULL)
 572     {
 573       edge_image=DestroyImage(edge_image);
 574       sharp_image=DestroyImage(sharp_image);
 575       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
 576     }
 577   (void) ResetMagickMemory(kernel,0,(size_t) width*sizeof(*kernel));
 578   for (i=0; i < (ssize_t) width; i+=2)
 579   {
 580     kernel[i]=(double *) AcquireQuantumMemory((size_t) (width-i),(width-i)*
 581       sizeof(**kernel));
 582     if (kernel[i] == (double *) NULL)
 583       break;
 584     normalize=0.0;
 585     j=(ssize_t) (width-i)/2;
 586     k=0;
 587     for (v=(-j); v <= j; v++)
 588     {
 589       for (u=(-j); u <= j; u++)
 590       {
 591         kernel[i][k]=(double) (-exp(-((double) u*u+v*v)/(2.0*MagickSigma*
 592           MagickSigma))/(2.0*MagickPI*MagickSigma*MagickSigma));
 593         normalize+=kernel[i][k];
 594         k++;
 595       }
 596     }
 597     if (fabs(normalize) <= MagickEpsilon)
 598       normalize=1.0;
 599     normalize=1.0/normalize;
 600     for (k=0; k < (j*j); k++)
 601       kernel[i][k]=normalize*kernel[i][k];
 602   }
 603   if (i < (ssize_t) width)
 604     {
 605       for (i-=2; i >= 0; i-=2)
 606         kernel[i]=(double *) RelinquishMagickMemory(kernel[i]);
 607       kernel=(double **) RelinquishMagickMemory(kernel);
 608       edge_image=DestroyImage(edge_image);
 609       sharp_image=DestroyImage(sharp_image);
 610       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
 611     }
 612   /*
 613     Adaptively sharpen image.
 614   */
 615   status=MagickTrue;
 616   progress=0;
 617   image_view=AcquireCacheView(image);
 618   edge_view=AcquireCacheView(edge_image);
 619   sharp_view=AcquireCacheView(sharp_image);
 620 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 621   #pragma omp parallel for schedule(static,4) shared(progress,status)
 622 #endif
 623   for (y=0; y < (ssize_t) sharp_image->rows; y++)
 624   {
 625     register const Quantum
 626       *restrict r;
 627
 628     register Quantum
 629       *restrict q;
 630
 631     register ssize_t
 632       x;
 633
 634     if (status == MagickFalse)
 635       continue;
 636     r=GetCacheViewVirtualPixels(edge_view,0,y,edge_image->columns,1,exception);
 637     q=QueueCacheViewAuthenticPixels(sharp_view,0,y,sharp_image->columns,1,
 638       exception);
 639     if ((r == (const Quantum *) NULL) || (q == (Quantum *) NULL))
 640       {
 641         status=MagickFalse;
 642         continue;
 643       }
 644     for (x=0; x < (ssize_t) sharp_image->columns; x++)
 645     {
 646       register const Quantum
 647         *restrict p;
 648
 649       register ssize_t
 650         i;
 651
 652       ssize_t
 653         center,
 654         j;
 655
 656       j=(ssize_t) ceil((double) width*QuantumScale*
 657         GetPixelIntensity(edge_image,r)-0.5);
 658       if (j < 0)
 659         j=0;
 660       else
 661         if (j > (ssize_t) width)
 662           j=(ssize_t) width;
 663       if ((j & 0x01) != 0)
 664         j--;
 665       p=GetCacheViewVirtualPixels(image_view,x-((ssize_t) (width-j)/2L),y-
 666         (ssize_t) ((width-j)/2L),width-j,width-j,exception);
 667       if (p == (const Quantum *) NULL)
 668         break;
 669       center=(ssize_t) GetPixelChannels(image)*(width-j)*
 670         ((width-j)/2L)+GetPixelChannels(image)*((width-j)/2);
 671       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
 672       {
 673         MagickRealType
 674           alpha,
 675           gamma,
 676           pixel;
 677
 678         PixelChannel
 679           channel;
 680
 681         PixelTrait
 682           sharp_traits,
 683           traits;
 684
 685         register const double
 686           *restrict k;
 687
 688         register const Quantum
 689           *restrict pixels;
 690
 691         register ssize_t
 692           u;
 693
 694         ssize_t
 695           v;
 696
 697         channel=GetPixelChannelMapChannel(image,i);
 698         traits=GetPixelChannelMapTraits(image,channel);
 699         sharp_traits=GetPixelChannelMapTraits(sharp_image,channel);
 700         if ((traits == UndefinedPixelTrait) ||
 701             (sharp_traits == UndefinedPixelTrait))
 702           continue;
 703         if ((sharp_traits & CopyPixelTrait) != 0)
 704           {
 705             SetPixelChannel(sharp_image,channel,p[center+i],q);
 706             continue;
 707           }
 708         k=kernel[j];
 709         pixels=p;
 710         pixel=bias;
 711         gamma=0.0;
 712         if ((sharp_traits & BlendPixelTrait) == 0)
 713           {
 714             /*
 715               No alpha blending.
 716             */
 717             for (v=0; v < (ssize_t) (width-j); v++)
 718             {
 719               for (u=0; u < (ssize_t) (width-j); u++)
 720               {
 721                 pixel+=(*k)*pixels[i];
 722                 gamma+=(*k);
 723                 k++;
 724                 pixels+=GetPixelChannels(image);
 725               }
 726             }
 727             gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
 728             SetPixelChannel(sharp_image,channel,ClampToQuantum(gamma*pixel),q);
 729             continue;
 730           }
 731         /*
 732           Alpha blending.
 733         */
 734         for (v=0; v < (ssize_t) (width-j); v++)
 735         {
 736           for (u=0; u < (ssize_t) (width-j); u++)
 737           {
 738             alpha=(MagickRealType) (QuantumScale*GetPixelAlpha(image,pixels));
 739             pixel+=(*k)*alpha*pixels[i];
 740             gamma+=(*k)*alpha;
 741             k++;
 742             pixels+=GetPixelChannels(image);
 743           }
 744         }
 745         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
 746         SetPixelChannel(sharp_image,channel,ClampToQuantum(gamma*pixel),q);
 747       }
 748       q+=GetPixelChannels(sharp_image);
 749       r+=GetPixelChannels(edge_image);
 750     }
 751     if (SyncCacheViewAuthenticPixels(sharp_view,exception) == MagickFalse)
 752       status=MagickFalse;
 753     if (image->progress_monitor != (MagickProgressMonitor) NULL)
 754       {
 755         MagickBooleanType
 756           proceed;
 757
 758 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 759   #pragma omp critical (MagickCore_AdaptiveSharpenImage)
 760 #endif
 761         proceed=SetImageProgress(image,AdaptiveSharpenImageTag,progress++,
 762           image->rows);
 763         if (proceed == MagickFalse)
 764           status=MagickFalse;
 765       }
 766   }
 767   sharp_image->type=image->type;
 768   sharp_view=DestroyCacheView(sharp_view);
 769   edge_view=DestroyCacheView(edge_view);
 770   image_view=DestroyCacheView(image_view);
 771   edge_image=DestroyImage(edge_image);
 772   for (i=0; i < (ssize_t) width;  i+=2)
 773     kernel[i]=(double *) RelinquishMagickMemory(kernel[i]);
 774   kernel=(double **) RelinquishMagickMemory(kernel);
 775   if (status == MagickFalse)
 776     sharp_image=DestroyImage(sharp_image);
 777   return(sharp_image);
 778 }
 779 \f
 780 /*
 781 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 782 %                                                                             %
 783 %                                                                             %
 784 %                                                                             %
 785 %     B l u r I m a g e                                                       %
 786 %                                                                             %
 787 %                                                                             %
 788 %                                                                             %
 789 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 790 %
 791 %  BlurImage() blurs an image.  We convolve the image with a Gaussian operator
 792 %  of the given radius and standard deviation (sigma).  For reasonable results,
 793 %  the radius should be larger than sigma.  Use a radius of 0 and BlurImage()
 794 %  selects a suitable radius for you.
 795 %
 796 %  BlurImage() differs from GaussianBlurImage() in that it uses a separable
 797 %  kernel which is faster but mathematically equivalent to the non-separable
 798 %  kernel.
 799 %
 800 %  The format of the BlurImage method is:
 801 %
 802 %      Image *BlurImage(const Image *image,const double radius,
 803 %        const double sigma,const double bias,ExceptionInfo *exception)
 804 %
 805 %  A description of each parameter follows:
 806 %
 807 %    o image: the image.
 808 %
 809 %    o radius: the radius of the Gaussian, in pixels, not counting the center
 810 %      pixel.
 811 %
 812 %    o sigma: the standard deviation of the Gaussian, in pixels.
 813 %
 814 %    o bias: the bias.
 815 %
 816 %    o exception: return any errors or warnings in this structure.
 817 %
 818 */
 819
 820 static double *GetBlurKernel(const size_t width,const double sigma)
 821 {
 822   double
 823     *kernel,
 824     normalize;
 825
 826   register ssize_t
 827     i;
 828
 829   ssize_t
 830     j,
 831     k;
 832
 833   /*
 834     Generate a 1-D convolution kernel.
 835   */
 836   (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
 837   kernel=(double *) AcquireQuantumMemory((size_t) width,sizeof(*kernel));
 838   if (kernel == (double *) NULL)
 839     return(0);
 840   normalize=0.0;
 841   j=(ssize_t) width/2;
 842   i=0;
 843   for (k=(-j); k <= j; k++)
 844   {
 845     kernel[i]=(double) (exp(-((double) k*k)/(2.0*MagickSigma*MagickSigma))/
 846       (MagickSQ2PI*MagickSigma));
 847     normalize+=kernel[i];
 848     i++;
 849   }
 850   for (i=0; i < (ssize_t) width; i++)
 851     kernel[i]/=normalize;
 852   return(kernel);
 853 }
 854
 855 MagickExport Image *BlurImage(const Image *image,const double radius,
 856   const double sigma,const double bias,ExceptionInfo *exception)
 857 {
 858 #define BlurImageTag  "Blur/Image"
 859
 860   CacheView
 861     *blur_view,
 862     *image_view;
 863
 864   double
 865     *kernel;
 866
 867   Image
 868     *blur_image;
 869
 870   MagickBooleanType
 871     status;
 872
 873   MagickOffsetType
 874     progress;
 875
 876   register ssize_t
 877     i;
 878
 879   size_t
 880     width;
 881
 882   ssize_t
 883     center,
 884     x,
 885     y;
 886
 887   /*
 888     Initialize blur image attributes.
 889   */
 890   assert(image != (Image *) NULL);
 891   assert(image->signature == MagickSignature);
 892   if (image->debug != MagickFalse)
 893     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
 894   assert(exception != (ExceptionInfo *) NULL);
 895   assert(exception->signature == MagickSignature);
 896   blur_image=CloneImage(image,0,0,MagickTrue,exception);
 897   if (blur_image == (Image *) NULL)
 898     return((Image *) NULL);
 899   if (fabs(sigma) <= MagickEpsilon)
 900     return(blur_image);
 901   if (SetImageStorageClass(blur_image,DirectClass,exception) == MagickFalse)
 902     {
 903       blur_image=DestroyImage(blur_image);
 904       return((Image *) NULL);
 905     }
 906   width=GetOptimalKernelWidth1D(radius,sigma);
 907   kernel=GetBlurKernel(width,sigma);
 908   if (kernel == (double *) NULL)
 909     {
 910       blur_image=DestroyImage(blur_image);
 911       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
 912     }
 913   if (image->debug != MagickFalse)
 914     {
 915       char
 916         format[MaxTextExtent],
 917         *message;
 918
 919       register const double
 920         *k;
 921
 922       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
 923         "  BlurImage with %.20g kernel:",(double) width);
 924       message=AcquireString("");
 925       k=kernel;
 926       for (i=0; i < (ssize_t) width; i++)
 927       {
 928         *message='\0';
 929         (void) FormatLocaleString(format,MaxTextExtent,"%.20g: ",(double) i);
 930         (void) ConcatenateString(&message,format);
 931         (void) FormatLocaleString(format,MaxTextExtent,"%g ",*k++);
 932         (void) ConcatenateString(&message,format);
 933         (void) LogMagickEvent(TransformEvent,GetMagickModule(),"%s",message);
 934       }
 935       message=DestroyString(message);
 936     }
 937   /*
 938     Blur rows.
 939   */
 940   status=MagickTrue;
 941   progress=0;
 942   center=(ssize_t) GetPixelChannels(image)*(width/2L);
 943   image_view=AcquireCacheView(image);
 944   blur_view=AcquireCacheView(blur_image);
 945 #if defined(MAGICKCORE_OPENMP_SUPPORT)
 946   #pragma omp parallel for schedule(static,4) shared(progress,status)
 947 #endif
 948   for (y=0; y < (ssize_t) image->rows; y++)
 949   {
 950     register const Quantum
 951       *restrict p;
 952
 953     register Quantum
 954       *restrict q;
 955
 956     register ssize_t
 957       x;
 958
 959     if (status == MagickFalse)
 960       continue;
 961     p=GetCacheViewVirtualPixels(image_view,-((ssize_t) width/2L),y,
 962       image->columns+width,1,exception);
 963     q=GetCacheViewAuthenticPixels(blur_view,0,y,blur_image->columns,1,
 964       exception);
 965     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
 966       {
 967         status=MagickFalse;
 968         continue;
 969       }
 970     for (x=0; x < (ssize_t) image->columns; x++)
 971     {
 972       register ssize_t
 973         i;
 974
 975       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
 976       {
 977         MagickRealType
 978           alpha,
 979           gamma,
 980           pixel;
 981
 982         PixelChannel
 983           channel;
 984
 985         PixelTrait
 986           blur_traits,
 987           traits;
 988
 989         register const double
 990           *restrict k;
 991
 992         register const Quantum
 993           *restrict pixels;
 994
 995         register ssize_t
 996           u;
 997
 998         channel=GetPixelChannelMapChannel(image,i);
 999         traits=GetPixelChannelMapTraits(image,channel);
1000         blur_traits=GetPixelChannelMapTraits(blur_image,channel);
1001         if ((traits == UndefinedPixelTrait) ||
1002             (blur_traits == UndefinedPixelTrait))
1003           continue;
1004         if ((blur_traits & CopyPixelTrait) != 0)
1005           {
1006             SetPixelChannel(blur_image,channel,p[center+i],q);
1007             continue;
1008           }
1009         k=kernel;
1010         pixels=p;
1011         pixel=0.0;
1012         if ((blur_traits & BlendPixelTrait) == 0)
1013           {
1014             /*
1015               No alpha blending.
1016             */
1017             for (u=0; u < (ssize_t) width; u++)
1018             {
1019               pixel+=(*k)*pixels[i];
1020               k++;
1021               pixels+=GetPixelChannels(image);
1022             }
1023             SetPixelChannel(blur_image,channel,ClampToQuantum(pixel),q);
1024             continue;
1025           }
1026         /*
1027           Alpha blending.
1028         */
1029         gamma=0.0;
1030         for (u=0; u < (ssize_t) width; u++)
1031         {
1032           alpha=(MagickRealType) (QuantumScale*GetPixelAlpha(image,pixels));
1033           pixel+=(*k)*alpha*pixels[i];
1034           gamma+=(*k)*alpha;
1035           k++;
1036           pixels+=GetPixelChannels(image);
1037         }
1038         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
1039         SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
1040       }
1041       p+=GetPixelChannels(image);
1042       q+=GetPixelChannels(blur_image);
1043     }
1044     if (SyncCacheViewAuthenticPixels(blur_view,exception) == MagickFalse)
1045       status=MagickFalse;
1046     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1047       {
1048         MagickBooleanType
1049           proceed;
1050
1051 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1052   #pragma omp critical (MagickCore_BlurImage)
1053 #endif
1054         proceed=SetImageProgress(image,BlurImageTag,progress++,blur_image->rows+
1055           blur_image->columns);
1056         if (proceed == MagickFalse)
1057           status=MagickFalse;
1058       }
1059   }
1060   blur_view=DestroyCacheView(blur_view);
1061   image_view=DestroyCacheView(image_view);
1062   /*
1063     Blur columns.
1064   */
1065   image_view=AcquireCacheView(blur_image);
1066   blur_view=AcquireCacheView(blur_image);
1067 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1068   #pragma omp parallel for schedule(static,4) shared(progress,status)
1069 #endif
1070   for (x=0; x < (ssize_t) blur_image->columns; x++)
1071   {
1072     register const Quantum
1073       *restrict p;
1074
1075     register Quantum
1076       *restrict q;
1077
1078     register ssize_t
1079       y;
1080
1081     if (status == MagickFalse)
1082       continue;
1083     p=GetCacheViewVirtualPixels(image_view,x,-((ssize_t) width/2L),1,
1084       blur_image->rows+width,exception);
1085     q=GetCacheViewAuthenticPixels(blur_view,x,0,1,blur_image->rows,exception);
1086     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
1087       {
1088         status=MagickFalse;
1089         continue;
1090       }
1091     for (y=0; y < (ssize_t) blur_image->rows; y++)
1092     {
1093       register ssize_t
1094         i;
1095
1096       for (i=0; i < (ssize_t) GetPixelChannels(blur_image); i++)
1097       {
1098         MagickRealType
1099           alpha,
1100           gamma,
1101           pixel;
1102
1103         PixelChannel
1104           channel;
1105
1106         PixelTrait
1107           blur_traits,
1108           traits;
1109
1110         register const double
1111           *restrict k;
1112
1113         register const Quantum
1114           *restrict pixels;
1115
1116         register ssize_t
1117           u;
1118
1119         channel=GetPixelChannelMapChannel(blur_image,i);
1120         traits=GetPixelChannelMapTraits(blur_image,channel);
1121         blur_traits=GetPixelChannelMapTraits(blur_image,channel);
1122         if ((traits == UndefinedPixelTrait) ||
1123             (blur_traits == UndefinedPixelTrait))
1124           continue;
1125         if ((blur_traits & CopyPixelTrait) != 0)
1126           {
1127             SetPixelChannel(blur_image,channel,p[center+i],q);
1128             continue;
1129           }
1130         k=kernel;
1131         pixels=p;
1132         pixel=0.0;
1133         if ((blur_traits & BlendPixelTrait) == 0)
1134           {
1135             /*
1136               No alpha blending.
1137             */
1138             for (u=0; u < (ssize_t) width; u++)
1139             {
1140               pixel+=(*k)*pixels[i];
1141               k++;
1142               pixels+=GetPixelChannels(blur_image);
1143             }
1144             SetPixelChannel(blur_image,channel,ClampToQuantum(pixel),q);
1145             continue;
1146           }
1147         /*
1148           Alpha blending.
1149         */
1150         gamma=0.0;
1151         for (u=0; u < (ssize_t) width; u++)
1152         {
1153           alpha=(MagickRealType) (QuantumScale*GetPixelAlpha(blur_image,
1154             pixels));
1155           pixel+=(*k)*alpha*pixels[i];
1156           gamma+=(*k)*alpha;
1157           k++;
1158           pixels+=GetPixelChannels(blur_image);
1159         }
1160         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
1161         SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
1162       }
1163       p+=GetPixelChannels(blur_image);
1164       q+=GetPixelChannels(blur_image);
1165     }
1166     if (SyncCacheViewAuthenticPixels(blur_view,exception) == MagickFalse)
1167       status=MagickFalse;
1168     if (blur_image->progress_monitor != (MagickProgressMonitor) NULL)
1169       {
1170         MagickBooleanType
1171           proceed;
1172
1173 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1174   #pragma omp critical (MagickCore_BlurImage)
1175 #endif
1176         proceed=SetImageProgress(blur_image,BlurImageTag,progress++,
1177           blur_image->rows+blur_image->columns);
1178         if (proceed == MagickFalse)
1179           status=MagickFalse;
1180       }
1181   }
1182   blur_view=DestroyCacheView(blur_view);
1183   image_view=DestroyCacheView(image_view);
1184   kernel=(double *) RelinquishMagickMemory(kernel);
1185   if (status == MagickFalse)
1186     blur_image=DestroyImage(blur_image);
1187   blur_image->type=image->type;
1188   return(blur_image);
1189 }
1190 \f
1191 /*
1192 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1193 %                                                                             %
1194 %                                                                             %
1195 %                                                                             %
1196 %     C o n v o l v e I m a g e                                               %
1197 %                                                                             %
1198 %                                                                             %
1199 %                                                                             %
1200 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1201 %
1202 %  ConvolveImage() applies a custom convolution kernel to the image.
1203 %
1204 %  The format of the ConvolveImage method is:
1205 %
1206 %      Image *ConvolveImage(const Image *image,const KernelInfo *kernel,
1207 %        ExceptionInfo *exception)
1208 %
1209 %  A description of each parameter follows:
1210 %
1211 %    o image: the image.
1212 %
1213 %    o kernel: the filtering kernel.
1214 %
1215 %    o exception: return any errors or warnings in this structure.
1216 %
1217 */
1218 MagickExport Image *ConvolveImage(const Image *image,
1219   const KernelInfo *kernel_info,ExceptionInfo *exception)
1220 {
1221 #define ConvolveImageTag  "Convolve/Image"
1222
1223   CacheView
1224     *convolve_view,
1225     *image_view;
1226
1227   Image
1228     *convolve_image;
1229
1230   MagickBooleanType
1231     status;
1232
1233   MagickOffsetType
1234     progress;
1235
1236   ssize_t
1237     center,
1238     y;
1239
1240   /*
1241     Initialize convolve image attributes.
1242   */
1243   assert(image != (Image *) NULL);
1244   assert(image->signature == MagickSignature);
1245   if (image->debug != MagickFalse)
1246     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1247   assert(exception != (ExceptionInfo *) NULL);
1248   assert(exception->signature == MagickSignature);
1249   if ((kernel_info->width % 2) == 0)
1250     ThrowImageException(OptionError,"KernelWidthMustBeAnOddNumber");
1251   convolve_image=CloneImage(image,image->columns,image->rows,MagickTrue,
1252     exception);
1253   if (convolve_image == (Image *) NULL)
1254     return((Image *) NULL);
1255   if (SetImageStorageClass(convolve_image,DirectClass,exception) == MagickFalse)
1256     {
1257       convolve_image=DestroyImage(convolve_image);
1258       return((Image *) NULL);
1259     }
1260   if (image->debug != MagickFalse)
1261     {
1262       char
1263         format[MaxTextExtent],
1264         *message;
1265
1266       register const MagickRealType
1267         *k;
1268
1269       register ssize_t
1270         u;
1271
1272       ssize_t
1273         v;
1274
1275       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1276         "  ConvolveImage with %.20gx%.20g kernel:",(double) kernel_info->width,
1277         (double) kernel_info->height);
1278       message=AcquireString("");
1279       k=kernel_info->values;
1280       for (v=0; v < (ssize_t) kernel_info->width; v++)
1281       {
1282         *message='\0';
1283         (void) FormatLocaleString(format,MaxTextExtent,"%.20g: ",(double) v);
1284         (void) ConcatenateString(&message,format);
1285         for (u=0; u < (ssize_t) kernel_info->height; u++)
1286         {
1287           (void) FormatLocaleString(format,MaxTextExtent,"%g ",*k++);
1288           (void) ConcatenateString(&message,format);
1289         }
1290         (void) LogMagickEvent(TransformEvent,GetMagickModule(),"%s",message);
1291       }
1292       message=DestroyString(message);
1293     }
1294   status=AccelerateConvolveImage(image,kernel_info,convolve_image,exception);
1295   if (status == MagickTrue)
1296     return(convolve_image);
1297   /*
1298     Convolve image.
1299   */
1300   center=(ssize_t) GetPixelChannels(image)*(image->columns+kernel_info->width)*
1301     (kernel_info->height/2L)+GetPixelChannels(image)*(kernel_info->width/2L);
1302   status=MagickTrue;
1303   progress=0;
1304   image_view=AcquireCacheView(image);
1305   convolve_view=AcquireCacheView(convolve_image);
1306 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1307   #pragma omp parallel for schedule(static,4) shared(progress,status)
1308 #endif
1309   for (y=0; y < (ssize_t) image->rows; y++)
1310   {
1311     register const Quantum
1312       *restrict p;
1313
1314     register Quantum
1315       *restrict q;
1316
1317     register ssize_t
1318       x;
1319
1320     if (status == MagickFalse)
1321       continue;
1322     p=GetCacheViewVirtualPixels(image_view,-((ssize_t) kernel_info->width/2L),y-
1323       (ssize_t) (kernel_info->height/2L),image->columns+kernel_info->width,
1324       kernel_info->height,exception);
1325     q=QueueCacheViewAuthenticPixels(convolve_view,0,y,convolve_image->columns,1,
1326       exception);
1327     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
1328       {
1329         status=MagickFalse;
1330         continue;
1331       }
1332     for (x=0; x < (ssize_t) image->columns; x++)
1333     {
1334       register ssize_t
1335         i;
1336
1337       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1338       {
1339         MagickRealType
1340           alpha,
1341           gamma,
1342           pixel;
1343
1344         PixelChannel
1345           channel;
1346
1347         PixelTrait
1348           convolve_traits,
1349           traits;
1350
1351         register const MagickRealType
1352           *restrict k;
1353
1354         register const Quantum
1355           *restrict pixels;
1356
1357         register ssize_t
1358           u;
1359
1360         ssize_t
1361           v;
1362
1363         channel=GetPixelChannelMapChannel(image,i);
1364         traits=GetPixelChannelMapTraits(image,channel);
1365         convolve_traits=GetPixelChannelMapTraits(convolve_image,channel);
1366         if ((traits == UndefinedPixelTrait) ||
1367             (convolve_traits == UndefinedPixelTrait))
1368           continue;
1369         if ((convolve_traits & CopyPixelTrait) != 0)
1370           {
1371             SetPixelChannel(convolve_image,channel,p[center+i],q);
1372             continue;
1373           }
1374         k=kernel_info->values;
1375         pixels=p;
1376         pixel=kernel_info->bias;
1377         if ((convolve_traits & BlendPixelTrait) == 0)
1378           {
1379             /*
1380               No alpha blending.
1381             */
1382             for (v=0; v < (ssize_t) kernel_info->height; v++)
1383             {
1384               for (u=0; u < (ssize_t) kernel_info->width; u++)
1385               {
1386                 pixel+=(*k)*pixels[i];
1387                 k++;
1388                 pixels+=GetPixelChannels(image);
1389               }
1390               pixels+=image->columns*GetPixelChannels(image);
1391             }
1392             SetPixelChannel(convolve_image,channel,ClampToQuantum(pixel),q);
1393             continue;
1394           }
1395         /*
1396           Alpha blending.
1397         */
1398         gamma=0.0;
1399         for (v=0; v < (ssize_t) kernel_info->height; v++)
1400         {
1401           for (u=0; u < (ssize_t) kernel_info->width; u++)
1402           {
1403             alpha=(MagickRealType) (QuantumScale*GetPixelAlpha(image,pixels));
1404             pixel+=(*k)*alpha*pixels[i];
1405             gamma+=(*k)*alpha;
1406             k++;
1407             pixels+=GetPixelChannels(image);
1408           }
1409           pixels+=image->columns*GetPixelChannels(image);
1410         }
1411         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
1412         SetPixelChannel(convolve_image,channel,ClampToQuantum(gamma*pixel),q);
1413       }
1414       p+=GetPixelChannels(image);
1415       q+=GetPixelChannels(convolve_image);
1416     }
1417     if (SyncCacheViewAuthenticPixels(convolve_view,exception) == MagickFalse)
1418       status=MagickFalse;
1419     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1420       {
1421         MagickBooleanType
1422           proceed;
1423
1424 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1425   #pragma omp critical (MagickCore_ConvolveImage)
1426 #endif
1427         proceed=SetImageProgress(image,ConvolveImageTag,progress++,image->rows);
1428         if (proceed == MagickFalse)
1429           status=MagickFalse;
1430       }
1431   }
1432   convolve_image->type=image->type;
1433   convolve_view=DestroyCacheView(convolve_view);
1434   image_view=DestroyCacheView(image_view);
1435   if (status == MagickFalse)
1436     convolve_image=DestroyImage(convolve_image);
1437   return(convolve_image);
1438 }
1439 \f
1440 /*
1441 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1442 %                                                                             %
1443 %                                                                             %
1444 %                                                                             %
1445 %     D e s p e c k l e I m a g e                                             %
1446 %                                                                             %
1447 %                                                                             %
1448 %                                                                             %
1449 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1450 %
1451 %  DespeckleImage() reduces the speckle noise in an image while perserving the
1452 %  edges of the original image.
1453 %
1454 %  The format of the DespeckleImage method is:
1455 %
1456 %      Image *DespeckleImage(const Image *image,ExceptionInfo *exception)
1457 %
1458 %  A description of each parameter follows:
1459 %
1460 %    o image: the image.
1461 %
1462 %    o exception: return any errors or warnings in this structure.
1463 %
1464 */
1465
1466 static void Hull(const ssize_t x_offset,const ssize_t y_offset,
1467   const size_t columns,const size_t rows,Quantum *f,Quantum *g,
1468   const int polarity)
1469 {
1470   MagickRealType
1471     v;
1472
1473   register Quantum
1474     *p,
1475     *q,
1476     *r,
1477     *s;
1478
1479   register ssize_t
1480     x;
1481
1482   ssize_t
1483     y;
1484
1485   assert(f != (Quantum *) NULL);
1486   assert(g != (Quantum *) NULL);
1487   p=f+(columns+2);
1488   q=g+(columns+2);
1489   r=p+(y_offset*((ssize_t) columns+2)+x_offset);
1490   for (y=0; y < (ssize_t) rows; y++)
1491   {
1492     p++;
1493     q++;
1494     r++;
1495     if (polarity > 0)
1496       for (x=(ssize_t) columns; x != 0; x--)
1497       {
1498         v=(MagickRealType) (*p);
1499         if ((MagickRealType) *r >= (v+(MagickRealType) ScaleCharToQuantum(2)))
1500           v+=ScaleCharToQuantum(1);
1501         *q=(Quantum) v;
1502         p++;
1503         q++;
1504         r++;
1505       }
1506     else
1507       for (x=(ssize_t) columns; x != 0; x--)
1508       {
1509         v=(MagickRealType) (*p);
1510         if ((MagickRealType) *r <= (v-(MagickRealType) ScaleCharToQuantum(2)))
1511           v-=(ssize_t) ScaleCharToQuantum(1);
1512         *q=(Quantum) v;
1513         p++;
1514         q++;
1515         r++;
1516       }
1517     p++;
1518     q++;
1519     r++;
1520   }
1521   p=f+(columns+2);
1522   q=g+(columns+2);
1523   r=q+(y_offset*((ssize_t) columns+2)+x_offset);
1524   s=q-(y_offset*((ssize_t) columns+2)+x_offset);
1525   for (y=0; y < (ssize_t) rows; y++)
1526   {
1527     p++;
1528     q++;
1529     r++;
1530     s++;
1531     if (polarity > 0)
1532       for (x=(ssize_t) columns; x != 0; x--)
1533       {
1534         v=(MagickRealType) (*q);
1535         if (((MagickRealType) *s >=
1536              (v+(MagickRealType) ScaleCharToQuantum(2))) &&
1537             ((MagickRealType) *r > v))
1538           v+=ScaleCharToQuantum(1);
1539         *p=(Quantum) v;
1540         p++;
1541         q++;
1542         r++;
1543         s++;
1544       }
1545     else
1546       for (x=(ssize_t) columns; x != 0; x--)
1547       {
1548         v=(MagickRealType) (*q);
1549         if (((MagickRealType) *s <=
1550              (v-(MagickRealType) ScaleCharToQuantum(2))) &&
1551             ((MagickRealType) *r < v))
1552           v-=(MagickRealType) ScaleCharToQuantum(1);
1553         *p=(Quantum) v;
1554         p++;
1555         q++;
1556         r++;
1557         s++;
1558       }
1559     p++;
1560     q++;
1561     r++;
1562     s++;
1563   }
1564 }
1565
1566 MagickExport Image *DespeckleImage(const Image *image,ExceptionInfo *exception)
1567 {
1568 #define DespeckleImageTag  "Despeckle/Image"
1569
1570   CacheView
1571     *despeckle_view,
1572     *image_view;
1573
1574   Image
1575     *despeckle_image;
1576
1577   MagickBooleanType
1578     status;
1579
1580   Quantum
1581     *restrict buffers,
1582     *restrict pixels;
1583
1584   register ssize_t
1585     i;
1586
1587   size_t
1588     length;
1589
1590   static const ssize_t
1591     X[4] = {0, 1, 1,-1},
1592     Y[4] = {1, 0, 1, 1};
1593
1594   /*
1595     Allocate despeckled image.
1596   */
1597   assert(image != (const Image *) NULL);
1598   assert(image->signature == MagickSignature);
1599   if (image->debug != MagickFalse)
1600     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1601   assert(exception != (ExceptionInfo *) NULL);
1602   assert(exception->signature == MagickSignature);
1603   despeckle_image=CloneImage(image,0,0,MagickTrue,exception);
1604   if (despeckle_image == (Image *) NULL)
1605     return((Image *) NULL);
1606   status=SetImageStorageClass(despeckle_image,DirectClass,exception);
1607   if (status == MagickFalse)
1608     {
1609       despeckle_image=DestroyImage(despeckle_image);
1610       return((Image *) NULL);
1611     }
1612   /*
1613     Allocate image buffers.
1614   */
1615   length=(size_t) ((image->columns+2)*(image->rows+2));
1616   pixels=(Quantum *) AcquireQuantumMemory(length,2*sizeof(*pixels));
1617   buffers=(Quantum *) AcquireQuantumMemory(length,2*sizeof(*pixels));
1618   if ((pixels == (Quantum *) NULL) || (buffers == (Quantum *) NULL))
1619     {
1620       if (buffers != (Quantum *) NULL)
1621         buffers=(Quantum *) RelinquishMagickMemory(buffers);
1622       if (pixels != (Quantum *) NULL)
1623         pixels=(Quantum *) RelinquishMagickMemory(pixels);
1624       despeckle_image=DestroyImage(despeckle_image);
1625       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1626     }
1627   /*
1628     Reduce speckle in the image.
1629   */
1630   status=MagickTrue;
1631   image_view=AcquireCacheView(image);
1632   despeckle_view=AcquireCacheView(despeckle_image);
1633   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1634   {
1635     PixelChannel
1636        channel;
1637
1638     PixelTrait
1639       despeckle_traits,
1640       traits;
1641
1642     register Quantum
1643       *buffer,
1644       *pixel;
1645
1646     register ssize_t
1647       k,
1648       x;
1649
1650     ssize_t
1651       j,
1652       y;
1653
1654     if (status == MagickFalse)
1655       continue;
1656     channel=GetPixelChannelMapChannel(image,i);
1657     traits=GetPixelChannelMapTraits(image,channel);
1658     despeckle_traits=GetPixelChannelMapTraits(despeckle_image,channel);
1659     if ((traits == UndefinedPixelTrait) ||
1660         (despeckle_traits == UndefinedPixelTrait))
1661       continue;
1662     if ((despeckle_traits & CopyPixelTrait) != 0)
1663       continue;
1664     pixel=pixels;
1665     (void) ResetMagickMemory(pixel,0,length*sizeof(*pixel));
1666     buffer=buffers;
1667     j=(ssize_t) image->columns+2;
1668     for (y=0; y < (ssize_t) image->rows; y++)
1669     {
1670       register const Quantum
1671         *restrict p;
1672
1673       p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
1674       if (p == (const Quantum *) NULL)
1675         {
1676           status=MagickFalse;
1677           continue;
1678         }
1679       j++;
1680       for (x=0; x < (ssize_t) image->columns; x++)
1681       {
1682         pixel[j++]=p[i];
1683         p+=GetPixelChannels(image);
1684       }
1685       j++;
1686     }
1687     (void) ResetMagickMemory(buffer,0,length*sizeof(*buffer));
1688     for (k=0; k < 4; k++)
1689     {
1690       Hull(X[k],Y[k],image->columns,image->rows,pixel,buffer,1);
1691       Hull(-X[k],-Y[k],image->columns,image->rows,pixel,buffer,1);
1692       Hull(-X[k],-Y[k],image->columns,image->rows,pixel,buffer,-1);
1693       Hull(X[k],Y[k],image->columns,image->rows,pixel,buffer,-1);
1694     }
1695     j=(ssize_t) image->columns+2;
1696     for (y=0; y < (ssize_t) image->rows; y++)
1697     {
1698       MagickBooleanType
1699         sync;
1700
1701       register Quantum
1702         *restrict q;
1703
1704       q=GetCacheViewAuthenticPixels(despeckle_view,0,y,despeckle_image->columns,
1705         1,exception);
1706       if (q == (Quantum *) NULL)
1707         {
1708           status=MagickFalse;
1709           continue;
1710         }
1711       j++;
1712       for (x=0; x < (ssize_t) image->columns; x++)
1713       {
1714         SetPixelChannel(despeckle_image,channel,pixel[j++],q);
1715         q+=GetPixelChannels(despeckle_image);
1716       }
1717       sync=SyncCacheViewAuthenticPixels(despeckle_view,exception);
1718       if (sync == MagickFalse)
1719         status=MagickFalse;
1720       j++;
1721     }
1722     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1723       {
1724         MagickBooleanType
1725           proceed;
1726
1727         proceed=SetImageProgress(image,DespeckleImageTag,(MagickOffsetType) i,
1728           GetPixelChannels(image));
1729         if (proceed == MagickFalse)
1730           status=MagickFalse;
1731       }
1732   }
1733   despeckle_view=DestroyCacheView(despeckle_view);
1734   image_view=DestroyCacheView(image_view);
1735   buffers=(Quantum *) RelinquishMagickMemory(buffers);
1736   pixels=(Quantum *) RelinquishMagickMemory(pixels);
1737   despeckle_image->type=image->type;
1738   if (status == MagickFalse)
1739     despeckle_image=DestroyImage(despeckle_image);
1740   return(despeckle_image);
1741 }
1742 \f
1743 /*
1744 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1745 %                                                                             %
1746 %                                                                             %
1747 %                                                                             %
1748 %     E d g e I m a g e                                                       %
1749 %                                                                             %
1750 %                                                                             %
1751 %                                                                             %
1752 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1753 %
1754 %  EdgeImage() finds edges in an image.  Radius defines the radius of the
1755 %  convolution filter.  Use a radius of 0 and EdgeImage() selects a suitable
1756 %  radius for you.
1757 %
1758 %  The format of the EdgeImage method is:
1759 %
1760 %      Image *EdgeImage(const Image *image,const double radius,
1761 %        const double sigma,ExceptionInfo *exception)
1762 %
1763 %  A description of each parameter follows:
1764 %
1765 %    o image: the image.
1766 %
1767 %    o radius: the radius of the pixel neighborhood.
1768 %
1769 %    o sigma: the standard deviation of the Gaussian, in pixels.
1770 %
1771 %    o exception: return any errors or warnings in this structure.
1772 %
1773 */
1774 MagickExport Image *EdgeImage(const Image *image,const double radius,
1775   const double sigma,ExceptionInfo *exception)
1776 {
1777   Image
1778     *edge_image;
1779
1780   KernelInfo
1781     *kernel_info;
1782
1783   register ssize_t
1784     i;
1785
1786   size_t
1787     width;
1788
1789   ssize_t
1790     j,
1791     u,
1792     v;
1793
1794   assert(image != (const Image *) NULL);
1795   assert(image->signature == MagickSignature);
1796   if (image->debug != MagickFalse)
1797     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1798   assert(exception != (ExceptionInfo *) NULL);
1799   assert(exception->signature == MagickSignature);
1800   width=GetOptimalKernelWidth1D(radius,sigma);
1801   kernel_info=AcquireKernelInfo((const char *) NULL);
1802   if (kernel_info == (KernelInfo *) NULL)
1803     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1804   kernel_info->width=width;
1805   kernel_info->height=width;
1806   kernel_info->values=(MagickRealType *) AcquireAlignedMemory(
1807     kernel_info->width,kernel_info->width*sizeof(*kernel_info->values));
1808   if (kernel_info->values == (MagickRealType *) NULL)
1809     {
1810       kernel_info=DestroyKernelInfo(kernel_info);
1811       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1812     }
1813   j=(ssize_t) kernel_info->width/2;
1814   i=0;
1815   for (v=(-j); v <= j; v++)
1816   {
1817     for (u=(-j); u <= j; u++)
1818     {
1819       kernel_info->values[i]=(-1.0);
1820       i++;
1821     }
1822   }
1823   kernel_info->values[i/2]=(double) (width*width-1.0);
1824   kernel_info->bias=image->bias;   /* FUTURE: User bias on a edge image? */
1825   edge_image=ConvolveImage(image,kernel_info,exception);
1826   kernel_info=DestroyKernelInfo(kernel_info);
1827   return(edge_image);
1828 }
1829 \f
1830 /*
1831 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1832 %                                                                             %
1833 %                                                                             %
1834 %                                                                             %
1835 %     E m b o s s I m a g e                                                   %
1836 %                                                                             %
1837 %                                                                             %
1838 %                                                                             %
1839 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1840 %
1841 %  EmbossImage() returns a grayscale image with a three-dimensional effect.
1842 %  We convolve the image with a Gaussian operator of the given radius and
1843 %  standard deviation (sigma).  For reasonable results, radius should be
1844 %  larger than sigma.  Use a radius of 0 and Emboss() selects a suitable
1845 %  radius for you.
1846 %
1847 %  The format of the EmbossImage method is:
1848 %
1849 %      Image *EmbossImage(const Image *image,const double radius,
1850 %        const double sigma,ExceptionInfo *exception)
1851 %
1852 %  A description of each parameter follows:
1853 %
1854 %    o image: the image.
1855 %
1856 %    o radius: the radius of the pixel neighborhood.
1857 %
1858 %    o sigma: the standard deviation of the Gaussian, in pixels.
1859 %
1860 %    o exception: return any errors or warnings in this structure.
1861 %
1862 */
1863 MagickExport Image *EmbossImage(const Image *image,const double radius,
1864   const double sigma,ExceptionInfo *exception)
1865 {
1866   Image
1867     *emboss_image;
1868
1869   KernelInfo
1870     *kernel_info;
1871
1872   register ssize_t
1873     i;
1874
1875   size_t
1876     width;
1877
1878   ssize_t
1879     j,
1880     k,
1881     u,
1882     v;
1883
1884   assert(image != (const Image *) NULL);
1885   assert(image->signature == MagickSignature);
1886   if (image->debug != MagickFalse)
1887     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1888   assert(exception != (ExceptionInfo *) NULL);
1889   assert(exception->signature == MagickSignature);
1890   width=GetOptimalKernelWidth1D(radius,sigma);
1891   kernel_info=AcquireKernelInfo((const char *) NULL);
1892   if (kernel_info == (KernelInfo *) NULL)
1893     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1894   kernel_info->width=width;
1895   kernel_info->height=width;
1896   kernel_info->values=(MagickRealType *) AcquireAlignedMemory(
1897     kernel_info->width,kernel_info->width*sizeof(*kernel_info->values));
1898   if (kernel_info->values == (MagickRealType *) NULL)
1899     {
1900       kernel_info=DestroyKernelInfo(kernel_info);
1901       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1902     }
1903   j=(ssize_t) kernel_info->width/2;
1904   k=j;
1905   i=0;
1906   for (v=(-j); v <= j; v++)
1907   {
1908     for (u=(-j); u <= j; u++)
1909     {
1910       kernel_info->values[i]=(double) (((u < 0) || (v < 0) ? -8.0 : 8.0)*
1911         exp(-((double) u*u+v*v)/(2.0*MagickSigma*MagickSigma))/
1912         (2.0*MagickPI*MagickSigma*MagickSigma));
1913       if (u != k)
1914         kernel_info->values[i]=0.0;
1915       i++;
1916     }
1917     k--;
1918   }
1919   kernel_info->bias=image->bias;  /* FUTURE: user bias on an edge image */
1920   emboss_image=ConvolveImage(image,kernel_info,exception);
1921   kernel_info=DestroyKernelInfo(kernel_info);
1922   if (emboss_image != (Image *) NULL)
1923     (void) EqualizeImage(emboss_image,exception);
1924   return(emboss_image);
1925 }
1926 \f
1927 /*
1928 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1929 %                                                                             %
1930 %                                                                             %
1931 %                                                                             %
1932 %     G a u s s i a n B l u r I m a g e                                       %
1933 %                                                                             %
1934 %                                                                             %
1935 %                                                                             %
1936 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1937 %
1938 %  GaussianBlurImage() blurs an image.  We convolve the image with a
1939 %  Gaussian operator of the given radius and standard deviation (sigma).
1940 %  For reasonable results, the radius should be larger than sigma.  Use a
1941 %  radius of 0 and GaussianBlurImage() selects a suitable radius for you
1942 %
1943 %  The format of the GaussianBlurImage method is:
1944 %
1945 %      Image *GaussianBlurImage(const Image *image,onst double radius,
1946 %        const double sigma,const double bias,ExceptionInfo *exception)
1947 %
1948 %  A description of each parameter follows:
1949 %
1950 %    o image: the image.
1951 %
1952 %    o radius: the radius of the Gaussian, in pixels, not counting the center
1953 %      pixel.
1954 %
1955 %    o sigma: the standard deviation of the Gaussian, in pixels.
1956 %
1957 %    o bias: the bias.
1958 %
1959 %    o exception: return any errors or warnings in this structure.
1960 %
1961 */
1962 MagickExport Image *GaussianBlurImage(const Image *image,const double radius,
1963   const double sigma,const double bias,ExceptionInfo *exception)
1964 {
1965   Image
1966     *blur_image;
1967
1968   KernelInfo
1969     *kernel_info;
1970
1971   register ssize_t
1972     i;
1973
1974   size_t
1975     width;
1976
1977   ssize_t
1978     j,
1979     u,
1980     v;
1981
1982   assert(image != (const Image *) NULL);
1983   assert(image->signature == MagickSignature);
1984   if (image->debug != MagickFalse)
1985     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1986   assert(exception != (ExceptionInfo *) NULL);
1987   assert(exception->signature == MagickSignature);
1988   width=GetOptimalKernelWidth2D(radius,sigma);
1989   kernel_info=AcquireKernelInfo((const char *) NULL);
1990   if (kernel_info == (KernelInfo *) NULL)
1991     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
1992   (void) ResetMagickMemory(kernel_info,0,sizeof(*kernel_info));
1993   kernel_info->width=width;
1994   kernel_info->height=width;
1995   kernel_info->bias=bias;  /* FUTURE: user bias on Gaussian Blur! non-sense */
1996   kernel_info->signature=MagickSignature;
1997   kernel_info->values=(MagickRealType *) AcquireAlignedMemory(
1998     kernel_info->width,kernel_info->width*sizeof(*kernel_info->values));
1999   if (kernel_info->values == (MagickRealType *) NULL)
2000     {
2001       kernel_info=DestroyKernelInfo(kernel_info);
2002       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2003     }
2004   j=(ssize_t) kernel_info->width/2;
2005   i=0;
2006   for (v=(-j); v <= j; v++)
2007   {
2008     for (u=(-j); u <= j; u++)
2009     {
2010       kernel_info->values[i]=(double) (exp(-((double) u*u+v*v)/(2.0*
2011         MagickSigma*MagickSigma))/(2.0*MagickPI*MagickSigma*MagickSigma));
2012       i++;
2013     }
2014   }
2015   blur_image=ConvolveImage(image,kernel_info,exception);
2016   kernel_info=DestroyKernelInfo(kernel_info);
2017   return(blur_image);
2018 }
2019 \f
2020 /*
2021 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2022 %                                                                             %
2023 %                                                                             %
2024 %                                                                             %
2025 %     M o t i o n B l u r I m a g e                                           %
2026 %                                                                             %
2027 %                                                                             %
2028 %                                                                             %
2029 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2030 %
2031 %  MotionBlurImage() simulates motion blur.  We convolve the image with a
2032 %  Gaussian operator of the given radius and standard deviation (sigma).
2033 %  For reasonable results, radius should be larger than sigma.  Use a
2034 %  radius of 0 and MotionBlurImage() selects a suitable radius for you.
2035 %  Angle gives the angle of the blurring motion.
2036 %
2037 %  Andrew Protano contributed this effect.
2038 %
2039 %  The format of the MotionBlurImage method is:
2040 %
2041 %    Image *MotionBlurImage(const Image *image,const double radius,
2042 %      const double sigma,const double angle,const double bias,
2043 %      ExceptionInfo *exception)
2044 %
2045 %  A description of each parameter follows:
2046 %
2047 %    o image: the image.
2048 %
2049 %    o radius: the radius of the Gaussian, in pixels, not counting
2050 %      the center pixel.
2051 %
2052 %    o sigma: the standard deviation of the Gaussian, in pixels.
2053 %
2054 %    o angle: Apply the effect along this angle.
2055 %
2056 %    o bias: the bias.
2057 %
2058 %    o exception: return any errors or warnings in this structure.
2059 %
2060 */
2061
2062 static double *GetMotionBlurKernel(const size_t width,const double sigma)
2063 {
2064   double
2065     *kernel,
2066     normalize;
2067
2068   register ssize_t
2069     i;
2070
2071   /*
2072    Generate a 1-D convolution kernel.
2073   */
2074   (void) LogMagickEvent(TraceEvent,GetMagickModule(),"...");
2075   kernel=(double *) AcquireQuantumMemory((size_t) width,sizeof(*kernel));
2076   if (kernel == (double *) NULL)
2077     return(kernel);
2078   normalize=0.0;
2079   for (i=0; i < (ssize_t) width; i++)
2080   {
2081     kernel[i]=(double) (exp((-((double) i*i)/(double) (2.0*MagickSigma*
2082       MagickSigma)))/(MagickSQ2PI*MagickSigma));
2083     normalize+=kernel[i];
2084   }
2085   for (i=0; i < (ssize_t) width; i++)
2086     kernel[i]/=normalize;
2087   return(kernel);
2088 }
2089
2090 MagickExport Image *MotionBlurImage(const Image *image,const double radius,
2091   const double sigma,const double angle,const double bias,
2092   ExceptionInfo *exception)
2093 {
2094   CacheView
2095     *blur_view,
2096     *image_view;
2097
2098   double
2099     *kernel;
2100
2101   Image
2102     *blur_image;
2103
2104   MagickBooleanType
2105     status;
2106
2107   MagickOffsetType
2108     progress;
2109
2110   OffsetInfo
2111     *offset;
2112
2113   PointInfo
2114     point;
2115
2116   register ssize_t
2117     i;
2118
2119   size_t
2120     width;
2121
2122   ssize_t
2123     y;
2124
2125   assert(image != (Image *) NULL);
2126   assert(image->signature == MagickSignature);
2127   if (image->debug != MagickFalse)
2128     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2129   assert(exception != (ExceptionInfo *) NULL);
2130   width=GetOptimalKernelWidth1D(radius,sigma);
2131   kernel=GetMotionBlurKernel(width,sigma);
2132   if (kernel == (double *) NULL)
2133     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2134   offset=(OffsetInfo *) AcquireQuantumMemory(width,sizeof(*offset));
2135   if (offset == (OffsetInfo *) NULL)
2136     {
2137       kernel=(double *) RelinquishMagickMemory(kernel);
2138       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2139     }
2140   blur_image=CloneImage(image,image->columns,image->rows,MagickTrue,exception);
2141   if (blur_image == (Image *) NULL)
2142     {
2143       kernel=(double *) RelinquishMagickMemory(kernel);
2144       offset=(OffsetInfo *) RelinquishMagickMemory(offset);
2145       return((Image *) NULL);
2146     }
2147   if (SetImageStorageClass(blur_image,DirectClass,exception) == MagickFalse)
2148     {
2149       kernel=(double *) RelinquishMagickMemory(kernel);
2150       offset=(OffsetInfo *) RelinquishMagickMemory(offset);
2151       blur_image=DestroyImage(blur_image);
2152       return((Image *) NULL);
2153     }
2154   point.x=(double) width*sin(DegreesToRadians(angle));
2155   point.y=(double) width*cos(DegreesToRadians(angle));
2156   for (i=0; i < (ssize_t) width; i++)
2157   {
2158     offset[i].x=(ssize_t) ceil((double) (i*point.y)/hypot(point.x,point.y)-0.5);
2159     offset[i].y=(ssize_t) ceil((double) (i*point.x)/hypot(point.x,point.y)-0.5);
2160   }
2161   /*
2162     Motion blur image.
2163   */
2164   status=MagickTrue;
2165   progress=0;
2166   image_view=AcquireCacheView(image);
2167   blur_view=AcquireCacheView(blur_image);
2168 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2169   #pragma omp parallel for schedule(static,4) shared(progress,status) omp_throttle(1)
2170 #endif
2171   for (y=0; y < (ssize_t) image->rows; y++)
2172   {
2173     register const Quantum
2174       *restrict p;
2175
2176     register Quantum
2177       *restrict q;
2178
2179     register ssize_t
2180       x;
2181
2182     if (status == MagickFalse)
2183       continue;
2184     p=GetCacheViewVirtualPixels(blur_view,0,y,image->columns,1,exception);
2185     q=GetCacheViewAuthenticPixels(blur_view,0,y,blur_image->columns,1,
2186       exception);
2187     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
2188       {
2189         status=MagickFalse;
2190         continue;
2191       }
2192     for (x=0; x < (ssize_t) image->columns; x++)
2193     {
2194       register ssize_t
2195         i;
2196
2197       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2198       {
2199         MagickRealType
2200           alpha,
2201           gamma,
2202           pixel;
2203
2204         PixelChannel
2205           channel;
2206
2207         PixelTrait
2208           blur_traits,
2209           traits;
2210
2211         register const Quantum
2212           *restrict r;
2213
2214         register double
2215           *restrict k;
2216
2217         register ssize_t
2218           j;
2219
2220         channel=GetPixelChannelMapChannel(image,i);
2221         traits=GetPixelChannelMapTraits(image,channel);
2222         blur_traits=GetPixelChannelMapTraits(blur_image,channel);
2223         if ((traits == UndefinedPixelTrait) ||
2224             (blur_traits == UndefinedPixelTrait))
2225           continue;
2226         if ((blur_traits & CopyPixelTrait) != 0)
2227           {
2228             SetPixelChannel(blur_image,channel,p[i],q);
2229             continue;
2230           }
2231         k=kernel;
2232         pixel=bias;
2233         if ((blur_traits & BlendPixelTrait) == 0)
2234           {
2235             for (j=0; j < (ssize_t) width; j++)
2236             {
2237               r=GetCacheViewVirtualPixels(image_view,x+offset[j].x,y+
2238                 offset[j].y,1,1,exception);
2239               if (r == (const Quantum *) NULL)
2240                 {
2241                   status=MagickFalse;
2242                   continue;
2243                 }
2244               pixel+=(*k)*r[i];
2245               k++;
2246             }
2247             SetPixelChannel(blur_image,channel,ClampToQuantum(pixel),q);
2248             continue;
2249           }
2250         alpha=0.0;
2251         gamma=0.0;
2252         for (j=0; j < (ssize_t) width; j++)
2253         {
2254           r=GetCacheViewVirtualPixels(image_view,x+offset[j].x,y+offset[j].y,1,
2255             1,exception);
2256           if (r == (const Quantum *) NULL)
2257             {
2258               status=MagickFalse;
2259               continue;
2260             }
2261           alpha=(MagickRealType) (QuantumScale*GetPixelAlpha(image,r));
2262           pixel+=(*k)*alpha*r[i];
2263           gamma+=(*k)*alpha;
2264           k++;
2265         }
2266         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
2267         SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
2268       }
2269       p+=GetPixelChannels(image);
2270       q+=GetPixelChannels(blur_image);
2271     }
2272     if (SyncCacheViewAuthenticPixels(blur_view,exception) == MagickFalse)
2273       status=MagickFalse;
2274     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2275       {
2276         MagickBooleanType
2277           proceed;
2278
2279 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2280   #pragma omp critical (MagickCore_MotionBlurImage)
2281 #endif
2282         proceed=SetImageProgress(image,BlurImageTag,progress++,image->rows);
2283         if (proceed == MagickFalse)
2284           status=MagickFalse;
2285       }
2286   }
2287   blur_view=DestroyCacheView(blur_view);
2288   image_view=DestroyCacheView(image_view);
2289   kernel=(double *) RelinquishMagickMemory(kernel);
2290   offset=(OffsetInfo *) RelinquishMagickMemory(offset);
2291   if (status == MagickFalse)
2292     blur_image=DestroyImage(blur_image);
2293   return(blur_image);
2294 }
2295 \f
2296 /*
2297 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2298 %                                                                             %
2299 %                                                                             %
2300 %                                                                             %
2301 %     P r e v i e w I m a g e                                                 %
2302 %                                                                             %
2303 %                                                                             %
2304 %                                                                             %
2305 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2306 %
2307 %  PreviewImage() tiles 9 thumbnails of the specified image with an image
2308 %  processing operation applied with varying parameters.  This may be helpful
2309 %  pin-pointing an appropriate parameter for a particular image processing
2310 %  operation.
2311 %
2312 %  The format of the PreviewImages method is:
2313 %
2314 %      Image *PreviewImages(const Image *image,const PreviewType preview,
2315 %        ExceptionInfo *exception)
2316 %
2317 %  A description of each parameter follows:
2318 %
2319 %    o image: the image.
2320 %
2321 %    o preview: the image processing operation.
2322 %
2323 %    o exception: return any errors or warnings in this structure.
2324 %
2325 */
2326 MagickExport Image *PreviewImage(const Image *image,const PreviewType preview,
2327   ExceptionInfo *exception)
2328 {
2329 #define NumberTiles  9
2330 #define PreviewImageTag  "Preview/Image"
2331 #define DefaultPreviewGeometry  "204x204+10+10"
2332
2333   char
2334     factor[MaxTextExtent],
2335     label[MaxTextExtent];
2336
2337   double
2338     degrees,
2339     gamma,
2340     percentage,
2341     radius,
2342     sigma,
2343     threshold;
2344
2345   Image
2346     *images,
2347     *montage_image,
2348     *preview_image,
2349     *thumbnail;
2350
2351   ImageInfo
2352     *preview_info;
2353
2354   MagickBooleanType
2355     proceed;
2356
2357   MontageInfo
2358     *montage_info;
2359
2360   QuantizeInfo
2361     quantize_info;
2362
2363   RectangleInfo
2364     geometry;
2365
2366   register ssize_t
2367     i,
2368     x;
2369
2370   size_t
2371     colors;
2372
2373   ssize_t
2374     y;
2375
2376   /*
2377     Open output image file.
2378   */
2379   assert(image != (Image *) NULL);
2380   assert(image->signature == MagickSignature);
2381   if (image->debug != MagickFalse)
2382     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2383   colors=2;
2384   degrees=0.0;
2385   gamma=(-0.2f);
2386   preview_info=AcquireImageInfo();
2387   SetGeometry(image,&geometry);
2388   (void) ParseMetaGeometry(DefaultPreviewGeometry,&geometry.x,&geometry.y,
2389     &geometry.width,&geometry.height);
2390   images=NewImageList();
2391   percentage=12.5;
2392   GetQuantizeInfo(&quantize_info);
2393   radius=0.0;
2394   sigma=1.0;
2395   threshold=0.0;
2396   x=0;
2397   y=0;
2398   for (i=0; i < NumberTiles; i++)
2399   {
2400     thumbnail=ThumbnailImage(image,geometry.width,geometry.height,exception);
2401     if (thumbnail == (Image *) NULL)
2402       break;
2403     (void) SetImageProgressMonitor(thumbnail,(MagickProgressMonitor) NULL,
2404       (void *) NULL);
2405     (void) SetImageProperty(thumbnail,"label",DefaultTileLabel,exception);
2406     if (i == (NumberTiles/2))
2407       {
2408         (void) QueryColorCompliance("#dfdfdf",AllCompliance,
2409           &thumbnail->matte_color,exception);
2410         AppendImageToList(&images,thumbnail);
2411         continue;
2412       }
2413     switch (preview)
2414     {
2415       case RotatePreview:
2416       {
2417         degrees+=45.0;
2418         preview_image=RotateImage(thumbnail,degrees,exception);
2419         (void) FormatLocaleString(label,MaxTextExtent,"rotate %g",degrees);
2420         break;
2421       }
2422       case ShearPreview:
2423       {
2424         degrees+=5.0;
2425         preview_image=ShearImage(thumbnail,degrees,degrees,exception);
2426         (void) FormatLocaleString(label,MaxTextExtent,"shear %gx%g",
2427           degrees,2.0*degrees);
2428         break;
2429       }
2430       case RollPreview:
2431       {
2432         x=(ssize_t) ((i+1)*thumbnail->columns)/NumberTiles;
2433         y=(ssize_t) ((i+1)*thumbnail->rows)/NumberTiles;
2434         preview_image=RollImage(thumbnail,x,y,exception);
2435         (void) FormatLocaleString(label,MaxTextExtent,"roll %+.20gx%+.20g",
2436           (double) x,(double) y);
2437         break;
2438       }
2439       case HuePreview:
2440       {
2441         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2442         if (preview_image == (Image *) NULL)
2443           break;
2444         (void) FormatLocaleString(factor,MaxTextExtent,"100,100,%g",
2445           2.0*percentage);
2446         (void) ModulateImage(preview_image,factor,exception);
2447         (void) FormatLocaleString(label,MaxTextExtent,"modulate %s",factor);
2448         break;
2449       }
2450       case SaturationPreview:
2451       {
2452         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2453         if (preview_image == (Image *) NULL)
2454           break;
2455         (void) FormatLocaleString(factor,MaxTextExtent,"100,%g",
2456           2.0*percentage);
2457         (void) ModulateImage(preview_image,factor,exception);
2458         (void) FormatLocaleString(label,MaxTextExtent,"modulate %s",factor);
2459         break;
2460       }
2461       case BrightnessPreview:
2462       {
2463         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2464         if (preview_image == (Image *) NULL)
2465           break;
2466         (void) FormatLocaleString(factor,MaxTextExtent,"%g",2.0*percentage);
2467         (void) ModulateImage(preview_image,factor,exception);
2468         (void) FormatLocaleString(label,MaxTextExtent,"modulate %s",factor);
2469         break;
2470       }
2471       case GammaPreview:
2472       default:
2473       {
2474         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2475         if (preview_image == (Image *) NULL)
2476           break;
2477         gamma+=0.4f;
2478         (void) GammaImage(preview_image,gamma,exception);
2479         (void) FormatLocaleString(label,MaxTextExtent,"gamma %g",gamma);
2480         break;
2481       }
2482       case SpiffPreview:
2483       {
2484         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2485         if (preview_image != (Image *) NULL)
2486           for (x=0; x < i; x++)
2487             (void) ContrastImage(preview_image,MagickTrue,exception);
2488         (void) FormatLocaleString(label,MaxTextExtent,"contrast (%.20g)",
2489           (double) i+1);
2490         break;
2491       }
2492       case DullPreview:
2493       {
2494         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2495         if (preview_image == (Image *) NULL)
2496           break;
2497         for (x=0; x < i; x++)
2498           (void) ContrastImage(preview_image,MagickFalse,exception);
2499         (void) FormatLocaleString(label,MaxTextExtent,"+contrast (%.20g)",
2500           (double) i+1);
2501         break;
2502       }
2503       case GrayscalePreview:
2504       {
2505         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2506         if (preview_image == (Image *) NULL)
2507           break;
2508         colors<<=1;
2509         quantize_info.number_colors=colors;
2510         quantize_info.colorspace=GRAYColorspace;
2511         (void) QuantizeImage(&quantize_info,preview_image,exception);
2512         (void) FormatLocaleString(label,MaxTextExtent,
2513           "-colorspace gray -colors %.20g",(double) colors);
2514         break;
2515       }
2516       case QuantizePreview:
2517       {
2518         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2519         if (preview_image == (Image *) NULL)
2520           break;
2521         colors<<=1;
2522         quantize_info.number_colors=colors;
2523         (void) QuantizeImage(&quantize_info,preview_image,exception);
2524         (void) FormatLocaleString(label,MaxTextExtent,"colors %.20g",(double)
2525           colors);
2526         break;
2527       }
2528       case DespecklePreview:
2529       {
2530         for (x=0; x < (i-1); x++)
2531         {
2532           preview_image=DespeckleImage(thumbnail,exception);
2533           if (preview_image == (Image *) NULL)
2534             break;
2535           thumbnail=DestroyImage(thumbnail);
2536           thumbnail=preview_image;
2537         }
2538         preview_image=DespeckleImage(thumbnail,exception);
2539         if (preview_image == (Image *) NULL)
2540           break;
2541         (void) FormatLocaleString(label,MaxTextExtent,"despeckle (%.20g)",
2542           (double) i+1);
2543         break;
2544       }
2545       case ReduceNoisePreview:
2546       {
2547         preview_image=StatisticImage(thumbnail,NonpeakStatistic,(size_t) radius,
2548           (size_t) radius,exception);
2549         (void) FormatLocaleString(label,MaxTextExtent,"noise %g",radius);
2550         break;
2551       }
2552       case AddNoisePreview:
2553       {
2554         switch ((int) i)
2555         {
2556           case 0:
2557           {
2558             (void) CopyMagickString(factor,"uniform",MaxTextExtent);
2559             break;
2560           }
2561           case 1:
2562           {
2563             (void) CopyMagickString(factor,"gaussian",MaxTextExtent);
2564             break;
2565           }
2566           case 2:
2567           {
2568             (void) CopyMagickString(factor,"multiplicative",MaxTextExtent);
2569             break;
2570           }
2571           case 3:
2572           {
2573             (void) CopyMagickString(factor,"impulse",MaxTextExtent);
2574             break;
2575           }
2576           case 4:
2577           {
2578             (void) CopyMagickString(factor,"laplacian",MaxTextExtent);
2579             break;
2580           }
2581           case 5:
2582           {
2583             (void) CopyMagickString(factor,"Poisson",MaxTextExtent);
2584             break;
2585           }
2586           default:
2587           {
2588             (void) CopyMagickString(thumbnail->magick,"NULL",MaxTextExtent);
2589             break;
2590           }
2591         }
2592         preview_image=StatisticImage(thumbnail,NonpeakStatistic,(size_t) i,
2593           (size_t) i,exception);
2594         (void) FormatLocaleString(label,MaxTextExtent,"+noise %s",factor);
2595         break;
2596       }
2597       case SharpenPreview:
2598       {
2599         /* FUTURE: user bias on sharpen! This is non-sensical! */
2600         preview_image=SharpenImage(thumbnail,radius,sigma,image->bias,
2601           exception);
2602         (void) FormatLocaleString(label,MaxTextExtent,"sharpen %gx%g",
2603           radius,sigma);
2604         break;
2605       }
2606       case BlurPreview:
2607       {
2608         /* FUTURE: user bias on blur! This is non-sensical! */
2609         preview_image=BlurImage(thumbnail,radius,sigma,image->bias,exception);
2610         (void) FormatLocaleString(label,MaxTextExtent,"blur %gx%g",radius,
2611           sigma);
2612         break;
2613       }
2614       case ThresholdPreview:
2615       {
2616         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2617         if (preview_image == (Image *) NULL)
2618           break;
2619         (void) BilevelImage(thumbnail,(double) (percentage*((MagickRealType)
2620           QuantumRange+1.0))/100.0,exception);
2621         (void) FormatLocaleString(label,MaxTextExtent,"threshold %g",
2622           (double) (percentage*((MagickRealType) QuantumRange+1.0))/100.0);
2623         break;
2624       }
2625       case EdgeDetectPreview:
2626       {
2627         preview_image=EdgeImage(thumbnail,radius,sigma,exception);
2628         (void) FormatLocaleString(label,MaxTextExtent,"edge %g",radius);
2629         break;
2630       }
2631       case SpreadPreview:
2632       {
2633         preview_image=SpreadImage(thumbnail,radius,thumbnail->interpolate,
2634           exception);
2635         (void) FormatLocaleString(label,MaxTextExtent,"spread %g",
2636           radius+0.5);
2637         break;
2638       }
2639       case SolarizePreview:
2640       {
2641         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2642         if (preview_image == (Image *) NULL)
2643           break;
2644         (void) SolarizeImage(preview_image,(double) QuantumRange*
2645           percentage/100.0,exception);
2646         (void) FormatLocaleString(label,MaxTextExtent,"solarize %g",
2647           (QuantumRange*percentage)/100.0);
2648         break;
2649       }
2650       case ShadePreview:
2651       {
2652         degrees+=10.0;
2653         preview_image=ShadeImage(thumbnail,MagickTrue,degrees,degrees,
2654           exception);
2655         (void) FormatLocaleString(label,MaxTextExtent,"shade %gx%g",
2656           degrees,degrees);
2657         break;
2658       }
2659       case RaisePreview:
2660       {
2661         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2662         if (preview_image == (Image *) NULL)
2663           break;
2664         geometry.width=(size_t) (2*i+2);
2665         geometry.height=(size_t) (2*i+2);
2666         geometry.x=i/2;
2667         geometry.y=i/2;
2668         (void) RaiseImage(preview_image,&geometry,MagickTrue,exception);
2669         (void) FormatLocaleString(label,MaxTextExtent,
2670           "raise %.20gx%.20g%+.20g%+.20g",(double) geometry.width,(double)
2671           geometry.height,(double) geometry.x,(double) geometry.y);
2672         break;
2673       }
2674       case SegmentPreview:
2675       {
2676         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2677         if (preview_image == (Image *) NULL)
2678           break;
2679         threshold+=0.4f;
2680         (void) SegmentImage(preview_image,RGBColorspace,MagickFalse,threshold,
2681           threshold,exception);
2682         (void) FormatLocaleString(label,MaxTextExtent,"segment %gx%g",
2683           threshold,threshold);
2684         break;
2685       }
2686       case SwirlPreview:
2687       {
2688         preview_image=SwirlImage(thumbnail,degrees,image->interpolate,
2689           exception);
2690         (void) FormatLocaleString(label,MaxTextExtent,"swirl %g",degrees);
2691         degrees+=45.0;
2692         break;
2693       }
2694       case ImplodePreview:
2695       {
2696         degrees+=0.1f;
2697         preview_image=ImplodeImage(thumbnail,degrees,image->interpolate,
2698           exception);
2699         (void) FormatLocaleString(label,MaxTextExtent,"implode %g",degrees);
2700         break;
2701       }
2702       case WavePreview:
2703       {
2704         degrees+=5.0f;
2705         preview_image=WaveImage(thumbnail,0.5*degrees,2.0*degrees,
2706           image->interpolate,exception);
2707         (void) FormatLocaleString(label,MaxTextExtent,"wave %gx%g",
2708           0.5*degrees,2.0*degrees);
2709         break;
2710       }
2711       case OilPaintPreview:
2712       {
2713         preview_image=OilPaintImage(thumbnail,(double) radius,(double) sigma,
2714           exception);
2715         (void) FormatLocaleString(label,MaxTextExtent,"charcoal %gx%g",
2716           radius,sigma);
2717         break;
2718       }
2719       case CharcoalDrawingPreview:
2720       {
2721         /* FUTURE: user bias on charcoal! This is non-sensical! */
2722         preview_image=CharcoalImage(thumbnail,(double) radius,(double) sigma,
2723           image->bias,exception);
2724         (void) FormatLocaleString(label,MaxTextExtent,"charcoal %gx%g",
2725           radius,sigma);
2726         break;
2727       }
2728       case JPEGPreview:
2729       {
2730         char
2731           filename[MaxTextExtent];
2732
2733         int
2734           file;
2735
2736         MagickBooleanType
2737           status;
2738
2739         preview_image=CloneImage(thumbnail,0,0,MagickTrue,exception);
2740         if (preview_image == (Image *) NULL)
2741           break;
2742         preview_info->quality=(size_t) percentage;
2743         (void) FormatLocaleString(factor,MaxTextExtent,"%.20g",(double)
2744           preview_info->quality);
2745         file=AcquireUniqueFileResource(filename);
2746         if (file != -1)
2747           file=close(file)-1;
2748         (void) FormatLocaleString(preview_image->filename,MaxTextExtent,
2749           "jpeg:%s",filename);
2750         status=WriteImage(preview_info,preview_image,exception);
2751         if (status != MagickFalse)
2752           {
2753             Image
2754               *quality_image;
2755
2756             (void) CopyMagickString(preview_info->filename,
2757               preview_image->filename,MaxTextExtent);
2758             quality_image=ReadImage(preview_info,exception);
2759             if (quality_image != (Image *) NULL)
2760               {
2761                 preview_image=DestroyImage(preview_image);
2762                 preview_image=quality_image;
2763               }
2764           }
2765         (void) RelinquishUniqueFileResource(preview_image->filename);
2766         if ((GetBlobSize(preview_image)/1024) >= 1024)
2767           (void) FormatLocaleString(label,MaxTextExtent,"quality %s\n%gmb ",
2768             factor,(double) ((MagickOffsetType) GetBlobSize(preview_image))/
2769             1024.0/1024.0);
2770         else
2771           if (GetBlobSize(preview_image) >= 1024)
2772             (void) FormatLocaleString(label,MaxTextExtent,
2773               "quality %s\n%gkb ",factor,(double) ((MagickOffsetType)
2774               GetBlobSize(preview_image))/1024.0);
2775           else
2776             (void) FormatLocaleString(label,MaxTextExtent,"quality %s\n%.20gb ",
2777               factor,(double) ((MagickOffsetType) GetBlobSize(thumbnail)));
2778         break;
2779       }
2780     }
2781     thumbnail=DestroyImage(thumbnail);
2782     percentage+=12.5;
2783     radius+=0.5;
2784     sigma+=0.25;
2785     if (preview_image == (Image *) NULL)
2786       break;
2787     (void) DeleteImageProperty(preview_image,"label");
2788     (void) SetImageProperty(preview_image,"label",label,exception);
2789     AppendImageToList(&images,preview_image);
2790     proceed=SetImageProgress(image,PreviewImageTag,(MagickOffsetType) i,
2791       NumberTiles);
2792     if (proceed == MagickFalse)
2793       break;
2794   }
2795   if (images == (Image *) NULL)
2796     {
2797       preview_info=DestroyImageInfo(preview_info);
2798       return((Image *) NULL);
2799     }
2800   /*
2801     Create the montage.
2802   */
2803   montage_info=CloneMontageInfo(preview_info,(MontageInfo *) NULL);
2804   (void) CopyMagickString(montage_info->filename,image->filename,MaxTextExtent);
2805   montage_info->shadow=MagickTrue;
2806   (void) CloneString(&montage_info->tile,"3x3");
2807   (void) CloneString(&montage_info->geometry,DefaultPreviewGeometry);
2808   (void) CloneString(&montage_info->frame,DefaultTileFrame);
2809   montage_image=MontageImages(images,montage_info,exception);
2810   montage_info=DestroyMontageInfo(montage_info);
2811   images=DestroyImageList(images);
2812   if (montage_image == (Image *) NULL)
2813     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2814   if (montage_image->montage != (char *) NULL)
2815     {
2816       /*
2817         Free image directory.
2818       */
2819       montage_image->montage=(char *) RelinquishMagickMemory(
2820         montage_image->montage);
2821       if (image->directory != (char *) NULL)
2822         montage_image->directory=(char *) RelinquishMagickMemory(
2823           montage_image->directory);
2824     }
2825   preview_info=DestroyImageInfo(preview_info);
2826   return(montage_image);
2827 }
2828 \f
2829 /*
2830 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2831 %                                                                             %
2832 %                                                                             %
2833 %                                                                             %
2834 %     R a d i a l B l u r I m a g e                                           %
2835 %                                                                             %
2836 %                                                                             %
2837 %                                                                             %
2838 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2839 %
2840 %  RadialBlurImage() applies a radial blur to the image.
2841 %
2842 %  Andrew Protano contributed this effect.
2843 %
2844 %  The format of the RadialBlurImage method is:
2845 %
2846 %    Image *RadialBlurImage(const Image *image,const double angle,
2847 %      const double blur,ExceptionInfo *exception)
2848 %
2849 %  A description of each parameter follows:
2850 %
2851 %    o image: the image.
2852 %
2853 %    o angle: the angle of the radial blur.
2854 %
2855 %    o blur: the blur.
2856 %
2857 %    o exception: return any errors or warnings in this structure.
2858 %
2859 */
2860 MagickExport Image *RadialBlurImage(const Image *image,const double angle,
2861   const double bias,ExceptionInfo *exception)
2862 {
2863   CacheView
2864     *blur_view,
2865     *image_view;
2866
2867   Image
2868     *blur_image;
2869
2870   MagickBooleanType
2871     status;
2872
2873   MagickOffsetType
2874     progress;
2875
2876   MagickRealType
2877     blur_radius,
2878     *cos_theta,
2879     offset,
2880     *sin_theta,
2881     theta;
2882
2883   PointInfo
2884     blur_center;
2885
2886   register ssize_t
2887     i;
2888
2889   size_t
2890     n;
2891
2892   ssize_t
2893     y;
2894
2895   /*
2896     Allocate blur image.
2897   */
2898   assert(image != (Image *) NULL);
2899   assert(image->signature == MagickSignature);
2900   if (image->debug != MagickFalse)
2901     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2902   assert(exception != (ExceptionInfo *) NULL);
2903   assert(exception->signature == MagickSignature);
2904   blur_image=CloneImage(image,image->columns,image->rows,MagickTrue,exception);
2905   if (blur_image == (Image *) NULL)
2906     return((Image *) NULL);
2907   if (SetImageStorageClass(blur_image,DirectClass,exception) == MagickFalse)
2908     {
2909       blur_image=DestroyImage(blur_image);
2910       return((Image *) NULL);
2911     }
2912   blur_center.x=(double) image->columns/2.0;
2913   blur_center.y=(double) image->rows/2.0;
2914   blur_radius=hypot(blur_center.x,blur_center.y);
2915   n=(size_t) fabs(4.0*DegreesToRadians(angle)*sqrt((double) blur_radius)+2UL);
2916   theta=DegreesToRadians(angle)/(MagickRealType) (n-1);
2917   cos_theta=(MagickRealType *) AcquireQuantumMemory((size_t) n,
2918     sizeof(*cos_theta));
2919   sin_theta=(MagickRealType *) AcquireQuantumMemory((size_t) n,
2920     sizeof(*sin_theta));
2921   if ((cos_theta == (MagickRealType *) NULL) ||
2922       (sin_theta == (MagickRealType *) NULL))
2923     {
2924       blur_image=DestroyImage(blur_image);
2925       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
2926     }
2927   offset=theta*(MagickRealType) (n-1)/2.0;
2928   for (i=0; i < (ssize_t) n; i++)
2929   {
2930     cos_theta[i]=cos((double) (theta*i-offset));
2931     sin_theta[i]=sin((double) (theta*i-offset));
2932   }
2933   /*
2934     Radial blur image.
2935   */
2936   status=MagickTrue;
2937   progress=0;
2938   image_view=AcquireCacheView(image);
2939   blur_view=AcquireCacheView(blur_image);
2940 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2941   #pragma omp parallel for schedule(static,4) shared(progress,status)
2942 #endif
2943   for (y=0; y < (ssize_t) image->rows; y++)
2944   {
2945     register const Quantum
2946       *restrict p;
2947
2948     register Quantum
2949       *restrict q;
2950
2951     register ssize_t
2952       x;
2953
2954     if (status == MagickFalse)
2955       continue;
2956     p=GetCacheViewVirtualPixels(blur_view,0,y,image->columns,1,exception);
2957     q=GetCacheViewAuthenticPixels(blur_view,0,y,blur_image->columns,1,
2958       exception);
2959     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
2960       {
2961         status=MagickFalse;
2962         continue;
2963       }
2964     for (x=0; x < (ssize_t) image->columns; x++)
2965     {
2966       MagickRealType
2967         radius;
2968
2969       PointInfo
2970         center;
2971
2972       register ssize_t
2973         i;
2974
2975       size_t
2976         step;
2977
2978       center.x=(double) x-blur_center.x;
2979       center.y=(double) y-blur_center.y;
2980       radius=hypot((double) center.x,center.y);
2981       if (radius == 0)
2982         step=1;
2983       else
2984         {
2985           step=(size_t) (blur_radius/radius);
2986           if (step == 0)
2987             step=1;
2988           else
2989             if (step >= n)
2990               step=n-1;
2991         }
2992       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2993       {
2994         MagickRealType
2995           gamma,
2996           pixel;
2997
2998         PixelChannel
2999           channel;
3000
3001         PixelTrait
3002           blur_traits,
3003           traits;
3004
3005         register const Quantum
3006           *restrict r;
3007
3008         register ssize_t
3009           j;
3010
3011         channel=GetPixelChannelMapChannel(image,i);
3012         traits=GetPixelChannelMapTraits(image,channel);
3013         blur_traits=GetPixelChannelMapTraits(blur_image,channel);
3014         if ((traits == UndefinedPixelTrait) ||
3015             (blur_traits == UndefinedPixelTrait))
3016           continue;
3017         if ((blur_traits & CopyPixelTrait) != 0)
3018           {
3019             SetPixelChannel(blur_image,channel,p[i],q);
3020             continue;
3021           }
3022         gamma=0.0;
3023         pixel=bias;
3024         if ((blur_traits & BlendPixelTrait) == 0)
3025           {
3026             for (j=0; j < (ssize_t) n; j+=(ssize_t) step)
3027             {
3028               r=GetCacheViewVirtualPixels(image_view, (ssize_t) (blur_center.x+
3029                 center.x*cos_theta[j]-center.y*sin_theta[j]+0.5),(ssize_t)
3030                 (blur_center.y+center.x*sin_theta[j]+center.y*cos_theta[j]+0.5),
3031                 1,1,exception);
3032               if (r == (const Quantum *) NULL)
3033                 {
3034                   status=MagickFalse;
3035                   continue;
3036                 }
3037               pixel+=r[i];
3038               gamma++;
3039             }
3040             gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
3041             SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
3042             continue;
3043           }
3044         for (j=0; j < (ssize_t) n; j+=(ssize_t) step)
3045         {
3046           r=GetCacheViewVirtualPixels(image_view, (ssize_t) (blur_center.x+
3047             center.x*cos_theta[j]-center.y*sin_theta[j]+0.5),(ssize_t)
3048             (blur_center.y+center.x*sin_theta[j]+center.y*cos_theta[j]+0.5),
3049             1,1,exception);
3050           if (r == (const Quantum *) NULL)
3051             {
3052               status=MagickFalse;
3053               continue;
3054             }
3055           pixel+=GetPixelAlpha(image,r)*r[i];
3056           gamma+=GetPixelAlpha(image,r);
3057         }
3058         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
3059         SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
3060       }
3061       p+=GetPixelChannels(image);
3062       q+=GetPixelChannels(blur_image);
3063     }
3064     if (SyncCacheViewAuthenticPixels(blur_view,exception) == MagickFalse)
3065       status=MagickFalse;
3066     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3067       {
3068         MagickBooleanType
3069           proceed;
3070
3071 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3072   #pragma omp critical (MagickCore_RadialBlurImage)
3073 #endif
3074         proceed=SetImageProgress(image,BlurImageTag,progress++,image->rows);
3075         if (proceed == MagickFalse)
3076           status=MagickFalse;
3077       }
3078   }
3079   blur_view=DestroyCacheView(blur_view);
3080   image_view=DestroyCacheView(image_view);
3081   cos_theta=(MagickRealType *) RelinquishMagickMemory(cos_theta);
3082   sin_theta=(MagickRealType *) RelinquishMagickMemory(sin_theta);
3083   if (status == MagickFalse)
3084     blur_image=DestroyImage(blur_image);
3085   return(blur_image);
3086 }
3087 \f
3088 /*
3089 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3090 %                                                                             %
3091 %                                                                             %
3092 %                                                                             %
3093 %     S e l e c t i v e B l u r I m a g e                                     %
3094 %                                                                             %
3095 %                                                                             %
3096 %                                                                             %
3097 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3098 %
3099 %  SelectiveBlurImage() selectively blur pixels within a contrast threshold.
3100 %  It is similar to the unsharpen mask that sharpens everything with contrast
3101 %  above a certain threshold.
3102 %
3103 %  The format of the SelectiveBlurImage method is:
3104 %
3105 %      Image *SelectiveBlurImage(const Image *image,const double radius,
3106 %        const double sigma,const double threshold,const double bias,
3107 %        ExceptionInfo *exception)
3108 %
3109 %  A description of each parameter follows:
3110 %
3111 %    o image: the image.
3112 %
3113 %    o radius: the radius of the Gaussian, in pixels, not counting the center
3114 %      pixel.
3115 %
3116 %    o sigma: the standard deviation of the Gaussian, in pixels.
3117 %
3118 %    o threshold: only pixels within this contrast threshold are included
3119 %      in the blur operation.
3120 %
3121 %    o bias: the bias.
3122 %
3123 %    o exception: return any errors or warnings in this structure.
3124 %
3125 */
3126 MagickExport Image *SelectiveBlurImage(const Image *image,const double radius,
3127   const double sigma,const double threshold,const double bias,
3128   ExceptionInfo *exception)
3129 {
3130 #define SelectiveBlurImageTag  "SelectiveBlur/Image"
3131
3132   CacheView
3133     *blur_view,
3134     *image_view;
3135
3136   double
3137     *kernel;
3138
3139   Image
3140     *blur_image;
3141
3142   MagickBooleanType
3143     status;
3144
3145   MagickOffsetType
3146     progress;
3147
3148   register ssize_t
3149     i;
3150
3151   size_t
3152     width;
3153
3154   ssize_t
3155     center,
3156     j,
3157     u,
3158     v,
3159     y;
3160
3161   /*
3162     Initialize blur image attributes.
3163   */
3164   assert(image != (Image *) NULL);
3165   assert(image->signature == MagickSignature);
3166   if (image->debug != MagickFalse)
3167     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3168   assert(exception != (ExceptionInfo *) NULL);
3169   assert(exception->signature == MagickSignature);
3170   width=GetOptimalKernelWidth1D(radius,sigma);
3171   kernel=(double *) AcquireQuantumMemory((size_t) width,width*sizeof(*kernel));
3172   if (kernel == (double *) NULL)
3173     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
3174   j=(ssize_t) width/2;
3175   i=0;
3176   for (v=(-j); v <= j; v++)
3177   {
3178     for (u=(-j); u <= j; u++)
3179       kernel[i++]=(double) (exp(-((double) u*u+v*v)/(2.0*MagickSigma*
3180         MagickSigma))/(2.0*MagickPI*MagickSigma*MagickSigma));
3181   }
3182   if (image->debug != MagickFalse)
3183     {
3184       char
3185         format[MaxTextExtent],
3186         *message;
3187
3188       register const double
3189         *k;
3190
3191       ssize_t
3192         u,
3193         v;
3194
3195       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
3196         "  SelectiveBlurImage with %.20gx%.20g kernel:",(double) width,(double)
3197         width);
3198       message=AcquireString("");
3199       k=kernel;
3200       for (v=0; v < (ssize_t) width; v++)
3201       {
3202         *message='\0';
3203         (void) FormatLocaleString(format,MaxTextExtent,"%.20g: ",(double) v);
3204         (void) ConcatenateString(&message,format);
3205         for (u=0; u < (ssize_t) width; u++)
3206         {
3207           (void) FormatLocaleString(format,MaxTextExtent,"%+f ",*k++);
3208           (void) ConcatenateString(&message,format);
3209         }
3210         (void) LogMagickEvent(TransformEvent,GetMagickModule(),"%s",message);
3211       }
3212       message=DestroyString(message);
3213     }
3214   blur_image=CloneImage(image,image->columns,image->rows,MagickTrue,exception);
3215   if (blur_image == (Image *) NULL)
3216     return((Image *) NULL);
3217   if (SetImageStorageClass(blur_image,DirectClass,exception) == MagickFalse)
3218     {
3219       blur_image=DestroyImage(blur_image);
3220       return((Image *) NULL);
3221     }
3222   /*
3223     Threshold blur image.
3224   */
3225   status=MagickTrue;
3226   progress=0;
3227   center=(ssize_t) (GetPixelChannels(image)*(image->columns+width)*(width/2L)+
3228     GetPixelChannels(image)*(width/2L));
3229   image_view=AcquireCacheView(image);
3230   blur_view=AcquireCacheView(blur_image);
3231 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3232   #pragma omp parallel for schedule(static,4) shared(progress,status)
3233 #endif
3234   for (y=0; y < (ssize_t) image->rows; y++)
3235   {
3236     double
3237       contrast;
3238
3239     MagickBooleanType
3240       sync;
3241
3242     register const Quantum
3243       *restrict p;
3244
3245     register Quantum
3246       *restrict q;
3247
3248     register ssize_t
3249       x;
3250
3251     if (status == MagickFalse)
3252       continue;
3253     p=GetCacheViewVirtualPixels(image_view,-((ssize_t) width/2L),y-(ssize_t)
3254       (width/2L),image->columns+width,width,exception);
3255     q=GetCacheViewAuthenticPixels(blur_view,0,y,blur_image->columns,1,
3256       exception);
3257     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
3258       {
3259         status=MagickFalse;
3260         continue;
3261       }
3262     for (x=0; x < (ssize_t) image->columns; x++)
3263     {
3264       register ssize_t
3265         i;
3266
3267       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
3268       {
3269         MagickRealType
3270           alpha,
3271           gamma,
3272           intensity,
3273           pixel;
3274
3275         PixelChannel
3276           channel;
3277
3278         PixelTrait
3279           blur_traits,
3280           traits;
3281
3282         register const double
3283           *restrict k;
3284
3285         register const Quantum
3286           *restrict pixels;
3287
3288         register ssize_t
3289           u;
3290
3291         ssize_t
3292           v;
3293
3294         channel=GetPixelChannelMapChannel(image,i);
3295         traits=GetPixelChannelMapTraits(image,channel);
3296         blur_traits=GetPixelChannelMapTraits(blur_image,channel);
3297         if ((traits == UndefinedPixelTrait) ||
3298             (blur_traits == UndefinedPixelTrait))
3299           continue;
3300         if ((blur_traits & CopyPixelTrait) != 0)
3301           {
3302             SetPixelChannel(blur_image,channel,p[center+i],q);
3303             continue;
3304           }
3305         k=kernel;
3306         pixel=bias;
3307         pixels=p;
3308         intensity=(MagickRealType) GetPixelIntensity(image,p+center);
3309         gamma=0.0;
3310         if ((blur_traits & BlendPixelTrait) == 0)
3311           {
3312             for (v=0; v < (ssize_t) width; v++)
3313             {
3314               for (u=0; u < (ssize_t) width; u++)
3315               {
3316                 contrast=GetPixelIntensity(image,pixels)-intensity;
3317                 if (fabs(contrast) < threshold)
3318                   {
3319                     pixel+=(*k)*pixels[i];
3320                     gamma+=(*k);
3321                   }
3322                 k++;
3323                 pixels+=GetPixelChannels(image);
3324               }
3325               pixels+=image->columns*GetPixelChannels(image);
3326             }
3327             if (fabs((double) gamma) < MagickEpsilon)
3328               {
3329                 SetPixelChannel(blur_image,channel,p[center+i],q);
3330                 continue;
3331               }
3332             gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
3333             SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
3334             continue;
3335           }
3336         for (v=0; v < (ssize_t) width; v++)
3337         {
3338           for (u=0; u < (ssize_t) width; u++)
3339           {
3340             contrast=GetPixelIntensity(image,pixels)-intensity;
3341             if (fabs(contrast) < threshold)
3342               {
3343                 alpha=(MagickRealType) (QuantumScale*
3344                   GetPixelAlpha(image,pixels));
3345                 pixel+=(*k)*alpha*pixels[i];
3346                 gamma+=(*k)*alpha;
3347               }
3348             k++;
3349             pixels+=GetPixelChannels(image);
3350           }
3351           pixels+=image->columns*GetPixelChannels(image);
3352         }
3353         if (fabs((double) gamma) < MagickEpsilon)
3354           {
3355             SetPixelChannel(blur_image,channel,p[center+i],q);
3356             continue;
3357           }
3358         gamma=1.0/(fabs((double) gamma) <= MagickEpsilon ? 1.0 : gamma);
3359         SetPixelChannel(blur_image,channel,ClampToQuantum(gamma*pixel),q);
3360       }
3361       p+=GetPixelChannels(image);
3362       q+=GetPixelChannels(blur_image);
3363     }
3364     sync=SyncCacheViewAuthenticPixels(blur_view,exception);
3365     if (sync == MagickFalse)
3366       status=MagickFalse;
3367     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3368       {
3369         MagickBooleanType
3370           proceed;
3371
3372 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3373   #pragma omp critical (MagickCore_SelectiveBlurImage)
3374 #endif
3375         proceed=SetImageProgress(image,SelectiveBlurImageTag,progress++,
3376           image->rows);
3377         if (proceed == MagickFalse)
3378           status=MagickFalse;
3379       }
3380   }
3381   blur_image->type=image->type;
3382   blur_view=DestroyCacheView(blur_view);
3383   image_view=DestroyCacheView(image_view);
3384   kernel=(double *) RelinquishMagickMemory(kernel);
3385   if (status == MagickFalse)
3386     blur_image=DestroyImage(blur_image);
3387   return(blur_image);
3388 }
3389 \f
3390 /*
3391 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3392 %                                                                             %
3393 %                                                                             %
3394 %                                                                             %
3395 %     S h a d e I m a g e                                                     %
3396 %                                                                             %
3397 %                                                                             %
3398 %                                                                             %
3399 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3400 %
3401 %  ShadeImage() shines a distant light on an image to create a
3402 %  three-dimensional effect. You control the positioning of the light with
3403 %  azimuth and elevation; azimuth is measured in degrees off the x axis
3404 %  and elevation is measured in pixels above the Z axis.
3405 %
3406 %  The format of the ShadeImage method is:
3407 %
3408 %      Image *ShadeImage(const Image *image,const MagickBooleanType gray,
3409 %        const double azimuth,const double elevation,ExceptionInfo *exception)
3410 %
3411 %  A description of each parameter follows:
3412 %
3413 %    o image: the image.
3414 %
3415 %    o gray: A value other than zero shades the intensity of each pixel.
3416 %
3417 %    o azimuth, elevation:  Define the light source direction.
3418 %
3419 %    o exception: return any errors or warnings in this structure.
3420 %
3421 */
3422 MagickExport Image *ShadeImage(const Image *image,const MagickBooleanType gray,
3423   const double azimuth,const double elevation,ExceptionInfo *exception)
3424 {
3425 #define ShadeImageTag  "Shade/Image"
3426
3427   CacheView
3428     *image_view,
3429     *shade_view;
3430
3431   Image
3432     *shade_image;
3433
3434   MagickBooleanType
3435     status;
3436
3437   MagickOffsetType
3438     progress;
3439
3440   PrimaryInfo
3441     light;
3442
3443   ssize_t
3444     y;
3445
3446   /*
3447     Initialize shaded image attributes.
3448   */
3449   assert(image != (const Image *) NULL);
3450   assert(image->signature == MagickSignature);
3451   if (image->debug != MagickFalse)
3452     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3453   assert(exception != (ExceptionInfo *) NULL);
3454   assert(exception->signature == MagickSignature);
3455   shade_image=CloneImage(image,image->columns,image->rows,MagickTrue,exception);
3456   if (shade_image == (Image *) NULL)
3457     return((Image *) NULL);
3458   if (SetImageStorageClass(shade_image,DirectClass,exception) == MagickFalse)
3459     {
3460       shade_image=DestroyImage(shade_image);
3461       return((Image *) NULL);
3462     }
3463   /*
3464     Compute the light vector.
3465   */
3466   light.x=(double) QuantumRange*cos(DegreesToRadians(azimuth))*
3467     cos(DegreesToRadians(elevation));
3468   light.y=(double) QuantumRange*sin(DegreesToRadians(azimuth))*
3469     cos(DegreesToRadians(elevation));
3470   light.z=(double) QuantumRange*sin(DegreesToRadians(elevation));
3471   /*
3472     Shade image.
3473   */
3474   status=MagickTrue;
3475   progress=0;
3476   image_view=AcquireCacheView(image);
3477   shade_view=AcquireCacheView(shade_image);
3478 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3479   #pragma omp parallel for schedule(static,4) shared(progress,status)
3480 #endif
3481   for (y=0; y < (ssize_t) image->rows; y++)
3482   {
3483     MagickRealType
3484       distance,
3485       normal_distance,
3486       shade;
3487
3488     PrimaryInfo
3489       normal;
3490
3491     register const Quantum
3492       *restrict center,
3493       *restrict p,
3494       *restrict post,
3495       *restrict pre;
3496
3497     register Quantum
3498       *restrict q;
3499
3500     register ssize_t
3501       x;
3502
3503     if (status == MagickFalse)
3504       continue;
3505     p=GetCacheViewVirtualPixels(image_view,-1,y-1,image->columns+2,3,exception);
3506     q=QueueCacheViewAuthenticPixels(shade_view,0,y,shade_image->columns,1,
3507       exception);
3508     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
3509       {
3510         status=MagickFalse;
3511         continue;
3512       }
3513     /*
3514       Shade this row of pixels.
3515     */
3516     normal.z=2.0*(double) QuantumRange;  /* constant Z of surface normal */
3517     pre=p+GetPixelChannels(image);
3518     center=pre+(image->columns+2)*GetPixelChannels(image);
3519     post=center+(image->columns+2)*GetPixelChannels(image);
3520     for (x=0; x < (ssize_t) image->columns; x++)
3521     {
3522       register ssize_t
3523         i;
3524
3525       /*
3526         Determine the surface normal and compute shading.
3527       */
3528       normal.x=(double) (GetPixelIntensity(image,pre-GetPixelChannels(image))+
3529         GetPixelIntensity(image,center-GetPixelChannels(image))+
3530         GetPixelIntensity(image,post-GetPixelChannels(image))-
3531         GetPixelIntensity(image,pre+GetPixelChannels(image))-
3532         GetPixelIntensity(image,center+GetPixelChannels(image))-
3533         GetPixelIntensity(image,post+GetPixelChannels(image)));
3534       normal.y=(double) (GetPixelIntensity(image,post-GetPixelChannels(image))+
3535         GetPixelIntensity(image,post)+GetPixelIntensity(image,post+
3536         GetPixelChannels(image))-GetPixelIntensity(image,pre-
3537         GetPixelChannels(image))-GetPixelIntensity(image,pre)-
3538         GetPixelIntensity(image,pre+GetPixelChannels(image)));
3539       if ((normal.x == 0.0) && (normal.y == 0.0))
3540         shade=light.z;
3541       else
3542         {
3543           shade=0.0;
3544           distance=normal.x*light.x+normal.y*light.y+normal.z*light.z;
3545           if (distance > MagickEpsilon)
3546             {
3547               normal_distance=
3548                 normal.x*normal.x+normal.y*normal.y+normal.z*normal.z;
3549               if (normal_distance > (MagickEpsilon*MagickEpsilon))
3550                 shade=distance/sqrt((double) normal_distance);
3551             }
3552         }
3553       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
3554       {
3555         PixelChannel
3556           channel;
3557
3558         PixelTrait
3559           shade_traits,
3560           traits;
3561
3562         channel=GetPixelChannelMapChannel(image,i);
3563         traits=GetPixelChannelMapTraits(image,channel);
3564         shade_traits=GetPixelChannelMapTraits(shade_image,channel);
3565         if ((traits == UndefinedPixelTrait) ||
3566             (shade_traits == UndefinedPixelTrait))
3567           continue;
3568         if ((shade_traits & CopyPixelTrait) != 0)
3569           {
3570             SetPixelChannel(shade_image,channel,center[i],q);
3571             continue;
3572           }
3573         if (gray != MagickFalse)
3574           {
3575             SetPixelChannel(shade_image,channel,ClampToQuantum(shade),q);
3576             continue;
3577           }
3578         SetPixelChannel(shade_image,channel,ClampToQuantum(QuantumScale*shade*
3579           center[i]),q);
3580       }
3581       pre+=GetPixelChannels(image);
3582       center+=GetPixelChannels(image);
3583       post+=GetPixelChannels(image);
3584       q+=GetPixelChannels(shade_image);
3585     }
3586     if (SyncCacheViewAuthenticPixels(shade_view,exception) == MagickFalse)
3587       status=MagickFalse;
3588     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3589       {
3590         MagickBooleanType
3591           proceed;
3592
3593 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3594   #pragma omp critical (MagickCore_ShadeImage)
3595 #endif
3596         proceed=SetImageProgress(image,ShadeImageTag,progress++,image->rows);
3597         if (proceed == MagickFalse)
3598           status=MagickFalse;
3599       }
3600   }
3601   shade_view=DestroyCacheView(shade_view);
3602   image_view=DestroyCacheView(image_view);
3603   if (status == MagickFalse)
3604     shade_image=DestroyImage(shade_image);
3605   return(shade_image);
3606 }
3607 \f
3608 /*
3609 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3610 %                                                                             %
3611 %                                                                             %
3612 %                                                                             %
3613 %     S h a r p e n I m a g e                                                 %
3614 %                                                                             %
3615 %                                                                             %
3616 %                                                                             %
3617 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3618 %
3619 %  SharpenImage() sharpens the image.  We convolve the image with a Gaussian
3620 %  operator of the given radius and standard deviation (sigma).  For
3621 %  reasonable results, radius should be larger than sigma.  Use a radius of 0
3622 %  and SharpenImage() selects a suitable radius for you.
3623 %
3624 %  Using a separable kernel would be faster, but the negative weights cancel
3625 %  out on the corners of the kernel producing often undesirable ringing in the
3626 %  filtered result; this can be avoided by using a 2D gaussian shaped image
3627 %  sharpening kernel instead.
3628 %
3629 %  The format of the SharpenImage method is:
3630 %
3631 %    Image *SharpenImage(const Image *image,const double radius,
3632 %      const double sigma,const double bias,ExceptionInfo *exception)
3633 %
3634 %  A description of each parameter follows:
3635 %
3636 %    o image: the image.
3637 %
3638 %    o radius: the radius of the Gaussian, in pixels, not counting the center
3639 %      pixel.
3640 %
3641 %    o sigma: the standard deviation of the Laplacian, in pixels.
3642 %
3643 %    o bias: bias.
3644 %
3645 %    o exception: return any errors or warnings in this structure.
3646 %
3647 */
3648 MagickExport Image *SharpenImage(const Image *image,const double radius,
3649   const double sigma,const double bias,ExceptionInfo *exception)
3650 {
3651   double
3652     normalize;
3653
3654   Image
3655     *sharp_image;
3656
3657   KernelInfo
3658     *kernel_info;
3659
3660   register ssize_t
3661     i;
3662
3663   size_t
3664     width;
3665
3666   ssize_t
3667     j,
3668     u,
3669     v;
3670
3671   assert(image != (const Image *) NULL);
3672   assert(image->signature == MagickSignature);
3673   if (image->debug != MagickFalse)
3674     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3675   assert(exception != (ExceptionInfo *) NULL);
3676   assert(exception->signature == MagickSignature);
3677   width=GetOptimalKernelWidth2D(radius,sigma);
3678   kernel_info=AcquireKernelInfo((const char *) NULL);
3679   if (kernel_info == (KernelInfo *) NULL)
3680     ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
3681   (void) ResetMagickMemory(kernel_info,0,sizeof(*kernel_info));
3682   kernel_info->width=width;
3683   kernel_info->height=width;
3684   kernel_info->bias=bias;   /* FUTURE: user bias - non-sensical! */
3685   kernel_info->signature=MagickSignature;
3686   kernel_info->values=(MagickRealType *) AcquireAlignedMemory(
3687     kernel_info->width,kernel_info->width*sizeof(*kernel_info->values));
3688   if (kernel_info->values == (MagickRealType *) NULL)
3689     {
3690       kernel_info=DestroyKernelInfo(kernel_info);
3691       ThrowImageException(ResourceLimitError,"MemoryAllocationFailed");
3692     }
3693   normalize=0.0;
3694   j=(ssize_t) kernel_info->width/2;
3695   i=0;
3696   for (v=(-j); v <= j; v++)
3697   {
3698     for (u=(-j); u <= j; u++)
3699     {
3700       kernel_info->values[i]=(double) (-exp(-((double) u*u+v*v)/(2.0*
3701         MagickSigma*MagickSigma))/(2.0*MagickPI*MagickSigma*MagickSigma));
3702       normalize+=kernel_info->values[i];
3703       i++;
3704     }
3705   }
3706   kernel_info->values[i/2]=(double) ((-2.0)*normalize);
3707   sharp_image=ConvolveImage(image,kernel_info,exception);
3708   kernel_info=DestroyKernelInfo(kernel_info);
3709   return(sharp_image);
3710 }
3711 \f
3712 /*
3713 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3714 %                                                                             %
3715 %                                                                             %
3716 %                                                                             %
3717 %     S p r e a d I m a g e                                                   %
3718 %                                                                             %
3719 %                                                                             %
3720 %                                                                             %
3721 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3722 %
3723 %  SpreadImage() is a special effects method that randomly displaces each
3724 %  pixel in a block defined by the radius parameter.
3725 %
3726 %  The format of the SpreadImage method is:
3727 %
3728 %      Image *SpreadImage(const Image *image,const double radius,
3729 %        const PixelInterpolateMethod method,ExceptionInfo *exception)
3730 %
3731 %  A description of each parameter follows:
3732 %
3733 %    o image: the image.
3734 %
3735 %    o radius:  choose a random pixel in a neighborhood of this extent.
3736 %
3737 %    o method:  the pixel interpolation method.
3738 %
3739 %    o exception: return any errors or warnings in this structure.
3740 %
3741 */
3742 MagickExport Image *SpreadImage(const Image *image,const double radius,
3743   const PixelInterpolateMethod method,ExceptionInfo *exception)
3744 {
3745 #define SpreadImageTag  "Spread/Image"
3746
3747   CacheView
3748     *image_view,
3749     *spread_view;
3750
3751   Image
3752     *spread_image;
3753
3754   MagickBooleanType
3755     status;
3756
3757   MagickOffsetType
3758     progress;
3759
3760   RandomInfo
3761     **restrict random_info;
3762
3763   size_t
3764     width;
3765
3766   ssize_t
3767     y;
3768
3769   /*
3770     Initialize spread image attributes.
3771   */
3772   assert(image != (Image *) NULL);
3773   assert(image->signature == MagickSignature);
3774   if (image->debug != MagickFalse)
3775     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3776   assert(exception != (ExceptionInfo *) NULL);
3777   assert(exception->signature == MagickSignature);
3778   spread_image=CloneImage(image,image->columns,image->rows,MagickTrue,
3779     exception);
3780   if (spread_image == (Image *) NULL)
3781     return((Image *) NULL);
3782   if (SetImageStorageClass(spread_image,DirectClass,exception) == MagickFalse)
3783     {
3784       spread_image=DestroyImage(spread_image);
3785       return((Image *) NULL);
3786     }
3787   /*
3788     Spread image.
3789   */
3790   status=MagickTrue;
3791   progress=0;
3792   width=GetOptimalKernelWidth1D(radius,0.5);
3793   random_info=AcquireRandomInfoThreadSet();
3794   image_view=AcquireCacheView(image);
3795   spread_view=AcquireCacheView(spread_image);
3796 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3797   #pragma omp parallel for schedule(static,4) shared(progress,status) omp_throttle(1)
3798 #endif
3799   for (y=0; y < (ssize_t) image->rows; y++)
3800   {
3801     const int
3802       id = GetOpenMPThreadId();
3803
3804     register const Quantum
3805       *restrict p;
3806
3807     register Quantum
3808       *restrict q;
3809
3810     register ssize_t
3811       x;
3812
3813     if (status == MagickFalse)
3814       continue;
3815     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
3816     q=QueueCacheViewAuthenticPixels(spread_view,0,y,spread_image->columns,1,
3817       exception);
3818     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
3819       {
3820         status=MagickFalse;
3821         continue;
3822       }
3823     for (x=0; x < (ssize_t) image->columns; x++)
3824     {
3825       PointInfo
3826         point;
3827
3828       point.x=GetPseudoRandomValue(random_info[id]);
3829       point.y=GetPseudoRandomValue(random_info[id]);
3830       status=InterpolatePixelChannels(image,image_view,spread_image,method,
3831         (double) x+width*point.x-0.5,(double) y+width*point.y-0.5,q,exception);
3832       q+=GetPixelChannels(spread_image);
3833     }
3834     if (SyncCacheViewAuthenticPixels(spread_view,exception) == MagickFalse)
3835       status=MagickFalse;
3836     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3837       {
3838         MagickBooleanType
3839           proceed;
3840
3841 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3842   #pragma omp critical (MagickCore_SpreadImage)
3843 #endif
3844         proceed=SetImageProgress(image,SpreadImageTag,progress++,image->rows);
3845         if (proceed == MagickFalse)
3846           status=MagickFalse;
3847       }
3848   }
3849   spread_view=DestroyCacheView(spread_view);
3850   image_view=DestroyCacheView(image_view);
3851   random_info=DestroyRandomInfoThreadSet(random_info);
3852   return(spread_image);
3853 }
3854 \f
3855 /*
3856 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3857 %                                                                             %
3858 %                                                                             %
3859 %                                                                             %
3860 %     U n s h a r p M a s k I m a g e                                         %
3861 %                                                                             %
3862 %                                                                             %
3863 %                                                                             %
3864 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3865 %
3866 %  UnsharpMaskImage() sharpens one or more image channels.  We convolve the
3867 %  image with a Gaussian operator of the given radius and standard deviation
3868 %  (sigma).  For reasonable results, radius should be larger than sigma.  Use a
3869 %  radius of 0 and UnsharpMaskImage() selects a suitable radius for you.
3870 %
3871 %  The format of the UnsharpMaskImage method is:
3872 %
3873 %    Image *UnsharpMaskImage(const Image *image,const double radius,
3874 %      const double sigma,const double amount,const double threshold,
3875 %      ExceptionInfo *exception)
3876 %
3877 %  A description of each parameter follows:
3878 %
3879 %    o image: the image.
3880 %
3881 %    o radius: the radius of the Gaussian, in pixels, not counting the center
3882 %      pixel.
3883 %
3884 %    o sigma: the standard deviation of the Gaussian, in pixels.
3885 %
3886 %    o amount: the percentage of the difference between the original and the
3887 %      blur image that is added back into the original.
3888 %
3889 %    o threshold: the threshold in pixels needed to apply the diffence amount.
3890 %
3891 %    o exception: return any errors or warnings in this structure.
3892 %
3893 */
3894 MagickExport Image *UnsharpMaskImage(const Image *image,const double radius,
3895   const double sigma,const double amount,const double threshold,
3896   ExceptionInfo *exception)
3897 {
3898 #define SharpenImageTag  "Sharpen/Image"
3899
3900   CacheView
3901     *image_view,
3902     *unsharp_view;
3903
3904   Image
3905     *unsharp_image;
3906
3907   MagickBooleanType
3908     status;
3909
3910   MagickOffsetType
3911     progress;
3912
3913   MagickRealType
3914     quantum_threshold;
3915
3916   ssize_t
3917     y;
3918
3919   assert(image != (const Image *) NULL);
3920   assert(image->signature == MagickSignature);
3921   if (image->debug != MagickFalse)
3922     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3923   assert(exception != (ExceptionInfo *) NULL);
3924
3925
3926   /* FUTURE:  use of bias on sharpen is non-sensical */
3927   unsharp_image=BlurImage(image,radius,sigma,image->bias,exception);
3928
3929   if (unsharp_image == (Image *) NULL)
3930     return((Image *) NULL);
3931   quantum_threshold=(MagickRealType) QuantumRange*threshold;
3932   /*
3933     Unsharp-mask image.
3934   */
3935   status=MagickTrue;
3936   progress=0;
3937   image_view=AcquireCacheView(image);
3938   unsharp_view=AcquireCacheView(unsharp_image);
3939 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3940   #pragma omp parallel for schedule(static,4) shared(progress,status)
3941 #endif
3942   for (y=0; y < (ssize_t) image->rows; y++)
3943   {
3944     register const Quantum
3945       *restrict p;
3946
3947     register Quantum
3948       *restrict q;
3949
3950     register ssize_t
3951       x;
3952
3953     if (status == MagickFalse)
3954       continue;
3955     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
3956     q=GetCacheViewAuthenticPixels(unsharp_view,0,y,unsharp_image->columns,1,
3957       exception);
3958     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
3959       {
3960         status=MagickFalse;
3961         continue;
3962       }
3963     for (x=0; x < (ssize_t) image->columns; x++)
3964     {
3965       register ssize_t
3966         i;
3967
3968       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
3969       {
3970         MagickRealType
3971           pixel;
3972
3973         PixelChannel
3974           channel;
3975
3976         PixelTrait
3977           traits,
3978           unsharp_traits;
3979
3980         channel=GetPixelChannelMapChannel(image,i);
3981         traits=GetPixelChannelMapTraits(image,channel);
3982         unsharp_traits=GetPixelChannelMapTraits(unsharp_image,channel);
3983         if ((traits == UndefinedPixelTrait) ||
3984             (unsharp_traits == UndefinedPixelTrait))
3985           continue;
3986         if ((unsharp_traits & CopyPixelTrait) != 0)
3987           {
3988             SetPixelChannel(unsharp_image,channel,p[i],q);
3989             continue;
3990           }
3991         pixel=p[i]-(MagickRealType) GetPixelChannel(unsharp_image,channel,q);
3992         if (fabs(2.0*pixel) < quantum_threshold)
3993           pixel=(MagickRealType) p[i];
3994         else
3995           pixel=(MagickRealType) p[i]+amount*pixel;
3996         SetPixelChannel(unsharp_image,channel,ClampToQuantum(pixel),q);
3997       }
3998       p+=GetPixelChannels(image);
3999       q+=GetPixelChannels(unsharp_image);
4000     }
4001     if (SyncCacheViewAuthenticPixels(unsharp_view,exception) == MagickFalse)
4002       status=MagickFalse;
4003     if (image->progress_monitor != (MagickProgressMonitor) NULL)
4004       {
4005         MagickBooleanType
4006           proceed;
4007
4008 #if defined(MAGICKCORE_OPENMP_SUPPORT)
4009   #pragma omp critical (MagickCore_UnsharpMaskImage)
4010 #endif
4011         proceed=SetImageProgress(image,SharpenImageTag,progress++,image->rows);
4012         if (proceed == MagickFalse)
4013           status=MagickFalse;
4014       }
4015   }
4016   unsharp_image->type=image->type;
4017   unsharp_view=DestroyCacheView(unsharp_view);
4018   image_view=DestroyCacheView(image_view);
4019   if (status == MagickFalse)
4020     unsharp_image=DestroyImage(unsharp_image);
4021   return(unsharp_image);
4022 }