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