]> granicus.if.org Git - imagemagick/blob - MagickCore/enhance.c
comment typos
[imagemagick] / MagickCore / enhance.c
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 %                                John Cristy                                  %
17 %                                 July 1992                                   %
18 %                                                                             %
19 %                                                                             %
20 %  Copyright 1999-2012 ImageMagick Studio LLC, a non-profit organization      %
21 %  dedicated to making software imaging solutions freely available.           %
22 %                                                                             %
23 %  You may not use this file except in compliance with the License.  You may  %
24 %  obtain a copy of the License at                                            %
25 %                                                                             %
26 %    http://www.imagemagick.org/script/license.php                            %
27 %                                                                             %
28 %  Unless required by applicable law or agreed to in writing, software        %
29 %  distributed under the License is distributed on an "AS IS" BASIS,          %
30 %  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.   %
31 %  See the License for the specific language governing permissions and        %
32 %  limitations under the License.                                             %
33 %                                                                             %
34 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
35 %
36 %
37 %
38 */
39 \f
40 /*
41   Include declarations.
42 */
43 #include "MagickCore/studio.h"
44 #include "MagickCore/artifact.h"
45 #include "MagickCore/cache.h"
46 #include "MagickCore/cache-view.h"
47 #include "MagickCore/color.h"
48 #include "MagickCore/color-private.h"
49 #include "MagickCore/colorspace.h"
50 #include "MagickCore/colorspace-private.h"
51 #include "MagickCore/composite-private.h"
52 #include "MagickCore/enhance.h"
53 #include "MagickCore/exception.h"
54 #include "MagickCore/exception-private.h"
55 #include "MagickCore/fx.h"
56 #include "MagickCore/gem.h"
57 #include "MagickCore/gem-private.h"
58 #include "MagickCore/geometry.h"
59 #include "MagickCore/histogram.h"
60 #include "MagickCore/image.h"
61 #include "MagickCore/image-private.h"
62 #include "MagickCore/memory_.h"
63 #include "MagickCore/monitor.h"
64 #include "MagickCore/monitor-private.h"
65 #include "MagickCore/option.h"
66 #include "MagickCore/pixel.h"
67 #include "MagickCore/pixel-accessor.h"
68 #include "MagickCore/quantum.h"
69 #include "MagickCore/quantum-private.h"
70 #include "MagickCore/resample.h"
71 #include "MagickCore/resample-private.h"
72 #include "MagickCore/resource_.h"
73 #include "MagickCore/statistic.h"
74 #include "MagickCore/string_.h"
75 #include "MagickCore/string-private.h"
76 #include "MagickCore/thread-private.h"
77 #include "MagickCore/threshold.h"
78 #include "MagickCore/token.h"
79 #include "MagickCore/xml-tree.h"
80 #include "MagickCore/xml-tree-private.h"
81 \f
82 /*
83 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
84 %                                                                             %
85 %                                                                             %
86 %                                                                             %
87 %     A u t o G a m m a I m a g e                                             %
88 %                                                                             %
89 %                                                                             %
90 %                                                                             %
91 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
92 %
93 %  AutoGammaImage() extract the 'mean' from the image and adjust the image
94 %  to try make set its gamma appropriatally.
95 %
96 %  The format of the AutoGammaImage method is:
97 %
98 %      MagickBooleanType AutoGammaImage(Image *image,ExceptionInfo *exception)
99 %
100 %  A description of each parameter follows:
101 %
102 %    o image: The image to auto-level
103 %
104 %    o exception: return any errors or warnings in this structure.
105 %
106 */
107 MagickExport MagickBooleanType AutoGammaImage(Image *image,
108   ExceptionInfo *exception)
109 {
110   double
111     gamma,
112     log_mean,
113     mean,
114     sans;
115
116   MagickStatusType
117     status;
118
119   register ssize_t
120     i;
121
122   log_mean=log(0.5);
123   if (image->channel_mask == DefaultChannels)
124     {
125       /*
126         Apply gamma correction equally across all given channels.
127       */
128       (void) GetImageMean(image,&mean,&sans,exception);
129       gamma=log(mean*QuantumScale)/log_mean;
130       return(LevelImage(image,0.0,(double) QuantumRange,gamma,exception));
131     }
132   /*
133     Auto-gamma each channel separately.
134   */
135   status=MagickTrue;
136   for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
137   {
138     ChannelType
139       channel_mask;
140
141     PixelChannel
142       channel;
143
144     PixelTrait
145       traits;
146
147     channel=GetPixelChannelChannel(image,i);
148     traits=GetPixelChannelTraits(image,channel);
149     if ((traits & UpdatePixelTrait) == 0)
150       continue;
151     channel_mask=SetImageChannelMask(image,(ChannelType) (1 << i));
152     status=GetImageMean(image,&mean,&sans,exception);
153     gamma=log(mean*QuantumScale)/log_mean;
154     status&=LevelImage(image,0.0,(double) QuantumRange,gamma,exception);
155     (void) SetImageChannelMask(image,channel_mask);
156     if (status == MagickFalse)
157       break;
158   }
159   return(status != 0 ? MagickTrue : MagickFalse);
160 }
161 \f
162 /*
163 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
164 %                                                                             %
165 %                                                                             %
166 %                                                                             %
167 %     A u t o L e v e l I m a g e                                             %
168 %                                                                             %
169 %                                                                             %
170 %                                                                             %
171 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
172 %
173 %  AutoLevelImage() adjusts the levels of a particular image channel by
174 %  scaling the minimum and maximum values to the full quantum range.
175 %
176 %  The format of the LevelImage method is:
177 %
178 %      MagickBooleanType AutoLevelImage(Image *image,ExceptionInfo *exception)
179 %
180 %  A description of each parameter follows:
181 %
182 %    o image: The image to auto-level
183 %
184 %    o exception: return any errors or warnings in this structure.
185 %
186 */
187 MagickExport MagickBooleanType AutoLevelImage(Image *image,
188   ExceptionInfo *exception)
189 {
190   return(MinMaxStretchImage(image,0.0,0.0,1.0,exception));
191 }
192 \f
193 /*
194 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
195 %                                                                             %
196 %                                                                             %
197 %                                                                             %
198 %     B r i g h t n e s s C o n t r a s t I m a g e                           %
199 %                                                                             %
200 %                                                                             %
201 %                                                                             %
202 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
203 %
204 %  BrightnessContrastImage() changes the brightness and/or contrast of an
205 %  image.  It converts the brightness and contrast parameters into slope and
206 %  intercept and calls a polynomical function to apply to the image.
207 %
208 %  The format of the BrightnessContrastImage method is:
209 %
210 %      MagickBooleanType BrightnessContrastImage(Image *image,
211 %        const double brightness,const double contrast,ExceptionInfo *exception)
212 %
213 %  A description of each parameter follows:
214 %
215 %    o image: the image.
216 %
217 %    o brightness: the brightness percent (-100 .. 100).
218 %
219 %    o contrast: the contrast percent (-100 .. 100).
220 %
221 %    o exception: return any errors or warnings in this structure.
222 %
223 */
224 MagickExport MagickBooleanType BrightnessContrastImage(Image *image,
225   const double brightness,const double contrast,ExceptionInfo *exception)
226 {
227 #define BrightnessContastImageTag  "BrightnessContast/Image"
228
229   double
230     alpha,
231     coefficients[2],
232     intercept,
233     slope;
234
235   MagickBooleanType
236     status;
237
238   /*
239     Compute slope and intercept.
240   */
241   assert(image != (Image *) NULL);
242   assert(image->signature == MagickSignature);
243   if (image->debug != MagickFalse)
244     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
245   alpha=contrast;
246   slope=tan((double) (MagickPI*(alpha/100.0+1.0)/4.0));
247   if (slope < 0.0)
248     slope=0.0;
249   intercept=brightness/100.0+((100-brightness)/200.0)*(1.0-slope);
250   coefficients[0]=slope;
251   coefficients[1]=intercept;
252   status=FunctionImage(image,PolynomialFunction,2,coefficients,exception);
253   return(status);
254 }
255 \f
256 /*
257 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
258 %                                                                             %
259 %                                                                             %
260 %                                                                             %
261 %     C l u t I m a g e                                                       %
262 %                                                                             %
263 %                                                                             %
264 %                                                                             %
265 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
266 %
267 %  ClutImage() replaces each color value in the given image, by using it as an
268 %  index to lookup a replacement color value in a Color Look UP Table in the
269 %  form of an image.  The values are extracted along a diagonal of the CLUT
270 %  image so either a horizontal or vertial gradient image can be used.
271 %
272 %  Typically this is used to either re-color a gray-scale image according to a
273 %  color gradient in the CLUT image, or to perform a freeform histogram
274 %  (level) adjustment according to the (typically gray-scale) gradient in the
275 %  CLUT image.
276 %
277 %  When the 'channel' mask includes the matte/alpha transparency channel but
278 %  one image has no such channel it is assumed that that image is a simple
279 %  gray-scale image that will effect the alpha channel values, either for
280 %  gray-scale coloring (with transparent or semi-transparent colors), or
281 %  a histogram adjustment of existing alpha channel values.   If both images
282 %  have matte channels, direct and normal indexing is applied, which is rarely
283 %  used.
284 %
285 %  The format of the ClutImage method is:
286 %
287 %      MagickBooleanType ClutImage(Image *image,Image *clut_image,
288 %        const PixelInterpolateMethod method,ExceptionInfo *exception)
289 %
290 %  A description of each parameter follows:
291 %
292 %    o image: the image, which is replaced by indexed CLUT values
293 %
294 %    o clut_image: the color lookup table image for replacement color values.
295 %
296 %    o method: the pixel interpolation method.
297 %
298 %    o exception: return any errors or warnings in this structure.
299 %
300 */
301 MagickExport MagickBooleanType ClutImage(Image *image,const Image *clut_image,
302   const PixelInterpolateMethod method,ExceptionInfo *exception)
303 {
304 #define ClutImageTag  "Clut/Image"
305
306   CacheView
307     *clut_view,
308     *image_view;
309
310   MagickBooleanType
311     status;
312
313   MagickOffsetType
314     progress;
315
316   PixelInfo
317     *clut_map;
318
319   register ssize_t
320     i;
321
322   ssize_t adjust,
323     y;
324
325   assert(image != (Image *) NULL);
326   assert(image->signature == MagickSignature);
327   if (image->debug != MagickFalse)
328     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
329   assert(clut_image != (Image *) NULL);
330   assert(clut_image->signature == MagickSignature);
331   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
332     return(MagickFalse);
333   if (IsGrayColorspace(image->colorspace) != MagickFalse)
334     (void) TransformImageColorspace(image,RGBColorspace,exception);
335   clut_map=(PixelInfo *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*clut_map));
336   if (clut_map == (PixelInfo *) NULL)
337     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
338       image->filename);
339   /*
340     Clut image.
341   */
342   status=MagickTrue;
343   progress=0;
344   adjust=(ssize_t) (clut_image->interpolate == IntegerInterpolatePixel ? 0 : 1);
345   clut_view=AcquireVirtualCacheView(clut_image,exception);
346 #if defined(MAGICKCORE_OPENMP_SUPPORT)
347   #pragma omp parallel for schedule(static,4) \
348     dynamic_number_threads(image,image->columns,1,1)
349 #endif
350   for (i=0; i <= (ssize_t) MaxMap; i++)
351   {
352     GetPixelInfo(clut_image,clut_map+i);
353     (void) InterpolatePixelInfo(clut_image,clut_view,method,
354       QuantumScale*i*(clut_image->columns-adjust),QuantumScale*i*
355       (clut_image->rows-adjust),clut_map+i,exception);
356   }
357   clut_view=DestroyCacheView(clut_view);
358   image_view=AcquireAuthenticCacheView(image,exception);
359 #if defined(MAGICKCORE_OPENMP_SUPPORT)
360   #pragma omp parallel for schedule(static,4) shared(progress,status) \
361     dynamic_number_threads(image,image->columns,image->rows,1)
362 #endif
363   for (y=0; y < (ssize_t) image->rows; y++)
364   {
365     PixelInfo
366       pixel;
367
368     register Quantum
369       *restrict q;
370
371     register ssize_t
372       x;
373
374     if (status == MagickFalse)
375       continue;
376     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
377     if (q == (Quantum *) NULL)
378       {
379         status=MagickFalse;
380         continue;
381       }
382     GetPixelInfo(image,&pixel);
383     for (x=0; x < (ssize_t) image->columns; x++)
384     {
385       if (GetPixelMask(image,q) != 0)
386         {
387           q+=GetPixelChannels(image);
388           continue;
389         }
390       GetPixelInfoPixel(image,q,&pixel);
391       pixel.red=clut_map[ScaleQuantumToMap(
392         ClampToQuantum(pixel.red))].red;
393       pixel.green=clut_map[ScaleQuantumToMap(
394         ClampToQuantum(pixel.green))].green;
395       pixel.blue=clut_map[ScaleQuantumToMap(
396         ClampToQuantum(pixel.blue))].blue;
397       pixel.black=clut_map[ScaleQuantumToMap(
398         ClampToQuantum(pixel.black))].black;
399       pixel.alpha=clut_map[ScaleQuantumToMap(
400         ClampToQuantum(pixel.alpha))].alpha;
401       SetPixelInfoPixel(image,&pixel,q);
402       q+=GetPixelChannels(image);
403     }
404     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
405       status=MagickFalse;
406     if (image->progress_monitor != (MagickProgressMonitor) NULL)
407       {
408         MagickBooleanType
409           proceed;
410
411 #if defined(MAGICKCORE_OPENMP_SUPPORT)
412         #pragma omp critical (MagickCore_ClutImage)
413 #endif
414         proceed=SetImageProgress(image,ClutImageTag,progress++,image->rows);
415         if (proceed == MagickFalse)
416           status=MagickFalse;
417       }
418   }
419   image_view=DestroyCacheView(image_view);
420   clut_map=(PixelInfo *) RelinquishMagickMemory(clut_map);
421   if ((clut_image->alpha_trait == BlendPixelTrait) &&
422       ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0))
423     (void) SetImageAlphaChannel(image,ActivateAlphaChannel,exception);
424   return(status);
425 }
426 \f
427 /*
428 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
429 %                                                                             %
430 %                                                                             %
431 %                                                                             %
432 %     C o l o r D e c i s i o n L i s t I m a g e                             %
433 %                                                                             %
434 %                                                                             %
435 %                                                                             %
436 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
437 %
438 %  ColorDecisionListImage() accepts a lightweight Color Correction Collection
439 %  (CCC) file which solely contains one or more color corrections and applies
440 %  the correction to the image.  Here is a sample CCC file:
441 %
442 %    <ColorCorrectionCollection xmlns="urn:ASC:CDL:v1.2">
443 %          <ColorCorrection id="cc03345">
444 %                <SOPNode>
445 %                     <Slope> 0.9 1.2 0.5 </Slope>
446 %                     <Offset> 0.4 -0.5 0.6 </Offset>
447 %                     <Power> 1.0 0.8 1.5 </Power>
448 %                </SOPNode>
449 %                <SATNode>
450 %                     <Saturation> 0.85 </Saturation>
451 %                </SATNode>
452 %          </ColorCorrection>
453 %    </ColorCorrectionCollection>
454 %
455 %  which includes the slop, offset, and power for each of the RGB channels
456 %  as well as the saturation.
457 %
458 %  The format of the ColorDecisionListImage method is:
459 %
460 %      MagickBooleanType ColorDecisionListImage(Image *image,
461 %        const char *color_correction_collection,ExceptionInfo *exception)
462 %
463 %  A description of each parameter follows:
464 %
465 %    o image: the image.
466 %
467 %    o color_correction_collection: the color correction collection in XML.
468 %
469 %    o exception: return any errors or warnings in this structure.
470 %
471 */
472 MagickExport MagickBooleanType ColorDecisionListImage(Image *image,
473   const char *color_correction_collection,ExceptionInfo *exception)
474 {
475 #define ColorDecisionListCorrectImageTag  "ColorDecisionList/Image"
476
477   typedef struct _Correction
478   {
479     double
480       slope,
481       offset,
482       power;
483   } Correction;
484
485   typedef struct _ColorCorrection
486   {
487     Correction
488       red,
489       green,
490       blue;
491
492     double
493       saturation;
494   } ColorCorrection;
495
496   CacheView
497     *image_view;
498
499   char
500     token[MaxTextExtent];
501
502   ColorCorrection
503     color_correction;
504
505   const char
506     *content,
507     *p;
508
509   MagickBooleanType
510     status;
511
512   MagickOffsetType
513     progress;
514
515   PixelInfo
516     *cdl_map;
517
518   register ssize_t
519     i;
520
521   ssize_t
522     y;
523
524   XMLTreeInfo
525     *cc,
526     *ccc,
527     *sat,
528     *sop;
529
530   /*
531     Allocate and initialize cdl maps.
532   */
533   assert(image != (Image *) NULL);
534   assert(image->signature == MagickSignature);
535   if (image->debug != MagickFalse)
536     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
537   if (color_correction_collection == (const char *) NULL)
538     return(MagickFalse);
539   ccc=NewXMLTree((const char *) color_correction_collection,exception);
540   if (ccc == (XMLTreeInfo *) NULL)
541     return(MagickFalse);
542   cc=GetXMLTreeChild(ccc,"ColorCorrection");
543   if (cc == (XMLTreeInfo *) NULL)
544     {
545       ccc=DestroyXMLTree(ccc);
546       return(MagickFalse);
547     }
548   color_correction.red.slope=1.0;
549   color_correction.red.offset=0.0;
550   color_correction.red.power=1.0;
551   color_correction.green.slope=1.0;
552   color_correction.green.offset=0.0;
553   color_correction.green.power=1.0;
554   color_correction.blue.slope=1.0;
555   color_correction.blue.offset=0.0;
556   color_correction.blue.power=1.0;
557   color_correction.saturation=0.0;
558   sop=GetXMLTreeChild(cc,"SOPNode");
559   if (sop != (XMLTreeInfo *) NULL)
560     {
561       XMLTreeInfo
562         *offset,
563         *power,
564         *slope;
565
566       slope=GetXMLTreeChild(sop,"Slope");
567       if (slope != (XMLTreeInfo *) NULL)
568         {
569           content=GetXMLTreeContent(slope);
570           p=(const char *) content;
571           for (i=0; (*p != '\0') && (i < 3); i++)
572           {
573             GetMagickToken(p,&p,token);
574             if (*token == ',')
575               GetMagickToken(p,&p,token);
576             switch (i)
577             {
578               case 0:
579               {
580                 color_correction.red.slope=StringToDouble(token,(char **) NULL);
581                 break;
582               }
583               case 1:
584               {
585                 color_correction.green.slope=StringToDouble(token,
586                   (char **) NULL);
587                 break;
588               }
589               case 2:
590               {
591                 color_correction.blue.slope=StringToDouble(token,
592                   (char **) NULL);
593                 break;
594               }
595             }
596           }
597         }
598       offset=GetXMLTreeChild(sop,"Offset");
599       if (offset != (XMLTreeInfo *) NULL)
600         {
601           content=GetXMLTreeContent(offset);
602           p=(const char *) content;
603           for (i=0; (*p != '\0') && (i < 3); i++)
604           {
605             GetMagickToken(p,&p,token);
606             if (*token == ',')
607               GetMagickToken(p,&p,token);
608             switch (i)
609             {
610               case 0:
611               {
612                 color_correction.red.offset=StringToDouble(token,
613                   (char **) NULL);
614                 break;
615               }
616               case 1:
617               {
618                 color_correction.green.offset=StringToDouble(token,
619                   (char **) NULL);
620                 break;
621               }
622               case 2:
623               {
624                 color_correction.blue.offset=StringToDouble(token,
625                   (char **) NULL);
626                 break;
627               }
628             }
629           }
630         }
631       power=GetXMLTreeChild(sop,"Power");
632       if (power != (XMLTreeInfo *) NULL)
633         {
634           content=GetXMLTreeContent(power);
635           p=(const char *) content;
636           for (i=0; (*p != '\0') && (i < 3); i++)
637           {
638             GetMagickToken(p,&p,token);
639             if (*token == ',')
640               GetMagickToken(p,&p,token);
641             switch (i)
642             {
643               case 0:
644               {
645                 color_correction.red.power=StringToDouble(token,(char **) NULL);
646                 break;
647               }
648               case 1:
649               {
650                 color_correction.green.power=StringToDouble(token,
651                   (char **) NULL);
652                 break;
653               }
654               case 2:
655               {
656                 color_correction.blue.power=StringToDouble(token,
657                   (char **) NULL);
658                 break;
659               }
660             }
661           }
662         }
663     }
664   sat=GetXMLTreeChild(cc,"SATNode");
665   if (sat != (XMLTreeInfo *) NULL)
666     {
667       XMLTreeInfo
668         *saturation;
669
670       saturation=GetXMLTreeChild(sat,"Saturation");
671       if (saturation != (XMLTreeInfo *) NULL)
672         {
673           content=GetXMLTreeContent(saturation);
674           p=(const char *) content;
675           GetMagickToken(p,&p,token);
676           color_correction.saturation=StringToDouble(token,(char **) NULL);
677         }
678     }
679   ccc=DestroyXMLTree(ccc);
680   if (image->debug != MagickFalse)
681     {
682       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
683         "  Color Correction Collection:");
684       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
685         "  color_correction.red.slope: %g",color_correction.red.slope);
686       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
687         "  color_correction.red.offset: %g",color_correction.red.offset);
688       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
689         "  color_correction.red.power: %g",color_correction.red.power);
690       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
691         "  color_correction.green.slope: %g",color_correction.green.slope);
692       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
693         "  color_correction.green.offset: %g",color_correction.green.offset);
694       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
695         "  color_correction.green.power: %g",color_correction.green.power);
696       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
697         "  color_correction.blue.slope: %g",color_correction.blue.slope);
698       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
699         "  color_correction.blue.offset: %g",color_correction.blue.offset);
700       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
701         "  color_correction.blue.power: %g",color_correction.blue.power);
702       (void) LogMagickEvent(TransformEvent,GetMagickModule(),
703         "  color_correction.saturation: %g",color_correction.saturation);
704     }
705   cdl_map=(PixelInfo *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*cdl_map));
706   if (cdl_map == (PixelInfo *) NULL)
707     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
708       image->filename);
709 #if defined(MAGICKCORE_OPENMP_SUPPORT)
710   #pragma omp parallel for schedule(static,4) \
711     dynamic_number_threads(image,image->columns,1,1)
712 #endif
713   for (i=0; i <= (ssize_t) MaxMap; i++)
714   {
715     cdl_map[i].red=(double) ScaleMapToQuantum((double)
716       (MaxMap*(pow(color_correction.red.slope*i/MaxMap+
717       color_correction.red.offset,color_correction.red.power))));
718     cdl_map[i].green=(double) ScaleMapToQuantum((double)
719       (MaxMap*(pow(color_correction.green.slope*i/MaxMap+
720       color_correction.green.offset,color_correction.green.power))));
721     cdl_map[i].blue=(double) ScaleMapToQuantum((double)
722       (MaxMap*(pow(color_correction.blue.slope*i/MaxMap+
723       color_correction.blue.offset,color_correction.blue.power))));
724   }
725   if (image->storage_class == PseudoClass)
726     for (i=0; i < (ssize_t) image->colors; i++)
727     {
728       /*
729         Apply transfer function to colormap.
730       */
731       double
732         luma;
733
734       luma=0.21267*image->colormap[i].red+0.71526*image->colormap[i].green+
735         0.07217*image->colormap[i].blue;
736       image->colormap[i].red=luma+color_correction.saturation*cdl_map[
737         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))].red-luma;
738       image->colormap[i].green=luma+color_correction.saturation*cdl_map[
739         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].green))].green-luma;
740       image->colormap[i].blue=luma+color_correction.saturation*cdl_map[
741         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].blue))].blue-luma;
742     }
743   /*
744     Apply transfer function to image.
745   */
746   status=MagickTrue;
747   progress=0;
748   image_view=AcquireAuthenticCacheView(image,exception);
749 #if defined(MAGICKCORE_OPENMP_SUPPORT)
750   #pragma omp parallel for schedule(static,4) shared(progress,status) \
751     dynamic_number_threads(image,image->columns,image->rows,1)
752 #endif
753   for (y=0; y < (ssize_t) image->rows; y++)
754   {
755     double
756       luma;
757
758     register Quantum
759       *restrict q;
760
761     register ssize_t
762       x;
763
764     if (status == MagickFalse)
765       continue;
766     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
767     if (q == (Quantum *) NULL)
768       {
769         status=MagickFalse;
770         continue;
771       }
772     for (x=0; x < (ssize_t) image->columns; x++)
773     {
774       luma=0.21267*GetPixelRed(image,q)+0.71526*GetPixelGreen(image,q)+0.07217*
775         GetPixelBlue(image,q);
776       SetPixelRed(image,ClampToQuantum(luma+color_correction.saturation*
777         (cdl_map[ScaleQuantumToMap(GetPixelRed(image,q))].red-luma)),q);
778       SetPixelGreen(image,ClampToQuantum(luma+color_correction.saturation*
779         (cdl_map[ScaleQuantumToMap(GetPixelGreen(image,q))].green-luma)),q);
780       SetPixelBlue(image,ClampToQuantum(luma+color_correction.saturation*
781         (cdl_map[ScaleQuantumToMap(GetPixelBlue(image,q))].blue-luma)),q);
782       q+=GetPixelChannels(image);
783     }
784     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
785       status=MagickFalse;
786     if (image->progress_monitor != (MagickProgressMonitor) NULL)
787       {
788         MagickBooleanType
789           proceed;
790
791 #if defined(MAGICKCORE_OPENMP_SUPPORT)
792         #pragma omp critical (MagickCore_ColorDecisionListImageChannel)
793 #endif
794         proceed=SetImageProgress(image,ColorDecisionListCorrectImageTag,
795           progress++,image->rows);
796         if (proceed == MagickFalse)
797           status=MagickFalse;
798       }
799   }
800   image_view=DestroyCacheView(image_view);
801   cdl_map=(PixelInfo *) RelinquishMagickMemory(cdl_map);
802   return(status);
803 }
804 \f
805 /*
806 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
807 %                                                                             %
808 %                                                                             %
809 %                                                                             %
810 %     C o n t r a s t I m a g e                                               %
811 %                                                                             %
812 %                                                                             %
813 %                                                                             %
814 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
815 %
816 %  ContrastImage() enhances the intensity differences between the lighter and
817 %  darker elements of the image.  Set sharpen to a MagickTrue to increase the
818 %  image contrast otherwise the contrast is reduced.
819 %
820 %  The format of the ContrastImage method is:
821 %
822 %      MagickBooleanType ContrastImage(Image *image,
823 %        const MagickBooleanType sharpen,ExceptionInfo *exception)
824 %
825 %  A description of each parameter follows:
826 %
827 %    o image: the image.
828 %
829 %    o sharpen: Increase or decrease image contrast.
830 %
831 %    o exception: return any errors or warnings in this structure.
832 %
833 */
834
835 static void Contrast(const int sign,double *red,double *green,double *blue)
836 {
837   double
838     brightness,
839     hue,
840     saturation;
841
842   /*
843     Enhance contrast: dark color become darker, light color become lighter.
844   */
845   assert(red != (double *) NULL);
846   assert(green != (double *) NULL);
847   assert(blue != (double *) NULL);
848   hue=0.0;
849   saturation=0.0;
850   brightness=0.0;
851   ConvertRGBToHSB(*red,*green,*blue,&hue,&saturation,&brightness);
852   brightness+=0.5*sign*(0.5*(sin((double) (MagickPI*(brightness-0.5)))+1.0)-
853     brightness);
854   if (brightness > 1.0)
855     brightness=1.0;
856   else
857     if (brightness < 0.0)
858       brightness=0.0;
859   ConvertHSBToRGB(hue,saturation,brightness,red,green,blue);
860 }
861
862 MagickExport MagickBooleanType ContrastImage(Image *image,
863   const MagickBooleanType sharpen,ExceptionInfo *exception)
864 {
865 #define ContrastImageTag  "Contrast/Image"
866
867   CacheView
868     *image_view;
869
870   int
871     sign;
872
873   MagickBooleanType
874     status;
875
876   MagickOffsetType
877     progress;
878
879   register ssize_t
880     i;
881
882   ssize_t
883     y;
884
885   assert(image != (Image *) NULL);
886   assert(image->signature == MagickSignature);
887   if (image->debug != MagickFalse)
888     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
889   sign=sharpen != MagickFalse ? 1 : -1;
890   if (image->storage_class == PseudoClass)
891     {
892       /*
893         Contrast enhance colormap.
894       */
895       for (i=0; i < (ssize_t) image->colors; i++)
896         Contrast(sign,&image->colormap[i].red,&image->colormap[i].green,
897           &image->colormap[i].blue);
898     }
899   /*
900     Contrast enhance image.
901   */
902   status=MagickTrue;
903   progress=0;
904   image_view=AcquireAuthenticCacheView(image,exception);
905 #if defined(MAGICKCORE_OPENMP_SUPPORT)
906   #pragma omp parallel for schedule(static,4) shared(progress,status) \
907     dynamic_number_threads(image,image->columns,image->rows,1)
908 #endif
909   for (y=0; y < (ssize_t) image->rows; y++)
910   {
911     double
912       blue,
913       green,
914       red;
915
916     register Quantum
917       *restrict q;
918
919     register ssize_t
920       x;
921
922     if (status == MagickFalse)
923       continue;
924     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
925     if (q == (Quantum *) NULL)
926       {
927         status=MagickFalse;
928         continue;
929       }
930     for (x=0; x < (ssize_t) image->columns; x++)
931     {
932       red=(double) GetPixelRed(image,q);
933       green=(double) GetPixelGreen(image,q);
934       blue=(double) GetPixelBlue(image,q);
935       Contrast(sign,&red,&green,&blue);
936       SetPixelRed(image,ClampToQuantum(red),q);
937       SetPixelGreen(image,ClampToQuantum(green),q);
938       SetPixelBlue(image,ClampToQuantum(blue),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 critical (MagickCore_ContrastImage)
950 #endif
951         proceed=SetImageProgress(image,ContrastImageTag,progress++,image->rows);
952         if (proceed == MagickFalse)
953           status=MagickFalse;
954       }
955   }
956   image_view=DestroyCacheView(image_view);
957   return(status);
958 }
959 \f
960 /*
961 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
962 %                                                                             %
963 %                                                                             %
964 %                                                                             %
965 %     C o n t r a s t S t r e t c h I m a g e                                 %
966 %                                                                             %
967 %                                                                             %
968 %                                                                             %
969 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
970 %
971 %  ContrastStretchImage() is a simple image enhancement technique that attempts
972 %  to improve the contrast in an image by 'stretching' the range of intensity
973 %  values it contains to span a desired range of values. It differs from the
974 %  more sophisticated histogram equalization in that it can only apply a
975 %  linear scaling function to the image pixel values.  As a result the
976 %  'enhancement' is less harsh.
977 %
978 %  The format of the ContrastStretchImage method is:
979 %
980 %      MagickBooleanType ContrastStretchImage(Image *image,
981 %        const char *levels,ExceptionInfo *exception)
982 %
983 %  A description of each parameter follows:
984 %
985 %    o image: the image.
986 %
987 %    o black_point: the black point.
988 %
989 %    o white_point: the white point.
990 %
991 %    o levels: Specify the levels where the black and white points have the
992 %      range of 0 to number-of-pixels (e.g. 1%, 10x90%, etc.).
993 %
994 %    o exception: return any errors or warnings in this structure.
995 %
996 */
997 MagickExport MagickBooleanType ContrastStretchImage(Image *image,
998   const double black_point,const double white_point,ExceptionInfo *exception)
999 {
1000 #define MaxRange(color)  ((double) ScaleQuantumToMap((Quantum) (color)))
1001 #define ContrastStretchImageTag  "ContrastStretch/Image"
1002
1003   CacheView
1004     *image_view;
1005
1006   MagickBooleanType
1007     status;
1008
1009   MagickOffsetType
1010     progress;
1011
1012   double
1013     *black,
1014     *histogram,
1015     *stretch_map,
1016     *white;
1017
1018   register ssize_t
1019     i;
1020
1021   size_t
1022     number_channels;
1023
1024   ssize_t
1025     y;
1026
1027   /*
1028     Allocate histogram and stretch map.
1029   */
1030   assert(image != (Image *) NULL);
1031   assert(image->signature == MagickSignature);
1032   if (image->debug != MagickFalse)
1033     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1034   black=(double *) AcquireQuantumMemory(GetPixelChannels(image),sizeof(*black));
1035   white=(double *) AcquireQuantumMemory(GetPixelChannels(image),sizeof(*white));
1036   histogram=(double *) AcquireQuantumMemory(MaxMap+1UL,GetPixelChannels(image)*
1037     sizeof(*histogram));
1038   stretch_map=(double *) AcquireQuantumMemory(MaxMap+1UL,
1039     GetPixelChannels(image)*sizeof(*stretch_map));
1040   if ((black == (double *) NULL) || (white == (double *) NULL) ||
1041       (histogram == (double *) NULL) || (stretch_map == (double *) NULL))
1042     {
1043       if (stretch_map != (double *) NULL)
1044         stretch_map=(double *) RelinquishMagickMemory(stretch_map);
1045       if (histogram != (double *) NULL)
1046         histogram=(double *) RelinquishMagickMemory(histogram);
1047       if (white != (double *) NULL)
1048         white=(double *) RelinquishMagickMemory(white);
1049       if (black != (double *) NULL)
1050         black=(double *) RelinquishMagickMemory(black);
1051       ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
1052         image->filename);
1053     }
1054   /*
1055     Form histogram.
1056   */
1057   status=MagickTrue;
1058   (void) ResetMagickMemory(histogram,0,(MaxMap+1)*GetPixelChannels(image)*
1059     sizeof(*histogram));
1060   image_view=AcquireVirtualCacheView(image,exception);
1061   for (y=0; y < (ssize_t) image->rows; y++)
1062   {
1063     register const Quantum
1064       *restrict p;
1065
1066     register ssize_t
1067       x;
1068
1069     if (status == MagickFalse)
1070       continue;
1071     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
1072     if (p == (const Quantum *) NULL)
1073       {
1074         status=MagickFalse;
1075         continue;
1076       }
1077     for (x=0; x < (ssize_t) image->columns; x++)
1078     {
1079       double
1080         pixel;
1081
1082       register ssize_t
1083         i;
1084
1085       pixel=GetPixelIntensity(image,p);
1086       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1087       {
1088         if (image->channel_mask != DefaultChannels)
1089           pixel=(double) p[i];
1090         histogram[GetPixelChannels(image)*ScaleQuantumToMap(
1091           ClampToQuantum(pixel))+i]++;
1092       }
1093       p+=GetPixelChannels(image);
1094     }
1095   }
1096   image_view=DestroyCacheView(image_view);
1097   /*
1098     Find the histogram boundaries by locating the black/white levels.
1099   */
1100   number_channels=GetPixelChannels(image);
1101 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1102   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1103     dynamic_number_threads(image,image->columns,1,1)
1104 #endif
1105   for (i=0; i < (ssize_t) number_channels; i++)
1106   {
1107     double
1108       intensity;
1109
1110     register ssize_t
1111       j;
1112
1113     black[i]=0.0;
1114     white[i]=MaxRange(QuantumRange);
1115     intensity=0.0;
1116     for (j=0; j <= (ssize_t) MaxMap; j++)
1117     {
1118       intensity+=histogram[GetPixelChannels(image)*j+i];
1119       if (intensity > black_point)
1120         break;
1121     }
1122     black[i]=(double) j;
1123     intensity=0.0;
1124     for (j=(ssize_t) MaxMap; j != 0; j--)
1125     {
1126       intensity+=histogram[GetPixelChannels(image)*j+i];
1127       if (intensity > ((double) image->columns*image->rows-white_point))
1128         break;
1129     }
1130     white[i]=(double) j;
1131   }
1132   histogram=(double *) RelinquishMagickMemory(histogram);
1133   /*
1134     Stretch the histogram to create the stretched image mapping.
1135   */
1136   (void) ResetMagickMemory(stretch_map,0,(MaxMap+1)*GetPixelChannels(image)*
1137     sizeof(*stretch_map));
1138   number_channels=GetPixelChannels(image);
1139 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1140   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1141     dynamic_number_threads(image,image->columns,1,1)
1142 #endif
1143   for (i=0; i < (ssize_t) number_channels; i++)
1144   {
1145     register ssize_t
1146       j;
1147
1148     for (j=0; j <= (ssize_t) MaxMap; j++)
1149     {
1150       if (j < (ssize_t) black[i])
1151         stretch_map[GetPixelChannels(image)*j+i]=0.0;
1152       else
1153         if (j > (ssize_t) white[i])
1154           stretch_map[GetPixelChannels(image)*j+i]=(double)
1155             QuantumRange;
1156         else
1157           if (black[i] != white[i])
1158             stretch_map[GetPixelChannels(image)*j+i]=(double)
1159               ScaleMapToQuantum((double) (MaxMap*(j-black[i])/
1160               (white[i]-black[i])));
1161     }
1162   }
1163   if (image->storage_class == PseudoClass)
1164     {
1165       register ssize_t
1166         j;
1167
1168       /*
1169         Stretch-contrast colormap.
1170       */
1171       for (j=0; j < (ssize_t) image->colors; j++)
1172       {
1173         if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1174           {
1175             i=GetPixelChannelChannel(image,RedPixelChannel);
1176             if (black[i] != white[i])
1177               image->colormap[j].red=stretch_map[GetPixelChannels(image)*
1178                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].red))]+i;
1179           }
1180         if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1181           {
1182             i=GetPixelChannelChannel(image,GreenPixelChannel);
1183             if (black[i] != white[i])
1184               image->colormap[j].green=stretch_map[GetPixelChannels(image)*
1185                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].green))]+i;
1186           }
1187         if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1188           {
1189             i=GetPixelChannelChannel(image,BluePixelChannel);
1190             if (black[i] != white[i])
1191               image->colormap[j].blue=stretch_map[GetPixelChannels(image)*
1192                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].blue))]+i;
1193           }
1194         if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
1195           {
1196             i=GetPixelChannelChannel(image,AlphaPixelChannel);
1197             if (black[i] != white[i])
1198               image->colormap[j].alpha=stretch_map[GetPixelChannels(image)*
1199                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].alpha))]+i;
1200           }
1201       }
1202     }
1203   /*
1204     Stretch-contrast image.
1205   */
1206   status=MagickTrue;
1207   progress=0;
1208   image_view=AcquireAuthenticCacheView(image,exception);
1209 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1210   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1211     dynamic_number_threads(image,image->columns,image->rows,1)
1212 #endif
1213   for (y=0; y < (ssize_t) image->rows; y++)
1214   {
1215     register Quantum
1216       *restrict q;
1217
1218     register ssize_t
1219       x;
1220
1221     if (status == MagickFalse)
1222       continue;
1223     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1224     if (q == (Quantum *) NULL)
1225       {
1226         status=MagickFalse;
1227         continue;
1228       }
1229     for (x=0; x < (ssize_t) image->columns; x++)
1230     {
1231       register ssize_t
1232         i;
1233
1234       if (GetPixelMask(image,q) != 0)
1235         {
1236           q+=GetPixelChannels(image);
1237           continue;
1238         }
1239       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1240       {
1241         PixelChannel
1242           channel;
1243
1244         PixelTrait
1245           traits;
1246
1247         channel=GetPixelChannelChannel(image,i);
1248         traits=GetPixelChannelTraits(image,channel);
1249         if (((traits & UpdatePixelTrait) == 0) || (black[i] == white[i]))
1250           continue;
1251         q[i]=ClampToQuantum(stretch_map[GetPixelChannels(image)*
1252           ScaleQuantumToMap(q[i])+i]);
1253       }
1254       q+=GetPixelChannels(image);
1255     }
1256     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1257       status=MagickFalse;
1258     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1259       {
1260         MagickBooleanType
1261           proceed;
1262
1263 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1264         #pragma omp critical (MagickCore_ContrastStretchImage)
1265 #endif
1266         proceed=SetImageProgress(image,ContrastStretchImageTag,progress++,
1267           image->rows);
1268         if (proceed == MagickFalse)
1269           status=MagickFalse;
1270       }
1271   }
1272   image_view=DestroyCacheView(image_view);
1273   stretch_map=(double *) RelinquishMagickMemory(stretch_map);
1274   white=(double *) RelinquishMagickMemory(white);
1275   black=(double *) RelinquishMagickMemory(black);
1276   return(status);
1277 }
1278 \f
1279 /*
1280 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1281 %                                                                             %
1282 %                                                                             %
1283 %                                                                             %
1284 %     E n h a n c e I m a g e                                                 %
1285 %                                                                             %
1286 %                                                                             %
1287 %                                                                             %
1288 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1289 %
1290 %  EnhanceImage() applies a digital filter that improves the quality of a
1291 %  noisy image.
1292 %
1293 %  The format of the EnhanceImage method is:
1294 %
1295 %      Image *EnhanceImage(const Image *image,ExceptionInfo *exception)
1296 %
1297 %  A description of each parameter follows:
1298 %
1299 %    o image: the image.
1300 %
1301 %    o exception: return any errors or warnings in this structure.
1302 %
1303 */
1304 MagickExport Image *EnhanceImage(const Image *image,ExceptionInfo *exception)
1305 {
1306 #define EnhancePixel(weight) \
1307   mean=((double) r[i]+GetPixelChannel(enhance_image,channel,q))/2.0; \
1308   distance=(double) r[i]-(double) GetPixelChannel( \
1309     enhance_image,channel,q); \
1310   distance_squared=QuantumScale*(2.0*((double) QuantumRange+1.0)+ \
1311     mean)*distance*distance; \
1312   if (distance_squared < ((double) QuantumRange*(double) \
1313       QuantumRange/25.0f)) \
1314     { \
1315       aggregate+=(weight)*r[i]; \
1316       total_weight+=(weight); \
1317     } \
1318   r+=GetPixelChannels(image);
1319 #define EnhanceImageTag  "Enhance/Image"
1320
1321   CacheView
1322     *enhance_view,
1323     *image_view;
1324
1325   Image
1326     *enhance_image;
1327
1328   MagickBooleanType
1329     status;
1330
1331   MagickOffsetType
1332     progress;
1333
1334   ssize_t
1335     y;
1336
1337   /*
1338     Initialize enhanced image attributes.
1339   */
1340   assert(image != (const Image *) NULL);
1341   assert(image->signature == MagickSignature);
1342   if (image->debug != MagickFalse)
1343     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1344   assert(exception != (ExceptionInfo *) NULL);
1345   assert(exception->signature == MagickSignature);
1346   enhance_image=CloneImage(image,image->columns,image->rows,MagickTrue,
1347     exception);
1348   if (enhance_image == (Image *) NULL)
1349     return((Image *) NULL);
1350   if (SetImageStorageClass(enhance_image,DirectClass,exception) == MagickFalse)
1351     {
1352       enhance_image=DestroyImage(enhance_image);
1353       return((Image *) NULL);
1354     }
1355   /*
1356     Enhance image.
1357   */
1358   status=MagickTrue;
1359   progress=0;
1360   image_view=AcquireVirtualCacheView(image,exception);
1361   enhance_view=AcquireAuthenticCacheView(enhance_image,exception);
1362 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1363   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1364     dynamic_number_threads(image,image->columns,image->rows,1)
1365 #endif
1366   for (y=0; y < (ssize_t) image->rows; y++)
1367   {
1368     register const Quantum
1369       *restrict p;
1370
1371     register Quantum
1372       *restrict q;
1373
1374     register ssize_t
1375       x;
1376
1377     ssize_t
1378       center;
1379
1380     if (status == MagickFalse)
1381       continue;
1382     p=GetCacheViewVirtualPixels(image_view,-2,y-2,image->columns+4,5,exception);
1383     q=QueueCacheViewAuthenticPixels(enhance_view,0,y,enhance_image->columns,1,
1384       exception);
1385     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
1386       {
1387         status=MagickFalse;
1388         continue;
1389       }
1390     center=(ssize_t) GetPixelChannels(image)*(2*(image->columns+4)+2);
1391     for (x=0; x < (ssize_t) image->columns; x++)
1392     {
1393       register ssize_t
1394         i;
1395
1396       if (GetPixelMask(image,p) != 0)
1397         {
1398           p+=GetPixelChannels(image);
1399           q+=GetPixelChannels(enhance_image);
1400           continue;
1401         }
1402       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1403       {
1404         double
1405           aggregate,
1406           distance,
1407           distance_squared,
1408           mean,
1409           total_weight;
1410
1411         PixelChannel
1412           channel;
1413
1414         PixelTrait
1415           enhance_traits,
1416           traits;
1417
1418         register const Quantum
1419           *restrict r;
1420
1421         channel=GetPixelChannelChannel(image,i);
1422         traits=GetPixelChannelTraits(image,channel);
1423         enhance_traits=GetPixelChannelTraits(enhance_image,channel);
1424         if ((traits == UndefinedPixelTrait) ||
1425             (enhance_traits == UndefinedPixelTrait))
1426           continue;
1427         SetPixelChannel(enhance_image,channel,p[center+i],q);
1428         if ((enhance_traits & CopyPixelTrait) != 0)
1429           continue;
1430         /*
1431           Compute weighted average of target pixel color components.
1432         */
1433         aggregate=0.0;
1434         total_weight=0.0;
1435         r=p;
1436         EnhancePixel(5.0); EnhancePixel(8.0); EnhancePixel(10.0);
1437           EnhancePixel(8.0); EnhancePixel(5.0);
1438         r=p+1*GetPixelChannels(image)*(image->columns+4);
1439         EnhancePixel(8.0); EnhancePixel(20.0); EnhancePixel(40.0);
1440           EnhancePixel(20.0); EnhancePixel(8.0);
1441         r=p+2*GetPixelChannels(image)*(image->columns+4);
1442         EnhancePixel(10.0); EnhancePixel(40.0); EnhancePixel(80.0);
1443           EnhancePixel(40.0); EnhancePixel(10.0);
1444         r=p+3*GetPixelChannels(image)*(image->columns+4);
1445         EnhancePixel(8.0); EnhancePixel(20.0); EnhancePixel(40.0);
1446           EnhancePixel(20.0); EnhancePixel(8.0);
1447         r=p+4*GetPixelChannels(image)*(image->columns+4);
1448         EnhancePixel(5.0); EnhancePixel(8.0); EnhancePixel(10.0);
1449           EnhancePixel(8.0); EnhancePixel(5.0);
1450         SetPixelChannel(enhance_image,channel,ClampToQuantum(aggregate/
1451           total_weight),q);
1452       }
1453       p+=GetPixelChannels(image);
1454       q+=GetPixelChannels(enhance_image);
1455     }
1456     if (SyncCacheViewAuthenticPixels(enhance_view,exception) == MagickFalse)
1457       status=MagickFalse;
1458     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1459       {
1460         MagickBooleanType
1461           proceed;
1462
1463 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1464         #pragma omp critical (MagickCore_EnhanceImage)
1465 #endif
1466         proceed=SetImageProgress(image,EnhanceImageTag,progress++,image->rows);
1467         if (proceed == MagickFalse)
1468           status=MagickFalse;
1469       }
1470   }
1471   enhance_view=DestroyCacheView(enhance_view);
1472   image_view=DestroyCacheView(image_view);
1473   return(enhance_image);
1474 }
1475 \f
1476 /*
1477 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1478 %                                                                             %
1479 %                                                                             %
1480 %                                                                             %
1481 %     E q u a l i z e I m a g e                                               %
1482 %                                                                             %
1483 %                                                                             %
1484 %                                                                             %
1485 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1486 %
1487 %  EqualizeImage() applies a histogram equalization to the image.
1488 %
1489 %  The format of the EqualizeImage method is:
1490 %
1491 %      MagickBooleanType EqualizeImage(Image *image,ExceptionInfo *exception)
1492 %
1493 %  A description of each parameter follows:
1494 %
1495 %    o image: the image.
1496 %
1497 %    o exception: return any errors or warnings in this structure.
1498 %
1499 */
1500 MagickExport MagickBooleanType EqualizeImage(Image *image,
1501   ExceptionInfo *exception)
1502 {
1503 #define EqualizeImageTag  "Equalize/Image"
1504
1505   CacheView
1506     *image_view;
1507
1508   MagickBooleanType
1509     status;
1510
1511   MagickOffsetType
1512     progress;
1513
1514   double
1515     black[CompositePixelChannel],
1516     *equalize_map,
1517     *histogram,
1518     *map,
1519     white[CompositePixelChannel];
1520
1521   register ssize_t
1522     i;
1523
1524   size_t
1525     number_channels;
1526
1527   ssize_t
1528     y;
1529
1530   /*
1531     Allocate and initialize histogram arrays.
1532   */
1533   assert(image != (Image *) NULL);
1534   assert(image->signature == MagickSignature);
1535   if (image->debug != MagickFalse)
1536     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1537   equalize_map=(double *) AcquireQuantumMemory(MaxMap+1UL,
1538     GetPixelChannels(image)*sizeof(*equalize_map));
1539   histogram=(double *) AcquireQuantumMemory(MaxMap+1UL,
1540     GetPixelChannels(image)*sizeof(*histogram));
1541   map=(double *) AcquireQuantumMemory(MaxMap+1UL,
1542     GetPixelChannels(image)*sizeof(*map));
1543   if ((equalize_map == (double *) NULL) ||
1544       (histogram == (double *) NULL) ||
1545       (map == (double *) NULL))
1546     {
1547       if (map != (double *) NULL)
1548         map=(double *) RelinquishMagickMemory(map);
1549       if (histogram != (double *) NULL)
1550         histogram=(double *) RelinquishMagickMemory(histogram);
1551       if (equalize_map != (double *) NULL)
1552         equalize_map=(double *) RelinquishMagickMemory(equalize_map);
1553       ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
1554         image->filename);
1555     }
1556   /*
1557     Form histogram.
1558   */
1559   status=MagickTrue;
1560   (void) ResetMagickMemory(histogram,0,(MaxMap+1)*GetPixelChannels(image)*
1561     sizeof(*histogram));
1562   image_view=AcquireVirtualCacheView(image,exception);
1563   for (y=0; y < (ssize_t) image->rows; y++)
1564   {
1565     register const Quantum
1566       *restrict p;
1567
1568     register ssize_t
1569       x;
1570
1571     if (status == MagickFalse)
1572       continue;
1573     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
1574     if (p == (const Quantum *) NULL)
1575       {
1576         status=MagickFalse;
1577         continue;
1578       }
1579     for (x=0; x < (ssize_t) image->columns; x++)
1580     {
1581       register ssize_t
1582         i;
1583
1584       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1585         histogram[GetPixelChannels(image)*ScaleQuantumToMap(p[i])+i]++;
1586       p+=GetPixelChannels(image);
1587     }
1588   }
1589   image_view=DestroyCacheView(image_view);
1590   /*
1591     Integrate the histogram to get the equalization map.
1592   */
1593   number_channels=GetPixelChannels(image);
1594 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1595   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1596     dynamic_number_threads(image,image->columns,1,1)
1597 #endif
1598   for (i=0; i < (ssize_t) number_channels; i++)
1599   {
1600     double
1601       intensity;
1602
1603     register ssize_t
1604       j;
1605
1606     intensity=0.0;
1607     for (j=0; j <= (ssize_t) MaxMap; j++)
1608     {
1609       intensity+=histogram[GetPixelChannels(image)*j+i];
1610       map[GetPixelChannels(image)*j+i]=intensity;
1611     }
1612   }
1613   (void) ResetMagickMemory(equalize_map,0,(MaxMap+1)*GetPixelChannels(image)*
1614     sizeof(*equalize_map));
1615   number_channels=GetPixelChannels(image);
1616 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1617   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1618     dynamic_number_threads(image,image->columns,1,1)
1619 #endif
1620   for (i=0; i < (ssize_t) number_channels; i++)
1621   {
1622     register ssize_t
1623       j;
1624
1625     black[i]=map[i];
1626     white[i]=map[GetPixelChannels(image)*MaxMap+i];
1627     if (black[i] != white[i])
1628       for (j=0; j <= (ssize_t) MaxMap; j++)
1629         equalize_map[GetPixelChannels(image)*j+i]=(double)
1630           ScaleMapToQuantum((double) ((MaxMap*(map[
1631           GetPixelChannels(image)*j+i]-black[i]))/(white[i]-black[i])));
1632   }
1633   histogram=(double *) RelinquishMagickMemory(histogram);
1634   map=(double *) RelinquishMagickMemory(map);
1635   if (image->storage_class == PseudoClass)
1636     {
1637       PixelChannel
1638         channel;
1639
1640       register ssize_t
1641         j;
1642
1643       /*
1644         Equalize colormap.
1645       */
1646       for (j=0; j < (ssize_t) image->colors; j++)
1647       {
1648         if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1649           {
1650             channel=GetPixelChannelChannel(image,RedPixelChannel);
1651             if (black[channel] != white[channel])
1652               image->colormap[j].red=equalize_map[GetPixelChannels(image)*
1653                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].red))]+
1654                 channel;
1655           }
1656         if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1657           {
1658             channel=GetPixelChannelChannel(image,GreenPixelChannel);
1659             if (black[channel] != white[channel])
1660               image->colormap[j].green=equalize_map[GetPixelChannels(image)*
1661                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].green))]+
1662                 channel;
1663           }
1664         if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1665           {
1666             channel=GetPixelChannelChannel(image,BluePixelChannel);
1667             if (black[channel] != white[channel])
1668               image->colormap[j].blue=equalize_map[GetPixelChannels(image)*
1669                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].blue))]+
1670                 channel;
1671           }
1672         if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
1673           {
1674             channel=GetPixelChannelChannel(image,AlphaPixelChannel);
1675             if (black[channel] != white[channel])
1676               image->colormap[j].alpha=equalize_map[GetPixelChannels(image)*
1677                 ScaleQuantumToMap(ClampToQuantum(image->colormap[j].alpha))]+
1678                 channel;
1679           }
1680       }
1681     }
1682   /*
1683     Equalize image.
1684   */
1685   progress=0;
1686   image_view=AcquireAuthenticCacheView(image,exception);
1687 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1688   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1689     dynamic_number_threads(image,image->columns,image->rows,1)
1690 #endif
1691   for (y=0; y < (ssize_t) image->rows; y++)
1692   {
1693     register Quantum
1694       *restrict q;
1695
1696     register ssize_t
1697       x;
1698
1699     if (status == MagickFalse)
1700       continue;
1701     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1702     if (q == (Quantum *) NULL)
1703       {
1704         status=MagickFalse;
1705         continue;
1706       }
1707     for (x=0; x < (ssize_t) image->columns; x++)
1708     {
1709       register ssize_t
1710         i;
1711
1712       if (GetPixelMask(image,q) != 0)
1713         {
1714           q+=GetPixelChannels(image);
1715           continue;
1716         }
1717       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1718       {
1719         PixelChannel
1720           channel;
1721
1722         PixelTrait
1723           traits;
1724
1725         channel=GetPixelChannelChannel(image,i);
1726         traits=GetPixelChannelTraits(image,channel);
1727         if (((traits & UpdatePixelTrait) == 0) || (black[i] == white[i]))
1728           continue;
1729         q[i]=ClampToQuantum(equalize_map[GetPixelChannels(image)*
1730           ScaleQuantumToMap(q[i])+i]);
1731       }
1732       q+=GetPixelChannels(image);
1733     }
1734     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1735       status=MagickFalse;
1736     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1737       {
1738         MagickBooleanType
1739           proceed;
1740
1741 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1742         #pragma omp critical (MagickCore_EqualizeImage)
1743 #endif
1744         proceed=SetImageProgress(image,EqualizeImageTag,progress++,image->rows);
1745         if (proceed == MagickFalse)
1746           status=MagickFalse;
1747       }
1748   }
1749   image_view=DestroyCacheView(image_view);
1750   equalize_map=(double *) RelinquishMagickMemory(equalize_map);
1751   return(status);
1752 }
1753 \f
1754 /*
1755 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1756 %                                                                             %
1757 %                                                                             %
1758 %                                                                             %
1759 %     G a m m a I m a g e                                                     %
1760 %                                                                             %
1761 %                                                                             %
1762 %                                                                             %
1763 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1764 %
1765 %  GammaImage() gamma-corrects a particular image channel.  The same
1766 %  image viewed on different devices will have perceptual differences in the
1767 %  way the image's intensities are represented on the screen.  Specify
1768 %  individual gamma levels for the red, green, and blue channels, or adjust
1769 %  all three with the gamma parameter.  Values typically range from 0.8 to 2.3.
1770 %
1771 %  You can also reduce the influence of a particular channel with a gamma
1772 %  value of 0.
1773 %
1774 %  The format of the GammaImage method is:
1775 %
1776 %      MagickBooleanType GammaImage(Image *image,const double gamma,
1777 %        ExceptionInfo *exception)
1778 %
1779 %  A description of each parameter follows:
1780 %
1781 %    o image: the image.
1782 %
1783 %    o level: the image gamma as a string (e.g. 1.6,1.2,1.0).
1784 %
1785 %    o gamma: the image gamma.
1786 %
1787 */
1788 MagickExport MagickBooleanType GammaImage(Image *image,const double gamma,
1789   ExceptionInfo *exception)
1790 {
1791 #define GammaCorrectImageTag  "GammaCorrect/Image"
1792
1793   CacheView
1794     *image_view;
1795
1796   MagickBooleanType
1797     status;
1798
1799   MagickOffsetType
1800     progress;
1801
1802   Quantum
1803     *gamma_map;
1804
1805   register ssize_t
1806     i;
1807
1808   ssize_t
1809     y;
1810
1811   /*
1812     Allocate and initialize gamma maps.
1813   */
1814   assert(image != (Image *) NULL);
1815   assert(image->signature == MagickSignature);
1816   if (image->debug != MagickFalse)
1817     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
1818   if (gamma == 1.0)
1819     return(MagickTrue);
1820   gamma_map=(Quantum *) AcquireQuantumMemory(MaxMap+1UL,sizeof(*gamma_map));
1821   if (gamma_map == (Quantum *) NULL)
1822     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
1823       image->filename);
1824   (void) ResetMagickMemory(gamma_map,0,(MaxMap+1)*sizeof(*gamma_map));
1825   if (gamma != 0.0)
1826 #if defined(MAGICKCORE_OPENMP_SUPPORT) && (MaxMap > 256)
1827     #pragma omp parallel for \
1828       dynamic_number_threads(image,image->columns,1,1)
1829 #endif
1830     for (i=0; i <= (ssize_t) MaxMap; i++)
1831       gamma_map[i]=ScaleMapToQuantum((double) (MaxMap*pow((double) i/
1832         MaxMap,1.0/gamma)));
1833   if (image->storage_class == PseudoClass)
1834     for (i=0; i < (ssize_t) image->colors; i++)
1835     {
1836       /*
1837         Gamma-correct colormap.
1838       */
1839       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
1840         image->colormap[i].red=(double) gamma_map[
1841           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))];
1842       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
1843         image->colormap[i].green=(double) gamma_map[
1844           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].green))];
1845       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
1846         image->colormap[i].blue=(double) gamma_map[
1847           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].blue))];
1848       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
1849         image->colormap[i].alpha=(double) gamma_map[
1850           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].alpha))];
1851     }
1852   /*
1853     Gamma-correct image.
1854   */
1855   status=MagickTrue;
1856   progress=0;
1857   image_view=AcquireAuthenticCacheView(image,exception);
1858 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1859   #pragma omp parallel for schedule(static,4) shared(progress,status) \
1860     dynamic_number_threads(image,image->columns,image->rows,1)
1861 #endif
1862   for (y=0; y < (ssize_t) image->rows; y++)
1863   {
1864     register Quantum
1865       *restrict q;
1866
1867     register ssize_t
1868       x;
1869
1870     if (status == MagickFalse)
1871       continue;
1872     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
1873     if (q == (Quantum *) NULL)
1874       {
1875         status=MagickFalse;
1876         continue;
1877       }
1878     for (x=0; x < (ssize_t) image->columns; x++)
1879     {
1880       register ssize_t
1881         i;
1882
1883       if (GetPixelMask(image,q) != 0)
1884         {
1885           q+=GetPixelChannels(image);
1886           continue;
1887         }
1888       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
1889       {
1890         PixelChannel
1891           channel;
1892
1893         PixelTrait
1894           traits;
1895
1896         channel=GetPixelChannelChannel(image,i);
1897         traits=GetPixelChannelTraits(image,channel);
1898         if ((traits & UpdatePixelTrait) == 0)
1899           continue;
1900         q[i]=gamma_map[ScaleQuantumToMap(q[i])];
1901       }
1902       q+=GetPixelChannels(image);
1903     }
1904     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
1905       status=MagickFalse;
1906     if (image->progress_monitor != (MagickProgressMonitor) NULL)
1907       {
1908         MagickBooleanType
1909           proceed;
1910
1911 #if defined(MAGICKCORE_OPENMP_SUPPORT)
1912         #pragma omp critical (MagickCore_GammaImage)
1913 #endif
1914         proceed=SetImageProgress(image,GammaCorrectImageTag,progress++,
1915           image->rows);
1916         if (proceed == MagickFalse)
1917           status=MagickFalse;
1918       }
1919   }
1920   image_view=DestroyCacheView(image_view);
1921   gamma_map=(Quantum *) RelinquishMagickMemory(gamma_map);
1922   if (image->gamma != 0.0)
1923     image->gamma*=gamma;
1924   return(status);
1925 }
1926 \f
1927 /*
1928 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1929 %                                                                             %
1930 %                                                                             %
1931 %                                                                             %
1932 %     H a l d C l u t I m a g e                                               %
1933 %                                                                             %
1934 %                                                                             %
1935 %                                                                             %
1936 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1937 %
1938 %  HaldClutImage() applies a Hald color lookup table to the image.  A Hald
1939 %  color lookup table is a 3-dimensional color cube mapped to 2 dimensions.
1940 %  Create it with the HALD coder.  You can apply any color transformation to
1941 %  the Hald image and then use this method to apply the transform to the
1942 %  image.
1943 %
1944 %  The format of the HaldClutImage method is:
1945 %
1946 %      MagickBooleanType HaldClutImage(Image *image,Image *hald_image,
1947 %        ExceptionInfo *exception)
1948 %
1949 %  A description of each parameter follows:
1950 %
1951 %    o image: the image, which is replaced by indexed CLUT values
1952 %
1953 %    o hald_image: the color lookup table image for replacement color values.
1954 %
1955 %    o exception: return any errors or warnings in this structure.
1956 %
1957 */
1958
1959 static inline size_t MagickMin(const size_t x,const size_t y)
1960 {
1961   if (x < y)
1962     return(x);
1963   return(y);
1964 }
1965
1966 MagickExport MagickBooleanType HaldClutImage(Image *image,
1967   const Image *hald_image,ExceptionInfo *exception)
1968 {
1969 #define HaldClutImageTag  "Clut/Image"
1970
1971   typedef struct _HaldInfo
1972   {
1973     double
1974       x,
1975       y,
1976       z;
1977   } HaldInfo;
1978
1979   CacheView
1980     *hald_view,
1981     *image_view;
1982
1983   double
1984     width;
1985
1986   MagickBooleanType
1987     status;
1988
1989   MagickOffsetType
1990     progress;
1991
1992   PixelInfo
1993     zero;
1994
1995   size_t
1996     cube_size,
1997     length,
1998     level;
1999
2000   ssize_t
2001     y;
2002
2003   assert(image != (Image *) NULL);
2004   assert(image->signature == MagickSignature);
2005   if (image->debug != MagickFalse)
2006     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2007   assert(hald_image != (Image *) NULL);
2008   assert(hald_image->signature == MagickSignature);
2009   if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
2010     return(MagickFalse);
2011   if (IsGrayColorspace(image->colorspace) != MagickFalse)
2012     (void) TransformImageColorspace(image,RGBColorspace,exception);
2013   if (image->alpha_trait != BlendPixelTrait)
2014     (void) SetImageAlphaChannel(image,OpaqueAlphaChannel,exception);
2015   /*
2016     Hald clut image.
2017   */
2018   status=MagickTrue;
2019   progress=0;
2020   length=MagickMin(hald_image->columns,hald_image->rows);
2021   for (level=2; (level*level*level) < length; level++) ;
2022   level*=level;
2023   cube_size=level*level;
2024   width=(double) hald_image->columns;
2025   GetPixelInfo(hald_image,&zero);
2026   hald_view=AcquireVirtualCacheView(hald_image,exception);
2027   image_view=AcquireAuthenticCacheView(image,exception);
2028 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2029   #pragma omp parallel for schedule(static,4) shared(progress,status) \
2030     dynamic_number_threads(image,image->columns,image->rows,1)
2031 #endif
2032   for (y=0; y < (ssize_t) image->rows; y++)
2033   {
2034     register Quantum
2035       *restrict q;
2036
2037     register ssize_t
2038       x;
2039
2040     if (status == MagickFalse)
2041       continue;
2042     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2043     if (q == (Quantum *) NULL)
2044       {
2045         status=MagickFalse;
2046         continue;
2047       }
2048     for (x=0; x < (ssize_t) image->columns; x++)
2049     {
2050       double
2051         offset;
2052
2053       HaldInfo
2054         point;
2055
2056       PixelInfo
2057         pixel,
2058         pixel1,
2059         pixel2,
2060         pixel3,
2061         pixel4;
2062
2063       point.x=QuantumScale*(level-1.0)*GetPixelRed(image,q);
2064       point.y=QuantumScale*(level-1.0)*GetPixelGreen(image,q);
2065       point.z=QuantumScale*(level-1.0)*GetPixelBlue(image,q);
2066       offset=point.x+level*floor(point.y)+cube_size*floor(point.z);
2067       point.x-=floor(point.x);
2068       point.y-=floor(point.y);
2069       point.z-=floor(point.z);
2070       pixel1=zero;
2071       (void) InterpolatePixelInfo(image,hald_view,image->interpolate,
2072         fmod(offset,width),floor(offset/width),&pixel1,exception);
2073       pixel2=zero;
2074       (void) InterpolatePixelInfo(image,hald_view,image->interpolate,
2075         fmod(offset+level,width),floor((offset+level)/width),&pixel2,exception);
2076       pixel3=zero;
2077       CompositePixelInfoAreaBlend(&pixel1,pixel1.alpha,&pixel2,pixel2.alpha,
2078         point.y,&pixel3);
2079       offset+=cube_size;
2080       (void) InterpolatePixelInfo(image,hald_view,image->interpolate,
2081         fmod(offset,width),floor(offset/width),&pixel1,exception);
2082       (void) InterpolatePixelInfo(image,hald_view,image->interpolate,
2083         fmod(offset+level,width),floor((offset+level)/width),&pixel2,exception);
2084       pixel4=zero;
2085       CompositePixelInfoAreaBlend(&pixel1,pixel1.alpha,&pixel2,pixel2.alpha,
2086         point.y,&pixel4);
2087       pixel=zero;
2088       CompositePixelInfoAreaBlend(&pixel3,pixel3.alpha,&pixel4,pixel4.alpha,
2089         point.z,&pixel);
2090       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2091         SetPixelRed(image,ClampToQuantum(pixel.red),q);
2092       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2093         SetPixelGreen(image,ClampToQuantum(pixel.green),q);
2094       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2095         SetPixelBlue(image,ClampToQuantum(pixel.blue),q);
2096       if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
2097           (image->colorspace == CMYKColorspace))
2098         SetPixelBlack(image,ClampToQuantum(pixel.black),q);
2099       if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
2100           (image->alpha_trait == BlendPixelTrait))
2101         SetPixelAlpha(image,ClampToQuantum(pixel.alpha),q);
2102       q+=GetPixelChannels(image);
2103     }
2104     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2105       status=MagickFalse;
2106     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2107       {
2108         MagickBooleanType
2109           proceed;
2110
2111 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2112         #pragma omp critical (MagickCore_HaldClutImage)
2113 #endif
2114         proceed=SetImageProgress(image,HaldClutImageTag,progress++,image->rows);
2115         if (proceed == MagickFalse)
2116           status=MagickFalse;
2117       }
2118   }
2119   hald_view=DestroyCacheView(hald_view);
2120   image_view=DestroyCacheView(image_view);
2121   return(status);
2122 }
2123 \f
2124 /*
2125 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2126 %                                                                             %
2127 %                                                                             %
2128 %                                                                             %
2129 %     L e v e l I m a g e                                                     %
2130 %                                                                             %
2131 %                                                                             %
2132 %                                                                             %
2133 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2134 %
2135 %  LevelImage() adjusts the levels of a particular image channel by
2136 %  scaling the colors falling between specified white and black points to
2137 %  the full available quantum range.
2138 %
2139 %  The parameters provided represent the black, and white points.  The black
2140 %  point specifies the darkest color in the image. Colors darker than the
2141 %  black point are set to zero.  White point specifies the lightest color in
2142 %  the image.  Colors brighter than the white point are set to the maximum
2143 %  quantum value.
2144 %
2145 %  If a '!' flag is given, map black and white colors to the given levels
2146 %  rather than mapping those levels to black and white.  See
2147 %  LevelizeImage() below.
2148 %
2149 %  Gamma specifies a gamma correction to apply to the image.
2150 %
2151 %  The format of the LevelImage method is:
2152 %
2153 %      MagickBooleanType LevelImage(Image *image,const double black_point,
2154 %        const double white_point,const double gamma,ExceptionInfo *exception)
2155 %
2156 %  A description of each parameter follows:
2157 %
2158 %    o image: the image.
2159 %
2160 %    o black_point: The level to map zero (black) to.
2161 %
2162 %    o white_point: The level to map QuantumRange (white) to.
2163 %
2164 %    o exception: return any errors or warnings in this structure.
2165 %
2166 */
2167
2168 static inline double LevelPixel(const double black_point,
2169   const double white_point,const double gamma,const double pixel)
2170 {
2171   double
2172     level_pixel,
2173     scale;
2174
2175   scale=(white_point != black_point) ? 1.0/(white_point-black_point) : 1.0;
2176   level_pixel=(double) QuantumRange*pow(scale*((double) pixel-
2177     black_point),1.0/gamma);
2178   return(level_pixel);
2179 }
2180
2181 MagickExport MagickBooleanType LevelImage(Image *image,const double black_point,
2182   const double white_point,const double gamma,ExceptionInfo *exception)
2183 {
2184 #define LevelImageTag  "Level/Image"
2185
2186   CacheView
2187     *image_view;
2188
2189   MagickBooleanType
2190     status;
2191
2192   MagickOffsetType
2193     progress;
2194
2195   register ssize_t
2196     i;
2197
2198   ssize_t
2199     y;
2200
2201   /*
2202     Allocate and initialize levels map.
2203   */
2204   assert(image != (Image *) NULL);
2205   assert(image->signature == MagickSignature);
2206   if (image->debug != MagickFalse)
2207     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2208   if (image->storage_class == PseudoClass)
2209     for (i=0; i < (ssize_t) image->colors; i++)
2210     {
2211       /*
2212         Level colormap.
2213       */
2214       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2215         image->colormap[i].red=(double) ClampToQuantum(LevelPixel(black_point,
2216           white_point,gamma,image->colormap[i].red));
2217       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2218         image->colormap[i].green=(double) ClampToQuantum(LevelPixel(black_point,
2219           white_point,gamma,image->colormap[i].green));
2220       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2221         image->colormap[i].blue=(double) ClampToQuantum(LevelPixel(black_point,
2222           white_point,gamma,image->colormap[i].blue));
2223       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
2224         image->colormap[i].alpha=(double) ClampToQuantum(LevelPixel(black_point,
2225           white_point,gamma,image->colormap[i].alpha));
2226       }
2227   /*
2228     Level image.
2229   */
2230   status=MagickTrue;
2231   progress=0;
2232   image_view=AcquireAuthenticCacheView(image,exception);
2233 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2234   #pragma omp parallel for schedule(static,4) shared(progress,status) \
2235     dynamic_number_threads(image,image->columns,image->rows,1)
2236 #endif
2237   for (y=0; y < (ssize_t) image->rows; y++)
2238   {
2239     register Quantum
2240       *restrict q;
2241
2242     register ssize_t
2243       x;
2244
2245     if (status == MagickFalse)
2246       continue;
2247     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2248     if (q == (Quantum *) NULL)
2249       {
2250         status=MagickFalse;
2251         continue;
2252       }
2253     for (x=0; x < (ssize_t) image->columns; x++)
2254     {
2255       register ssize_t
2256         i;
2257
2258       if (GetPixelMask(image,q) != 0)
2259         {
2260           q+=GetPixelChannels(image);
2261           continue;
2262         }
2263       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2264       {
2265         PixelChannel
2266           channel;
2267
2268         PixelTrait
2269           traits;
2270
2271         channel=GetPixelChannelChannel(image,i);
2272         traits=GetPixelChannelTraits(image,channel);
2273         if ((traits & UpdatePixelTrait) == 0)
2274           continue;
2275         q[i]=ClampToQuantum(LevelPixel(black_point,white_point,gamma,
2276           (double) q[i]));
2277       }
2278       q+=GetPixelChannels(image);
2279     }
2280     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2281       status=MagickFalse;
2282     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2283       {
2284         MagickBooleanType
2285           proceed;
2286
2287 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2288         #pragma omp critical (MagickCore_LevelImage)
2289 #endif
2290         proceed=SetImageProgress(image,LevelImageTag,progress++,image->rows);
2291         if (proceed == MagickFalse)
2292           status=MagickFalse;
2293       }
2294   }
2295   image_view=DestroyCacheView(image_view);
2296   if (status != MagickFalse)
2297     (void) ClampImage(image,exception);
2298   return(status);
2299 }
2300 \f
2301 /*
2302 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2303 %                                                                             %
2304 %                                                                             %
2305 %                                                                             %
2306 %     L e v e l i z e I m a g e                                               %
2307 %                                                                             %
2308 %                                                                             %
2309 %                                                                             %
2310 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2311 %
2312 %  LevelizeImage() applies the reversed LevelImage() operation to just
2313 %  the specific channels specified.  It compresses the full range of color
2314 %  values, so that they lie between the given black and white points. Gamma is
2315 %  applied before the values are mapped.
2316 %
2317 %  LevelizeImage() can be called with by using a +level command line
2318 %  API option, or using a '!' on a -level or LevelImage() geometry string.
2319 %
2320 %  It can be used to de-contrast a greyscale image to the exact levels
2321 %  specified.  Or by using specific levels for each channel of an image you
2322 %  can convert a gray-scale image to any linear color gradient, according to
2323 %  those levels.
2324 %
2325 %  The format of the LevelizeImage method is:
2326 %
2327 %      MagickBooleanType LevelizeImage(Image *image,const double black_point,
2328 %        const double white_point,const double gamma,ExceptionInfo *exception)
2329 %
2330 %  A description of each parameter follows:
2331 %
2332 %    o image: the image.
2333 %
2334 %    o black_point: The level to map zero (black) to.
2335 %
2336 %    o white_point: The level to map QuantumRange (white) to.
2337 %
2338 %    o gamma: adjust gamma by this factor before mapping values.
2339 %
2340 %    o exception: return any errors or warnings in this structure.
2341 %
2342 */
2343 MagickExport MagickBooleanType LevelizeImage(Image *image,
2344   const double black_point,const double white_point,const double gamma,
2345   ExceptionInfo *exception)
2346 {
2347 #define LevelizeImageTag  "Levelize/Image"
2348 #define LevelizeValue(x) (ClampToQuantum((pow((double) (QuantumScale*(x)), \
2349   1.0/gamma))*(white_point-black_point)+black_point))
2350
2351   CacheView
2352     *image_view;
2353
2354   MagickBooleanType
2355     status;
2356
2357   MagickOffsetType
2358     progress;
2359
2360   register ssize_t
2361     i;
2362
2363   ssize_t
2364     y;
2365
2366   /*
2367     Allocate and initialize levels map.
2368   */
2369   assert(image != (Image *) NULL);
2370   assert(image->signature == MagickSignature);
2371   if (image->debug != MagickFalse)
2372     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2373   if (IsGrayColorspace(image->colorspace) != MagickFalse)
2374     (void) SetImageColorspace(image,RGBColorspace,exception);
2375   if (image->storage_class == PseudoClass)
2376     for (i=0; i < (ssize_t) image->colors; i++)
2377     {
2378       /*
2379         Level colormap.
2380       */
2381       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2382         image->colormap[i].red=(double) LevelizeValue(
2383           image->colormap[i].red);
2384       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2385         image->colormap[i].green=(double) LevelizeValue(
2386           image->colormap[i].green);
2387       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2388         image->colormap[i].blue=(double) LevelizeValue(
2389           image->colormap[i].blue);
2390       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
2391         image->colormap[i].alpha=(double) LevelizeValue(
2392           image->colormap[i].alpha);
2393     }
2394   /*
2395     Level image.
2396   */
2397   status=MagickTrue;
2398   progress=0;
2399   image_view=AcquireAuthenticCacheView(image,exception);
2400 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2401   #pragma omp parallel for schedule(static,4) shared(progress,status) \
2402     dynamic_number_threads(image,image->columns,image->rows,1)
2403 #endif
2404   for (y=0; y < (ssize_t) image->rows; y++)
2405   {
2406     register Quantum
2407       *restrict q;
2408
2409     register ssize_t
2410       x;
2411
2412     if (status == MagickFalse)
2413       continue;
2414     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2415     if (q == (Quantum *) NULL)
2416       {
2417         status=MagickFalse;
2418         continue;
2419       }
2420     for (x=0; x < (ssize_t) image->columns; x++)
2421     {
2422       register ssize_t
2423         i;
2424
2425       if (GetPixelMask(image,q) != 0)
2426         {
2427           q+=GetPixelChannels(image);
2428           continue;
2429         }
2430       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
2431       {
2432         PixelChannel
2433           channel;
2434
2435         PixelTrait
2436           traits;
2437
2438         channel=GetPixelChannelChannel(image,i);
2439         traits=GetPixelChannelTraits(image,channel);
2440         if ((traits & UpdatePixelTrait) == 0)
2441           continue;
2442         q[i]=LevelizeValue(q[i]);
2443       }
2444       q+=GetPixelChannels(image);
2445     }
2446     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
2447       status=MagickFalse;
2448     if (image->progress_monitor != (MagickProgressMonitor) NULL)
2449       {
2450         MagickBooleanType
2451           proceed;
2452
2453 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2454         #pragma omp critical (MagickCore_LevelizeImage)
2455 #endif
2456         proceed=SetImageProgress(image,LevelizeImageTag,progress++,image->rows);
2457         if (proceed == MagickFalse)
2458           status=MagickFalse;
2459       }
2460   }
2461   image_view=DestroyCacheView(image_view);
2462   return(status);
2463 }
2464 \f
2465 /*
2466 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2467 %                                                                             %
2468 %                                                                             %
2469 %                                                                             %
2470 %     L e v e l I m a g e C o l o r s                                         %
2471 %                                                                             %
2472 %                                                                             %
2473 %                                                                             %
2474 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2475 %
2476 %  LevelImageColors() maps the given color to "black" and "white" values,
2477 %  linearly spreading out the colors, and level values on a channel by channel
2478 %  bases, as per LevelImage().  The given colors allows you to specify
2479 %  different level ranges for each of the color channels separately.
2480 %
2481 %  If the boolean 'invert' is set true the image values will modifyed in the
2482 %  reverse direction. That is any existing "black" and "white" colors in the
2483 %  image will become the color values given, with all other values compressed
2484 %  appropriatally.  This effectivally maps a greyscale gradient into the given
2485 %  color gradient.
2486 %
2487 %  The format of the LevelImageColors method is:
2488 %
2489 %    MagickBooleanType LevelImageColors(Image *image,
2490 %      const PixelInfo *black_color,const PixelInfo *white_color,
2491 %      const MagickBooleanType invert,ExceptionInfo *exception)
2492 %
2493 %  A description of each parameter follows:
2494 %
2495 %    o image: the image.
2496 %
2497 %    o black_color: The color to map black to/from
2498 %
2499 %    o white_point: The color to map white to/from
2500 %
2501 %    o invert: if true map the colors (levelize), rather than from (level)
2502 %
2503 %    o exception: return any errors or warnings in this structure.
2504 %
2505 */
2506 MagickExport MagickBooleanType LevelImageColors(Image *image,
2507   const PixelInfo *black_color,const PixelInfo *white_color,
2508   const MagickBooleanType invert,ExceptionInfo *exception)
2509 {
2510   ChannelType
2511     channel_mask;
2512
2513   MagickStatusType
2514     status;
2515
2516   /*
2517     Allocate and initialize levels map.
2518   */
2519   assert(image != (Image *) NULL);
2520   assert(image->signature == MagickSignature);
2521   if (image->debug != MagickFalse)
2522     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2523   status=MagickFalse;
2524   if (invert == MagickFalse)
2525     {
2526       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2527         {
2528           channel_mask=SetImageChannelMask(image,RedChannel);
2529           status|=LevelImage(image,black_color->red,white_color->red,1.0,
2530             exception);
2531           (void) SetImageChannelMask(image,channel_mask);
2532         }
2533       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2534         {
2535           channel_mask=SetImageChannelMask(image,GreenChannel);
2536           status|=LevelImage(image,black_color->green,white_color->green,1.0,
2537             exception);
2538           (void) SetImageChannelMask(image,channel_mask);
2539         }
2540       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2541         {
2542           channel_mask=SetImageChannelMask(image,BlueChannel);
2543           status|=LevelImage(image,black_color->blue,white_color->blue,1.0,
2544             exception);
2545           (void) SetImageChannelMask(image,channel_mask);
2546         }
2547       if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
2548           (image->colorspace == CMYKColorspace))
2549         {
2550           channel_mask=SetImageChannelMask(image,BlackChannel);
2551           status|=LevelImage(image,black_color->black,white_color->black,1.0,
2552             exception);
2553           (void) SetImageChannelMask(image,channel_mask);
2554         }
2555       if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
2556           (image->alpha_trait == BlendPixelTrait))
2557         {
2558           channel_mask=SetImageChannelMask(image,AlphaChannel);
2559           status|=LevelImage(image,black_color->alpha,white_color->alpha,1.0,
2560             exception);
2561           (void) SetImageChannelMask(image,channel_mask);
2562         }
2563     }
2564   else
2565     {
2566       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
2567         {
2568           channel_mask=SetImageChannelMask(image,RedChannel);
2569           status|=LevelizeImage(image,black_color->red,white_color->red,1.0,
2570             exception);
2571           (void) SetImageChannelMask(image,channel_mask);
2572         }
2573       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
2574         {
2575           channel_mask=SetImageChannelMask(image,GreenChannel);
2576           status|=LevelizeImage(image,black_color->green,white_color->green,1.0,
2577             exception);
2578           (void) SetImageChannelMask(image,channel_mask);
2579         }
2580       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
2581         {
2582           channel_mask=SetImageChannelMask(image,BlueChannel);
2583           status|=LevelizeImage(image,black_color->blue,white_color->blue,1.0,
2584             exception);
2585           (void) SetImageChannelMask(image,channel_mask);
2586         }
2587       if (((GetPixelBlackTraits(image) & UpdatePixelTrait) != 0) &&
2588           (image->colorspace == CMYKColorspace))
2589         {
2590           channel_mask=SetImageChannelMask(image,BlackChannel);
2591           status|=LevelizeImage(image,black_color->black,white_color->black,1.0,
2592             exception);
2593           (void) SetImageChannelMask(image,channel_mask);
2594         }
2595       if (((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0) &&
2596           (image->alpha_trait == BlendPixelTrait))
2597         {
2598           channel_mask=SetImageChannelMask(image,AlphaChannel);
2599           status|=LevelizeImage(image,black_color->alpha,white_color->alpha,1.0,
2600             exception);
2601           (void) SetImageChannelMask(image,channel_mask);
2602         }
2603     }
2604   return(status == 0 ? MagickFalse : MagickTrue);
2605 }
2606 \f
2607 /*
2608 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2609 %                                                                             %
2610 %                                                                             %
2611 %                                                                             %
2612 %     L i n e a r S t r e t c h I m a g e                                     %
2613 %                                                                             %
2614 %                                                                             %
2615 %                                                                             %
2616 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2617 %
2618 %  LinearStretchImage() discards any pixels below the black point and above
2619 %  the white point and levels the remaining pixels.
2620 %
2621 %  The format of the LinearStretchImage method is:
2622 %
2623 %      MagickBooleanType LinearStretchImage(Image *image,
2624 %        const double black_point,const double white_point,
2625 %        ExceptionInfo *exception)
2626 %
2627 %  A description of each parameter follows:
2628 %
2629 %    o image: the image.
2630 %
2631 %    o black_point: the black point.
2632 %
2633 %    o white_point: the white point.
2634 %
2635 %    o exception: return any errors or warnings in this structure.
2636 %
2637 */
2638 MagickExport MagickBooleanType LinearStretchImage(Image *image,
2639   const double black_point,const double white_point,ExceptionInfo *exception)
2640 {
2641 #define LinearStretchImageTag  "LinearStretch/Image"
2642
2643   CacheView
2644     *image_view;
2645
2646   MagickBooleanType
2647     status;
2648
2649   double
2650     *histogram,
2651     intensity;
2652
2653   ssize_t
2654     black,
2655     white,
2656     y;
2657
2658   /*
2659     Allocate histogram and linear map.
2660   */
2661   assert(image != (Image *) NULL);
2662   assert(image->signature == MagickSignature);
2663   histogram=(double *) AcquireQuantumMemory(MaxMap+1UL,
2664     sizeof(*histogram));
2665   if (histogram == (double *) NULL)
2666     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
2667       image->filename);
2668   /*
2669     Form histogram.
2670   */
2671   (void) ResetMagickMemory(histogram,0,(MaxMap+1)*sizeof(*histogram));
2672   image_view=AcquireVirtualCacheView(image,exception);
2673   for (y=0; y < (ssize_t) image->rows; y++)
2674   {
2675     register const Quantum
2676       *restrict p;
2677
2678     register ssize_t
2679       x;
2680
2681     p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
2682     if (p == (const Quantum *) NULL)
2683       break;
2684     for (x=0; x < (ssize_t) image->columns; x++)
2685     {
2686       double
2687         intensity;
2688
2689       intensity=GetPixelIntensity(image,p);
2690       histogram[ScaleQuantumToMap(ClampToQuantum(intensity))]++;
2691       p+=GetPixelChannels(image);
2692     }
2693   }
2694   image_view=DestroyCacheView(image_view);
2695   /*
2696     Find the histogram boundaries by locating the black and white point levels.
2697   */
2698   intensity=0.0;
2699   for (black=0; black < (ssize_t) MaxMap; black++)
2700   {
2701     intensity+=histogram[black];
2702     if (intensity >= black_point)
2703       break;
2704   }
2705   intensity=0.0;
2706   for (white=(ssize_t) MaxMap; white != 0; white--)
2707   {
2708     intensity+=histogram[white];
2709     if (intensity >= white_point)
2710       break;
2711   }
2712   histogram=(double *) RelinquishMagickMemory(histogram);
2713   status=LevelImage(image,(double) black,(double) white,1.0,exception);
2714   return(status);
2715 }
2716 \f
2717 /*
2718 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2719 %                                                                             %
2720 %                                                                             %
2721 %                                                                             %
2722 %     M o d u l a t e I m a g e                                               %
2723 %                                                                             %
2724 %                                                                             %
2725 %                                                                             %
2726 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
2727 %
2728 %  ModulateImage() lets you control the brightness, saturation, and hue
2729 %  of an image.  Modulate represents the brightness, saturation, and hue
2730 %  as one parameter (e.g. 90,150,100).  If the image colorspace is HSL, the
2731 %  modulation is lightness, saturation, and hue.  For HWB, use blackness,
2732 %  whiteness, and hue. And for HCL, use chrome, luma, and hue.
2733 %
2734 %  The format of the ModulateImage method is:
2735 %
2736 %      MagickBooleanType ModulateImage(Image *image,const char *modulate,
2737 %        ExceptionInfo *exception)
2738 %
2739 %  A description of each parameter follows:
2740 %
2741 %    o image: the image.
2742 %
2743 %    o modulate: Define the percent change in brightness, saturation, and hue.
2744 %
2745 %    o exception: return any errors or warnings in this structure.
2746 %
2747 */
2748
2749 static inline void ModulateHCL(const double percent_hue,
2750   const double percent_chroma,const double percent_luma,double *red,
2751   double *green,double *blue)
2752 {
2753   double
2754     hue,
2755     luma,
2756     chroma;
2757
2758   /*
2759     Increase or decrease color luma, chroma, or hue.
2760   */
2761   ConvertRGBToHCL(*red,*green,*blue,&hue,&chroma,&luma);
2762   hue+=0.5*(0.01*percent_hue-1.0);
2763   while (hue < 0.0)
2764     hue+=1.0;
2765   while (hue > 1.0)
2766     hue-=1.0;
2767   chroma*=0.01*percent_chroma;
2768   luma*=0.01*percent_luma;
2769   ConvertHCLToRGB(hue,chroma,luma,red,green,blue);
2770 }
2771
2772 static inline void ModulateHSB(const double percent_hue,
2773   const double percent_saturation,const double percent_brightness,double *red,
2774   double *green,double *blue)
2775 {
2776   double
2777     brightness,
2778     hue,
2779     saturation;
2780
2781   /*
2782     Increase or decrease color brightness, saturation, or hue.
2783   */
2784   ConvertRGBToHSB(*red,*green,*blue,&hue,&saturation,&brightness);
2785   hue+=0.5*(0.01*percent_hue-1.0);
2786   while (hue < 0.0)
2787     hue+=1.0;
2788   while (hue > 1.0)
2789     hue-=1.0;
2790   saturation*=0.01*percent_saturation;
2791   brightness*=0.01*percent_brightness;
2792   ConvertHSBToRGB(hue,saturation,brightness,red,green,blue);
2793 }
2794
2795 static inline void ModulateHSL(const double percent_hue,
2796   const double percent_saturation,const double percent_lightness,double *red,
2797   double *green,double *blue)
2798 {
2799   double
2800     hue,
2801     lightness,
2802     saturation;
2803
2804   /*
2805     Increase or decrease color lightness, saturation, or hue.
2806   */
2807   ConvertRGBToHSL(*red,*green,*blue,&hue,&saturation,&lightness);
2808   hue+=0.5*(0.01*percent_hue-1.0);
2809   while (hue < 0.0)
2810     hue+=1.0;
2811   while (hue > 1.0)
2812     hue-=1.0;
2813   saturation*=0.01*percent_saturation;
2814   lightness*=0.01*percent_lightness;
2815   ConvertHSLToRGB(hue,saturation,lightness,red,green,blue);
2816 }
2817
2818 static inline void ModulateHWB(const double percent_hue,
2819   const double percent_whiteness,const double percent_blackness,double *red,
2820   double *green,double *blue)
2821 {
2822   double
2823     blackness,
2824     hue,
2825     whiteness;
2826
2827   /*
2828     Increase or decrease color blackness, whiteness, or hue.
2829   */
2830   ConvertRGBToHWB(*red,*green,*blue,&hue,&whiteness,&blackness);
2831   hue+=0.5*(0.01*percent_hue-1.0);
2832   while (hue < 0.0)
2833     hue+=1.0;
2834   while (hue > 1.0)
2835     hue-=1.0;
2836   blackness*=0.01*percent_blackness;
2837   whiteness*=0.01*percent_whiteness;
2838   ConvertHWBToRGB(hue,whiteness,blackness,red,green,blue);
2839 }
2840
2841 MagickExport MagickBooleanType ModulateImage(Image *image,const char *modulate,
2842   ExceptionInfo *exception)
2843 {
2844 #define ModulateImageTag  "Modulate/Image"
2845
2846   CacheView
2847     *image_view;
2848
2849   ColorspaceType
2850     colorspace;
2851
2852   const char
2853     *artifact;
2854
2855   double
2856     percent_brightness,
2857     percent_hue,
2858     percent_saturation;
2859
2860   GeometryInfo
2861     geometry_info;
2862
2863   MagickBooleanType
2864     status;
2865
2866   MagickOffsetType
2867     progress;
2868
2869   MagickStatusType
2870     flags;
2871
2872   register ssize_t
2873     i;
2874
2875   ssize_t
2876     y;
2877
2878   /*
2879     Initialize modulate table.
2880   */
2881   assert(image != (Image *) NULL);
2882   assert(image->signature == MagickSignature);
2883   if (image->debug != MagickFalse)
2884     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
2885   if (modulate == (char *) NULL)
2886     return(MagickFalse);
2887   if (IssRGBCompatibleColorspace(image->colorspace) == MagickFalse)
2888     (void) TransformImageColorspace(image,sRGBColorspace,exception);
2889   flags=ParseGeometry(modulate,&geometry_info);
2890   percent_brightness=geometry_info.rho;
2891   percent_saturation=geometry_info.sigma;
2892   if ((flags & SigmaValue) == 0)
2893     percent_saturation=100.0;
2894   percent_hue=geometry_info.xi;
2895   if ((flags & XiValue) == 0)
2896     percent_hue=100.0;
2897   colorspace=UndefinedColorspace;
2898   artifact=GetImageArtifact(image,"modulate:colorspace");
2899   if (artifact != (const char *) NULL)
2900     colorspace=(ColorspaceType) ParseCommandOption(MagickColorspaceOptions,
2901       MagickFalse,artifact);
2902   if (image->storage_class == PseudoClass)
2903     for (i=0; i < (ssize_t) image->colors; i++)
2904     {
2905       double
2906         blue,
2907         green,
2908         red;
2909
2910       /*
2911         Modulate colormap.
2912       */
2913       red=image->colormap[i].red;
2914       green=image->colormap[i].green;
2915       blue=image->colormap[i].blue;
2916       switch (colorspace)
2917       {
2918         case HCLColorspace:
2919         {
2920           ModulateHCL(percent_hue,percent_saturation,percent_brightness,
2921             &red,&green,&blue);
2922           break;
2923         }
2924         case HSBColorspace:
2925         {
2926           ModulateHSB(percent_hue,percent_saturation,percent_brightness,
2927             &red,&green,&blue);
2928           break;
2929         }
2930         case HSLColorspace:
2931         default:
2932         {
2933           ModulateHSL(percent_hue,percent_saturation,percent_brightness,
2934             &red,&green,&blue);
2935           break;
2936         }
2937         case HWBColorspace:
2938         {
2939           ModulateHWB(percent_hue,percent_saturation,percent_brightness,
2940             &red,&green,&blue);
2941           break;
2942         }
2943       }
2944     }
2945   /*
2946     Modulate image.
2947   */
2948   status=MagickTrue;
2949   progress=0;
2950   image_view=AcquireAuthenticCacheView(image,exception);
2951 #if defined(MAGICKCORE_OPENMP_SUPPORT)
2952   #pragma omp parallel for schedule(static,4) shared(progress,status) \
2953     dynamic_number_threads(image,image->columns,image->rows,1)
2954 #endif
2955   for (y=0; y < (ssize_t) image->rows; y++)
2956   {
2957     register Quantum
2958       *restrict q;
2959
2960     register ssize_t
2961       x;
2962
2963     if (status == MagickFalse)
2964       continue;
2965     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
2966     if (q == (Quantum *) NULL)
2967       {
2968         status=MagickFalse;
2969         continue;
2970       }
2971     for (x=0; x < (ssize_t) image->columns; x++)
2972     {
2973       double
2974         blue,
2975         green,
2976         red;
2977
2978       red=(double) GetPixelRed(image,q);
2979       green=(double) GetPixelGreen(image,q);
2980       blue=(double) GetPixelBlue(image,q);
2981       switch (colorspace)
2982       {
2983         case HCLColorspace:
2984         {
2985           ModulateHCL(percent_hue,percent_saturation,percent_brightness,
2986             &red,&green,&blue);
2987           break;
2988         }
2989         case HSBColorspace:
2990         {
2991           ModulateHSB(percent_hue,percent_saturation,percent_brightness,
2992             &red,&green,&blue);
2993           break;
2994         }
2995         case HSLColorspace:
2996         default:
2997         {
2998           ModulateHSL(percent_hue,percent_saturation,percent_brightness,
2999             &red,&green,&blue);
3000           break;
3001         }
3002         case HWBColorspace:
3003         {
3004           ModulateHWB(percent_hue,percent_saturation,percent_brightness,
3005             &red,&green,&blue);
3006           break;
3007         }
3008       }
3009       SetPixelRed(image,ClampToQuantum(red),q);
3010       SetPixelGreen(image,ClampToQuantum(green),q);
3011       SetPixelBlue(image,ClampToQuantum(blue),q);
3012       q+=GetPixelChannels(image);
3013     }
3014     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
3015       status=MagickFalse;
3016     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3017       {
3018         MagickBooleanType
3019           proceed;
3020
3021 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3022         #pragma omp critical (MagickCore_ModulateImage)
3023 #endif
3024         proceed=SetImageProgress(image,ModulateImageTag,progress++,image->rows);
3025         if (proceed == MagickFalse)
3026           status=MagickFalse;
3027       }
3028   }
3029   image_view=DestroyCacheView(image_view);
3030   return(status);
3031 }
3032 \f
3033 /*
3034 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3035 %                                                                             %
3036 %                                                                             %
3037 %                                                                             %
3038 %     N e g a t e I m a g e                                                   %
3039 %                                                                             %
3040 %                                                                             %
3041 %                                                                             %
3042 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3043 %
3044 %  NegateImage() negates the colors in the reference image.  The grayscale
3045 %  option means that only grayscale values within the image are negated.
3046 %
3047 %  The format of the NegateImage method is:
3048 %
3049 %      MagickBooleanType NegateImage(Image *image,
3050 %        const MagickBooleanType grayscale,ExceptionInfo *exception)
3051 %
3052 %  A description of each parameter follows:
3053 %
3054 %    o image: the image.
3055 %
3056 %    o grayscale: If MagickTrue, only negate grayscale pixels within the image.
3057 %
3058 %    o exception: return any errors or warnings in this structure.
3059 %
3060 */
3061 MagickExport MagickBooleanType NegateImage(Image *image,
3062   const MagickBooleanType grayscale,ExceptionInfo *exception)
3063 {
3064 #define NegateImageTag  "Negate/Image"
3065
3066   CacheView
3067     *image_view;
3068
3069   MagickBooleanType
3070     status;
3071
3072   MagickOffsetType
3073     progress;
3074
3075   register ssize_t
3076     i;
3077
3078   ssize_t
3079     y;
3080
3081   assert(image != (Image *) NULL);
3082   assert(image->signature == MagickSignature);
3083   if (image->debug != MagickFalse)
3084     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3085   if (image->storage_class == PseudoClass)
3086     for (i=0; i < (ssize_t) image->colors; i++)
3087     {
3088       /*
3089         Negate colormap.
3090       */
3091       if (grayscale != MagickFalse)
3092         if ((image->colormap[i].red != image->colormap[i].green) ||
3093           (image->colormap[i].green != image->colormap[i].blue))
3094           continue;
3095       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
3096         image->colormap[i].red=QuantumRange-image->colormap[i].red;
3097       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
3098         image->colormap[i].green=QuantumRange-image->colormap[i].green;
3099       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
3100         image->colormap[i].blue=QuantumRange-image->colormap[i].blue;
3101     }
3102   /*
3103     Negate image.
3104   */
3105   status=MagickTrue;
3106   progress=0;
3107   image_view=AcquireAuthenticCacheView(image,exception);
3108   if (grayscale != MagickFalse)
3109     {
3110 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3111       #pragma omp parallel for schedule(static) shared(progress,status) \
3112         dynamic_number_threads(image,image->columns,image->rows,1)
3113 #endif
3114       for (y=0; y < (ssize_t) image->rows; y++)
3115       {
3116         MagickBooleanType
3117           sync;
3118
3119         register Quantum
3120           *restrict q;
3121
3122         register ssize_t
3123           x;
3124
3125         if (status == MagickFalse)
3126           continue;
3127         q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,
3128           exception);
3129         if (q == (Quantum *) NULL)
3130           {
3131             status=MagickFalse;
3132             continue;
3133           }
3134         for (x=0; x < (ssize_t) image->columns; x++)
3135         {
3136           register ssize_t
3137             i;
3138
3139           if ((GetPixelMask(image,q) != 0) ||
3140               (IsPixelGray(image,q) != MagickFalse))
3141             {
3142               q+=GetPixelChannels(image);
3143               continue;
3144             }
3145           for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
3146           {
3147             PixelChannel
3148               channel;
3149
3150             PixelTrait
3151               traits;
3152
3153             channel=GetPixelChannelChannel(image,i);
3154             traits=GetPixelChannelTraits(image,channel);
3155             if ((traits & UpdatePixelTrait) == 0)
3156               continue;
3157             q[i]=QuantumRange-q[i];
3158           }
3159           q+=GetPixelChannels(image);
3160         }
3161         sync=SyncCacheViewAuthenticPixels(image_view,exception);
3162         if (sync == MagickFalse)
3163           status=MagickFalse;
3164         if (image->progress_monitor != (MagickProgressMonitor) NULL)
3165           {
3166             MagickBooleanType
3167               proceed;
3168
3169 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3170             #pragma omp critical (MagickCore_NegateImage)
3171 #endif
3172             proceed=SetImageProgress(image,NegateImageTag,progress++,
3173               image->rows);
3174             if (proceed == MagickFalse)
3175               status=MagickFalse;
3176           }
3177       }
3178       image_view=DestroyCacheView(image_view);
3179       return(MagickTrue);
3180     }
3181   /*
3182     Negate image.
3183   */
3184 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3185   #pragma omp parallel for schedule(static) shared(progress,status) \
3186     dynamic_number_threads(image,image->columns,image->rows,1)
3187 #endif
3188   for (y=0; y < (ssize_t) image->rows; y++)
3189   {
3190     register Quantum
3191       *restrict q;
3192
3193     register ssize_t
3194       x;
3195
3196     if (status == MagickFalse)
3197       continue;
3198     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
3199     if (q == (Quantum *) NULL)
3200       {
3201         status=MagickFalse;
3202         continue;
3203       }
3204     for (x=0; x < (ssize_t) image->columns; x++)
3205     {
3206       register ssize_t
3207         i;
3208
3209       if (GetPixelMask(image,q) != 0)
3210         {
3211           q+=GetPixelChannels(image);
3212           continue;
3213         }
3214       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
3215       {
3216         PixelChannel
3217           channel;
3218
3219         PixelTrait
3220           traits;
3221
3222         channel=GetPixelChannelChannel(image,i);
3223         traits=GetPixelChannelTraits(image,channel);
3224         if ((traits & UpdatePixelTrait) == 0)
3225           continue;
3226         q[i]=QuantumRange-q[i];
3227       }
3228       q+=GetPixelChannels(image);
3229     }
3230     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
3231       status=MagickFalse;
3232     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3233       {
3234         MagickBooleanType
3235           proceed;
3236
3237 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3238         #pragma omp critical (MagickCore_NegateImage)
3239 #endif
3240         proceed=SetImageProgress(image,NegateImageTag,progress++,image->rows);
3241         if (proceed == MagickFalse)
3242           status=MagickFalse;
3243       }
3244   }
3245   image_view=DestroyCacheView(image_view);
3246   return(status);
3247 }
3248 \f
3249 /*
3250 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3251 %                                                                             %
3252 %                                                                             %
3253 %                                                                             %
3254 %     N o r m a l i z e I m a g e                                             %
3255 %                                                                             %
3256 %                                                                             %
3257 %                                                                             %
3258 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3259 %
3260 %  The NormalizeImage() method enhances the contrast of a color image by
3261 %  mapping the darkest 2 percent of all pixel to black and the brightest
3262 %  1 percent to white.
3263 %
3264 %  The format of the NormalizeImage method is:
3265 %
3266 %      MagickBooleanType NormalizeImage(Image *image,ExceptionInfo *exception)
3267 %
3268 %  A description of each parameter follows:
3269 %
3270 %    o image: the image.
3271 %
3272 %    o exception: return any errors or warnings in this structure.
3273 %
3274 */
3275 MagickExport MagickBooleanType NormalizeImage(Image *image,
3276   ExceptionInfo *exception)
3277 {
3278   double
3279     black_point,
3280     white_point;
3281
3282   black_point=(double) image->columns*image->rows*0.0015;
3283   white_point=(double) image->columns*image->rows*0.9995;
3284   return(ContrastStretchImage(image,black_point,white_point,exception));
3285 }
3286 \f
3287 /*
3288 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3289 %                                                                             %
3290 %                                                                             %
3291 %                                                                             %
3292 %     S i g m o i d a l C o n t r a s t I m a g e                             %
3293 %                                                                             %
3294 %                                                                             %
3295 %                                                                             %
3296 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3297 %
3298 %  SigmoidalContrastImage() adjusts the contrast of an image with a non-linear
3299 %  sigmoidal contrast algorithm.  Increase the contrast of the image using a
3300 %  sigmoidal transfer function without saturating highlights or shadows.
3301 %  Contrast indicates how much to increase the contrast (0 is none; 3 is
3302 %  typical; 20 is pushing it); mid-point indicates where midtones fall in the
3303 %  resultant image (0 is white; 50% is middle-gray; 100% is black).  Set
3304 %  sharpen to MagickTrue to increase the image contrast otherwise the contrast
3305 %  is reduced.
3306 %
3307 %  The format of the SigmoidalContrastImage method is:
3308 %
3309 %      MagickBooleanType SigmoidalContrastImage(Image *image,
3310 %        const MagickBooleanType sharpen,const char *levels,
3311 %        ExceptionInfo *exception)
3312 %
3313 %  A description of each parameter follows:
3314 %
3315 %    o image: the image.
3316 %
3317 %    o sharpen: Increase or decrease image contrast.
3318 %
3319 %    o alpha: strength of the contrast, the larger the number the more
3320 %      'threshold-like' it becomes.
3321 %
3322 %    o beta: midpoint of the function as a color value 0 to QuantumRange.
3323 %
3324 %    o exception: return any errors or warnings in this structure.
3325 %
3326 */
3327 MagickExport MagickBooleanType SigmoidalContrastImage(Image *image,
3328   const MagickBooleanType sharpen,const double contrast,const double midpoint,
3329   ExceptionInfo *exception)
3330 {
3331 #define SigmoidalContrastImageTag  "SigmoidalContrast/Image"
3332
3333   CacheView
3334     *image_view;
3335
3336   MagickBooleanType
3337     status;
3338
3339   MagickOffsetType
3340     progress;
3341
3342   Quantum
3343     *sigmoidal_map;
3344
3345   register ssize_t
3346     i;
3347
3348   ssize_t
3349     y;
3350
3351   /*
3352     Sigmoidal function Sig with inflexion point moved to b and "slope
3353     constant" set to a.
3354     The first version, based on the hyperbolic tangent tanh, when
3355     combined with the scaling step, is an exact arithmetic clone of the
3356     the sigmoid function based on the logistic curve. The equivalence is
3357     based on the identity
3358     1/(1+exp(-t)) = (1+tanh(t/2))/2
3359     (http://de.wikipedia.org/wiki/Sigmoidfunktion) and the fact that the
3360     scaled sigmoidal derivation is invariant under affine transformations
3361     of the ordinate.
3362     The tanh version is almost certainly more accurate and cheaper.
3363     The 0.5 factor in the argument is to clone the legacy ImageMagick
3364     behavior.
3365     The reason for making the define depend on atanh even though it only
3366     uses tanh has to do with the construction of the inverse of the scaled
3367     sigmoidal.
3368   */
3369 #if defined(MAGICKCORE_HAVE_ATANH)
3370 #define Sig(a,b,x) ( tanh((0.5*(a))*((x)-(b))) )
3371 #else
3372 #define Sig(a,b,x) ( 1.0/(1.0+exp((a)*((b)-(x)))) )
3373 #endif
3374   /*
3375     Scaled sigmoidal formula:
3376     ( Sig(a,b,x)-Sig(a,b,0) ) / ( Sig(a,b,1) - Sig(a,b,0) )
3377     See http://osdir.com/ml/video.image-magick.devel/2005-04/msg00006.html
3378     and http://www.cs.dartmouth.edu/farid/downloads/tutorials/fip.pdf.
3379     The limit of ScaledSig as a->0 is the identity, but a=0 gives a
3380     division by zero. This is fixed below by hardwiring the identity when a
3381     is small. This would appear to be safe because the series expansion of
3382     the sigmoidal function around x=b is 1/2-a*(b-x)/4+... so that s(1)-s(0)
3383     is about a/4.
3384   */
3385 #define ScaledSig(a,b,x) ( \
3386   (Sig((a),(b),(x))-Sig((a),(b),0.0)) / (Sig((a),(b),1.0)-Sig((a),(b),0.0)) )
3387   /*
3388     Inverse of ScaledSig, used for +sigmoidal-contrast:
3389   */
3390 #if defined(MAGICKCORE_HAVE_ATANH)
3391 #define InverseScaledSig(a,b,x) ( (b) +                             \
3392   atanh( (Sig((a),(b),1.0)-Sig((a),(b),0.0))*(x)+Sig((a),(b),0.0) ) \
3393   / (0.5*(a)) )
3394 #else
3395 #define InverseScaledSig(a,b,x) ( (b) -                                     \
3396   log( 1.0/((Sig((a),(b),1.0)-Sig((a),(b),0.0))*(x)+Sig((a),(b),0.0))-1.0 ) \
3397   / (a) )
3398 #endif
3399
3400   /*
3401     Allocate and initialize sigmoidal maps.
3402   */
3403   assert(image != (Image *) NULL);
3404   assert(image->signature == MagickSignature);
3405   if (image->debug != MagickFalse)
3406     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
3407   sigmoidal_map=(Quantum *) AcquireQuantumMemory(MaxMap+1UL,
3408     sizeof(*sigmoidal_map));
3409   if (sigmoidal_map == (Quantum *) NULL)
3410     ThrowBinaryException(ResourceLimitError,"MemoryAllocationFailed",
3411       image->filename);
3412   (void) ResetMagickMemory(sigmoidal_map,0,(MaxMap+1)*sizeof(*sigmoidal_map));
3413   if (contrast<4.0*MagickEpsilon)
3414     for (i=0; i <= (ssize_t) MaxMap; i++)
3415       sigmoidal_map[i]=ScaleMapToQuantum((double) i);
3416   else if (sharpen != MagickFalse)
3417     for (i=0; i <= (ssize_t) MaxMap; i++)
3418       sigmoidal_map[i]=ScaleMapToQuantum( (double) (MaxMap*
3419         ScaledSig(contrast,QuantumScale*midpoint,(double) i/MaxMap)));
3420   else
3421     for (i=0; i <= (ssize_t) MaxMap; i++)
3422       sigmoidal_map[i]=ScaleMapToQuantum((double) (MaxMap*
3423         InverseScaledSig(contrast,QuantumScale*midpoint,(double) i/MaxMap)));
3424   if (image->storage_class == PseudoClass)
3425     for (i=0; i < (ssize_t) image->colors; i++)
3426     {
3427       /*
3428         Sigmoidal-contrast enhance colormap.
3429       */
3430       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
3431         image->colormap[i].red=(double) ClampToQuantum((double) sigmoidal_map[
3432           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))]);
3433       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
3434         image->colormap[i].green=(double) ClampToQuantum((double) sigmoidal_map[
3435           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].green))]);
3436       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
3437         image->colormap[i].blue=(double) ClampToQuantum((double) sigmoidal_map[
3438           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].blue))]);
3439       if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
3440         image->colormap[i].alpha=(double) ClampToQuantum((double) sigmoidal_map[
3441           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].alpha))]);
3442     }
3443   /*
3444     Sigmoidal-contrast enhance image.
3445   */
3446   status=MagickTrue;
3447   progress=0;
3448   image_view=AcquireAuthenticCacheView(image,exception);
3449 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3450   #pragma omp parallel for schedule(static,4) shared(progress,status) \
3451     dynamic_number_threads(image,image->columns,image->rows,1)
3452 #endif
3453   for (y=0; y < (ssize_t) image->rows; y++)
3454   {
3455     register Quantum
3456       *restrict q;
3457
3458     register ssize_t
3459       x;
3460
3461     if (status == MagickFalse)
3462       continue;
3463     q=GetCacheViewAuthenticPixels(image_view,0,y,image->columns,1,exception);
3464     if (q == (Quantum *) NULL)
3465       {
3466         status=MagickFalse;
3467         continue;
3468       }
3469     for (x=0; x < (ssize_t) image->columns; x++)
3470     {
3471       register ssize_t
3472         i;
3473
3474       if (GetPixelMask(image,q) != 0)
3475         {
3476           q+=GetPixelChannels(image);
3477           continue;
3478         }
3479       for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
3480       {
3481         PixelChannel
3482           channel;
3483
3484         PixelTrait
3485           traits;
3486
3487         channel=GetPixelChannelChannel(image,i);
3488         traits=GetPixelChannelTraits(image,channel);
3489         if ((traits & UpdatePixelTrait) == 0)
3490           continue;
3491         q[i]=ClampToQuantum((double) sigmoidal_map[ScaleQuantumToMap(q[i])]);
3492       }
3493       q+=GetPixelChannels(image);
3494     }
3495     if (SyncCacheViewAuthenticPixels(image_view,exception) == MagickFalse)
3496       status=MagickFalse;
3497     if (image->progress_monitor != (MagickProgressMonitor) NULL)
3498       {
3499         MagickBooleanType
3500           proceed;
3501
3502 #if defined(MAGICKCORE_OPENMP_SUPPORT)
3503         #pragma omp critical (MagickCore_SigmoidalContrastImage)
3504 #endif
3505         proceed=SetImageProgress(image,SigmoidalContrastImageTag,progress++,
3506           image->rows);
3507         if (proceed == MagickFalse)
3508           status=MagickFalse;
3509       }
3510   }
3511   image_view=DestroyCacheView(image_view);
3512   sigmoidal_map=(Quantum *) RelinquishMagickMemory(sigmoidal_map);
3513   return(status);
3514 }