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