]> granicus.if.org Git - imagemagick/blobdiff - MagickCore/enhance.c
(no commit message)
[imagemagick] / MagickCore / enhance.c
index 0311039de668196d265ac33a7f5a23eac83d9065..44477517b2ba2477b79d1419fb0cf05b5009221a 100644 (file)
@@ -718,15 +718,15 @@ MagickExport MagickBooleanType ColorDecisionListImage(Image *image,
     for (i=0; i < (ssize_t) image->colors; i++)
     {
       /*
-       Apply transfer function to colormap.
+        Apply transfer function to colormap.
       */
       double
-       luma;
+        luma;
 
       luma=0.21267*image->colormap[i].red+0.71526*image->colormap[i].green+
-       0.07217*image->colormap[i].blue;
+        0.07217*image->colormap[i].blue;
       image->colormap[i].red=luma+color_correction.saturation*cdl_map[
-       ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))].red-luma;
+        ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))].red-luma;
       image->colormap[i].green=luma+color_correction.saturation*cdl_map[
         ScaleQuantumToMap(ClampToQuantum(image->colormap[i].green))].green-luma;
       image->colormap[i].blue=luma+color_correction.saturation*cdl_map[
@@ -1806,10 +1806,10 @@ MagickExport MagickBooleanType GammaImage(Image *image,const double gamma,
     for (i=0; i < (ssize_t) image->colors; i++)
     {
       /*
-       Gamma-correct colormap.
+        Gamma-correct colormap.
       */
       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
-       image->colormap[i].red=(double) gamma_map[
+        image->colormap[i].red=(double) gamma_map[
           ScaleQuantumToMap(ClampToQuantum(image->colormap[i].red))];
       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
         image->colormap[i].green=(double) gamma_map[
@@ -3058,18 +3058,18 @@ MagickExport MagickBooleanType NegateImage(Image *image,
     for (i=0; i < (ssize_t) image->colors; i++)
     {
       /*
-       Negate colormap.
+        Negate colormap.
       */
       if (grayscale != MagickFalse)
-       if ((image->colormap[i].red != image->colormap[i].green) ||
-         (image->colormap[i].green != image->colormap[i].blue))
-         continue;
+        if ((image->colormap[i].red != image->colormap[i].green) ||
+            (image->colormap[i].green != image->colormap[i].blue))
+          continue;
       if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
-       image->colormap[i].red=QuantumRange-image->colormap[i].red;
+        image->colormap[i].red=QuantumRange-image->colormap[i].red;
       if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
-       image->colormap[i].green=QuantumRange-image->colormap[i].green;
+        image->colormap[i].green=QuantumRange-image->colormap[i].green;
       if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
-       image->colormap[i].blue=QuantumRange-image->colormap[i].blue;
+        image->colormap[i].blue=QuantumRange-image->colormap[i].blue;
     }
   /*
     Negate image.
@@ -3292,11 +3292,93 @@ MagickExport MagickBooleanType NormalizeImage(Image *image,
 %    o exception: return any errors or warnings in this structure.
 %
 */
+
+/*
+  Sigmoidal function with inflexion point moved to b and "slope constant" set
+  to a.
+  The first version, based on the hyperbolic tangent tanh, when combined with
+  the scaling step, is an exact arithmetic clone of the the sigmoid function
+  based on the logistic curve. The equivalence is based on the identity
+
+  1/(1+exp(-t)) = (1+tanh(t/2))/2
+
+  (http://de.wikipedia.org/wiki/Sigmoidfunktion) and the fact that the scaled
+  sigmoidal derivation is invariant under affine transformations of the
+  ordinate.
+  The tanh version is almost certainly more accurate and cheaper.
+  The 0.5 factor in its argument is to clone the legacy ImageMagick behavior.
+  The reason for making the define depend on atanh even though it only uses
+  tanh has to do with the construction of the inverse of the scaled sigmoidal.
+*/
+#if defined(MAGICKCORE_HAVE_ATANH)
+#define Sigmoidal(a,b,x) ( tanh((0.5*(a))*((x)-(b))) )
+#else
+#define Sigmoidal(a,b,x) ( 1.0/(1.0+exp((a)*((b)-(x)))) )
+#endif
+/*
+  Scaled sigmoidal function:
+
+    ( Sigmoidal(a,b,x) - Sigmoidal(a,b,0) ) /
+    ( Sigmoidal(a,b,1) - Sigmoidal(a,b,0) )
+
+  See http://osdir.com/ml/video.image-magick.devel/2005-04/msg00006.html and
+  http://www.cs.dartmouth.edu/farid/downloads/tutorials/fip.pdf.  The limit
+  of ScaledSigmoidal as a->0 is the identity, but a=0 gives a division by
+  zero. This is fixed above by exiting immediately when contrast is small,
+  leaving the image (or colormap) unmodified. This appears to be safe because
+  the series expansion of the logistic sigmoidal function around x=b is
+
+  1/2-a*(b-x)/4+...
+
+  so that the key denominator s(1)-s(0) is about a/4 (a/2 with tanh).
+*/
+#define ScaledSigmoidal(a,b,x) (                    \
+  (Sigmoidal((a),(b),(x))-Sigmoidal((a),(b),0.0)) / \
+  (Sigmoidal((a),(b),1.0)-Sigmoidal((a),(b),0.0)) )
+/*
+  Inverse of ScaledSigmoidal, used for +sigmoidal-contrast.  Because b
+  may be 0 or 1, the argument of the hyperbolic tangent (resp. logistic
+  sigmoidal) may be outside of the interval (-1,1) (resp. (0,1)), even
+  when creating a LUT from in gamut values, hence the branching.  In
+  addition, HDRI may have out of gamut values.
+  InverseScaledSigmoidal is not a two-side inverse of ScaledSigmoidal:
+  It is only a right inverse. This is unavoidable.
+*/
+static inline double InverseScaledSigmoidal(const double a,const double b,
+  const double x)
+{
+  const double sig0=Sigmoidal(a,b,0.0);
+  const double argument=(Sigmoidal(a,b,1.0)-sig0)*x+sig0;
+  const double clamped=
+    (
+#if defined(MAGICKCORE_HAVE_ATANH)
+      argument < -1+MagickEpsilon
+      ?
+      -1+MagickEpsilon
+      :
+      ( argument > 1-MagickEpsilon ? 1-MagickEpsilon : argument )
+    );
+  return(b+(2.0/a)*atanh(clamped));
+#else
+      argument < MagickEpsilon
+      ?
+      MagickEpsilon
+      :
+      ( argument > 1-MagickEpsilon ? 1-MagickEpsilon : argument )
+    );
+  return(b+(-1.0/a)*log(1.0/clamped+-1.0));
+#endif
+}
+
 MagickExport MagickBooleanType SigmoidalContrastImage(Image *image,
   const MagickBooleanType sharpen,const double contrast,const double midpoint,
   ExceptionInfo *exception)
 {
 #define SigmoidalContrastImageTag  "SigmoidalContrast/Image"
+#define ScaledSig(x) ( ClampToQuantum(QuantumRange* \
+  ScaledSigmoidal(contrast,QuantumScale*midpoint,QuantumScale*(x))) )
+#define InverseScaledSig(x) ( ClampToQuantum(QuantumRange* \
+  InverseScaledSigmoidal(contrast,QuantumScale*midpoint,QuantumScale*(x))) )
 
   CacheView
     *image_view;
@@ -3311,115 +3393,53 @@ MagickExport MagickBooleanType SigmoidalContrastImage(Image *image,
     y;
 
   /*
-    Side effect: clamps values unless contrast<MagickEpsilon, in which
-    case nothing is done.
-  */
-  if (contrast<MagickEpsilon)
-    return(MagickTrue);
-
-  /*
-    Sigmoidal function with inflexion point moved to b and "slope
-    constant" set to a.
-    The first version, based on the hyperbolic tangent tanh, when
-    combined with the scaling step, is an exact arithmetic clone of the
-    the sigmoid function based on the logistic curve. The equivalence is
-    based on the identity
-    1/(1+exp(-t)) = (1+tanh(t/2))/2
-    (http://de.wikipedia.org/wiki/Sigmoidfunktion) and the fact that the
-    scaled sigmoidal derivation is invariant under affine transformations
-    of the ordinate.
-    The tanh version is almost certainly more accurate and cheaper.
-    The 0.5 factor in its argument is to clone the legacy ImageMagick
-    behavior.  The reason for making the define depend on atanh even
-    though it only uses tanh has to do with the construction of the
-    inverse of the scaled sigmoidal.
-  */
-#if defined(MAGICKCORE_HAVE_ATANH)
-#define Sigmoidal(a,b,x) ( tanh((0.5*(a))*((x)-(b))) )
-#else
-#define Sigmoidal(a,b,x) ( 1.0/(1.0+exp((a)*((b)-(x)))) )
-#endif
-  /*
-    Scaled sigmoidal formula:
-    ( Sigmoidal(a,b,x) - Sigmoidal(a,b,0) ) /
-    ( Sigmoidal(a,b,1) - Sigmoidal(a,b,0) )
-    See http://osdir.com/ml/video.image-magick.devel/2005-04/msg00006.html
-    and http://www.cs.dartmouth.edu/farid/downloads/tutorials/fip.pdf.
-    The limit of ScaledSigmoidal as a->0 is the identity, but a=0 gives a
-    division by zero. This is fixed above by exiting immediately when
-    contrast is small. This would appear to be safe because the series
-    expansion of the logistic sigmoidal function around x=b is
-    1/2-a*(b-x)/4+... so that s(1)-s(0) is about a/4. (With tanh, it's a/2.)
-  */
-#define ScaledSigmoidal(a,b,x) (                    \
-  (Sigmoidal((a),(b),(x))-Sigmoidal((a),(b),0.0)) / \
-  (Sigmoidal((a),(b),1.0)-Sigmoidal((a),(b),0.0)) )
-  /*
-    Inverse of ScaledSigmoidal, used for +sigmoidal-contrast:
-  */
-#if defined(MAGICKCORE_HAVE_ATANH)
-#define InverseScaledSigmoidal(a,b,x) ( (b) + (2.0/(a)) * atanh( \
-  (Sigmoidal((a),(b),1.0)-Sigmoidal((a),(b),0.0))*(x)+Sigmoidal((a),(b),0.0) ) )
-#else
-#define InverseScaledSigmoidal(a,b,x) ( (b) + (-1.0/(a)) * log( 1.0 /          \
-  ((Sigmoidal((a),(b),1.0)-Sigmoidal((a),(b),0.0))*(x)+Sigmoidal((a),(b),0.0)) \
-  + -1.0 ) )
-#endif
-  /*
-    Convenience macros. No clamping needed at the end because of monotonicity.
+    Convenience macros.
   */
