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