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