-#define ScaledSig(x) ( (Quantum) (QuantumRange*ScaledSigmoidal(contrast, \
-  QuantumScale*midpoint,QuantumScale*ClampToQuantum((x)))) )  
-#define InverseScaledSig(x) ( (Quantum) (QuantumRange*InverseScaledSigmoidal( \
-  contrast,QuantumScale*midpoint,QuantumScale*ClampToQuantum((x)))) )  
-
   assert(image != (Image *) NULL);
   assert(image->signature == MagickSignature);
   if (image->debug != MagickFalse)
     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
   /*
-    Sigmoidal-contrast enhance colormap:
+    Side effect: may clamp values unless contrast<MagickEpsilon, in which
+    case nothing is done.
+  */
+  if (contrast < MagickEpsilon)
+    return(MagickTrue);
+  /*
+    Sigmoidal-contrast enhance colormap.
   */
   if (image->storage_class == PseudoClass)
     {
       register ssize_t
-       i;
+        i;
 
       if (sharpen != MagickFalse)
-       {
-         for (i=0; i < (ssize_t) image->colors; i++)
-         {
-           if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].red=ScaledSig(image->colormap[i].red);
-           if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].green=ScaledSig(image->colormap[i].green);
-           if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].blue=ScaledSig(image->colormap[i].blue);
-           if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].alpha=ScaledSig(image->colormap[i].alpha);
-         }
-       }
+        for (i=0; i < (ssize_t) image->colors; i++)
+        {
+          if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].red=ScaledSig(image->colormap[i].red);
+          if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].green=ScaledSig(image->colormap[i].green);
+          if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].blue=ScaledSig(image->colormap[i].blue);
+          if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].alpha=ScaledSig(image->colormap[i].alpha);
+        }
       else
