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