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