-       {
-         for (i=0; i < (ssize_t) image->colors; i++)
-         {
-           if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].red=
-               InverseScaledSig(image->colormap[i].red);
-           if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].green=
-               InverseScaledSig(image->colormap[i].green);
-           if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].blue=
-               InverseScaledSig(image->colormap[i].blue);
-           if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
-             image->colormap[i].alpha=
-               InverseScaledSig(image->colormap[i].alpha);
-         }
-       }
+        for (i=0; i < (ssize_t) image->colors; i++)
+        {
+          if ((GetPixelRedTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].red=InverseScaledSig(image->colormap[i].red);
+          if ((GetPixelGreenTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].green=InverseScaledSig(image->colormap[i].green);
+          if ((GetPixelBlueTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].blue=InverseScaledSig(image->colormap[i].blue);
+          if ((GetPixelAlphaTraits(image) & UpdatePixelTrait) != 0)
+            image->colormap[i].alpha=InverseScaledSig(image->colormap[i].alpha);
+        }
     }
   /*
-    Sigmoidal-contrast enhance image:
+    Sigmoidal-contrast enhance image.
   */
   status=MagickTrue;
   progress=0;
@@ -3447,7 +3467,7 @@ MagickExport MagickBooleanType SigmoidalContrastImage(Image *image,
     for (x=0; x < (ssize_t) image->columns; x++)
     {
       register ssize_t
-       i;
+        i;
 
       if (GetPixelMask(image,q) != 0)
         {
@@ -3466,10 +3486,10 @@ MagickExport MagickBooleanType SigmoidalContrastImage(Image *image,
         traits=GetPixelChannelTraits(image,channel);
         if ((traits & UpdatePixelTrait) == 0)
           continue;
-       if (sharpen != MagickFalse)
-         q[i]=ScaledSig(q[i]);
-       else
-         q[i]=InverseScaledSig(q[i]);
+        if (sharpen != MagickFalse)
+          q[i]=ScaledSig(q[i]);
+        else
+          q[i]=InverseScaledSig(q[i]);
       }
       q+=GetPixelChannels(image);
     }