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

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