• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 /*
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 %                                                                             %
4 %                                                                             %
5 %                                                                             %
6 %              EEEEE  N   N  H   H   AAA   N   N   CCCC  EEEEE                %
7 %              E      NN  N  H   H  A   A  NN  N  C      E                    %
8 %              EEE    N N N  HHHHH  AAAAA  N N N  C      EEE                  %
9 %              E      N  NN  H   H  A   A  N  NN  C      E                    %
10 %              EEEEE  N   N  H   H  A   A  N   N   CCCC  EEEEE                %
11 %                                                                             %
12 %                                                                             %
13 %                    MagickCore Image Enhancement Methods                     %
14 %                                                                             %
15 %                              Software Design                                %
16 %                                   Cristy                                    %
17 %                                 July 1992                                   %
18 %                                                                             %
19 %                                                                             %
20 %  Copyright 1999-2019 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 %    https://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/artifact.h"
46 #include "MagickCore/attribute.h"
47 #include "MagickCore/cache.h"
48 #include "MagickCore/cache-private.h"
49 #include "MagickCore/cache-view.h"
50 #include "MagickCore/channel.h"
51 #include "MagickCore/color.h"
52 #include "MagickCore/color-private.h"
53 #include "MagickCore/colorspace.h"
54 #include "MagickCore/colorspace-private.h"
55 #include "MagickCore/composite-private.h"
56 #include "MagickCore/enhance.h"
57 #include "MagickCore/exception.h"
58 #include "MagickCore/exception-private.h"
59 #include "MagickCore/fx.h"
60 #include "MagickCore/gem.h"
61 #include "MagickCore/gem-private.h"
62 #include "MagickCore/geometry.h"
63 #include "MagickCore/histogram.h"
64 #include "MagickCore/image.h"
65 #include "MagickCore/image-private.h"
66 #include "MagickCore/memory_.h"
67 #include "MagickCore/monitor.h"
68 #include "MagickCore/monitor-private.h"
69 #include "MagickCore/option.h"
70 #include "MagickCore/pixel.h"
71 #include "MagickCore/pixel-accessor.h"
72 #include "MagickCore/quantum.h"
73 #include "MagickCore/quantum-private.h"
74 #include "MagickCore/resample.h"
75 #include "MagickCore/resample-private.h"
76 #include "MagickCore/resource_.h"
77 #include "MagickCore/statistic.h"
78 #include "MagickCore/string_.h"
79 #include "MagickCore/string-private.h"
80 #include "MagickCore/thread-private.h"
81 #include "MagickCore/threshold.h"
82 #include "MagickCore/token.h"
83 #include "MagickCore/xml-tree.h"
84 #include "MagickCore/xml-tree-private.h"
85 
86 /*
87 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
88 %                                                                             %
89 %                                                                             %
90 %                                                                             %
91 %     A u t o G a m m a I m a g e                                             %
92 %                                                                             %
93 %                                                                             %
94 %                                                                             %
95 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
96 %
97 %  AutoGammaImage() extract the 'mean' from the image and adjust the image
98 %  to try make set its gamma appropriatally.
99 %
100 %  The format of the AutoGammaImage method is:
101 %
102 %      MagickBooleanType AutoGammaImage(Image *image,ExceptionInfo *exception)
103 %
104 %  A description of each parameter follows:
105 %
106 %    o image: The image to auto-level
107 %
108 %    o exception: return any errors or warnings in this structure.
109 %
110 */
AutoGammaImage(Image * image,ExceptionInfo * exception)111 MagickExport MagickBooleanType AutoGammaImage(Image *image,
112   ExceptionInfo *exception)
113 {
114   double
115     gamma,
116     log_mean,
117     mean,
118     sans;
119 
120   MagickStatusType
121     status;
122 
123   register ssize_t
124     i;
125 
126   log_mean=log(0.5);
127   if (image->channel_mask == DefaultChannels)
128     {
129       /*
130         Apply gamma correction equally across all given channels.
131       */
132       (void) GetImageMean(image,&mean,&sans,exception);
133       gamma=log(mean*QuantumScale)/log_mean;
134       return(LevelImage(image,0.0,(double) QuantumRange,gamma,exception));
135     }
136   /*
137     Auto-gamma each channel separately.
138   */
139   status=MagickTrue;
140   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
141   {
142     ChannelType
143       channel_mask;
144 
145     PixelChannel channel = GetPixelChannelChannel(image,i);
146     PixelTrait traits = GetPixelChannelTraits(image,channel);
147     if ((traits & UpdatePixelTrait) == 0)
148       continue;
149     channel_mask=SetImageChannelMask(image,(ChannelType) (1UL << i));
150     status=GetImageMean(image,&mean,&sans,exception);
151     gamma=log(mean*QuantumScale)/log_mean;
152     status&=LevelImage(image,0.0,(double) QuantumRange,gamma,exception);
153     (void) SetImageChannelMask(image,channel_mask);
154     if (status == MagickFalse)
155       break;
156   }
157   return(status != 0 ? MagickTrue : MagickFalse);
158 }
159 
160 /*
161 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
162 %                                                                             %
163 %                                                                             %
164 %                                                                             %
165 %     A u t o L e v e l I m a g e                                             %
166 %                                                                             %
167 %                                                                             %
168 %                                                                             %
169 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
170 %
171 %  AutoLevelImage() adjusts the levels of a particular image channel by
172 %  scaling the minimum and maximum values to the full quantum range.
173 %
174 %  The format of the LevelImage method is:
175 %
176 %      MagickBooleanType AutoLevelImage(Image *image,ExceptionInfo *exception)
177 %
178 %  A description of each parameter follows:
179 %
180 %    o image: The image to auto-level
181 %
182 %    o exception: return any errors or warnings in this structure.
183 %
184 */
AutoLevelImage(Image * image,ExceptionInfo * exception)185 MagickExport MagickBooleanType AutoLevelImage(Image *image,
186   ExceptionInfo *exception)
187 {
188   return(MinMaxStretchImage(image,0.0,0.0,1.0,exception));
189 }
190 
191 /*
192 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
193 %                                                                             %
194 %                                                                             %
195 %                                                                             %
196 %     B r i g h t n e s s C o n t r a s t I m a g e                           %
197 %                                                                             %
198 %                                                                             %
199 %                                                                             %
200 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
201 %
202 %  BrightnessContrastImage() changes the brightness and/or contrast of an
203 %  image.  It converts the brightness and contrast parameters into slope and
204 %  intercept and calls a polynomical function to apply to the image.
205 %
206 %  The format of the BrightnessContrastImage method is:
207 %
208 %      MagickBooleanType BrightnessContrastImage(Image *image,
209 %        const double brightness,const double contrast,ExceptionInfo *exception)
210 %
211 %  A description of each parameter follows:
212 %
213 %    o image: the image.
214 %
215 %    o brightness: the brightness percent (-100 .. 100).
216 %
217 %    o contrast: the contrast percent (-100 .. 100).
218 %
219 %    o exception: return any errors or warnings in this structure.
220 %
221 */
BrightnessContrastImage(Image * image,const double brightness,const double contrast,ExceptionInfo * exception)222 MagickExport MagickBooleanType BrightnessContrastImage(Image *image,
223   const double brightness,const double contrast,ExceptionInfo *exception)
224 {
225 #define BrightnessContastImageTag  "BrightnessContast/Image"
226 
227   double
228     alpha,
229     coefficients[2],
230     intercept,
231     slope;
232 
233   MagickBooleanType
234     status;
235 
236   /*
237     Compute slope and intercept.
238   */
239   assert(image != (Image *) NULL);
240   assert(image->signature == MagickCoreSignature);
241   if (image->debug != MagickFalse)
242     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
243   alpha=contrast;
244   slope=tan((double) (MagickPI*(alpha/100.0+1.0)/4.0));
245   if (slope < 0.0)
246     slope=0.0;
247   intercept=brightness/100.0+((100-brightness)/200.0)*(1.0-slope);
248   coefficients[0]=slope;
249   coefficients[1]=intercept;
250   status=FunctionImage(image,PolynomialFunction,2,coefficients,exception);
251   return(status);
252 }
253 
254 /*
255 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
256 %                                                                             %
257 %                                                                             %
258 %                                                                             %
259 %     C L A H E I m a g e                                                     %
260 %                                                                             %
261 %                                                                             %
262 %                                                                             %
263 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
264 %
265 %  CLAHEImage() is a variant of adaptive histogram equalization in which the
266 %  contrast amplification is limited, so as to reduce this problem of noise
267 %  amplification.
268 %
269 %  Adapted from implementation by Karel Zuiderveld, karel@cv.ruu.nl in
270 %  "Graphics Gems IV", Academic Press, 1994.
271 %
272 %  The format of the CLAHEImage method is:
273 %
274 %      MagickBooleanType CLAHEImage(Image *image,const size_t width,
275 %        const size_t height,const size_t number_bins,const double clip_limit,
276 %        ExceptionInfo *exception)
277 %
278 %  A description of each parameter follows:
279 %
280 %    o image: the image.
281 %
282 %    o width: the width of the tile divisions to use in horizontal direction.
283 %
284 %    o height: the height of the tile divisions to use in vertical direction.
285 %
286 %    o number_bins: number of bins for histogram ("dynamic range").
287 %
288 %    o clip_limit: contrast limit for localised changes in contrast. A limit
289 %      less than 1 results in standard non-contrast limited AHE.
290 %
291 %    o exception: return any errors or warnings in this structure.
292 %
293 */
294 
295 typedef struct _RangeInfo
296 {
297   unsigned short
298     min,
299     max;
300 } RangeInfo;
301 
ClipCLAHEHistogram(const double clip_limit,const size_t number_bins,size_t * histogram)302 static void ClipCLAHEHistogram(const double clip_limit,const size_t number_bins,
303   size_t *histogram)
304 {
305 #define NumberCLAHEGrays  (65536)
306 
307   register ssize_t
308     i;
309 
310   size_t
311     cumulative_excess,
312     previous_excess,
313     step;
314 
315   ssize_t
316     excess;
317 
318   /*
319     Compute total number of excess pixels.
320   */
321   cumulative_excess=0;
322   for (i=0; i < (ssize_t) number_bins; i++)
323   {
324     excess=(ssize_t) histogram[i]-(ssize_t) clip_limit;
325     if (excess > 0)
326       cumulative_excess+=excess;
327   }
328   /*
329     Clip histogram and redistribute excess pixels across all bins.
330   */
331   step=cumulative_excess/number_bins;
332   excess=(ssize_t) (clip_limit-step);
333   for (i=0; i < (ssize_t) number_bins; i++)
334   {
335     if ((double) histogram[i] > clip_limit)
336       histogram[i]=(size_t) clip_limit;
337     else
338       if ((ssize_t) histogram[i] > excess)
339         {
340           cumulative_excess-=histogram[i]-excess;
341           histogram[i]=(size_t) clip_limit;
342         }
343       else
344         {
345           cumulative_excess-=step;
346           histogram[i]+=step;
347         }
348   }
349   /*
350     Redistribute remaining excess.
351   */
352   do
353   {
354     register size_t
355       *p;
356 
357     size_t
358       *q;
359 
360     previous_excess=cumulative_excess;
361     p=histogram;
362     q=histogram+number_bins;
363     while ((cumulative_excess != 0) && (p < q))
364     {
365       step=number_bins/cumulative_excess;
366       if (step < 1)
367         step=1;
368       for (p=histogram; (p < q) && (cumulative_excess != 0); p+=step)
369         if ((double) *p < clip_limit)
370           {
371             (*p)++;
372             cumulative_excess--;
373           }
374       p++;
375     }
376   } while ((cumulative_excess != 0) && (cumulative_excess < previous_excess));
377 }
378 
GenerateCLAHEHistogram(const RectangleInfo * clahe_info,const RectangleInfo * tile_info,const size_t number_bins,const unsigned short * lut,const unsigned short * pixels,size_t * histogram)379 static void GenerateCLAHEHistogram(const RectangleInfo *clahe_info,
380   const RectangleInfo *tile_info,const size_t number_bins,
381   const unsigned short *lut,const unsigned short *pixels,size_t *histogram)
382 {
383   register const unsigned short
384     *p;
385 
386   register ssize_t
387     i;
388 
389   /*
390     Classify the pixels into a gray histogram.
391   */
392   for (i=0; i < (ssize_t) number_bins; i++)
393     histogram[i]=0L;
394   p=pixels;
395   for (i=0; i < (ssize_t) tile_info->height; i++)
396   {
397     const unsigned short
398       *q;
399 
400     q=p+tile_info->width;
401     while (p < q)
402       histogram[lut[*p++]]++;
403     q+=clahe_info->width;
404     p=q-tile_info->width;
405   }
406 }
407 
InterpolateCLAHE(const RectangleInfo * clahe_info,const size_t * Q12,const size_t * Q22,const size_t * Q11,const size_t * Q21,const RectangleInfo * tile,const unsigned short * lut,unsigned short * pixels)408 static void InterpolateCLAHE(const RectangleInfo *clahe_info,const size_t *Q12,
409   const size_t *Q22,const size_t *Q11,const size_t *Q21,
410   const RectangleInfo *tile,const unsigned short *lut,unsigned short *pixels)
411 {
412   ssize_t
413     y;
414 
415   unsigned short
416     intensity;
417 
418   /*
419     Bilinear interpolate four tiles to eliminate boundary artifacts.
420   */
421   for (y=(ssize_t) tile->height; y > 0; y--)
422   {
423     register ssize_t
424       x;
425 
426     for (x=(ssize_t) tile->width; x > 0; x--)
427     {
428       intensity=lut[*pixels];
429       *pixels++=(unsigned short ) (PerceptibleReciprocal((double) tile->width*
430         tile->height)*(y*(x*Q12[intensity]+(tile->width-x)*Q22[intensity])+
431         (tile->height-y)*(x*Q11[intensity]+(tile->width-x)*Q21[intensity])));
432     }
433     pixels+=(clahe_info->width-tile->width);
434   }
435 }
436 
GenerateCLAHELut(const RangeInfo * range_info,const size_t number_bins,unsigned short * lut)437 static void GenerateCLAHELut(const RangeInfo *range_info,
438   const size_t number_bins,unsigned short *lut)
439 {
440   ssize_t
441     i;
442 
443   unsigned short
444     delta;
445 
446   /*
447     Scale input image [intensity min,max] to [0,number_bins-1].
448   */
449   delta=(unsigned short) ((range_info->max-range_info->min)/number_bins+1);
450   for (i=(ssize_t) range_info->min; i <= (ssize_t) range_info->max; i++)
451     lut[i]=(unsigned short) ((i-range_info->min)/delta);
452 }
453 
MapCLAHEHistogram(const RangeInfo * range_info,const size_t number_bins,const size_t number_pixels,size_t * histogram)454 static void MapCLAHEHistogram(const RangeInfo *range_info,
455   const size_t number_bins,const size_t number_pixels,size_t *histogram)
456 {
457   double
458     scale,
459     sum;
460 
461   register ssize_t
462     i;
463 
464   /*
465     Rescale histogram to range [min-intensity .. max-intensity].
466   */
467   scale=(double) (range_info->max-range_info->min)/number_pixels;
468   sum=0.0;
469   for (i=0; i < (ssize_t) number_bins; i++)
470   {
471     sum+=histogram[i];
472     histogram[i]=(size_t) (range_info->min+scale*sum);
473     if (histogram[i] > range_info->max)
474       histogram[i]=range_info->max;
475   }
476 }
477 
CLAHE(const RectangleInfo * clahe_info,const RectangleInfo * tile_info,const RangeInfo * range_info,const size_t number_bins,const double clip_limit,unsigned short * pixels)478 static MagickBooleanType CLAHE(const RectangleInfo *clahe_info,
479   const RectangleInfo *tile_info,const RangeInfo *range_info,
480   const size_t number_bins,const double clip_limit,unsigned short *pixels)
481 {
482   MemoryInfo
483     *tile_cache;
484 
485   register unsigned short
486     *p;
487 
488   size_t
489     limit,
490     *tiles;
491 
492   ssize_t
493     y;
494 
495   unsigned short
496     lut[NumberCLAHEGrays];
497 
498   /*
499     Constrast limited adapted histogram equalization.
500   */
501   if (clip_limit == 1.0)
502     return(MagickTrue);
503   tile_cache=AcquireVirtualMemory((size_t) clahe_info->x*clahe_info->y,
504     number_bins*sizeof(*tiles));
505   if (tile_cache == (MemoryInfo *) NULL)
506     return(MagickFalse);
507   tiles=(size_t *) GetVirtualMemoryBlob(tile_cache);
508   limit=(size_t) (clip_limit*(tile_info->width*tile_info->height)/number_bins);
509   if (limit < 1UL)
510     limit=1UL;
511   /*
512     Generate greylevel mappings for each tile.
513   */
514   GenerateCLAHELut(range_info,number_bins,lut);
515   p=pixels;
516   for (y=0; y < (ssize_t) clahe_info->y; y++)
517   {
518     register ssize_t
519       x;
520 
521     for (x=0; x < (ssize_t) clahe_info->x; x++)
522     {
523       size_t
524         *histogram;
525 
526       histogram=tiles+(number_bins*(y*clahe_info->x+x));
527       GenerateCLAHEHistogram(clahe_info,tile_info,number_bins,lut,p,histogram);
528       ClipCLAHEHistogram((double) limit,number_bins,histogram);
529       MapCLAHEHistogram(range_info,number_bins,tile_info->width*
530         tile_info->height,histogram);
531       p+=tile_info->width;
532     }
533     p+=clahe_info->width*(tile_info->height-1);
534   }
535   /*
536     Interpolate greylevel mappings to get CLAHE image.
537   */
538   p=pixels;
539   for (y=0; y <= (ssize_t) clahe_info->y; y++)
540   {
541     OffsetInfo
542       offset;
543 
544     RectangleInfo
545       tile;
546 
547     register ssize_t
548       x;
549 
550     tile.height=tile_info->height;
551     tile.y=y-1;
552     offset.y=tile.y+1;
553     if (y == 0)
554       {
555         /*
556           Top row.
557         */
558         tile.height=tile_info->height >> 1;
559         tile.y=0;
560         offset.y=0;
561       }
562     else
563       if (y == (ssize_t) clahe_info->y)
564         {
565           /*
566             Bottom row.
567           */
568           tile.height=(tile_info->height+1) >> 1;
569           tile.y=clahe_info->y-1;
570           offset.y=tile.y;
571         }
572     for (x=0; x <= (ssize_t) clahe_info->x; x++)
573     {
574       tile.width=tile_info->width;
575       tile.x=x-1;
576       offset.x=tile.x+1;
577       if (x == 0)
578         {
579           /*
580             Left column.
581           */
582           tile.width=tile_info->width >> 1;
583           tile.x=0;
584           offset.x=0;
585         }
586       else
587         if (x == (ssize_t) clahe_info->x)
588           {
589             /*
590               Right column.
591             */
592             tile.width=(tile_info->width+1) >> 1;
593             tile.x=clahe_info->x-1;
594             offset.x=tile.x;
595           }
596       InterpolateCLAHE(clahe_info,
597         tiles+(number_bins*(tile.y*clahe_info->x+tile.x)),     /* Q12 */
598         tiles+(number_bins*(tile.y*clahe_info->x+offset.x)),   /* Q22 */
599         tiles+(number_bins*(offset.y*clahe_info->x+tile.x)),   /* Q11 */
600         tiles+(number_bins*(offset.y*clahe_info->x+offset.x)), /* Q21 */
601         &tile,lut,p);
602       p+=tile.width;
603     }
604     p+=clahe_info->width*(tile.height-1);
605   }
606   tile_cache=RelinquishVirtualMemory(tile_cache);
607   return(MagickTrue);
608 }
609 
CLAHEImage(Image * image,const size_t width,const size_t height,const size_t number_bins,const double clip_limit,ExceptionInfo * exception)610 MagickExport MagickBooleanType CLAHEImage(Image *image,const size_t width,
611   const size_t height,const size_t number_bins,const double clip_limit,
612   ExceptionInfo *exception)
613 {
614 #define CLAHEImageTag  "CLAHE/Image"
615 
616   CacheView
617     *image_view;
618 
619   ColorspaceType
620     colorspace;
621 
622   MagickBooleanType
623     status;
624 
625   MagickOffsetType
626     progress;
627 
628   MemoryInfo
629     *pixel_cache;
630 
631   RangeInfo
632     range_info;
633 
634   RectangleInfo
635     clahe_info,
636     tile_info;
637 
638   size_t
639     n;
640 
641   ssize_t
642     y;
643 
644   unsigned short
645     *pixels;
646 
647   /*
648     Configure CLAHE parameters.
649   */
650   assert(image != (Image *) NULL);
651   assert(image->signature == MagickCoreSignature);
652   if (image->debug != MagickFalse)
653     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
654   range_info.min=0;
655   range_info.max=NumberCLAHEGrays-1;
656   tile_info.width=width;
657   if (tile_info.width == 0)
658     tile_info.width=image->columns >> 3;
659   tile_info.height=height;
660   if (tile_info.height == 0)
661     tile_info.height=image->rows >> 3;
662   tile_info.x=0;
663   if ((image->columns % tile_info.width) != 0)
664     tile_info.x=(ssize_t) tile_info.width-(image->columns % tile_info.width);
665   tile_info.y=0;
666   if ((image->rows % tile_info.height) != 0)
667     tile_info.y=(ssize_t) tile_info.height-(image->rows % tile_info.height);
668   clahe_info.width=image->columns+tile_info.x;
669   clahe_info.height=image->rows+tile_info.y;
670   clahe_info.x=(ssize_t) clahe_info.width/tile_info.width;
671   clahe_info.y=(ssize_t) clahe_info.height/tile_info.height;
672   pixel_cache=AcquireVirtualMemory(clahe_info.width,clahe_info.height*
673     sizeof(*pixels));
674   if (pixel_cache == (MemoryInfo *) NULL)
675     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
676       image->filename);
677   pixels=(unsigned short *) GetVirtualMemoryBlob(pixel_cache);
678   colorspace=image->colorspace;
679   if (TransformImageColorspace(image,LabColorspace,exception) == MagickFalse)
680     {
681       pixel_cache=RelinquishVirtualMemory(pixel_cache);
682       return(MagickFalse);
683     }
684   /*
685     Initialize CLAHE pixels.
686   */
687   image_view=AcquireVirtualCacheView(image,exception);
688   progress=0;
689   status=MagickTrue;
690   n=0;
691   for (y=0; y < (ssize_t) clahe_info.height; y++)
692   {
693     register const Quantum
694       *magick_restrict p;
695 
696     register ssize_t
697       x;
698 
699     if (status == MagickFalse)
700       continue;
701     p=GetCacheViewVirtualPixels(image_view,-(tile_info.x >> 1),y-
702       (tile_info.y >> 1),clahe_info.width,1,exception);
703     if (p == (const Quantum *) NULL)
704       {
705         status=MagickFalse;
706         continue;
707       }
708     for (x=0; x < (ssize_t) clahe_info.width; x++)
709     {
710       pixels[n++]=ScaleQuantumToShort(p[0]);
711       p+=GetPixelChannels(image);
712     }
713     if (image->progress_monitor != (MagickProgressMonitor) NULL)
714       {
715         MagickBooleanType
716           proceed;
717 
718 #if defined(MAGICKCORE_OPENMP_SUPPORT)
719         #pragma omp atomic
720 #endif
721         progress++;
722         proceed=SetImageProgress(image,CLAHEImageTag,progress,2*
723           GetPixelChannels(image));
724         if (proceed == MagickFalse)
725           status=MagickFalse;
726       }
727   }
728   image_view=DestroyCacheView(image_view);
729   status=CLAHE(&clahe_info,&tile_info,&range_info,number_bins == 0 ?
730     (size_t) 128 : MagickMin(number_bins,256),clip_limit,pixels);
731   if (status == MagickFalse)
732     (void) ThrowMagickException(exception,GetMagickModule(),
733       ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
734   /*
735     Push CLAHE pixels to CLAHE image.
736   */
737   image_view=AcquireAuthenticCacheView(image,exception);
738   n=clahe_info.width*(tile_info.y >> 1);
739   for (y=0; y < (ssize_t) image->rows; y++)
740   {
741     register Quantum
742       *magick_restrict q;
743 
744     register ssize_t
745       x;
746 
747     if (status == MagickFalse)
748       continue;
749     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
750     if (q == (Quantum *) NULL)
751       {
752         status=MagickFalse;
753         continue;
754       }
755     n+=tile_info.x >> 1;
756     for (x=0; x < (ssize_t) image->columns; x++)
757     {
758       q[0]=ScaleShortToQuantum(pixels[n++]);
759       q+=GetPixelChannels(image);
760     }
761     n+=(clahe_info.width-image->columns-(tile_info.x >> 1));
762     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
763       status=MagickFalse;
764     if (image->progress_monitor != (MagickProgressMonitor) NULL)
765       {
766         MagickBooleanType
767           proceed;
768 
769 #if defined(MAGICKCORE_OPENMP_SUPPORT)
770         #pragma omp atomic
771 #endif
772         progress++;
773         proceed=SetImageProgress(image,CLAHEImageTag,progress,2*
774           GetPixelChannels(image));
775         if (proceed == MagickFalse)
776           status=MagickFalse;
777       }
778   }
779   image_view=DestroyCacheView(image_view);
780   pixel_cache=RelinquishVirtualMemory(pixel_cache);
781   if (TransformImageColorspace(image,colorspace,exception) == MagickFalse)
782     status=MagickFalse;
783   return(status);
784 }
785 
786 /*
787 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
788 %                                                                             %
789 %                                                                             %
790 %                                                                             %
791 %     C l u t I m a g e                                                       %
792 %                                                                             %
793 %                                                                             %
794 %                                                                             %
795 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
796 %
797 %  ClutImage() replaces each color value in the given image, by using it as an
798 %  index to lookup a replacement color value in a Color Look UP Table in the
799 %  form of an image.  The values are extracted along a diagonal of the CLUT
800 %  image so either a horizontal or vertial gradient image can be used.
801 %
802 %  Typically this is used to either re-color a gray-scale image according to a
803 %  color gradient in the CLUT image, or to perform a freeform histogram
804 %  (level) adjustment according to the (typically gray-scale) gradient in the
805 %  CLUT image.
806 %
807 %  When the 'channel' mask includes the matte/alpha transparency channel but
808 %  one image has no such channel it is assumed that that image is a simple
809 %  gray-scale image that will effect the alpha channel values, either for
810 %  gray-scale coloring (with transparent or semi-transparent colors), or
811 %  a histogram adjustment of existing alpha channel values.   If both images
812 %  have matte channels, direct and normal indexing is applied, which is rarely
813 %  used.
814 %
815 %  The format of the ClutImage method is:
816 %
817 %      MagickBooleanType ClutImage(Image *image,Image *clut_image,
818 %        const PixelInterpolateMethod method,ExceptionInfo *exception)
819 %
820 %  A description of each parameter follows:
821 %
822 %    o image: the image, which is replaced by indexed CLUT values
823 %
824 %    o clut_image: the color lookup table image for replacement color values.
825 %
826 %    o method: the pixel interpolation method.
827 %
828 %    o exception: return any errors or warnings in this structure.
829 %
830 */
ClutImage(Image * image,const Image * clut_image,const PixelInterpolateMethod method,ExceptionInfo * exception)831 MagickExport MagickBooleanType ClutImage(Image *image,const Image *clut_image,
832   const PixelInterpolateMethod method,ExceptionInfo *exception)
833 {
834 #define ClutImageTag  "Clut/Image"
835 
836   CacheView
837     *clut_view,
838     *image_view;
839 
840   MagickBooleanType
841     status;
842 
843   MagickOffsetType
844     progress;
845 
846   PixelInfo
847     *clut_map;
848 
849   register ssize_t
850     i;
851 
852   ssize_t adjust,
853     y;
854 
855   assert(image != (Image *) NULL);
856   assert(image->signature == MagickCoreSignature);
857   if (image->debug != MagickFalse)
858     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
859   assert(clut_image != (Image *) NULL);
860   assert(clut_image->signature == MagickCoreSignature);
861   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
862     return(MagickFalse);
863   if ((IsGrayColorspace(image->colorspace) != MagickFalse) &&
864       (IsGrayColorspace(clut_image->colorspace) == MagickFalse))
865     (void) SetImageColorspace(image,sRGBColorspace,exception);
866   clut_map=(PixelInfo *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*clut_map));
867   if (clut_map == (PixelInfo *) NULL)
868     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
869       image->filename);
870   /*
871     Clut image.
872   */
873   status=MagickTrue;
874   progress=0;
875   adjust=(ssize_t) (clut_image->interpolate == IntegerInterpolatePixel ? 0 : 1);
876   clut_view=AcquireVirtualCacheView(clut_image,exception);
877   for (i=0; i <= (ssize_t) MaxMap; i++)
878   {
879     GetPixelInfo(clut_image,clut_map+i);
880     status=InterpolatePixelInfo(clut_image,clut_view,method,
881       (double) i*(clut_image->columns-adjust)/MaxMap,(double) i*
882       (clut_image->rows-adjust)/MaxMap,clut_map+i,exception);
883     if (status == MagickFalse)
884       break;
885   }
886   clut_view=DestroyCacheView(clut_view);
887   image_view=AcquireAuthenticCacheView(image,exception);
888 #if defined(MAGICKCORE_OPENMP_SUPPORT)
889   #pragma omp parallel for schedule(static) shared(progress,status) \
890     magick_number_threads(image,image,image->rows,1)
891 #endif
892   for (y=0; y < (ssize_t) image->rows; y++)
893   {
894     PixelInfo
895       pixel;
896 
897     register Quantum
898       *magick_restrict q;
899 
900     register ssize_t
901       x;
902 
903     if (status == MagickFalse)
904       continue;
905     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
906     if (q == (Quantum *) NULL)
907       {
908         status=MagickFalse;
909         continue;
910       }
911     GetPixelInfo(image,&pixel);
912     for (x=0; x < (ssize_t) image->columns; x++)
913     {
914       PixelTrait
915         traits;
916 
917       GetPixelInfoPixel(image,q,&pixel);
918       traits=GetPixelChannelTraits(image,RedPixelChannel);
919       if ((traits & UpdatePixelTrait) != 0)
920         pixel.red=clut_map[ScaleQuantumToMap(ClampToQuantum(
921           pixel.red))].red;
922       traits=GetPixelChannelTraits(image,GreenPixelChannel);
923       if ((traits & UpdatePixelTrait) != 0)
924         pixel.green=clut_map[ScaleQuantumToMap(ClampToQuantum(
925           pixel.green))].green;
926       traits=GetPixelChannelTraits(image,BluePixelChannel);
927       if ((traits & UpdatePixelTrait) != 0)
928         pixel.blue=clut_map[ScaleQuantumToMap(ClampToQuantum(
929           pixel.blue))].blue;
930       traits=GetPixelChannelTraits(image,BlackPixelChannel);
931       if ((traits & UpdatePixelTrait) != 0)
932         pixel.black=clut_map[ScaleQuantumToMap(ClampToQuantum(
933           pixel.black))].black;
934       traits=GetPixelChannelTraits(image,AlphaPixelChannel);
935       if ((traits & UpdatePixelTrait) != 0)
936         pixel.alpha=clut_map[ScaleQuantumToMap(ClampToQuantum(
937           pixel.alpha))].alpha;
938       SetPixelViaPixelInfo(image,&pixel,q);
939       q+=GetPixelChannels(image);
940     }
941     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
942       status=MagickFalse;
943     if (image->progress_monitor != (MagickProgressMonitor) NULL)
944       {
945         MagickBooleanType
946           proceed;
947 
948 #if defined(MAGICKCORE_OPENMP_SUPPORT)
949         #pragma omp atomic
950 #endif
951         progress++;
952         proceed=SetImageProgress(image,ClutImageTag,progress,image->rows);
953         if (proceed == MagickFalse)
954           status=MagickFalse;
955       }
956   }
957   image_view=DestroyCacheView(image_view);
958   clut_map=(PixelInfo *) RelinquishMagickMemory(clut_map);
959   if ((clut_image->alpha_trait != UndefinedPixelTrait) &&
960       ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0))
961     (void) SetImageAlphaChannel(image,ActivateAlphaChannel,exception);
962   return(status);
963 }
964 
965 /*
966 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
967 %                                                                             %
968 %                                                                             %
969 %                                                                             %
970 %     C o l o r D e c i s i o n L i s t I m a g e                             %
971 %                                                                             %
972 %                                                                             %
973 %                                                                             %
974 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
975 %
976 %  ColorDecisionListImage() accepts a lightweight Color Correction Collection
977 %  (CCC) file which solely contains one or more color corrections and applies
978 %  the correction to the image.  Here is a sample CCC file:
979 %
980 %    <ColorCorrectionCollection xmlns="urn:ASC:CDL:v1.2">
981 %          <ColorCorrection id="cc03345">
982 %                <SOPNode>
983 %                     <Slope> 0.9 1.2 0.5 </Slope>
984 %                     <Offset> 0.4 -0.5 0.6 </Offset>
985 %                     <Power> 1.0 0.8 1.5 </Power>
986 %                </SOPNode>
987 %                <SATNode>
988 %                     <Saturation> 0.85 </Saturation>
989 %                </SATNode>
990 %          </ColorCorrection>
991 %    </ColorCorrectionCollection>
992 %
993 %  which includes the slop, offset, and power for each of the RGB channels
994 %  as well as the saturation.
995 %
996 %  The format of the ColorDecisionListImage method is:
997 %
998 %      MagickBooleanType ColorDecisionListImage(Image *image,
999 %        const char *color_correction_collection,ExceptionInfo *exception)
1000 %
1001 %  A description of each parameter follows:
1002 %
1003 %    o image: the image.
1004 %
1005 %    o color_correction_collection: the color correction collection in XML.
1006 %
1007 %    o exception: return any errors or warnings in this structure.
1008 %
1009 */
ColorDecisionListImage(Image * image,const char * color_correction_collection,ExceptionInfo * exception)1010 MagickExport MagickBooleanType ColorDecisionListImage(Image *image,
1011   const char *color_correction_collection,ExceptionInfo *exception)
1012 {
1013 #define ColorDecisionListCorrectImageTag  "ColorDecisionList/Image"
1014 
1015   typedef struct _Correction
1016   {
1017     double
1018       slope,
1019       offset,
1020       power;
1021   } Correction;
1022 
1023   typedef struct _ColorCorrection
1024   {
1025     Correction
1026       red,
1027       green,
1028       blue;
1029 
1030     double
1031       saturation;
1032   } ColorCorrection;
1033 
1034   CacheView
1035     *image_view;
1036 
1037   char
1038     token[MagickPathExtent];
1039 
1040   ColorCorrection
1041     color_correction;
1042 
1043   const char
1044     *content,
1045     *p;
1046 
1047   MagickBooleanType
1048     status;
1049 
1050   MagickOffsetType
1051     progress;
1052 
1053   PixelInfo
1054     *cdl_map;
1055 
1056   register ssize_t
1057     i;
1058 
1059   ssize_t
1060     y;
1061 
1062   XMLTreeInfo
1063     *cc,
1064     *ccc,
1065     *sat,
1066     *sop;
1067 
1068   /*
1069     Allocate and initialize cdl maps.
1070   */
1071   assert(image != (Image *) NULL);
1072   assert(image->signature == MagickCoreSignature);
1073   if (image->debug != MagickFalse)
1074     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1075   if (color_correction_collection == (const char *) NULL)
1076     return(MagickFalse);
1077   ccc=NewXMLTree((const char *) color_correction_collection,exception);
1078   if (ccc == (XMLTreeInfo *) NULL)
1079     return(MagickFalse);
1080   cc=GetXMLTreeChild(ccc,"ColorCorrection");
1081   if (cc == (XMLTreeInfo *) NULL)
1082     {
1083       ccc=DestroyXMLTree(ccc);
1084       return(MagickFalse);
1085     }
1086   color_correction.red.slope=1.0;
1087   color_correction.red.offset=0.0;
1088   color_correction.red.power=1.0;
1089   color_correction.green.slope=1.0;
1090   color_correction.green.offset=0.0;
1091   color_correction.green.power=1.0;
1092   color_correction.blue.slope=1.0;
1093   color_correction.blue.offset=0.0;
1094   color_correction.blue.power=1.0;
1095   color_correction.saturation=0.0;
1096   sop=GetXMLTreeChild(cc,"SOPNode");
1097   if (sop != (XMLTreeInfo *) NULL)
1098     {
1099       XMLTreeInfo
1100         *offset,
1101         *power,
1102         *slope;
1103 
1104       slope=GetXMLTreeChild(sop,"Slope");
1105       if (slope != (XMLTreeInfo *) NULL)
1106         {
1107           content=GetXMLTreeContent(slope);
1108           p=(const char *) content;
1109           for (i=0; (*p != '\0') && (i < 3); i++)
1110           {
1111             GetNextToken(p,&p,MagickPathExtent,token);
1112             if (*token == ',')
1113               GetNextToken(p,&p,MagickPathExtent,token);
1114             switch (i)
1115             {
1116               case 0:
1117               {
1118                 color_correction.red.slope=StringToDouble(token,(char **) NULL);
1119                 break;
1120               }
1121               case 1:
1122               {
1123                 color_correction.green.slope=StringToDouble(token,
1124                   (char **) NULL);
1125                 break;
1126               }
1127               case 2:
1128               {
1129                 color_correction.blue.slope=StringToDouble(token,
1130                   (char **) NULL);
1131                 break;
1132               }
1133             }
1134           }
1135         }
1136       offset=GetXMLTreeChild(sop,"Offset");
1137       if (offset != (XMLTreeInfo *) NULL)
1138         {
1139           content=GetXMLTreeContent(offset);
1140           p=(const char *) content;
1141           for (i=0; (*p != '\0') && (i < 3); i++)
1142           {
1143             GetNextToken(p,&p,MagickPathExtent,token);
1144             if (*token == ',')
1145               GetNextToken(p,&p,MagickPathExtent,token);
1146             switch (i)
1147             {
1148               case 0:
1149               {
1150                 color_correction.red.offset=StringToDouble(token,
1151                   (char **) NULL);
1152                 break;
1153               }
1154               case 1:
1155               {
1156                 color_correction.green.offset=StringToDouble(token,
1157                   (char **) NULL);
1158                 break;
1159               }
1160               case 2:
1161               {
1162                 color_correction.blue.offset=StringToDouble(token,
1163                   (char **) NULL);
1164                 break;
1165               }
1166             }
1167           }
1168         }
1169       power=GetXMLTreeChild(sop,"Power");
1170       if (power != (XMLTreeInfo *) NULL)
1171         {
1172           content=GetXMLTreeContent(power);
1173           p=(const char *) content;
1174           for (i=0; (*p != '\0') && (i < 3); i++)
1175           {
1176             GetNextToken(p,&p,MagickPathExtent,token);
1177             if (*token == ',')
1178               GetNextToken(p,&p,MagickPathExtent,token);
1179             switch (i)
1180             {
1181               case 0:
1182               {
1183                 color_correction.red.power=StringToDouble(token,(char **) NULL);
1184                 break;
1185               }
1186               case 1:
1187               {
1188                 color_correction.green.power=StringToDouble(token,
1189                   (char **) NULL);
1190                 break;
1191               }
1192               case 2:
1193               {
1194                 color_correction.blue.power=StringToDouble(token,
1195                   (char **) NULL);
1196                 break;
1197               }
1198             }
1199           }
1200         }
1201     }
1202   sat=GetXMLTreeChild(cc,"SATNode");
1203   if (sat != (XMLTreeInfo *) NULL)
1204     {
1205       XMLTreeInfo
1206         *saturation;
1207 
1208       saturation=GetXMLTreeChild(sat,"Saturation");
1209       if (saturation != (XMLTreeInfo *) NULL)
1210         {
1211           content=GetXMLTreeContent(saturation);
1212           p=(const char *) content;
1213           GetNextToken(p,&p,MagickPathExtent,token);
1214           color_correction.saturation=StringToDouble(token,(char **) NULL);
1215         }
1216     }
1217   ccc=DestroyXMLTree(ccc);
1218   if (image->debug != MagickFalse)
1219     {
1220       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1221         "  Color Correction Collection:");
1222       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1223         "  color_correction.red.slope: %g",color_correction.red.slope);
1224       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1225         "  color_correction.red.offset: %g",color_correction.red.offset);
1226       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1227         "  color_correction.red.power: %g",color_correction.red.power);
1228       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1229         "  color_correction.green.slope: %g",color_correction.green.slope);
1230       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1231         "  color_correction.green.offset: %g",color_correction.green.offset);
1232       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1233         "  color_correction.green.power: %g",color_correction.green.power);
1234       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1235         "  color_correction.blue.slope: %g",color_correction.blue.slope);
1236       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1237         "  color_correction.blue.offset: %g",color_correction.blue.offset);
1238       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1239         "  color_correction.blue.power: %g",color_correction.blue.power);
1240       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
1241         "  color_correction.saturation: %g",color_correction.saturation);
1242     }
1243   cdl_map=(PixelInfo *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*cdl_map));
1244   if (cdl_map == (PixelInfo *) NULL)
1245     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
1246       image->filename);
1247   for (i=0; i <= (ssize_t) MaxMap; i++)
1248   {
1249     cdl_map[i].red=(double) ScaleMapToQuantum((double)
1250       (MaxMap*(pow(color_correction.red.slope*i/MaxMap+
1251       color_correction.red.offset,color_correction.red.power))));
1252     cdl_map[i].green=(double) ScaleMapToQuantum((double)
1253       (MaxMap*(pow(color_correction.green.slope*i/MaxMap+
1254       color_correction.green.offset,color_correction.green.power))));
1255     cdl_map[i].blue=(double) ScaleMapToQuantum((double)
1256       (MaxMap*(pow(color_correction.blue.slope*i/MaxMap+
1257       color_correction.blue.offset,color_correction.blue.power))));
1258   }
1259   if (image->storage_class == PseudoClass)
1260     for (i=0; i < (ssize_t) image->colors; i++)
1261     {
1262       /*
1263         Apply transfer function to colormap.
1264       */
1265       double
1266         luma;
1267 
1268       luma=0.21267f*image->colormap[i].red+0.71526*image->colormap[i].green+
1269         0.07217f*image->colormap[i].blue;
1270       image->colormap[i].red=luma+color_correction.saturation*cdl_map[
1271         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))].red-luma;
1272       image->colormap[i].green=luma+color_correction.saturation*cdl_map[
1273         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].green))].green-luma;
1274       image->colormap[i].blue=luma+color_correction.saturation*cdl_map[
1275         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].blue))].blue-luma;
1276     }
1277   /*
1278     Apply transfer function to image.
1279   */
1280   status=MagickTrue;
1281   progress=0;
1282   image_view=AcquireAuthenticCacheView(image,exception);
1283 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1284   #pragma omp parallel for schedule(static) shared(progress,status) \
1285     magick_number_threads(image,image,image->rows,1)
1286 #endif
1287   for (y=0; y < (ssize_t) image->rows; y++)
1288   {
1289     double
1290       luma;
1291 
1292     register Quantum
1293       *magick_restrict q;
1294 
1295     register ssize_t
1296       x;
1297 
1298     if (status == MagickFalse)
1299       continue;
1300     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1301     if (q == (Quantum *) NULL)
1302       {
1303         status=MagickFalse;
1304         continue;
1305       }
1306     for (x=0; x < (ssize_t) image->columns; x++)
1307     {
1308       luma=0.21267f*GetPixelRed(image,q)+0.71526*GetPixelGreen(image,q)+
1309         0.07217f*GetPixelBlue(image,q);
1310       SetPixelRed(image,ClampToQuantum(luma+color_correction.saturation*
1311         (cdl_map[ScaleQuantumToMap(GetPixelRed(image,q))].red-luma)),q);
1312       SetPixelGreen(image,ClampToQuantum(luma+color_correction.saturation*
1313         (cdl_map[ScaleQuantumToMap(GetPixelGreen(image,q))].green-luma)),q);
1314       SetPixelBlue(image,ClampToQuantum(luma+color_correction.saturation*
1315         (cdl_map[ScaleQuantumToMap(GetPixelBlue(image,q))].blue-luma)),q);
1316       q+=GetPixelChannels(image);
1317     }
1318     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1319       status=MagickFalse;
1320     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1321       {
1322         MagickBooleanType
1323           proceed;
1324 
1325 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1326         #pragma omp atomic
1327 #endif
1328         progress++;
1329         proceed=SetImageProgress(image,ColorDecisionListCorrectImageTag,
1330           progress,image->rows);
1331         if (proceed == MagickFalse)
1332           status=MagickFalse;
1333       }
1334   }
1335   image_view=DestroyCacheView(image_view);
1336   cdl_map=(PixelInfo *) RelinquishMagickMemory(cdl_map);
1337   return(status);
1338 }
1339 
1340 /*
1341 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1342 %                                                                             %
1343 %                                                                             %
1344 %                                                                             %
1345 %     C o n t r a s t I m a g e                                               %
1346 %                                                                             %
1347 %                                                                             %
1348 %                                                                             %
1349 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1350 %
1351 %  ContrastImage() enhances the intensity differences between the lighter and
1352 %  darker elements of the image.  Set sharpen to a MagickTrue to increase the
1353 %  image contrast otherwise the contrast is reduced.
1354 %
1355 %  The format of the ContrastImage method is:
1356 %
1357 %      MagickBooleanType ContrastImage(Image *image,
1358 %        const MagickBooleanType sharpen,ExceptionInfo *exception)
1359 %
1360 %  A description of each parameter follows:
1361 %
1362 %    o image: the image.
1363 %
1364 %    o sharpen: Increase or decrease image contrast.
1365 %
1366 %    o exception: return any errors or warnings in this structure.
1367 %
1368 */
1369 
Contrast(const int sign,double * red,double * green,double * blue)1370 static void Contrast(const int sign,double *red,double *green,double *blue)
1371 {
1372   double
1373     brightness,
1374     hue,
1375     saturation;
1376 
1377   /*
1378     Enhance contrast: dark color become darker, light color become lighter.
1379   */
1380   assert(red != (double *) NULL);
1381   assert(green != (double *) NULL);
1382   assert(blue != (double *) NULL);
1383   hue=0.0;
1384   saturation=0.0;
1385   brightness=0.0;
1386   ConvertRGBToHSB(*red,*green,*blue,&hue,&saturation,&brightness);
1387   brightness+=0.5*sign*(0.5*(sin((double) (MagickPI*(brightness-0.5)))+1.0)-
1388     brightness);
1389   if (brightness > 1.0)
1390     brightness=1.0;
1391   else
1392     if (brightness < 0.0)
1393       brightness=0.0;
1394   ConvertHSBToRGB(hue,saturation,brightness,red,green,blue);
1395 }
1396 
ContrastImage(Image * image,const MagickBooleanType sharpen,ExceptionInfo * exception)1397 MagickExport MagickBooleanType ContrastImage(Image *image,
1398   const MagickBooleanType sharpen,ExceptionInfo *exception)
1399 {
1400 #define ContrastImageTag  "Contrast/Image"
1401 
1402   CacheView
1403     *image_view;
1404 
1405   int
1406     sign;
1407 
1408   MagickBooleanType
1409     status;
1410 
1411   MagickOffsetType
1412     progress;
1413 
1414   register ssize_t
1415     i;
1416 
1417   ssize_t
1418     y;
1419 
1420   assert(image != (Image *) NULL);
1421   assert(image->signature == MagickCoreSignature);
1422 #if defined(MAGICKCORE_OPENCL_SUPPORT)
1423   if (AccelerateContrastImage(image,sharpen,exception) != MagickFalse)
1424     return(MagickTrue);
1425 #endif
1426   if (image->debug != MagickFalse)
1427     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1428   sign=sharpen != MagickFalse ? 1 : -1;
1429   if (image->storage_class == PseudoClass)
1430     {
1431       /*
1432         Contrast enhance colormap.
1433       */
1434       for (i=0; i < (ssize_t) image->colors; i++)
1435       {
1436         double
1437           blue,
1438           green,
1439           red;
1440 
1441         red=(double) image->colormap[i].red;
1442         green=(double) image->colormap[i].green;
1443         blue=(double) image->colormap[i].blue;
1444         Contrast(sign,&red,&green,&blue);
1445         image->colormap[i].red=(MagickRealType) red;
1446         image->colormap[i].green=(MagickRealType) green;
1447         image->colormap[i].blue=(MagickRealType) blue;
1448       }
1449     }
1450   /*
1451     Contrast enhance image.
1452   */
1453   status=MagickTrue;
1454   progress=0;
1455   image_view=AcquireAuthenticCacheView(image,exception);
1456 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1457   #pragma omp parallel for schedule(static) shared(progress,status) \
1458     magick_number_threads(image,image,image->rows,1)
1459 #endif
1460   for (y=0; y < (ssize_t) image->rows; y++)
1461   {
1462     double
1463       blue,
1464       green,
1465       red;
1466 
1467     register Quantum
1468       *magick_restrict q;
1469 
1470     register ssize_t
1471       x;
1472 
1473     if (status == MagickFalse)
1474       continue;
1475     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1476     if (q == (Quantum *) NULL)
1477       {
1478         status=MagickFalse;
1479         continue;
1480       }
1481     for (x=0; x < (ssize_t) image->columns; x++)
1482     {
1483       red=(double) GetPixelRed(image,q);
1484       green=(double) GetPixelGreen(image,q);
1485       blue=(double) GetPixelBlue(image,q);
1486       Contrast(sign,&red,&green,&blue);
1487       SetPixelRed(image,ClampToQuantum(red),q);
1488       SetPixelGreen(image,ClampToQuantum(green),q);
1489       SetPixelBlue(image,ClampToQuantum(blue),q);
1490       q+=GetPixelChannels(image);
1491     }
1492     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1493       status=MagickFalse;
1494     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1495       {
1496         MagickBooleanType
1497           proceed;
1498 
1499 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1500         #pragma omp atomic
1501 #endif
1502         progress++;
1503         proceed=SetImageProgress(image,ContrastImageTag,progress,image->rows);
1504         if (proceed == MagickFalse)
1505           status=MagickFalse;
1506       }
1507   }
1508   image_view=DestroyCacheView(image_view);
1509   return(status);
1510 }
1511 
1512 /*
1513 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1514 %                                                                             %
1515 %                                                                             %
1516 %                                                                             %
1517 %     C o n t r a s t S t r e t c h I m a g e                                 %
1518 %                                                                             %
1519 %                                                                             %
1520 %                                                                             %
1521 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1522 %
1523 %  ContrastStretchImage() is a simple image enhancement technique that attempts
1524 %  to improve the contrast in an image by 'stretching' the range of intensity
1525 %  values it contains to span a desired range of values. It differs from the
1526 %  more sophisticated histogram equalization in that it can only apply a
1527 %  linear scaling function to the image pixel values.  As a result the
1528 %  'enhancement' is less harsh.
1529 %
1530 %  The format of the ContrastStretchImage method is:
1531 %
1532 %      MagickBooleanType ContrastStretchImage(Image *image,
1533 %        const char *levels,ExceptionInfo *exception)
1534 %
1535 %  A description of each parameter follows:
1536 %
1537 %    o image: the image.
1538 %
1539 %    o black_point: the black point.
1540 %
1541 %    o white_point: the white point.
1542 %
1543 %    o levels: Specify the levels where the black and white points have the
1544 %      range of 0 to number-of-pixels (e.g. 1%, 10x90%, etc.).
1545 %
1546 %    o exception: return any errors or warnings in this structure.
1547 %
1548 */
ContrastStretchImage(Image * image,const double black_point,const double white_point,ExceptionInfo * exception)1549 MagickExport MagickBooleanType ContrastStretchImage(Image *image,
1550   const double black_point,const double white_point,ExceptionInfo *exception)
1551 {
1552 #define MaxRange(color)  ((double) ScaleQuantumToMap((Quantum) (color)))
1553 #define ContrastStretchImageTag  "ContrastStretch/Image"
1554 
1555   CacheView
1556     *image_view;
1557 
1558   double
1559     *black,
1560     *histogram,
1561     *stretch_map,
1562     *white;
1563 
1564   MagickBooleanType
1565     status;
1566 
1567   MagickOffsetType
1568     progress;
1569 
1570   register ssize_t
1571     i;
1572 
1573   ssize_t
1574     y;
1575 
1576   /*
1577     Allocate histogram and stretch map.
1578   */
1579   assert(image != (Image *) NULL);
1580   assert(image->signature == MagickCoreSignature);
1581   if (image->debug != MagickFalse)
1582     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1583   if (SetImageGray(image,exception) != MagickFalse)
1584     (void) SetImageColorspace(image,GRAYColorspace,exception);
1585   black=(double *) AcquireQuantumMemory(MaxPixelChannels,sizeof(*black));
1586   white=(double *) AcquireQuantumMemory(MaxPixelChannels,sizeof(*white));
1587   histogram=(double *) AcquireQuantumMemory(MaxMap+1UL,MaxPixelChannels*
1588     sizeof(*histogram));
1589   stretch_map=(double *) AcquireQuantumMemory(MaxMap+1UL,MaxPixelChannels*
1590     sizeof(*stretch_map));
1591   if ((black == (double *) NULL) || (white == (double *) NULL) ||
1592       (histogram == (double *) NULL) || (stretch_map == (double *) NULL))
1593     {
1594       if (stretch_map != (double *) NULL)
1595         stretch_map=(double *) RelinquishMagickMemory(stretch_map);
1596       if (histogram != (double *) NULL)
1597         histogram=(double *) RelinquishMagickMemory(histogram);
1598       if (white != (double *) NULL)
1599         white=(double *) RelinquishMagickMemory(white);
1600       if (black != (double *) NULL)
1601         black=(double *) RelinquishMagickMemory(black);
1602       ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
1603         image->filename);
1604     }
1605   /*
1606     Form histogram.
1607   */
1608   status=MagickTrue;
1609   (void) memset(histogram,0,(MaxMap+1)*GetPixelChannels(image)*
1610     sizeof(*histogram));
1611   image_view=AcquireVirtualCacheView(image,exception);
1612   for (y=0; y < (ssize_t) image->rows; y++)
1613   {
1614     register const Quantum
1615       *magick_restrict p;
1616 
1617     register ssize_t
1618       x;
1619 
1620     if (status == MagickFalse)
1621       continue;
1622     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
1623     if (p == (const Quantum *) NULL)
1624       {
1625         status=MagickFalse;
1626         continue;
1627       }
1628     for (x=0; x < (ssize_t) image->columns; x++)
1629     {
1630       double
1631         pixel;
1632 
1633       pixel=GetPixelIntensity(image,p);
1634       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1635       {
1636         if (image->channel_mask != DefaultChannels)
1637           pixel=(double) p[i];
1638         histogram[GetPixelChannels(image)*ScaleQuantumToMap(
1639           ClampToQuantum(pixel))+i]++;
1640       }
1641       p+=GetPixelChannels(image);
1642     }
1643   }
1644   image_view=DestroyCacheView(image_view);
1645   /*
1646     Find the histogram boundaries by locating the black/white levels.
1647   */
1648   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1649   {
1650     double
1651       intensity;
1652 
1653     register ssize_t
1654       j;
1655 
1656     black[i]=0.0;
1657     white[i]=MaxRange(QuantumRange);
1658     intensity=0.0;
1659     for (j=0; j <= (ssize_t) MaxMap; j++)
1660     {
1661       intensity+=histogram[GetPixelChannels(image)*j+i];
1662       if (intensity > black_point)
1663         break;
1664     }
1665     black[i]=(double) j;
1666     intensity=0.0;
1667     for (j=(ssize_t) MaxMap; j != 0; j--)
1668     {
1669       intensity+=histogram[GetPixelChannels(image)*j+i];
1670       if (intensity > ((double) image->columns*image->rows-white_point))
1671         break;
1672     }
1673     white[i]=(double) j;
1674   }
1675   histogram=(double *) RelinquishMagickMemory(histogram);
1676   /*
1677     Stretch the histogram to create the stretched image mapping.
1678   */
1679   (void) memset(stretch_map,0,(MaxMap+1)*GetPixelChannels(image)*
1680     sizeof(*stretch_map));
1681   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1682   {
1683     register ssize_t
1684       j;
1685 
1686     for (j=0; j <= (ssize_t) MaxMap; j++)
1687     {
1688       double
1689         gamma;
1690 
1691       gamma=PerceptibleReciprocal(white[i]-black[i]);
1692       if (j < (ssize_t) black[i])
1693         stretch_map[GetPixelChannels(image)*j+i]=0.0;
1694       else
1695         if (j > (ssize_t) white[i])
1696           stretch_map[GetPixelChannels(image)*j+i]=(double) QuantumRange;
1697         else
1698           if (black[i] != white[i])
1699             stretch_map[GetPixelChannels(image)*j+i]=(double) ScaleMapToQuantum(
1700               (double) (MaxMap*gamma*(j-black[i])));
1701     }
1702   }
1703   if (image->storage_class == PseudoClass)
1704     {
1705       register ssize_t
1706         j;
1707 
1708       /*
1709         Stretch-contrast colormap.
1710       */
1711       for (j=0; j < (ssize_t) image->colors; j++)
1712       {
1713         if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1714           {
1715             i=GetPixelChannelOffset(image,RedPixelChannel);
1716             image->colormap[j].red=stretch_map[GetPixelChannels(image)*
1717               ScaleQuantumToMap(ClampToQuantum(image->colormap[j].red))+i];
1718           }
1719         if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1720           {
1721             i=GetPixelChannelOffset(image,GreenPixelChannel);
1722             image->colormap[j].green=stretch_map[GetPixelChannels(image)*
1723               ScaleQuantumToMap(ClampToQuantum(image->colormap[j].green))+i];
1724           }
1725         if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1726           {
1727             i=GetPixelChannelOffset(image,BluePixelChannel);
1728             image->colormap[j].blue=stretch_map[GetPixelChannels(image)*
1729               ScaleQuantumToMap(ClampToQuantum(image->colormap[j].blue))+i];
1730           }
1731         if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
1732           {
1733             i=GetPixelChannelOffset(image,AlphaPixelChannel);
1734             image->colormap[j].alpha=stretch_map[GetPixelChannels(image)*
1735               ScaleQuantumToMap(ClampToQuantum(image->colormap[j].alpha))+i];
1736           }
1737       }
1738     }
1739   /*
1740     Stretch-contrast image.
1741   */
1742   status=MagickTrue;
1743   progress=0;
1744   image_view=AcquireAuthenticCacheView(image,exception);
1745 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1746   #pragma omp parallel for schedule(static) shared(progress,status) \
1747     magick_number_threads(image,image,image->rows,1)
1748 #endif
1749   for (y=0; y < (ssize_t) image->rows; y++)
1750   {
1751     register Quantum
1752       *magick_restrict q;
1753 
1754     register ssize_t
1755       x;
1756 
1757     if (status == MagickFalse)
1758       continue;
1759     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1760     if (q == (Quantum *) NULL)
1761       {
1762         status=MagickFalse;
1763         continue;
1764       }
1765     for (x=0; x < (ssize_t) image->columns; x++)
1766     {
1767       register ssize_t
1768         j;
1769 
1770       for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
1771       {
1772         PixelChannel channel = GetPixelChannelChannel(image,j);
1773         PixelTrait traits = GetPixelChannelTraits(image,channel);
1774         if ((traits & UpdatePixelTrait) == 0)
1775           continue;
1776         if (black[j] == white[j])
1777           continue;
1778         q[j]=ClampToQuantum(stretch_map[GetPixelChannels(image)*
1779           ScaleQuantumToMap(q[j])+j]);
1780       }
1781       q+=GetPixelChannels(image);
1782     }
1783     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1784       status=MagickFalse;
1785     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1786       {
1787         MagickBooleanType
1788           proceed;
1789 
1790 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1791         #pragma omp atomic
1792 #endif
1793         progress++;
1794         proceed=SetImageProgress(image,ContrastStretchImageTag,progress,
1795           image->rows);
1796         if (proceed == MagickFalse)
1797           status=MagickFalse;
1798       }
1799   }
1800   image_view=DestroyCacheView(image_view);
1801   stretch_map=(double *) RelinquishMagickMemory(stretch_map);
1802   white=(double *) RelinquishMagickMemory(white);
1803   black=(double *) RelinquishMagickMemory(black);
1804   return(status);
1805 }
1806 
1807 /*
1808 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1809 %                                                                             %
1810 %                                                                             %
1811 %                                                                             %
1812 %     E n h a n c e I m a g e                                                 %
1813 %                                                                             %
1814 %                                                                             %
1815 %                                                                             %
1816 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1817 %
1818 %  EnhanceImage() applies a digital filter that improves the quality of a
1819 %  noisy image.
1820 %
1821 %  The format of the EnhanceImage method is:
1822 %
1823 %      Image *EnhanceImage(const Image *image,ExceptionInfo *exception)
1824 %
1825 %  A description of each parameter follows:
1826 %
1827 %    o image: the image.
1828 %
1829 %    o exception: return any errors or warnings in this structure.
1830 %
1831 */
EnhanceImage(const Image * image,ExceptionInfo * exception)1832 MagickExport Image *EnhanceImage(const Image *image,ExceptionInfo *exception)
1833 {
1834 #define EnhanceImageTag  "Enhance/Image"
1835 #define EnhancePixel(weight) \
1836   mean=QuantumScale*((double) GetPixelRed(image,r)+pixel.red)/2.0; \
1837   distance=QuantumScale*((double) GetPixelRed(image,r)-pixel.red); \
1838   distance_squared=(4.0+mean)*distance*distance; \
1839   mean=QuantumScale*((double) GetPixelGreen(image,r)+pixel.green)/2.0; \
1840   distance=QuantumScale*((double) GetPixelGreen(image,r)-pixel.green); \
1841   distance_squared+=(7.0-mean)*distance*distance; \
1842   mean=QuantumScale*((double) GetPixelBlue(image,r)+pixel.blue)/2.0; \
1843   distance=QuantumScale*((double) GetPixelBlue(image,r)-pixel.blue); \
1844   distance_squared+=(5.0-mean)*distance*distance; \
1845   mean=QuantumScale*((double) GetPixelBlack(image,r)+pixel.black)/2.0; \
1846   distance=QuantumScale*((double) GetPixelBlack(image,r)-pixel.black); \
1847   distance_squared+=(5.0-mean)*distance*distance; \
1848   mean=QuantumScale*((double) GetPixelAlpha(image,r)+pixel.alpha)/2.0; \
1849   distance=QuantumScale*((double) GetPixelAlpha(image,r)-pixel.alpha); \
1850   distance_squared+=(5.0-mean)*distance*distance; \
1851   if (distance_squared < 0.069) \
1852     { \
1853       aggregate.red+=(weight)*GetPixelRed(image,r); \
1854       aggregate.green+=(weight)*GetPixelGreen(image,r); \
1855       aggregate.blue+=(weight)*GetPixelBlue(image,r); \
1856       aggregate.black+=(weight)*GetPixelBlack(image,r); \
1857       aggregate.alpha+=(weight)*GetPixelAlpha(image,r); \
1858       total_weight+=(weight); \
1859     } \
1860   r+=GetPixelChannels(image);
1861 
1862   CacheView
1863     *enhance_view,
1864     *image_view;
1865 
1866   Image
1867     *enhance_image;
1868 
1869   MagickBooleanType
1870     status;
1871 
1872   MagickOffsetType
1873     progress;
1874 
1875   ssize_t
1876     y;
1877 
1878   /*
1879     Initialize enhanced image attributes.
1880   */
1881   assert(image != (const Image *) NULL);
1882   assert(image->signature == MagickCoreSignature);
1883   if (image->debug != MagickFalse)
1884     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1885   assert(exception != (ExceptionInfo *) NULL);
1886   assert(exception->signature == MagickCoreSignature);
1887   enhance_image=CloneImage(image,0,0,MagickTrue,
1888     exception);
1889   if (enhance_image == (Image *) NULL)
1890     return((Image *) NULL);
1891   if (SetImageStorageClass(enhance_image,DirectClass,exception) == MagickFalse)
1892     {
1893       enhance_image=DestroyImage(enhance_image);
1894       return((Image *) NULL);
1895     }
1896   /*
1897     Enhance image.
1898   */
1899   status=MagickTrue;
1900   progress=0;
1901   image_view=AcquireVirtualCacheView(image,exception);
1902   enhance_view=AcquireAuthenticCacheView(enhance_image,exception);
1903 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1904   #pragma omp parallel for schedule(static) shared(progress,status) \
1905     magick_number_threads(image,enhance_image,image->rows,1)
1906 #endif
1907   for (y=0; y < (ssize_t) image->rows; y++)
1908   {
1909     PixelInfo
1910       pixel;
1911 
1912     register const Quantum
1913       *magick_restrict p;
1914 
1915     register Quantum
1916       *magick_restrict q;
1917 
1918     register ssize_t
1919       x;
1920 
1921     ssize_t
1922       center;
1923 
1924     if (status == MagickFalse)
1925       continue;
1926     p=GetCacheViewVirtualPixels(image_view,-2,y-2,image->columns+4,5,exception);
1927     q=QueueCacheViewAuthenticPixels(enhance_view,0,y,enhance_image->columns,1,
1928       exception);
1929     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
1930       {
1931         status=MagickFalse;
1932         continue;
1933       }
1934     center=(ssize_t) GetPixelChannels(image)*(2*(image->columns+4)+2);
1935     GetPixelInfo(image,&pixel);
1936     for (x=0; x < (ssize_t) image->columns; x++)
1937     {
1938       double
1939         distance,
1940         distance_squared,
1941         mean,
1942         total_weight;
1943 
1944       PixelInfo
1945         aggregate;
1946 
1947       register const Quantum
1948         *magick_restrict r;
1949 
1950       GetPixelInfo(image,&aggregate);
1951       total_weight=0.0;
1952       GetPixelInfoPixel(image,p+center,&pixel);
1953       r=p;
1954       EnhancePixel(5.0); EnhancePixel(8.0); EnhancePixel(10.0);
1955         EnhancePixel(8.0); EnhancePixel(5.0);
1956       r=p+GetPixelChannels(image)*(image->columns+4);
1957       EnhancePixel(8.0); EnhancePixel(20.0); EnhancePixel(40.0);
1958         EnhancePixel(20.0); EnhancePixel(8.0);
1959       r=p+2*GetPixelChannels(image)*(image->columns+4);
1960       EnhancePixel(10.0); EnhancePixel(40.0); EnhancePixel(80.0);
1961         EnhancePixel(40.0); EnhancePixel(10.0);
1962       r=p+3*GetPixelChannels(image)*(image->columns+4);
1963       EnhancePixel(8.0); EnhancePixel(20.0); EnhancePixel(40.0);
1964         EnhancePixel(20.0); EnhancePixel(8.0);
1965       r=p+4*GetPixelChannels(image)*(image->columns+4);
1966       EnhancePixel(5.0); EnhancePixel(8.0); EnhancePixel(10.0);
1967         EnhancePixel(8.0); EnhancePixel(5.0);
1968       if (total_weight > MagickEpsilon)
1969         {
1970           pixel.red=((aggregate.red+total_weight/2.0)/total_weight);
1971           pixel.green=((aggregate.green+total_weight/2.0)/total_weight);
1972           pixel.blue=((aggregate.blue+total_weight/2.0)/total_weight);
1973           pixel.black=((aggregate.black+total_weight/2.0)/total_weight);
1974           pixel.alpha=((aggregate.alpha+total_weight/2.0)/total_weight);
1975         }
1976       SetPixelViaPixelInfo(image,&pixel,q);
1977       p+=GetPixelChannels(image);
1978       q+=GetPixelChannels(enhance_image);
1979     }
1980     if (SyncCacheViewAuthenticPixels(enhance_view,exception) == MagickFalse)
1981       status=MagickFalse;
1982     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1983       {
1984         MagickBooleanType
1985           proceed;
1986 
1987 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1988         #pragma omp atomic
1989 #endif
1990         progress++;
1991         proceed=SetImageProgress(image,EnhanceImageTag,progress,image->rows);
1992         if (proceed == MagickFalse)
1993           status=MagickFalse;
1994       }
1995   }
1996   enhance_view=DestroyCacheView(enhance_view);
1997   image_view=DestroyCacheView(image_view);
1998   if (status == MagickFalse)
1999     enhance_image=DestroyImage(enhance_image);
2000   return(enhance_image);
2001 }
2002 
2003 /*
2004 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2005 %                                                                             %
2006 %                                                                             %
2007 %                                                                             %
2008 %     E q u a l i z e I m a g e                                               %
2009 %                                                                             %
2010 %                                                                             %
2011 %                                                                             %
2012 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2013 %
2014 %  EqualizeImage() applies a histogram equalization to the image.
2015 %
2016 %  The format of the EqualizeImage method is:
2017 %
2018 %      MagickBooleanType EqualizeImage(Image *image,ExceptionInfo *exception)
2019 %
2020 %  A description of each parameter follows:
2021 %
2022 %    o image: the image.
2023 %
2024 %    o exception: return any errors or warnings in this structure.
2025 %
2026 */
EqualizeImage(Image * image,ExceptionInfo * exception)2027 MagickExport MagickBooleanType EqualizeImage(Image *image,
2028   ExceptionInfo *exception)
2029 {
2030 #define EqualizeImageTag  "Equalize/Image"
2031 
2032   CacheView
2033     *image_view;
2034 
2035   double
2036     black[CompositePixelChannel+1],
2037     *equalize_map,
2038     *histogram,
2039     *map,
2040     white[CompositePixelChannel+1];
2041 
2042   MagickBooleanType
2043     status;
2044 
2045   MagickOffsetType
2046     progress;
2047 
2048   register ssize_t
2049     i;
2050 
2051   ssize_t
2052     y;
2053 
2054   /*
2055     Allocate and initialize histogram arrays.
2056   */
2057   assert(image != (Image *) NULL);
2058   assert(image->signature == MagickCoreSignature);
2059 #if defined(MAGICKCORE_OPENCL_SUPPORT)
2060   if (AccelerateEqualizeImage(image,exception) != MagickFalse)
2061     return(MagickTrue);
2062 #endif
2063   if (image->debug != MagickFalse)
2064     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2065   equalize_map=(double *) AcquireQuantumMemory(MaxMap+1UL,MaxPixelChannels*
2066     sizeof(*equalize_map));
2067   histogram=(double *) AcquireQuantumMemory(MaxMap+1UL,MaxPixelChannels*
2068     sizeof(*histogram));
2069   map=(double *) AcquireQuantumMemory(MaxMap+1UL,MaxPixelChannels*sizeof(*map));
2070   if ((equalize_map == (double *) NULL) || (histogram == (double *) NULL) ||
2071       (map == (double *) NULL))
2072     {
2073       if (map != (double *) NULL)
2074         map=(double *) RelinquishMagickMemory(map);
2075       if (histogram != (double *) NULL)
2076         histogram=(double *) RelinquishMagickMemory(histogram);
2077       if (equalize_map != (double *) NULL)
2078         equalize_map=(double *) RelinquishMagickMemory(equalize_map);
2079       ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
2080         image->filename);
2081     }
2082   /*
2083     Form histogram.
2084   */
2085   status=MagickTrue;
2086   (void) memset(histogram,0,(MaxMap+1)*GetPixelChannels(image)*
2087     sizeof(*histogram));
2088   image_view=AcquireVirtualCacheView(image,exception);
2089   for (y=0; y < (ssize_t) image->rows; y++)
2090   {
2091     register const Quantum
2092       *magick_restrict p;
2093 
2094     register ssize_t
2095       x;
2096 
2097     if (status == MagickFalse)
2098       continue;
2099     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
2100     if (p == (const Quantum *) NULL)
2101       {
2102         status=MagickFalse;
2103         continue;
2104       }
2105     for (x=0; x < (ssize_t) image->columns; x++)
2106     {
2107       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2108       {
2109         double
2110           intensity;
2111 
2112         intensity=(double) p[i];
2113         if ((image->channel_mask & SyncChannels) != 0)
2114           intensity=GetPixelIntensity(image,p);
2115         histogram[GetPixelChannels(image)*ScaleQuantumToMap(
2116           ClampToQuantum(intensity))+i]++;
2117       }
2118       p+=GetPixelChannels(image);
2119     }
2120   }
2121   image_view=DestroyCacheView(image_view);
2122   /*
2123     Integrate the histogram to get the equalization map.
2124   */
2125   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2126   {
2127     double
2128       intensity;
2129 
2130     register ssize_t
2131       j;
2132 
2133     intensity=0.0;
2134     for (j=0; j <= (ssize_t) MaxMap; j++)
2135     {
2136       intensity+=histogram[GetPixelChannels(image)*j+i];
2137       map[GetPixelChannels(image)*j+i]=intensity;
2138     }
2139   }
2140   (void) memset(equalize_map,0,(MaxMap+1)*GetPixelChannels(image)*
2141     sizeof(*equalize_map));
2142   (void) memset(black,0,sizeof(*black));
2143   (void) memset(white,0,sizeof(*white));
2144   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2145   {
2146     register ssize_t
2147       j;
2148 
2149     black[i]=map[i];
2150     white[i]=map[GetPixelChannels(image)*MaxMap+i];
2151     if (black[i] != white[i])
2152       for (j=0; j <= (ssize_t) MaxMap; j++)
2153         equalize_map[GetPixelChannels(image)*j+i]=(double)
2154           ScaleMapToQuantum((double) ((MaxMap*(map[
2155           GetPixelChannels(image)*j+i]-black[i]))/(white[i]-black[i])));
2156   }
2157   histogram=(double *) RelinquishMagickMemory(histogram);
2158   map=(double *) RelinquishMagickMemory(map);
2159   if (image->storage_class == PseudoClass)
2160     {
2161       register ssize_t
2162         j;
2163 
2164       /*
2165         Equalize colormap.
2166       */
2167       for (j=0; j < (ssize_t) image->colors; j++)
2168       {
2169         if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2170           {
2171             PixelChannel channel = GetPixelChannelChannel(image,
2172               RedPixelChannel);
2173             if (black[channel] != white[channel])
2174               image->colormap[j].red=equalize_map[GetPixelChannels(image)*
2175                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].red))+
2176                 channel];
2177           }
2178         if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2179           {
2180             PixelChannel channel = GetPixelChannelChannel(image,
2181               GreenPixelChannel);
2182             if (black[channel] != white[channel])
2183               image->colormap[j].green=equalize_map[GetPixelChannels(image)*
2184                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].green))+
2185                 channel];
2186           }
2187         if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2188           {
2189             PixelChannel channel = GetPixelChannelChannel(image,
2190               BluePixelChannel);
2191             if (black[channel] != white[channel])
2192               image->colormap[j].blue=equalize_map[GetPixelChannels(image)*
2193                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].blue))+
2194                 channel];
2195           }
2196         if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
2197           {
2198             PixelChannel channel = GetPixelChannelChannel(image,
2199               AlphaPixelChannel);
2200             if (black[channel] != white[channel])
2201               image->colormap[j].alpha=equalize_map[GetPixelChannels(image)*
2202                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].alpha))+
2203                 channel];
2204           }
2205       }
2206     }
2207   /*
2208     Equalize image.
2209   */
2210   progress=0;
2211   image_view=AcquireAuthenticCacheView(image,exception);
2212 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2213   #pragma omp parallel for schedule(static) shared(progress,status) \
2214     magick_number_threads(image,image,image->rows,1)
2215 #endif
2216   for (y=0; y < (ssize_t) image->rows; y++)
2217   {
2218     register Quantum
2219       *magick_restrict q;
2220 
2221     register ssize_t
2222       x;
2223 
2224     if (status == MagickFalse)
2225       continue;
2226     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2227     if (q == (Quantum *) NULL)
2228       {
2229         status=MagickFalse;
2230         continue;
2231       }
2232     for (x=0; x < (ssize_t) image->columns; x++)
2233     {
2234       register ssize_t
2235         j;
2236 
2237       for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
2238       {
2239         PixelChannel channel = GetPixelChannelChannel(image,j);
2240         PixelTrait traits = GetPixelChannelTraits(image,channel);
2241         if (((traits & UpdatePixelTrait) == 0) || (black[j] == white[j]))
2242           continue;
2243         q[j]=ClampToQuantum(equalize_map[GetPixelChannels(image)*
2244           ScaleQuantumToMap(q[j])+j]);
2245       }
2246       q+=GetPixelChannels(image);
2247     }
2248     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2249       status=MagickFalse;
2250     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2251       {
2252         MagickBooleanType
2253           proceed;
2254 
2255 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2256         #pragma omp atomic
2257 #endif
2258         progress++;
2259         proceed=SetImageProgress(image,EqualizeImageTag,progress,image->rows);
2260         if (proceed == MagickFalse)
2261           status=MagickFalse;
2262       }
2263   }
2264   image_view=DestroyCacheView(image_view);
2265   equalize_map=(double *) RelinquishMagickMemory(equalize_map);
2266   return(status);
2267 }
2268 
2269 /*
2270 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2271 %                                                                             %
2272 %                                                                             %
2273 %                                                                             %
2274 %     G a m m a I m a g e                                                     %
2275 %                                                                             %
2276 %                                                                             %
2277 %                                                                             %
2278 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2279 %
2280 %  GammaImage() gamma-corrects a particular image channel.  The same
2281 %  image viewed on different devices will have perceptual differences in the
2282 %  way the image's intensities are represented on the screen.  Specify
2283 %  individual gamma levels for the red, green, and blue channels, or adjust
2284 %  all three with the gamma parameter.  Values typically range from 0.8 to 2.3.
2285 %
2286 %  You can also reduce the influence of a particular channel with a gamma
2287 %  value of 0.
2288 %
2289 %  The format of the GammaImage method is:
2290 %
2291 %      MagickBooleanType GammaImage(Image *image,const double gamma,
2292 %        ExceptionInfo *exception)
2293 %
2294 %  A description of each parameter follows:
2295 %
2296 %    o image: the image.
2297 %
2298 %    o level: the image gamma as a string (e.g. 1.6,1.2,1.0).
2299 %
2300 %    o gamma: the image gamma.
2301 %
2302 */
2303 
gamma_pow(const double value,const double gamma)2304 static inline double gamma_pow(const double value,const double gamma)
2305 {
2306   return(value < 0.0 ? value : pow(value,gamma));
2307 }
2308 
GammaImage(Image * image,const double gamma,ExceptionInfo * exception)2309 MagickExport MagickBooleanType GammaImage(Image *image,const double gamma,
2310   ExceptionInfo *exception)
2311 {
2312 #define GammaCorrectImageTag  "GammaCorrect/Image"
2313 
2314   CacheView
2315     *image_view;
2316 
2317   MagickBooleanType
2318     status;
2319 
2320   MagickOffsetType
2321     progress;
2322 
2323   Quantum
2324     *gamma_map;
2325 
2326   register ssize_t
2327     i;
2328 
2329   ssize_t
2330     y;
2331 
2332   /*
2333     Allocate and initialize gamma maps.
2334   */
2335   assert(image != (Image *) NULL);
2336   assert(image->signature == MagickCoreSignature);
2337   if (image->debug != MagickFalse)
2338     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2339   if (gamma == 1.0)
2340     return(MagickTrue);
2341   gamma_map=(Quantum *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*gamma_map));
2342   if (gamma_map == (Quantum *) NULL)
2343     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
2344       image->filename);
2345   (void) memset(gamma_map,0,(MaxMap+1)*sizeof(*gamma_map));
2346   if (gamma != 0.0)
2347     for (i=0; i <= (ssize_t) MaxMap; i++)
2348       gamma_map[i]=ScaleMapToQuantum((double) (MaxMap*pow((double) i/
2349         MaxMap,1.0/gamma)));
2350   if (image->storage_class == PseudoClass)
2351     for (i=0; i < (ssize_t) image->colors; i++)
2352     {
2353       /*
2354         Gamma-correct colormap.
2355       */
2356       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2357         image->colormap[i].red=(double) gamma_map[ScaleQuantumToMap(
2358           ClampToQuantum(image->colormap[i].red))];
2359       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2360         image->colormap[i].green=(double) gamma_map[ScaleQuantumToMap(
2361           ClampToQuantum(image->colormap[i].green))];
2362       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2363         image->colormap[i].blue=(double) gamma_map[ScaleQuantumToMap(
2364           ClampToQuantum(image->colormap[i].blue))];
2365       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
2366         image->colormap[i].alpha=(double) gamma_map[ScaleQuantumToMap(
2367           ClampToQuantum(image->colormap[i].alpha))];
2368     }
2369   /*
2370     Gamma-correct image.
2371   */
2372   status=MagickTrue;
2373   progress=0;
2374   image_view=AcquireAuthenticCacheView(image,exception);
2375 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2376   #pragma omp parallel for schedule(static) shared(progress,status) \
2377     magick_number_threads(image,image,image->rows,1)
2378 #endif
2379   for (y=0; y < (ssize_t) image->rows; y++)
2380   {
2381     register Quantum
2382       *magick_restrict q;
2383 
2384     register ssize_t
2385       x;
2386 
2387     if (status == MagickFalse)
2388       continue;
2389     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2390     if (q == (Quantum *) NULL)
2391       {
2392         status=MagickFalse;
2393         continue;
2394       }
2395     for (x=0; x < (ssize_t) image->columns; x++)
2396     {
2397       register ssize_t
2398         j;
2399 
2400       for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
2401       {
2402         PixelChannel channel = GetPixelChannelChannel(image,j);
2403         PixelTrait traits = GetPixelChannelTraits(image,channel);
2404         if ((traits & UpdatePixelTrait) == 0)
2405           continue;
2406         q[j]=gamma_map[ScaleQuantumToMap(ClampToQuantum((MagickRealType)
2407           q[j]))];
2408       }
2409       q+=GetPixelChannels(image);
2410     }
2411     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2412       status=MagickFalse;
2413     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2414       {
2415         MagickBooleanType
2416           proceed;
2417 
2418 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2419         #pragma omp atomic
2420 #endif
2421         progress++;
2422         proceed=SetImageProgress(image,GammaCorrectImageTag,progress, image->rows);
2423         if (proceed == MagickFalse)
2424           status=MagickFalse;
2425       }
2426   }
2427   image_view=DestroyCacheView(image_view);
2428   gamma_map=(Quantum *) RelinquishMagickMemory(gamma_map);
2429   if (image->gamma != 0.0)
2430     image->gamma*=gamma;
2431   return(status);
2432 }
2433 
2434 /*
2435 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2436 %                                                                             %
2437 %                                                                             %
2438 %                                                                             %
2439 %     G r a y s c a l e I m a g e                                             %
2440 %                                                                             %
2441 %                                                                             %
2442 %                                                                             %
2443 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2444 %
2445 %  GrayscaleImage() converts the image to grayscale.
2446 %
2447 %  The format of the GrayscaleImage method is:
2448 %
2449 %      MagickBooleanType GrayscaleImage(Image *image,
2450 %        const PixelIntensityMethod method ,ExceptionInfo *exception)
2451 %
2452 %  A description of each parameter follows:
2453 %
2454 %    o image: the image.
2455 %
2456 %    o method: the pixel intensity method.
2457 %
2458 %    o exception: return any errors or warnings in this structure.
2459 %
2460 */
GrayscaleImage(Image * image,const PixelIntensityMethod method,ExceptionInfo * exception)2461 MagickExport MagickBooleanType GrayscaleImage(Image *image,
2462   const PixelIntensityMethod method,ExceptionInfo *exception)
2463 {
2464 #define GrayscaleImageTag  "Grayscale/Image"
2465 
2466   CacheView
2467     *image_view;
2468 
2469   MagickBooleanType
2470     status;
2471 
2472   MagickOffsetType
2473     progress;
2474 
2475   ssize_t
2476     y;
2477 
2478   assert(image != (Image *) NULL);
2479   assert(image->signature == MagickCoreSignature);
2480   if (image->debug != MagickFalse)
2481     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2482   if (image->storage_class == PseudoClass)
2483     {
2484       if (SyncImage(image,exception) == MagickFalse)
2485         return(MagickFalse);
2486       if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
2487         return(MagickFalse);
2488     }
2489 #if defined(MAGICKCORE_OPENCL_SUPPORT)
2490   if (AccelerateGrayscaleImage(image,method,exception) != MagickFalse)
2491     {
2492       image->intensity=method;
2493       image->type=GrayscaleType;
2494       if ((method == Rec601LuminancePixelIntensityMethod) ||
2495           (method == Rec709LuminancePixelIntensityMethod))
2496         return(SetImageColorspace(image,LinearGRAYColorspace,exception));
2497       return(SetImageColorspace(image,GRAYColorspace,exception));
2498     }
2499 #endif
2500   /*
2501     Grayscale image.
2502   */
2503   status=MagickTrue;
2504   progress=0;
2505   image_view=AcquireAuthenticCacheView(image,exception);
2506 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2507   #pragma omp parallel for schedule(static) shared(progress,status) \
2508     magick_number_threads(image,image,image->rows,1)
2509 #endif
2510   for (y=0; y < (ssize_t) image->rows; y++)
2511   {
2512     register Quantum
2513       *magick_restrict q;
2514 
2515     register ssize_t
2516       x;
2517 
2518     if (status == MagickFalse)
2519       continue;
2520     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2521     if (q == (Quantum *) NULL)
2522       {
2523         status=MagickFalse;
2524         continue;
2525       }
2526     for (x=0; x < (ssize_t) image->columns; x++)
2527     {
2528       MagickRealType
2529         blue,
2530         green,
2531         red,
2532         intensity;
2533 
2534       red=(MagickRealType) GetPixelRed(image,q);
2535       green=(MagickRealType) GetPixelGreen(image,q);
2536       blue=(MagickRealType) GetPixelBlue(image,q);
2537       intensity=0.0;
2538       switch (method)
2539       {
2540         case AveragePixelIntensityMethod:
2541         {
2542           intensity=(red+green+blue)/3.0;
2543           break;
2544         }
2545         case BrightnessPixelIntensityMethod:
2546         {
2547           intensity=MagickMax(MagickMax(red,green),blue);
2548           break;
2549         }
2550         case LightnessPixelIntensityMethod:
2551         {
2552           intensity=(MagickMin(MagickMin(red,green),blue)+
2553             MagickMax(MagickMax(red,green),blue))/2.0;
2554           break;
2555         }
2556         case MSPixelIntensityMethod:
2557         {
2558           intensity=(MagickRealType) (((double) red*red+green*green+
2559             blue*blue)/3.0);
2560           break;
2561         }
2562         case Rec601LumaPixelIntensityMethod:
2563         {
2564           if (image->colorspace == RGBColorspace)
2565             {
2566               red=EncodePixelGamma(red);
2567               green=EncodePixelGamma(green);
2568               blue=EncodePixelGamma(blue);
2569             }
2570           intensity=0.298839*red+0.586811*green+0.114350*blue;
2571           break;
2572         }
2573         case Rec601LuminancePixelIntensityMethod:
2574         {
2575           if (image->colorspace == sRGBColorspace)
2576             {
2577               red=DecodePixelGamma(red);
2578               green=DecodePixelGamma(green);
2579               blue=DecodePixelGamma(blue);
2580             }
2581           intensity=0.298839*red+0.586811*green+0.114350*blue;
2582           break;
2583         }
2584         case Rec709LumaPixelIntensityMethod:
2585         default:
2586         {
2587           if (image->colorspace == RGBColorspace)
2588             {
2589               red=EncodePixelGamma(red);
2590               green=EncodePixelGamma(green);
2591               blue=EncodePixelGamma(blue);
2592             }
2593           intensity=0.212656*red+0.715158*green+0.072186*blue;
2594           break;
2595         }
2596         case Rec709LuminancePixelIntensityMethod:
2597         {
2598           if (image->colorspace == sRGBColorspace)
2599             {
2600               red=DecodePixelGamma(red);
2601               green=DecodePixelGamma(green);
2602               blue=DecodePixelGamma(blue);
2603             }
2604           intensity=0.212656*red+0.715158*green+0.072186*blue;
2605           break;
2606         }
2607         case RMSPixelIntensityMethod:
2608         {
2609           intensity=(MagickRealType) (sqrt((double) red*red+green*green+
2610             blue*blue)/sqrt(3.0));
2611           break;
2612         }
2613       }
2614       SetPixelGray(image,ClampToQuantum(intensity),q);
2615       q+=GetPixelChannels(image);
2616     }
2617     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2618       status=MagickFalse;
2619     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2620       {
2621         MagickBooleanType
2622           proceed;
2623 
2624 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2625         #pragma omp atomic
2626 #endif
2627         progress++;
2628         proceed=SetImageProgress(image,GrayscaleImageTag,progress,image->rows);
2629         if (proceed == MagickFalse)
2630           status=MagickFalse;
2631       }
2632   }
2633   image_view=DestroyCacheView(image_view);
2634   image->intensity=method;
2635   image->type=GrayscaleType;
2636   if ((method == Rec601LuminancePixelIntensityMethod) ||
2637       (method == Rec709LuminancePixelIntensityMethod))
2638     return(SetImageColorspace(image,LinearGRAYColorspace,exception));
2639   return(SetImageColorspace(image,GRAYColorspace,exception));
2640 }
2641 
2642 /*
2643 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2644 %                                                                             %
2645 %                                                                             %
2646 %                                                                             %
2647 %     H a l d C l u t I m a g e                                               %
2648 %                                                                             %
2649 %                                                                             %
2650 %                                                                             %
2651 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2652 %
2653 %  HaldClutImage() applies a Hald color lookup table to the image.  A Hald
2654 %  color lookup table is a 3-dimensional color cube mapped to 2 dimensions.
2655 %  Create it with the HALD coder.  You can apply any color transformation to
2656 %  the Hald image and then use this method to apply the transform to the
2657 %  image.
2658 %
2659 %  The format of the HaldClutImage method is:
2660 %
2661 %      MagickBooleanType HaldClutImage(Image *image,Image *hald_image,
2662 %        ExceptionInfo *exception)
2663 %
2664 %  A description of each parameter follows:
2665 %
2666 %    o image: the image, which is replaced by indexed CLUT values
2667 %
2668 %    o hald_image: the color lookup table image for replacement color values.
2669 %
2670 %    o exception: return any errors or warnings in this structure.
2671 %
2672 */
HaldClutImage(Image * image,const Image * hald_image,ExceptionInfo * exception)2673 MagickExport MagickBooleanType HaldClutImage(Image *image,
2674   const Image *hald_image,ExceptionInfo *exception)
2675 {
2676 #define HaldClutImageTag  "Clut/Image"
2677 
2678   typedef struct _HaldInfo
2679   {
2680     double
2681       x,
2682       y,
2683       z;
2684   } HaldInfo;
2685 
2686   CacheView
2687     *hald_view,
2688     *image_view;
2689 
2690   double
2691     width;
2692 
2693   MagickBooleanType
2694     status;
2695 
2696   MagickOffsetType
2697     progress;
2698 
2699   PixelInfo
2700     zero;
2701 
2702   size_t
2703     cube_size,
2704     length,
2705     level;
2706 
2707   ssize_t
2708     y;
2709 
2710   assert(image != (Image *) NULL);
2711   assert(image->signature == MagickCoreSignature);
2712   if (image->debug != MagickFalse)
2713     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2714   assert(hald_image != (Image *) NULL);
2715   assert(hald_image->signature == MagickCoreSignature);
2716   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
2717     return(MagickFalse);
2718   if (image->alpha_trait == UndefinedPixelTrait)
2719     (void) SetImageAlphaChannel(image,OpaqueAlphaChannel,exception);
2720   /*
2721     Hald clut image.
2722   */
2723   status=MagickTrue;
2724   progress=0;
2725   length=(size_t) MagickMin((MagickRealType) hald_image->columns,
2726     (MagickRealType) hald_image->rows);
2727   for (level=2; (level*level*level) < length; level++) ;
2728   level*=level;
2729   cube_size=level*level;
2730   width=(double) hald_image->columns;
2731   GetPixelInfo(hald_image,&zero);
2732   hald_view=AcquireVirtualCacheView(hald_image,exception);
2733   image_view=AcquireAuthenticCacheView(image,exception);
2734 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2735   #pragma omp parallel for schedule(static) shared(progress,status) \
2736     magick_number_threads(image,image,image->rows,1)
2737 #endif
2738   for (y=0; y < (ssize_t) image->rows; y++)
2739   {
2740     register Quantum
2741       *magick_restrict q;
2742 
2743     register ssize_t
2744       x;
2745 
2746     if (status == MagickFalse)
2747       continue;
2748     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2749     if (q == (Quantum *) NULL)
2750       {
2751         status=MagickFalse;
2752         continue;
2753       }
2754     for (x=0; x < (ssize_t) image->columns; x++)
2755     {
2756       double
2757         offset;
2758 
2759       HaldInfo
2760         point;
2761 
2762       PixelInfo
2763         pixel,
2764         pixel1,
2765         pixel2,
2766         pixel3,
2767         pixel4;
2768 
2769       point.x=QuantumScale*(level-1.0)*GetPixelRed(image,q);
2770       point.y=QuantumScale*(level-1.0)*GetPixelGreen(image,q);
2771       point.z=QuantumScale*(level-1.0)*GetPixelBlue(image,q);
2772       offset=point.x+level*floor(point.y)+cube_size*floor(point.z);
2773       point.x-=floor(point.x);
2774       point.y-=floor(point.y);
2775       point.z-=floor(point.z);
2776       pixel1=zero;
2777       status=InterpolatePixelInfo(hald_image,hald_view,hald_image->interpolate,
2778         fmod(offset,width),floor(offset/width),&pixel1,exception);
2779       if (status == MagickFalse)
2780         break;
2781       pixel2=zero;
2782       status=InterpolatePixelInfo(hald_image,hald_view,hald_image->interpolate,
2783         fmod(offset+level,width),floor((offset+level)/width),&pixel2,exception);
2784       if (status == MagickFalse)
2785         break;
2786       pixel3=zero;
2787       CompositePixelInfoAreaBlend(&pixel1,pixel1.alpha,&pixel2,pixel2.alpha,
2788         point.y,&pixel3);
2789       offset+=cube_size;
2790       status=InterpolatePixelInfo(hald_image,hald_view,hald_image->interpolate,
2791         fmod(offset,width),floor(offset/width),&pixel1,exception);
2792       if (status == MagickFalse)
2793         break;
2794       status=InterpolatePixelInfo(hald_image,hald_view,hald_image->interpolate,
2795         fmod(offset+level,width),floor((offset+level)/width),&pixel2,exception);
2796       if (status == MagickFalse)
2797         break;
2798       pixel4=zero;
2799       CompositePixelInfoAreaBlend(&pixel1,pixel1.alpha,&pixel2,pixel2.alpha,
2800         point.y,&pixel4);
2801       pixel=zero;
2802       CompositePixelInfoAreaBlend(&pixel3,pixel3.alpha,&pixel4,pixel4.alpha,
2803         point.z,&pixel);
2804       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2805         SetPixelRed(image,ClampToQuantum(pixel.red),q);
2806       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2807         SetPixelGreen(image,ClampToQuantum(pixel.green),q);
2808       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2809         SetPixelBlue(image,ClampToQuantum(pixel.blue),q);
2810       if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
2811           (image->colorspace == CMYKColorspace))
2812         SetPixelBlack(image,ClampToQuantum(pixel.black),q);
2813       if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
2814           (image->alpha_trait != UndefinedPixelTrait))
2815         SetPixelAlpha(image,ClampToQuantum(pixel.alpha),q);
2816       q+=GetPixelChannels(image);
2817     }
2818     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2819       status=MagickFalse;
2820     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2821       {
2822         MagickBooleanType
2823           proceed;
2824 
2825 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2826         #pragma omp atomic
2827 #endif
2828         progress++;
2829         proceed=SetImageProgress(image,HaldClutImageTag,progress,image->rows);
2830         if (proceed == MagickFalse)
2831           status=MagickFalse;
2832       }
2833   }
2834   hald_view=DestroyCacheView(hald_view);
2835   image_view=DestroyCacheView(image_view);
2836   return(status);
2837 }
2838 
2839 /*
2840 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2841 %                                                                             %
2842 %                                                                             %
2843 %                                                                             %
2844 %     L e v e l I m a g e                                                     %
2845 %                                                                             %
2846 %                                                                             %
2847 %                                                                             %
2848 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2849 %
2850 %  LevelImage() adjusts the levels of a particular image channel by
2851 %  scaling the colors falling between specified white and black points to
2852 %  the full available quantum range.
2853 %
2854 %  The parameters provided represent the black, and white points.  The black
2855 %  point specifies the darkest color in the image. Colors darker than the
2856 %  black point are set to zero.  White point specifies the lightest color in
2857 %  the image.  Colors brighter than the white point are set to the maximum
2858 %  quantum value.
2859 %
2860 %  If a '!' flag is given, map black and white colors to the given levels
2861 %  rather than mapping those levels to black and white.  See
2862 %  LevelizeImage() below.
2863 %
2864 %  Gamma specifies a gamma correction to apply to the image.
2865 %
2866 %  The format of the LevelImage method is:
2867 %
2868 %      MagickBooleanType LevelImage(Image *image,const double black_point,
2869 %        const double white_point,const double gamma,ExceptionInfo *exception)
2870 %
2871 %  A description of each parameter follows:
2872 %
2873 %    o image: the image.
2874 %
2875 %    o black_point: The level to map zero (black) to.
2876 %
2877 %    o white_point: The level to map QuantumRange (white) to.
2878 %
2879 %    o exception: return any errors or warnings in this structure.
2880 %
2881 */
2882 
LevelPixel(const double black_point,const double white_point,const double gamma,const double pixel)2883 static inline double LevelPixel(const double black_point,
2884   const double white_point,const double gamma,const double pixel)
2885 {
2886   double
2887     level_pixel,
2888     scale;
2889 
2890   scale=PerceptibleReciprocal(white_point-black_point);
2891   level_pixel=QuantumRange*gamma_pow(scale*((double) pixel-black_point),
2892     1.0/gamma);
2893   return(level_pixel);
2894 }
2895 
LevelImage(Image * image,const double black_point,const double white_point,const double gamma,ExceptionInfo * exception)2896 MagickExport MagickBooleanType LevelImage(Image *image,const double black_point,
2897   const double white_point,const double gamma,ExceptionInfo *exception)
2898 {
2899 #define LevelImageTag  "Level/Image"
2900 
2901   CacheView
2902     *image_view;
2903 
2904   MagickBooleanType
2905     status;
2906 
2907   MagickOffsetType
2908     progress;
2909 
2910   register ssize_t
2911     i;
2912 
2913   ssize_t
2914     y;
2915 
2916   /*
2917     Allocate and initialize levels map.
2918   */
2919   assert(image != (Image *) NULL);
2920   assert(image->signature == MagickCoreSignature);
2921   if (image->debug != MagickFalse)
2922     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2923   if (image->storage_class == PseudoClass)
2924     for (i=0; i < (ssize_t) image->colors; i++)
2925     {
2926       /*
2927         Level colormap.
2928       */
2929       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2930         image->colormap[i].red=(double) ClampToQuantum(LevelPixel(black_point,
2931           white_point,gamma,image->colormap[i].red));
2932       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2933         image->colormap[i].green=(double) ClampToQuantum(LevelPixel(black_point,
2934           white_point,gamma,image->colormap[i].green));
2935       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2936         image->colormap[i].blue=(double) ClampToQuantum(LevelPixel(black_point,
2937           white_point,gamma,image->colormap[i].blue));
2938       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
2939         image->colormap[i].alpha=(double) ClampToQuantum(LevelPixel(black_point,
2940           white_point,gamma,image->colormap[i].alpha));
2941     }
2942   /*
2943     Level image.
2944   */
2945   status=MagickTrue;
2946   progress=0;
2947   image_view=AcquireAuthenticCacheView(image,exception);
2948 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2949   #pragma omp parallel for schedule(static) shared(progress,status) \
2950     magick_number_threads(image,image,image->rows,1)
2951 #endif
2952   for (y=0; y < (ssize_t) image->rows; y++)
2953   {
2954     register Quantum
2955       *magick_restrict q;
2956 
2957     register ssize_t
2958       x;
2959 
2960     if (status == MagickFalse)
2961       continue;
2962     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2963     if (q == (Quantum *) NULL)
2964       {
2965         status=MagickFalse;
2966         continue;
2967       }
2968     for (x=0; x < (ssize_t) image->columns; x++)
2969     {
2970       register ssize_t
2971         j;
2972 
2973       for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
2974       {
2975         PixelChannel channel = GetPixelChannelChannel(image,j);
2976         PixelTrait traits = GetPixelChannelTraits(image,channel);
2977         if ((traits & UpdatePixelTrait) == 0)
2978           continue;
2979         q[j]=ClampToQuantum(LevelPixel(black_point,white_point,gamma,
2980           (double) q[j]));
2981       }
2982       q+=GetPixelChannels(image);
2983     }
2984     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2985       status=MagickFalse;
2986     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2987       {
2988         MagickBooleanType
2989           proceed;
2990 
2991 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2992         #pragma omp atomic
2993 #endif
2994         progress++;
2995         proceed=SetImageProgress(image,LevelImageTag,progress,image->rows);
2996         if (proceed == MagickFalse)
2997           status=MagickFalse;
2998       }
2999   }
3000   image_view=DestroyCacheView(image_view);
3001   (void) ClampImage(image,exception);
3002   return(status);
3003 }
3004 
3005 /*
3006 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3007 %                                                                             %
3008 %                                                                             %
3009 %                                                                             %
3010 %     L e v e l i z e I m a g e                                               %
3011 %                                                                             %
3012 %                                                                             %
3013 %                                                                             %
3014 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3015 %
3016 %  LevelizeImage() applies the reversed LevelImage() operation to just
3017 %  the specific channels specified.  It compresses the full range of color
3018 %  values, so that they lie between the given black and white points. Gamma is
3019 %  applied before the values are mapped.
3020 %
3021 %  LevelizeImage() can be called with by using a +level command line
3022 %  API option, or using a '!' on a -level or LevelImage() geometry string.
3023 %
3024 %  It can be used to de-contrast a greyscale image to the exact levels
3025 %  specified.  Or by using specific levels for each channel of an image you
3026 %  can convert a gray-scale image to any linear color gradient, according to
3027 %  those levels.
3028 %
3029 %  The format of the LevelizeImage method is:
3030 %
3031 %      MagickBooleanType LevelizeImage(Image *image,const double black_point,
3032 %        const double white_point,const double gamma,ExceptionInfo *exception)
3033 %
3034 %  A description of each parameter follows:
3035 %
3036 %    o image: the image.
3037 %
3038 %    o black_point: The level to map zero (black) to.
3039 %
3040 %    o white_point: The level to map QuantumRange (white) to.
3041 %
3042 %    o gamma: adjust gamma by this factor before mapping values.
3043 %
3044 %    o exception: return any errors or warnings in this structure.
3045 %
3046 */
LevelizeImage(Image * image,const double black_point,const double white_point,const double gamma,ExceptionInfo * exception)3047 MagickExport MagickBooleanType LevelizeImage(Image *image,
3048   const double black_point,const double white_point,const double gamma,
3049   ExceptionInfo *exception)
3050 {
3051 #define LevelizeImageTag  "Levelize/Image"
3052 #define LevelizeValue(x) ClampToQuantum(((MagickRealType) gamma_pow((double) \
3053   (QuantumScale*(x)),gamma))*(white_point-black_point)+black_point)
3054 
3055   CacheView
3056     *image_view;
3057 
3058   MagickBooleanType
3059     status;
3060 
3061   MagickOffsetType
3062     progress;
3063 
3064   register ssize_t
3065     i;
3066 
3067   ssize_t
3068     y;
3069 
3070   /*
3071     Allocate and initialize levels map.
3072   */
3073   assert(image != (Image *) NULL);
3074   assert(image->signature == MagickCoreSignature);
3075   if (image->debug != MagickFalse)
3076     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3077   if (image->storage_class == PseudoClass)
3078     for (i=0; i < (ssize_t) image->colors; i++)
3079     {
3080       /*
3081         Level colormap.
3082       */
3083       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
3084         image->colormap[i].red=(double) LevelizeValue(image->colormap[i].red);
3085       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
3086         image->colormap[i].green=(double) LevelizeValue(
3087           image->colormap[i].green);
3088       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
3089         image->colormap[i].blue=(double) LevelizeValue(image->colormap[i].blue);
3090       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
3091         image->colormap[i].alpha=(double) LevelizeValue(
3092           image->colormap[i].alpha);
3093     }
3094   /*
3095     Level image.
3096   */
3097   status=MagickTrue;
3098   progress=0;
3099   image_view=AcquireAuthenticCacheView(image,exception);
3100 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3101   #pragma omp parallel for schedule(static) shared(progress,status) \
3102     magick_number_threads(image,image,image->rows,1)
3103 #endif
3104   for (y=0; y < (ssize_t) image->rows; y++)
3105   {
3106     register Quantum
3107       *magick_restrict q;
3108 
3109     register ssize_t
3110       x;
3111 
3112     if (status == MagickFalse)
3113       continue;
3114     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
3115     if (q == (Quantum *) NULL)
3116       {
3117         status=MagickFalse;
3118         continue;
3119       }
3120     for (x=0; x < (ssize_t) image->columns; x++)
3121     {
3122       register ssize_t
3123         j;
3124 
3125       for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
3126       {
3127         PixelChannel channel = GetPixelChannelChannel(image,j);
3128         PixelTrait traits = GetPixelChannelTraits(image,channel);
3129         if ((traits & UpdatePixelTrait) == 0)
3130           continue;
3131         q[j]=LevelizeValue(q[j]);
3132       }
3133       q+=GetPixelChannels(image);
3134     }
3135     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
3136       status=MagickFalse;
3137     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3138       {
3139         MagickBooleanType
3140           proceed;
3141 
3142 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3143         #pragma omp atomic
3144 #endif
3145         progress++;
3146         proceed=SetImageProgress(image,LevelizeImageTag,progress,image->rows);
3147         if (proceed == MagickFalse)
3148           status=MagickFalse;
3149       }
3150   }
3151   image_view=DestroyCacheView(image_view);
3152   return(status);
3153 }
3154 
3155 /*
3156 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3157 %                                                                             %
3158 %                                                                             %
3159 %                                                                             %
3160 %     L e v e l I m a g e C o l o r s                                         %
3161 %                                                                             %
3162 %                                                                             %
3163 %                                                                             %
3164 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3165 %
3166 %  LevelImageColors() maps the given color to "black" and "white" values,
3167 %  linearly spreading out the colors, and level values on a channel by channel
3168 %  bases, as per LevelImage().  The given colors allows you to specify
3169 %  different level ranges for each of the color channels separately.
3170 %
3171 %  If the boolean 'invert' is set true the image values will modifyed in the
3172 %  reverse direction. That is any existing "black" and "white" colors in the
3173 %  image will become the color values given, with all other values compressed
3174 %  appropriatally.  This effectivally maps a greyscale gradient into the given
3175 %  color gradient.
3176 %
3177 %  The format of the LevelImageColors method is:
3178 %
3179 %    MagickBooleanType LevelImageColors(Image *image,
3180 %      const PixelInfo *black_color,const PixelInfo *white_color,
3181 %      const MagickBooleanType invert,ExceptionInfo *exception)
3182 %
3183 %  A description of each parameter follows:
3184 %
3185 %    o image: the image.
3186 %
3187 %    o black_color: The color to map black to/from
3188 %
3189 %    o white_point: The color to map white to/from
3190 %
3191 %    o invert: if true map the colors (levelize), rather than from (level)
3192 %
3193 %    o exception: return any errors or warnings in this structure.
3194 %
3195 */
LevelImageColors(Image * image,const PixelInfo * black_color,const PixelInfo * white_color,const MagickBooleanType invert,ExceptionInfo * exception)3196 MagickExport MagickBooleanType LevelImageColors(Image *image,
3197   const PixelInfo *black_color,const PixelInfo *white_color,
3198   const MagickBooleanType invert,ExceptionInfo *exception)
3199 {
3200   ChannelType
3201     channel_mask;
3202 
3203   MagickStatusType
3204     status;
3205 
3206   /*
3207     Allocate and initialize levels map.
3208   */
3209   assert(image != (Image *) NULL);
3210   assert(image->signature == MagickCoreSignature);
3211   if (image->debug != MagickFalse)
3212     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3213   if ((IsGrayColorspace(image->colorspace) != MagickFalse) &&
3214       ((IsGrayColorspace(black_color->colorspace) == MagickFalse) ||
3215        (IsGrayColorspace(white_color->colorspace) == MagickFalse)))
3216     (void) SetImageColorspace(image,sRGBColorspace,exception);
3217   status=MagickTrue;
3218   if (invert == MagickFalse)
3219     {
3220       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
3221         {
3222           channel_mask=SetImageChannelMask(image,RedChannel);
3223           status&=LevelImage(image,black_color->red,white_color->red,1.0,
3224             exception);
3225           (void) SetImageChannelMask(image,channel_mask);
3226         }
3227       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
3228         {
3229           channel_mask=SetImageChannelMask(image,GreenChannel);
3230           status&=LevelImage(image,black_color->green,white_color->green,1.0,
3231             exception);
3232           (void) SetImageChannelMask(image,channel_mask);
3233         }
3234       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
3235         {
3236           channel_mask=SetImageChannelMask(image,BlueChannel);
3237           status&=LevelImage(image,black_color->blue,white_color->blue,1.0,
3238             exception);
3239           (void) SetImageChannelMask(image,channel_mask);
3240         }
3241       if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
3242           (image->colorspace == CMYKColorspace))
3243         {
3244           channel_mask=SetImageChannelMask(image,BlackChannel);
3245           status&=LevelImage(image,black_color->black,white_color->black,1.0,
3246             exception);
3247           (void) SetImageChannelMask(image,channel_mask);
3248         }
3249       if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
3250           (image->alpha_trait != UndefinedPixelTrait))
3251         {
3252           channel_mask=SetImageChannelMask(image,AlphaChannel);
3253           status&=LevelImage(image,black_color->alpha,white_color->alpha,1.0,
3254             exception);
3255           (void) SetImageChannelMask(image,channel_mask);
3256         }
3257     }
3258   else
3259     {
3260       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
3261         {
3262           channel_mask=SetImageChannelMask(image,RedChannel);
3263           status&=LevelizeImage(image,black_color->red,white_color->red,1.0,
3264             exception);
3265           (void) SetImageChannelMask(image,channel_mask);
3266         }
3267       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
3268         {
3269           channel_mask=SetImageChannelMask(image,GreenChannel);
3270           status&=LevelizeImage(image,black_color->green,white_color->green,1.0,
3271             exception);
3272           (void) SetImageChannelMask(image,channel_mask);
3273         }
3274       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
3275         {
3276           channel_mask=SetImageChannelMask(image,BlueChannel);
3277           status&=LevelizeImage(image,black_color->blue,white_color->blue,1.0,
3278             exception);
3279           (void) SetImageChannelMask(image,channel_mask);
3280         }
3281       if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
3282           (image->colorspace == CMYKColorspace))
3283         {
3284           channel_mask=SetImageChannelMask(image,BlackChannel);
3285           status&=LevelizeImage(image,black_color->black,white_color->black,1.0,
3286             exception);
3287           (void) SetImageChannelMask(image,channel_mask);
3288         }
3289       if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
3290           (image->alpha_trait != UndefinedPixelTrait))
3291         {
3292           channel_mask=SetImageChannelMask(image,AlphaChannel);
3293           status&=LevelizeImage(image,black_color->alpha,white_color->alpha,1.0,
3294             exception);
3295           (void) SetImageChannelMask(image,channel_mask);
3296         }
3297     }
3298   return(status != 0 ? MagickTrue : MagickFalse);
3299 }
3300 
3301 /*
3302 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3303 %                                                                             %
3304 %                                                                             %
3305 %                                                                             %
3306 %     L i n e a r S t r e t c h I m a g e                                     %
3307 %                                                                             %
3308 %                                                                             %
3309 %                                                                             %
3310 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3311 %
3312 %  LinearStretchImage() discards any pixels below the black point and above
3313 %  the white point and levels the remaining pixels.
3314 %
3315 %  The format of the LinearStretchImage method is:
3316 %
3317 %      MagickBooleanType LinearStretchImage(Image *image,
3318 %        const double black_point,const double white_point,
3319 %        ExceptionInfo *exception)
3320 %
3321 %  A description of each parameter follows:
3322 %
3323 %    o image: the image.
3324 %
3325 %    o black_point: the black point.
3326 %
3327 %    o white_point: the white point.
3328 %
3329 %    o exception: return any errors or warnings in this structure.
3330 %
3331 */
LinearStretchImage(Image * image,const double black_point,const double white_point,ExceptionInfo * exception)3332 MagickExport MagickBooleanType LinearStretchImage(Image *image,
3333   const double black_point,const double white_point,ExceptionInfo *exception)
3334 {
3335 #define LinearStretchImageTag  "LinearStretch/Image"
3336 
3337   CacheView
3338     *image_view;
3339 
3340   double
3341     *histogram,
3342     intensity;
3343 
3344   MagickBooleanType
3345     status;
3346 
3347   ssize_t
3348     black,
3349     white,
3350     y;
3351 
3352   /*
3353     Allocate histogram and linear map.
3354   */
3355   assert(image != (Image *) NULL);
3356   assert(image->signature == MagickCoreSignature);
3357   histogram=(double *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*histogram));
3358   if (histogram == (double *) NULL)
3359     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
3360       image->filename);
3361   /*
3362     Form histogram.
3363   */
3364   (void) memset(histogram,0,(MaxMap+1)*sizeof(*histogram));
3365   image_view=AcquireVirtualCacheView(image,exception);
3366   for (y=0; y < (ssize_t) image->rows; y++)
3367   {
3368     register const Quantum
3369       *magick_restrict p;
3370 
3371     register ssize_t
3372       x;
3373 
3374     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
3375     if (p == (const Quantum *) NULL)
3376       break;
3377     for (x=0; x < (ssize_t) image->columns; x++)
3378     {
3379       intensity=GetPixelIntensity(image,p);
3380       histogram[ScaleQuantumToMap(ClampToQuantum(intensity))]++;
3381       p+=GetPixelChannels(image);
3382     }
3383   }
3384   image_view=DestroyCacheView(image_view);
3385   /*
3386     Find the histogram boundaries by locating the black and white point levels.
3387   */
3388   intensity=0.0;
3389   for (black=0; black < (ssize_t) MaxMap; black++)
3390   {
3391     intensity+=histogram[black];
3392     if (intensity >= black_point)
3393       break;
3394   }
3395   intensity=0.0;
3396   for (white=(ssize_t) MaxMap; white != 0; white--)
3397   {
3398     intensity+=histogram[white];
3399     if (intensity >= white_point)
3400       break;
3401   }
3402   histogram=(double *) RelinquishMagickMemory(histogram);
3403   status=LevelImage(image,(double) ScaleMapToQuantum((MagickRealType) black),
3404     (double) ScaleMapToQuantum((MagickRealType) white),1.0,exception);
3405   return(status);
3406 }
3407 
3408 
3409 /*
3410 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3411 %                                                                             %
3412 %                                                                             %
3413 %                                                                             %
3414 %     M o d u l a t e I m a g e                                               %
3415 %                                                                             %
3416 %                                                                             %
3417 %                                                                             %
3418 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3419 %
3420 %  ModulateImage() lets you control the brightness, saturation, and hue
3421 %  of an image.  Modulate represents the brightness, saturation, and hue
3422 %  as one parameter (e.g. 90,150,100).  If the image colorspace is HSL, the
3423 %  modulation is lightness, saturation, and hue.  For HWB, use blackness,
3424 %  whiteness, and hue. And for HCL, use chrome, luma, and hue.
3425 %
3426 %  The format of the ModulateImage method is:
3427 %
3428 %      MagickBooleanType ModulateImage(Image *image,const char *modulate,
3429 %        ExceptionInfo *exception)
3430 %
3431 %  A description of each parameter follows:
3432 %
3433 %    o image: the image.
3434 %
3435 %    o modulate: Define the percent change in brightness, saturation, and hue.
3436 %
3437 %    o exception: return any errors or warnings in this structure.
3438 %
3439 */
3440 
ModulateHCL(const double percent_hue,const double percent_chroma,const double percent_luma,double * red,double * green,double * blue)3441 static inline void ModulateHCL(const double percent_hue,
3442   const double percent_chroma,const double percent_luma,double *red,
3443   double *green,double *blue)
3444 {
3445   double
3446     hue,
3447     luma,
3448     chroma;
3449 
3450   /*
3451     Increase or decrease color luma, chroma, or hue.
3452   */
3453   ConvertRGBToHCL(*red,*green,*blue,&hue,&chroma,&luma);
3454   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3455   chroma*=0.01*percent_chroma;
3456   luma*=0.01*percent_luma;
3457   ConvertHCLToRGB(hue,chroma,luma,red,green,blue);
3458 }
3459 
ModulateHCLp(const double percent_hue,const double percent_chroma,const double percent_luma,double * red,double * green,double * blue)3460 static inline void ModulateHCLp(const double percent_hue,
3461   const double percent_chroma,const double percent_luma,double *red,
3462   double *green,double *blue)
3463 {
3464   double
3465     hue,
3466     luma,
3467     chroma;
3468 
3469   /*
3470     Increase or decrease color luma, chroma, or hue.
3471   */
3472   ConvertRGBToHCLp(*red,*green,*blue,&hue,&chroma,&luma);
3473   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3474   chroma*=0.01*percent_chroma;
3475   luma*=0.01*percent_luma;
3476   ConvertHCLpToRGB(hue,chroma,luma,red,green,blue);
3477 }
3478 
ModulateHSB(const double percent_hue,const double percent_saturation,const double percent_brightness,double * red,double * green,double * blue)3479 static inline void ModulateHSB(const double percent_hue,
3480   const double percent_saturation,const double percent_brightness,double *red,
3481   double *green,double *blue)
3482 {
3483   double
3484     brightness,
3485     hue,
3486     saturation;
3487 
3488   /*
3489     Increase or decrease color brightness, saturation, or hue.
3490   */
3491   ConvertRGBToHSB(*red,*green,*blue,&hue,&saturation,&brightness);
3492   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3493   saturation*=0.01*percent_saturation;
3494   brightness*=0.01*percent_brightness;
3495   ConvertHSBToRGB(hue,saturation,brightness,red,green,blue);
3496 }
3497 
ModulateHSI(const double percent_hue,const double percent_saturation,const double percent_intensity,double * red,double * green,double * blue)3498 static inline void ModulateHSI(const double percent_hue,
3499   const double percent_saturation,const double percent_intensity,double *red,
3500   double *green,double *blue)
3501 {
3502   double
3503     intensity,
3504     hue,
3505     saturation;
3506 
3507   /*
3508     Increase or decrease color intensity, saturation, or hue.
3509   */
3510   ConvertRGBToHSI(*red,*green,*blue,&hue,&saturation,&intensity);
3511   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3512   saturation*=0.01*percent_saturation;
3513   intensity*=0.01*percent_intensity;
3514   ConvertHSIToRGB(hue,saturation,intensity,red,green,blue);
3515 }
3516 
ModulateHSL(const double percent_hue,const double percent_saturation,const double percent_lightness,double * red,double * green,double * blue)3517 static inline void ModulateHSL(const double percent_hue,
3518   const double percent_saturation,const double percent_lightness,double *red,
3519   double *green,double *blue)
3520 {
3521   double
3522     hue,
3523     lightness,
3524     saturation;
3525 
3526   /*
3527     Increase or decrease color lightness, saturation, or hue.
3528   */
3529   ConvertRGBToHSL(*red,*green,*blue,&hue,&saturation,&lightness);
3530   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3531   saturation*=0.01*percent_saturation;
3532   lightness*=0.01*percent_lightness;
3533   ConvertHSLToRGB(hue,saturation,lightness,red,green,blue);
3534 }
3535 
ModulateHSV(const double percent_hue,const double percent_saturation,const double percent_value,double * red,double * green,double * blue)3536 static inline void ModulateHSV(const double percent_hue,
3537   const double percent_saturation,const double percent_value,double *red,
3538   double *green,double *blue)
3539 {
3540   double
3541     hue,
3542     saturation,
3543     value;
3544 
3545   /*
3546     Increase or decrease color value, saturation, or hue.
3547   */
3548   ConvertRGBToHSV(*red,*green,*blue,&hue,&saturation,&value);
3549   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3550   saturation*=0.01*percent_saturation;
3551   value*=0.01*percent_value;
3552   ConvertHSVToRGB(hue,saturation,value,red,green,blue);
3553 }
3554 
ModulateHWB(const double percent_hue,const double percent_whiteness,const double percent_blackness,double * red,double * green,double * blue)3555 static inline void ModulateHWB(const double percent_hue,
3556   const double percent_whiteness,const double percent_blackness,double *red,
3557   double *green,double *blue)
3558 {
3559   double
3560     blackness,
3561     hue,
3562     whiteness;
3563 
3564   /*
3565     Increase or decrease color blackness, whiteness, or hue.
3566   */
3567   ConvertRGBToHWB(*red,*green,*blue,&hue,&whiteness,&blackness);
3568   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3569   blackness*=0.01*percent_blackness;
3570   whiteness*=0.01*percent_whiteness;
3571   ConvertHWBToRGB(hue,whiteness,blackness,red,green,blue);
3572 }
3573 
ModulateLCHab(const double percent_luma,const double percent_chroma,const double percent_hue,double * red,double * green,double * blue)3574 static inline void ModulateLCHab(const double percent_luma,
3575   const double percent_chroma,const double percent_hue,double *red,
3576   double *green,double *blue)
3577 {
3578   double
3579     hue,
3580     luma,
3581     chroma;
3582 
3583   /*
3584     Increase or decrease color luma, chroma, or hue.
3585   */
3586   ConvertRGBToLCHab(*red,*green,*blue,&luma,&chroma,&hue);
3587   luma*=0.01*percent_luma;
3588   chroma*=0.01*percent_chroma;
3589   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3590   ConvertLCHabToRGB(luma,chroma,hue,red,green,blue);
3591 }
3592 
ModulateLCHuv(const double percent_luma,const double percent_chroma,const double percent_hue,double * red,double * green,double * blue)3593 static inline void ModulateLCHuv(const double percent_luma,
3594   const double percent_chroma,const double percent_hue,double *red,
3595   double *green,double *blue)
3596 {
3597   double
3598     hue,
3599     luma,
3600     chroma;
3601 
3602   /*
3603     Increase or decrease color luma, chroma, or hue.
3604   */
3605   ConvertRGBToLCHuv(*red,*green,*blue,&luma,&chroma,&hue);
3606   luma*=0.01*percent_luma;
3607   chroma*=0.01*percent_chroma;
3608   hue+=fmod((percent_hue-100.0),200.0)/200.0;
3609   ConvertLCHuvToRGB(luma,chroma,hue,red,green,blue);
3610 }
3611 
ModulateImage(Image * image,const char * modulate,ExceptionInfo * exception)3612 MagickExport MagickBooleanType ModulateImage(Image *image,const char *modulate,
3613   ExceptionInfo *exception)
3614 {
3615 #define ModulateImageTag  "Modulate/Image"
3616 
3617   CacheView
3618     *image_view;
3619 
3620   ColorspaceType
3621     colorspace;
3622 
3623   const char
3624     *artifact;
3625 
3626   double
3627     percent_brightness,
3628     percent_hue,
3629     percent_saturation;
3630 
3631   GeometryInfo
3632     geometry_info;
3633 
3634   MagickBooleanType
3635     status;
3636 
3637   MagickOffsetType
3638     progress;
3639 
3640   MagickStatusType
3641     flags;
3642 
3643   register ssize_t
3644     i;
3645 
3646   ssize_t
3647     y;
3648 
3649   /*
3650     Initialize modulate table.
3651   */
3652   assert(image != (Image *) NULL);
3653   assert(image->signature == MagickCoreSignature);
3654   if (image->debug != MagickFalse)
3655     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3656   if (modulate == (char *) NULL)
3657     return(MagickFalse);
3658   if (IssRGBCompatibleColorspace(image->colorspace) == MagickFalse)
3659     (void) SetImageColorspace(image,sRGBColorspace,exception);
3660   flags=ParseGeometry(modulate,&geometry_info);
3661   percent_brightness=geometry_info.rho;
3662   percent_saturation=geometry_info.sigma;
3663   if ((flags & SigmaValue) == 0)
3664     percent_saturation=100.0;
3665   percent_hue=geometry_info.xi;
3666   if ((flags & XiValue) == 0)
3667     percent_hue=100.0;
3668   colorspace=UndefinedColorspace;
3669   artifact=GetImageArtifact(image,"modulate:colorspace");
3670   if (artifact != (const char *) NULL)
3671     colorspace=(ColorspaceType) ParseCommandOption(MagickColorspaceOptions,
3672       MagickFalse,artifact);
3673   if (image->storage_class == PseudoClass)
3674     for (i=0; i < (ssize_t) image->colors; i++)
3675     {
3676       double
3677         blue,
3678         green,
3679         red;
3680 
3681       /*
3682         Modulate image colormap.
3683       */
3684       red=(double) image->colormap[i].red;
3685       green=(double) image->colormap[i].green;
3686       blue=(double) image->colormap[i].blue;
3687       switch (colorspace)
3688       {
3689         case HCLColorspace:
3690         {
3691           ModulateHCL(percent_hue,percent_saturation,percent_brightness,
3692             &red,&green,&blue);
3693           break;
3694         }
3695         case HCLpColorspace:
3696         {
3697           ModulateHCLp(percent_hue,percent_saturation,percent_brightness,
3698             &red,&green,&blue);
3699           break;
3700         }
3701         case HSBColorspace:
3702         {
3703           ModulateHSB(percent_hue,percent_saturation,percent_brightness,
3704             &red,&green,&blue);
3705           break;
3706         }
3707         case HSIColorspace:
3708         {
3709           ModulateHSI(percent_hue,percent_saturation,percent_brightness,
3710             &red,&green,&blue);
3711           break;
3712         }
3713         case HSLColorspace:
3714         default:
3715         {
3716           ModulateHSL(percent_hue,percent_saturation,percent_brightness,
3717             &red,&green,&blue);
3718           break;
3719         }
3720         case HSVColorspace:
3721         {
3722           ModulateHSV(percent_hue,percent_saturation,percent_brightness,
3723             &red,&green,&blue);
3724           break;
3725         }
3726         case HWBColorspace:
3727         {
3728           ModulateHWB(percent_hue,percent_saturation,percent_brightness,
3729             &red,&green,&blue);
3730           break;
3731         }
3732         case LCHColorspace:
3733         case LCHabColorspace:
3734         {
3735           ModulateLCHab(percent_brightness,percent_saturation,percent_hue,
3736             &red,&green,&blue);
3737           break;
3738         }
3739         case LCHuvColorspace:
3740         {
3741           ModulateLCHuv(percent_brightness,percent_saturation,percent_hue,
3742             &red,&green,&blue);
3743           break;
3744         }
3745       }
3746       image->colormap[i].red=red;
3747       image->colormap[i].green=green;
3748       image->colormap[i].blue=blue;
3749     }
3750   /*
3751     Modulate image.
3752   */
3753 #if defined(MAGICKCORE_OPENCL_SUPPORT)
3754   if (AccelerateModulateImage(image,percent_brightness,percent_hue,
3755         percent_saturation,colorspace,exception) != MagickFalse)
3756     return(MagickTrue);
3757 #endif
3758   status=MagickTrue;
3759   progress=0;
3760   image_view=AcquireAuthenticCacheView(image,exception);
3761 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3762   #pragma omp parallel for schedule(static) shared(progress,status) \
3763     magick_number_threads(image,image,image->rows,1)
3764 #endif
3765   for (y=0; y < (ssize_t) image->rows; y++)
3766   {
3767     register Quantum
3768       *magick_restrict q;
3769 
3770     register ssize_t
3771       x;
3772 
3773     if (status == MagickFalse)
3774       continue;
3775     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
3776     if (q == (Quantum *) NULL)
3777       {
3778         status=MagickFalse;
3779         continue;
3780       }
3781     for (x=0; x < (ssize_t) image->columns; x++)
3782     {
3783       double
3784         blue,
3785         green,
3786         red;
3787 
3788       red=(double) GetPixelRed(image,q);
3789       green=(double) GetPixelGreen(image,q);
3790       blue=(double) GetPixelBlue(image,q);
3791       switch (colorspace)
3792       {
3793         case HCLColorspace:
3794         {
3795           ModulateHCL(percent_hue,percent_saturation,percent_brightness,
3796             &red,&green,&blue);
3797           break;
3798         }
3799         case HCLpColorspace:
3800         {
3801           ModulateHCLp(percent_hue,percent_saturation,percent_brightness,
3802             &red,&green,&blue);
3803           break;
3804         }
3805         case HSBColorspace:
3806         {
3807           ModulateHSB(percent_hue,percent_saturation,percent_brightness,
3808             &red,&green,&blue);
3809           break;
3810         }
3811         case HSLColorspace:
3812         default:
3813         {
3814           ModulateHSL(percent_hue,percent_saturation,percent_brightness,
3815             &red,&green,&blue);
3816           break;
3817         }
3818         case HSVColorspace:
3819         {
3820           ModulateHSV(percent_hue,percent_saturation,percent_brightness,
3821             &red,&green,&blue);
3822           break;
3823         }
3824         case HWBColorspace:
3825         {
3826           ModulateHWB(percent_hue,percent_saturation,percent_brightness,
3827             &red,&green,&blue);
3828           break;
3829         }
3830         case LCHabColorspace:
3831         {
3832           ModulateLCHab(percent_brightness,percent_saturation,percent_hue,
3833             &red,&green,&blue);
3834           break;
3835         }
3836         case LCHColorspace:
3837         case LCHuvColorspace:
3838         {
3839           ModulateLCHuv(percent_brightness,percent_saturation,percent_hue,
3840             &red,&green,&blue);
3841           break;
3842         }
3843       }
3844       SetPixelRed(image,ClampToQuantum(red),q);
3845       SetPixelGreen(image,ClampToQuantum(green),q);
3846       SetPixelBlue(image,ClampToQuantum(blue),q);
3847       q+=GetPixelChannels(image);
3848     }
3849     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
3850       status=MagickFalse;
3851     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3852       {
3853         MagickBooleanType
3854           proceed;
3855 
3856 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3857         #pragma omp atomic
3858 #endif
3859         progress++;
3860         proceed=SetImageProgress(image,ModulateImageTag,progress,image->rows);
3861         if (proceed == MagickFalse)
3862           status=MagickFalse;
3863       }
3864   }
3865   image_view=DestroyCacheView(image_view);
3866   return(status);
3867 }
3868 
3869 /*
3870 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3871 %                                                                             %
3872 %                                                                             %
3873 %                                                                             %
3874 %     N e g a t e I m a g e                                                   %
3875 %                                                                             %
3876 %                                                                             %
3877 %                                                                             %
3878 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3879 %
3880 %  NegateImage() negates the colors in the reference image.  The grayscale
3881 %  option means that only grayscale values within the image are negated.
3882 %
3883 %  The format of the NegateImage method is:
3884 %
3885 %      MagickBooleanType NegateImage(Image *image,
3886 %        const MagickBooleanType grayscale,ExceptionInfo *exception)
3887 %
3888 %  A description of each parameter follows:
3889 %
3890 %    o image: the image.
3891 %
3892 %    o grayscale: If MagickTrue, only negate grayscale pixels within the image.
3893 %
3894 %    o exception: return any errors or warnings in this structure.
3895 %
3896 */
NegateImage(Image * image,const MagickBooleanType grayscale,ExceptionInfo * exception)3897 MagickExport MagickBooleanType NegateImage(Image *image,
3898   const MagickBooleanType grayscale,ExceptionInfo *exception)
3899 {
3900 #define NegateImageTag  "Negate/Image"
3901 
3902   CacheView
3903     *image_view;
3904 
3905   MagickBooleanType
3906     status;
3907 
3908   MagickOffsetType
3909     progress;
3910 
3911   register ssize_t
3912     i;
3913 
3914   ssize_t
3915     y;
3916 
3917   assert(image != (Image *) NULL);
3918   assert(image->signature == MagickCoreSignature);
3919   if (image->debug != MagickFalse)
3920     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3921   if (image->storage_class == PseudoClass)
3922     for (i=0; i < (ssize_t) image->colors; i++)
3923     {
3924       /*
3925         Negate colormap.
3926       */
3927       if( grayscale != MagickFalse )
3928         if ((image->colormap[i].red != image->colormap[i].green) ||
3929             (image->colormap[i].green != image->colormap[i].blue))
3930           continue;
3931       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
3932         image->colormap[i].red=QuantumRange-image->colormap[i].red;
3933       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
3934         image->colormap[i].green=QuantumRange-image->colormap[i].green;
3935       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
3936         image->colormap[i].blue=QuantumRange-image->colormap[i].blue;
3937     }
3938   /*
3939     Negate image.
3940   */
3941   status=MagickTrue;
3942   progress=0;
3943   image_view=AcquireAuthenticCacheView(image,exception);
3944   if( grayscale != MagickFalse )
3945     {
3946       for (y=0; y < (ssize_t) image->rows; y++)
3947       {
3948         MagickBooleanType
3949           sync;
3950 
3951         register Quantum
3952           *magick_restrict q;
3953 
3954         register ssize_t
3955           x;
3956 
3957         if (status == MagickFalse)
3958           continue;
3959         q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,
3960           exception);
3961         if (q == (Quantum *) NULL)
3962           {
3963             status=MagickFalse;
3964             continue;
3965           }
3966         for (x=0; x < (ssize_t) image->columns; x++)
3967         {
3968           register ssize_t
3969             j;
3970 
3971           if (IsPixelGray(image,q) != MagickFalse)
3972             {
3973               q+=GetPixelChannels(image);
3974               continue;
3975             }
3976           for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
3977           {
3978             PixelChannel channel = GetPixelChannelChannel(image,j);
3979             PixelTrait traits = GetPixelChannelTraits(image,channel);
3980             if ((traits & UpdatePixelTrait) == 0)
3981               continue;
3982             q[j]=QuantumRange-q[j];
3983           }
3984           q+=GetPixelChannels(image);
3985         }
3986         sync=SyncCacheViewAuthenticPixels(image_view,exception);
3987         if (sync == MagickFalse)
3988           status=MagickFalse;
3989         if (image->progress_monitor != (MagickProgressMonitor) NULL)
3990           {
3991             MagickBooleanType
3992               proceed;
3993 
3994 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3995             #pragma omp atomic
3996 #endif
3997             progress++;
3998             proceed=SetImageProgress(image,NegateImageTag,progress,image->rows);
3999             if (proceed == MagickFalse)
4000               status=MagickFalse;
4001           }
4002       }
4003       image_view=DestroyCacheView(image_view);
4004       return(MagickTrue);
4005     }
4006   /*
4007     Negate image.
4008   */
4009 #if defined(MAGICKCORE_OPENMP_SUPPORT)
4010   #pragma omp parallel for schedule(static) shared(progress,status) \
4011     magick_number_threads(image,image,image->rows,1)
4012 #endif
4013   for (y=0; y < (ssize_t) image->rows; y++)
4014   {
4015     register Quantum
4016       *magick_restrict q;
4017 
4018     register ssize_t
4019       x;
4020 
4021     if (status == MagickFalse)
4022       continue;
4023     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
4024     if (q == (Quantum *) NULL)
4025       {
4026         status=MagickFalse;
4027         continue;
4028       }
4029     for (x=0; x < (ssize_t) image->columns; x++)
4030     {
4031       register ssize_t
4032         j;
4033 
4034       for (j=0; j < (ssize_t) GetPixelChannels(image); j++)
4035       {
4036         PixelChannel channel = GetPixelChannelChannel(image,j);
4037         PixelTrait traits = GetPixelChannelTraits(image,channel);
4038         if ((traits & UpdatePixelTrait) == 0)
4039           continue;
4040         q[j]=QuantumRange-q[j];
4041       }
4042       q+=GetPixelChannels(image);
4043     }
4044     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
4045       status=MagickFalse;
4046     if (image->progress_monitor != (MagickProgressMonitor) NULL)
4047       {
4048         MagickBooleanType
4049           proceed;
4050 
4051 #if defined(MAGICKCORE_OPENMP_SUPPORT)
4052         #pragma omp atomic
4053 #endif
4054         progress++;
4055         proceed=SetImageProgress(image,NegateImageTag,progress,image->rows);
4056         if (proceed == MagickFalse)
4057           status=MagickFalse;
4058       }
4059   }
4060   image_view=DestroyCacheView(image_view);
4061   return(status);
4062 }
4063 
4064 /*
4065 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4066 %                                                                             %
4067 %                                                                             %
4068 %                                                                             %
4069 %     N o r m a l i z e I m a g e                                             %
4070 %                                                                             %
4071 %                                                                             %
4072 %                                                                             %
4073 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4074 %
4075 %  The NormalizeImage() method enhances the contrast of a color image by
4076 %  mapping the darkest 2 percent of all pixel to black and the brightest
4077 %  1 percent to white.
4078 %
4079 %  The format of the NormalizeImage method is:
4080 %
4081 %      MagickBooleanType NormalizeImage(Image *image,ExceptionInfo *exception)
4082 %
4083 %  A description of each parameter follows:
4084 %
4085 %    o image: the image.
4086 %
4087 %    o exception: return any errors or warnings in this structure.
4088 %
4089 */
NormalizeImage(Image * image,ExceptionInfo * exception)4090 MagickExport MagickBooleanType NormalizeImage(Image *image,
4091   ExceptionInfo *exception)
4092 {
4093   double
4094     black_point,
4095     white_point;
4096 
4097   black_point=(double) image->columns*image->rows*0.0015;
4098   white_point=(double) image->columns*image->rows*0.9995;
4099   return(ContrastStretchImage(image,black_point,white_point,exception));
4100 }
4101 
4102 /*
4103 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4104 %                                                                             %
4105 %                                                                             %
4106 %                                                                             %
4107 %     S i g m o i d a l C o n t r a s t I m a g e                             %
4108 %                                                                             %
4109 %                                                                             %
4110 %                                                                             %
4111 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
4112 %
4113 %  SigmoidalContrastImage() adjusts the contrast of an image with a non-linear
4114 %  sigmoidal contrast algorithm.  Increase the contrast of the image using a
4115 %  sigmoidal transfer function without saturating highlights or shadows.
4116 %  Contrast indicates how much to increase the contrast (0 is none; 3 is
4117 %  typical; 20 is pushing it); mid-point indicates where midtones fall in the
4118 %  resultant image (0 is white; 50% is middle-gray; 100% is black).  Set
4119 %  sharpen to MagickTrue to increase the image contrast otherwise the contrast
4120 %  is reduced.
4121 %
4122 %  The format of the SigmoidalContrastImage method is:
4123 %
4124 %      MagickBooleanType SigmoidalContrastImage(Image *image,
4125 %        const MagickBooleanType sharpen,const char *levels,
4126 %        ExceptionInfo *exception)
4127 %
4128 %  A description of each parameter follows:
4129 %
4130 %    o image: the image.
4131 %
4132 %    o sharpen: Increase or decrease image contrast.
4133 %
4134 %    o contrast: strength of the contrast, the larger the number the more
4135 %      'threshold-like' it becomes.
4136 %
4137 %    o midpoint: midpoint of the function as a color value 0 to QuantumRange.
4138 %
4139 %    o exception: return any errors or warnings in this structure.
4140 %
4141 */
4142 
4143 /*
4144   ImageMagick 6 has a version of this function which uses LUTs.
4145 */
4146 
4147 /*
4148   Sigmoidal function Sigmoidal with inflexion point moved to b and "slope
4149   constant" set to a.
4150 
4151   The first version, based on the hyperbolic tangent tanh, when combined with
4152   the scaling step, is an exact arithmetic clone of the the sigmoid function
4153   based on the logistic curve. The equivalence is based on the identity
4154 
4155     1/(1+exp(-t)) = (1+tanh(t/2))/2
4156 
4157   (http://de.wikipedia.org/wiki/Sigmoidfunktion) and the fact that the
4158   scaled sigmoidal derivation is invariant under affine transformations of
4159   the ordinate.
4160 
4161   The tanh version is almost certainly more accurate and cheaper.  The 0.5
4162   factor in the argument is to clone the legacy ImageMagick behavior. The
4163   reason for making the define depend on atanh even though it only uses tanh
4164   has to do with the construction of the inverse of the scaled sigmoidal.
4165 */
4166 #if defined(MAGICKCORE_HAVE_ATANH)
4167 #define Sigmoidal(a,b,x) ( tanh((0.5*(a))*((x)-(b))) )
4168 #else
4169 #define Sigmoidal(a,b,x) ( 1.0/(1.0+exp((a)*((b)-(x)))) )
4170 #endif
4171 /*
4172   Scaled sigmoidal function:
4173 
4174     ( Sigmoidal(a,b,x) - Sigmoidal(a,b,0) ) /
4175     ( Sigmoidal(a,b,1) - Sigmoidal(a,b,0) )
4176 
4177   See http://osdir.com/ml/video.image-magick.devel/2005-04/msg00006.html and
4178   http://www.cs.dartmouth.edu/farid/downloads/tutorials/fip.pdf.  The limit
4179   of ScaledSigmoidal as a->0 is the identity, but a=0 gives a division by
4180   zero. This is fixed below by exiting immediately when contrast is small,
4181   leaving the image (or colormap) unmodified. This appears to be safe because
4182   the series expansion of the logistic sigmoidal function around x=b is
4183 
4184   1/2-a*(b-x)/4+...
4185 
4186   so that the key denominator s(1)-s(0) is about a/4 (a/2 with tanh).
4187 */
4188 #define ScaledSigmoidal(a,b,x) (                    \
4189   (Sigmoidal((a),(b),(x))-Sigmoidal((a),(b),0.0)) / \
4190   (Sigmoidal((a),(b),1.0)-Sigmoidal((a),(b),0.0)) )
4191 /*
4192   Inverse of ScaledSigmoidal, used for +sigmoidal-contrast.  Because b
4193   may be 0 or 1, the argument of the hyperbolic tangent (resp. logistic
4194   sigmoidal) may be outside of the interval (-1,1) (resp. (0,1)), even
4195   when creating a LUT from in gamut values, hence the branching.  In
4196   addition, HDRI may have out of gamut values.
4197   InverseScaledSigmoidal is not a two-sided inverse of ScaledSigmoidal:
4198   It is only a right inverse. This is unavoidable.
4199 */
InverseScaledSigmoidal(const double a,const double b,const double x)4200 static inline double InverseScaledSigmoidal(const double a,const double b,
4201   const double x)
4202 {
4203   const double sig0=Sigmoidal(a,b,0.0);
4204   const double sig1=Sigmoidal(a,b,1.0);
4205   const double argument=(sig1-sig0)*x+sig0;
4206   const double clamped=
4207     (
4208 #if defined(MAGICKCORE_HAVE_ATANH)
4209       argument < -1+MagickEpsilon
4210       ?
4211       -1+MagickEpsilon
4212       :
4213       ( argument > 1-MagickEpsilon ? 1-MagickEpsilon : argument )
4214     );
4215   return(b+(2.0/a)*atanh(clamped));
4216 #else
4217       argument < MagickEpsilon
4218       ?
4219       MagickEpsilon
4220       :
4221       ( argument > 1-MagickEpsilon ? 1-MagickEpsilon : argument )
4222     );
4223   return(b-log(1.0/clamped-1.0)/a);
4224 #endif
4225 }
4226 
SigmoidalContrastImage(Image * image,const MagickBooleanType sharpen,const double contrast,const double midpoint,ExceptionInfo * exception)4227 MagickExport MagickBooleanType SigmoidalContrastImage(Image *image,
4228   const MagickBooleanType sharpen,const double contrast,const double midpoint,
4229   ExceptionInfo *exception)
4230 {
4231 #define SigmoidalContrastImageTag  "SigmoidalContrast/Image"
4232 #define ScaledSig(x) ( ClampToQuantum(QuantumRange* \
4233   ScaledSigmoidal(contrast,QuantumScale*midpoint,QuantumScale*(x))) )
4234 #define InverseScaledSig(x) ( ClampToQuantum(QuantumRange* \
4235   InverseScaledSigmoidal(contrast,QuantumScale*midpoint,QuantumScale*(x))) )
4236 
4237   CacheView
4238     *image_view;
4239 
4240   MagickBooleanType
4241     status;
4242 
4243   MagickOffsetType
4244     progress;
4245 
4246   ssize_t
4247     y;
4248 
4249   /*
4250     Convenience macros.
4251   */
4252   assert(image != (Image *) NULL);
4253   assert(image->signature == MagickCoreSignature);
4254   if (image->debug != MagickFalse)
4255     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
4256   /*
4257     Side effect: may clamp values unless contrast<MagickEpsilon, in which
4258     case nothing is done.
4259   */
4260   if (contrast < MagickEpsilon)
4261     return(MagickTrue);
4262   /*
4263     Sigmoidal-contrast enhance colormap.
4264   */
4265   if (image->storage_class == PseudoClass)
4266     {
4267       register ssize_t
4268         i;
4269 
4270       if( sharpen != MagickFalse )
4271         for (i=0; i < (ssize_t) image->colors; i++)
4272         {
4273           if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
4274             image->colormap[i].red=(MagickRealType) ScaledSig(
4275               image->colormap[i].red);
4276           if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
4277             image->colormap[i].green=(MagickRealType) ScaledSig(
4278               image->colormap[i].green);
4279           if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
4280             image->colormap[i].blue=(MagickRealType) ScaledSig(
4281               image->colormap[i].blue);
4282           if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
4283             image->colormap[i].alpha=(MagickRealType) ScaledSig(
4284               image->colormap[i].alpha);
4285         }
4286       else
4287         for (i=0; i < (ssize_t) image->colors; i++)
4288         {
4289           if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
4290             image->colormap[i].red=(MagickRealType) InverseScaledSig(
4291               image->colormap[i].red);
4292           if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
4293             image->colormap[i].green=(MagickRealType) InverseScaledSig(
4294               image->colormap[i].green);
4295           if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
4296             image->colormap[i].blue=(MagickRealType) InverseScaledSig(
4297               image->colormap[i].blue);
4298           if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
4299             image->colormap[i].alpha=(MagickRealType) InverseScaledSig(
4300               image->colormap[i].alpha);
4301         }
4302     }
4303   /*
4304     Sigmoidal-contrast enhance image.
4305   */
4306   status=MagickTrue;
4307   progress=0;
4308   image_view=AcquireAuthenticCacheView(image,exception);
4309 #if defined(MAGICKCORE_OPENMP_SUPPORT)
4310   #pragma omp parallel for schedule(static) shared(progress,status) \
4311     magick_number_threads(image,image,image->rows,1)
4312 #endif
4313   for (y=0; y < (ssize_t) image->rows; y++)
4314   {
4315     register Quantum
4316       *magick_restrict q;
4317 
4318     register ssize_t
4319       x;
4320 
4321     if (status == MagickFalse)
4322       continue;
4323     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
4324     if (q == (Quantum *) NULL)
4325       {
4326         status=MagickFalse;
4327         continue;
4328       }
4329     for (x=0; x < (ssize_t) image->columns; x++)
4330     {
4331       register ssize_t
4332         i;
4333 
4334       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
4335       {
4336         PixelChannel channel = GetPixelChannelChannel(image,i);
4337         PixelTrait traits = GetPixelChannelTraits(image,channel);
4338         if ((traits & UpdatePixelTrait) == 0)
4339           continue;
4340         if( sharpen != MagickFalse )
4341           q[i]=ScaledSig(q[i]);
4342         else
4343           q[i]=InverseScaledSig(q[i]);
4344       }
4345       q+=GetPixelChannels(image);
4346     }
4347     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
4348       status=MagickFalse;
4349     if (image->progress_monitor != (MagickProgressMonitor) NULL)
4350       {
4351         MagickBooleanType
4352           proceed;
4353 
4354 #if defined(MAGICKCORE_OPENMP_SUPPORT)
4355         #pragma omp atomic
4356 #endif
4357         progress++;
4358         proceed=SetImageProgress(image,SigmoidalContrastImageTag,progress,
4359           image->rows);
4360         if (proceed == MagickFalse)
4361           status=MagickFalse;
4362       }
4363   }
4364   image_view=DestroyCacheView(image_view);
4365   return(status);
4366 }
4367