]> granicus.if.org Git - imagemagick/blobdiff - MagickCore/compare.c
(no commit message)
[imagemagick] / MagickCore / compare.c
index cce6749492a6a84d1977ef5ac3d5e5c5ba71a4a6..2ceffd8bf27a97698029f88d26dc1691a01a810d 100644 (file)
 %                      MagickCore Image Comparison Methods                    %
 %                                                                             %
 %                              Software Design                                %
-%                                   Cristy                                    %
+%                                  Cristy                                     %
 %                               December 2003                                 %
 %                                                                             %
 %                                                                             %
-%  Copyright 1999-2014 ImageMagick Studio LLC, a non-profit organization      %
+%  Copyright 1999-2015 ImageMagick Studio LLC, a non-profit organization      %
 %  dedicated to making software imaging solutions freely available.           %
 %                                                                             %
 %  You may not use this file except in compliance with the License.  You may  %
 %    o exception: return any errors or warnings in this structure.
 %
 */
+
+static size_t GetImageChannels(const Image *image)
+{
+  register ssize_t
+    i;
+
+  size_t
+    channels;
+
+  channels=0;
+  for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
+  {
+    PixelChannel channel=GetPixelChannelChannel(image,i);
+    PixelTrait traits=GetPixelChannelTraits(image,channel);
+    if ((traits & UpdatePixelTrait) != 0)
+      channels++;
+  }
+  return(channels == 0 ? 1 : channels);
+}
+
+static inline MagickBooleanType ValidateImageMorphology(
+  const Image *restrict image,const Image *restrict reconstruct_image)
+{
+  /*
+    Does the image match the reconstructed image morphology?
+  */
+  return(MagickTrue);
+}
+
 MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
   const MetricType metric,double *distortion,ExceptionInfo *exception)
 {
@@ -112,6 +141,9 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
     *image_view,
     *reconstruct_view;
 
+  double
+    fuzz;
+
   const char
     *artifact;
 
@@ -126,6 +158,13 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
     highlight,
     lowlight;
 
+  RectangleInfo
+    geometry;
+
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
@@ -140,19 +179,22 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
   if (image->debug != MagickFalse)
     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
   if (metric != PerceptualHashErrorMetric)
-    if ((reconstruct_image->columns != image->columns) ||
-        (reconstruct_image->rows != image->rows))
-      ThrowImageException(ImageError,"ImageSizeDiffers");
+    if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
+      ThrowImageException(ImageError,"ImageMorphologyDiffers");
   status=GetImageDistortion(image,reconstruct_image,metric,distortion,
     exception);
   if (status == MagickFalse)
     return((Image *) NULL);
-  difference_image=CloneImage(image,0,0,MagickTrue,exception);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  SetGeometry(image,&geometry);
+  geometry.width=columns;
+  geometry.height=rows;
+  difference_image=ExtentImage(image,&geometry,exception);
   if (difference_image == (Image *) NULL)
     return((Image *) NULL);
   (void) SetImageAlphaChannel(difference_image,OpaqueAlphaChannel,exception);
-  highlight_image=CloneImage(image,image->columns,image->rows,MagickTrue,
-    exception);
+  highlight_image=CloneImage(image,columns,rows,MagickTrue,exception);
   if (highlight_image == (Image *) NULL)
     {
       difference_image=DestroyImage(difference_image);
@@ -178,14 +220,15 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
     Generate difference image.
   */
   status=MagickTrue;
+  fuzz=GetFuzzyColorDistance(image,reconstruct_image);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
   highlight_view=AcquireAuthenticCacheView(highlight_image,exception);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp parallel for schedule(static,4) shared(status) \
-    magick_threads(image,highlight_image,image->rows,1)
+    magick_threads(image,highlight_image,rows,1)
 #endif
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     MagickBooleanType
       sync;
@@ -202,18 +245,16 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
-    r=QueueCacheViewAuthenticPixels(highlight_view,0,y,highlight_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
+    r=QueueCacheViewAuthenticPixels(highlight_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL) ||
         (r == (Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -227,7 +268,7 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
 
       if (GetPixelReadMask(image,p) == 0)
         {
-          SetPixelInfoPixel(highlight_image,&lowlight,r);
+          SetPixelViaPixelInfo(highlight_image,&lowlight,r);
           p+=GetPixelChannels(image);
           q+=GetPixelChannels(reconstruct_image);
           r+=GetPixelChannels(highlight_image);
@@ -250,13 +291,16 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
             ((reconstruct_traits & UpdatePixelTrait) == 0))
           continue;
         distance=Sa*p[i]-Da*GetPixelChannel(reconstruct_image,channel,q);
-        if (fabs((double) distance) >= MagickEpsilon)
-          difference=MagickTrue;
+        if ((distance*distance) > fuzz)
+          {
+            difference=MagickTrue;
+            break;
+          }
       }
       if (difference == MagickFalse)
-        SetPixelInfoPixel(highlight_image,&lowlight,r);
+        SetPixelViaPixelInfo(highlight_image,&lowlight,r);
       else
-        SetPixelInfoPixel(highlight_image,&highlight,r);
+        SetPixelViaPixelInfo(highlight_image,&highlight,r);
       p+=GetPixelChannels(image);
       q+=GetPixelChannels(reconstruct_image);
       r+=GetPixelChannels(highlight_image);
@@ -310,13 +354,6 @@ MagickExport Image *CompareImages(Image *image,const Image *reconstruct_image,
 %
 */
 
-static inline double MagickMax(const double x,const double y)
-{
-  if (x > y)
-    return(x);
-  return(y);
-}
-
 static MagickBooleanType GetAbsoluteDistortion(const Image *image,
   const Image *reconstruct_image,double *distortion,ExceptionInfo *exception)
 {
@@ -330,6 +367,10 @@ static MagickBooleanType GetAbsoluteDistortion(const Image *image,
   MagickBooleanType
     status;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
@@ -337,23 +378,16 @@ static MagickBooleanType GetAbsoluteDistortion(const Image *image,
     Compute the absolute difference in pixels between two images.
   */
   status=MagickTrue;
-  if (image->fuzz == 0.0)
-    fuzz=MagickMax(reconstruct_image->fuzz,MagickSQ1_2)*
-      MagickMax(reconstruct_image->fuzz,MagickSQ1_2);
-  else
-    if (reconstruct_image->fuzz == 0.0)
-      fuzz=MagickMax(image->fuzz,MagickSQ1_2)*
-        MagickMax(image->fuzz,MagickSQ1_2);
-    else
-      fuzz=MagickMax(image->fuzz,MagickSQ1_2)*
-        MagickMax(reconstruct_image->fuzz,MagickSQ1_2);
+  fuzz=GetFuzzyColorDistance(image,reconstruct_image);
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp parallel for schedule(static,4) shared(status) \
-    magick_threads(image,image,image->rows,1)
+    magick_threads(image,image,rows,1)
 #endif
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     double
       channel_distortion[MaxPixelChannels+1];
@@ -368,16 +402,15 @@ static MagickBooleanType GetAbsoluteDistortion(const Image *image,
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
     (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -414,15 +447,12 @@ static MagickBooleanType GetAbsoluteDistortion(const Image *image,
         distance=Sa*p[i]-Da*GetPixelChannel(reconstruct_image,channel,q);
         if ((distance*distance) > fuzz)
           {
+            channel_distortion[i]++;
             difference=MagickTrue;
-            break;
           }
       }
       if (difference != MagickFalse)
-        {
-          channel_distortion[i]++;
-          channel_distortion[CompositePixelChannel]++;
-        }
+        channel_distortion[CompositePixelChannel]++;
       p+=GetPixelChannels(image);
       q+=GetPixelChannels(reconstruct_image);
     }
@@ -437,25 +467,6 @@ static MagickBooleanType GetAbsoluteDistortion(const Image *image,
   return(status);
 }
 
-static size_t GetImageChannels(const Image *image)
-{
-  register ssize_t
-    i;
-
-  size_t
-    channels;
-
-  channels=0;
-  for (i=0; i < (ssize_t) GetPixelChannels(image); i++)
-  {
-    PixelChannel channel=GetPixelChannelChannel(image,i);
-    PixelTrait traits=GetPixelChannelTraits(image,channel);
-    if ((traits & UpdatePixelTrait) != 0)
-      channels++;
-  }
-  return(channels);
-}
-
 static MagickBooleanType GetFuzzDistortion(const Image *image,
   const Image *reconstruct_image,double *distortion,ExceptionInfo *exception)
 {
@@ -469,17 +480,23 @@ static MagickBooleanType GetFuzzDistortion(const Image *image,
   register ssize_t
     i;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
   status=MagickTrue;
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp parallel for schedule(static,4) shared(status) \
-    magick_threads(image,image,image->rows,1)
+    magick_threads(image,image,rows,1)
 #endif
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     double
       channel_distortion[MaxPixelChannels+1];
@@ -494,16 +511,15 @@ static MagickBooleanType GetFuzzDistortion(const Image *image,
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
     (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -550,7 +566,7 @@ static MagickBooleanType GetFuzzDistortion(const Image *image,
   reconstruct_view=DestroyCacheView(reconstruct_view);
   image_view=DestroyCacheView(image_view);
   for (i=0; i <= MaxPixelChannels; i++)
-    distortion[i]/=((double) image->columns*image->rows);
+    distortion[i]/=((double) columns*rows);
   distortion[CompositePixelChannel]/=(double) GetImageChannels(image);
   distortion[CompositePixelChannel]=sqrt(distortion[CompositePixelChannel]);
   return(status);
@@ -569,17 +585,23 @@ static MagickBooleanType GetMeanAbsoluteDistortion(const Image *image,
   register ssize_t
     i;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
   status=MagickTrue;
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp parallel for schedule(static,4) shared(status) \
-    magick_threads(image,image,image->rows,1)
+    magick_threads(image,image,rows,1)
 #endif
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     double
       channel_distortion[MaxPixelChannels+1];
@@ -594,16 +616,15 @@ static MagickBooleanType GetMeanAbsoluteDistortion(const Image *image,
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
     (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -650,7 +671,7 @@ static MagickBooleanType GetMeanAbsoluteDistortion(const Image *image,
   reconstruct_view=DestroyCacheView(reconstruct_view);
   image_view=DestroyCacheView(image_view);
   for (i=0; i <= MaxPixelChannels; i++)
-    distortion[i]/=((double) image->columns*image->rows);
+    distortion[i]/=((double) columns*rows);
   distortion[CompositePixelChannel]/=(double) GetImageChannels(image);
   return(status);
 }
@@ -670,6 +691,10 @@ static MagickBooleanType GetMeanErrorPerPixel(Image *image,
     maximum_error,
     mean_error;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
@@ -677,9 +702,11 @@ static MagickBooleanType GetMeanErrorPerPixel(Image *image,
   area=0.0;
   maximum_error=0.0;
   mean_error=0.0;
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     register const Quantum
       *restrict p,
@@ -688,15 +715,14 @@ static MagickBooleanType GetMeanErrorPerPixel(Image *image,
     register ssize_t
       x;
 
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL))
       {
         status=MagickFalse;
         break;
       }
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -726,8 +752,7 @@ static MagickBooleanType GetMeanErrorPerPixel(Image *image,
             (reconstruct_traits == UndefinedPixelTrait) ||
             ((reconstruct_traits & UpdatePixelTrait) == 0))
           continue;
-        distance=fabs((double) (Sa*p[i]-Da*GetPixelChannel(reconstruct_image,
-          channel,q)));
+        distance=fabs(Sa*p[i]-Da*GetPixelChannel(reconstruct_image,channel,q));
         distortion[i]+=distance;
         distortion[CompositePixelChannel]+=distance;
         mean_error+=distance*distance;
@@ -760,17 +785,23 @@ static MagickBooleanType GetMeanSquaredDistortion(const Image *image,
   register ssize_t
     i;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
   status=MagickTrue;
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp parallel for schedule(static,4) shared(status) \
-    magick_threads(image,image,image->rows,1)
+    magick_threads(image,image,rows,1)
 #endif
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     double
       channel_distortion[MaxPixelChannels+1];
@@ -785,16 +816,15 @@ static MagickBooleanType GetMeanSquaredDistortion(const Image *image,
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
     (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -841,7 +871,7 @@ static MagickBooleanType GetMeanSquaredDistortion(const Image *image,
   reconstruct_view=DestroyCacheView(reconstruct_view);
   image_view=DestroyCacheView(image_view);
   for (i=0; i <= MaxPixelChannels; i++)
-    distortion[i]/=((double) image->columns*image->rows);
+    distortion[i]/=((double) columns*rows);
   distortion[CompositePixelChannel]/=GetImageChannels(image);
   return(status);
 }
@@ -872,6 +902,10 @@ static MagickBooleanType GetNormalizedCrossCorrelationDistortion(
   register ssize_t
     i;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
@@ -895,10 +929,12 @@ static MagickBooleanType GetNormalizedCrossCorrelationDistortion(
   progress=0;
   for (i=0; i <= MaxPixelChannels; i++)
     distortion[i]=0.0;
-  area=1.0/((double) image->columns*image->rows);
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
+  area=1.0/((double) columns*rows);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     register const Quantum
       *restrict p,
@@ -909,15 +945,14 @@ static MagickBooleanType GetNormalizedCrossCorrelationDistortion(
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -956,8 +991,7 @@ static MagickBooleanType GetNormalizedCrossCorrelationDistortion(
         MagickBooleanType
           proceed;
 
-        proceed=SetImageProgress(image,SimilarityImageTag,progress++,
-          image->rows);
+        proceed=SetImageProgress(image,SimilarityImageTag,progress++,rows);
         if (proceed == MagickFalse)
           status=MagickFalse;
       }
@@ -1002,17 +1036,23 @@ static MagickBooleanType GetPeakAbsoluteDistortion(const Image *image,
   MagickBooleanType
     status;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
   status=MagickTrue;
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp parallel for schedule(static,4) shared(status) \
-    magick_threads(image,image,image->rows,1)
+    magick_threads(image,image,rows,1)
 #endif
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     double
       channel_distortion[MaxPixelChannels+1];
@@ -1027,16 +1067,15 @@ static MagickBooleanType GetPeakAbsoluteDistortion(const Image *image,
 
     if (status == MagickFalse)
       continue;
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (const Quantum *) NULL))
       {
         status=MagickFalse;
         continue;
       }
     (void) ResetMagickMemory(channel_distortion,0,sizeof(channel_distortion));
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       double
         Da,
@@ -1088,6 +1127,15 @@ static MagickBooleanType GetPeakAbsoluteDistortion(const Image *image,
   return(status);
 }
 
+static inline double MagickLog10(const double x)
+{
+#define Log10Epsilon  (1.0e-11)
+
+ if (fabs(x) < Log10Epsilon)
+   return(log10(Log10Epsilon));
+ return(log10(fabs(x)));
+}
+
 static MagickBooleanType GetPeakSignalToNoiseRatio(const Image *image,
   const Image *reconstruct_image,double *distortion,ExceptionInfo *exception)
 {
@@ -1099,156 +1147,100 @@ static MagickBooleanType GetPeakSignalToNoiseRatio(const Image *image,
 
   status=GetMeanSquaredDistortion(image,reconstruct_image,distortion,exception);
   for (i=0; i <= MaxPixelChannels; i++)
-    distortion[i]=20.0*log10((double) 1.0/sqrt(distortion[i]));
+    distortion[i]=20.0*MagickLog10((double) 1.0/sqrt(distortion[i]));
   return(status);
 }
 
-static Image *PerceptualHashImage(const Image *image,
-  const ColorspaceType colorspace,ExceptionInfo *exception)
-{
-  Image
-    *clone_image,
-    *phash_image;
-
-  MagickBooleanType
-    status;
-
-  /*
-    Transform colorspace then blur perceptual hash image.
-  */
-  clone_image=CloneImage(image,0,0,MagickTrue,exception);
-  if (clone_image == (Image *) NULL)
-    return((Image *) NULL);
-  status=TransformImageColorspace(clone_image,colorspace,exception);
-  if (status == MagickFalse)
-    {
-      clone_image=DestroyImage(clone_image);
-      return((Image *) NULL);
-    }
-  phash_image=BlurImage(clone_image,0.0,1.0,exception);
-  clone_image=DestroyImage(clone_image);
-  if (phash_image != (Image *) NULL)
-    phash_image->depth=8;
-  return(phash_image);
-}
-
 static MagickBooleanType GetPerceptualHashDistortion(const Image *image,
   const Image *reconstruct_image,double *distortion,ExceptionInfo *exception)
 {
-  ChannelMoments
-    *image_moments,
-    *reconstruct_moments;
-
-  Image
-    *phash_image;
+  ChannelPerceptualHash
+    *image_phash,
+    *reconstruct_phash;
 
-  register ssize_t
-    i;
+  ssize_t
+    channel;
 
   /*
-    Compute perceptual hash in the native image colorspace.
+    Compute perceptual hash in the sRGB colorspace.
   */
-  phash_image=PerceptualHashImage(image,image->colorspace,exception);
-  if (phash_image == (Image *) NULL)
-    return(MagickFalse);
-  image_moments=GetImageMoments(phash_image,exception);
-  phash_image=DestroyImage(phash_image);
-  if (image_moments == (ChannelMoments *) NULL)
+  image_phash=GetImagePerceptualHash(image,exception);
+  if (image_phash == (ChannelPerceptualHash *) NULL)
     return(MagickFalse);
-  phash_image=PerceptualHashImage(reconstruct_image,
-    reconstruct_image->colorspace,exception);
-  if (phash_image == (Image *) NULL)
-    {
-      image_moments=(ChannelMoments *) RelinquishMagickMemory(image_moments);
-      return(MagickFalse);
-    }
-  reconstruct_moments=GetImageMoments(phash_image,exception);
-  phash_image=DestroyImage(phash_image);
-  if (reconstruct_moments == (ChannelMoments *) NULL)
+  reconstruct_phash=GetImagePerceptualHash(reconstruct_image,exception);
+  if (reconstruct_phash == (ChannelPerceptualHash *) NULL)
     {
-      image_moments=(ChannelMoments *) RelinquishMagickMemory(image_moments);
+      image_phash=(ChannelPerceptualHash *) RelinquishMagickMemory(image_phash);
       return(MagickFalse);
     }
-  for (i=0; i < 7; i++)
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+  #pragma omp parallel for schedule(static,4)
+#endif
+  for (channel=0; channel < MaxPixelChannels; channel++)
   {
-   ssize_t
-      channel;
-      
-    /*
-      Compute sum of moment differences squared.
-    */
-    for (channel=0; channel < MaxPixelChannels; channel++)
+    double
+      difference;
+
+    register ssize_t
+      i;
+
+    difference=0.0;
+    for (i=0; i < MaximumNumberOfImageMoments; i++)
     {
       double
-        difference;
+        alpha,
+        beta;
 
-      if ((fabs(image_moments[channel].ellipse_intensity) < MagickEpsilon) ||
-          (fabs(reconstruct_moments[channel].ellipse_intensity) < MagickEpsilon))
-        continue;
-      difference=log10(fabs(reconstruct_moments[channel].I[i]))-
-         log10(fabs(image_moments[channel].I[i]));
-      distortion[channel]+=difference*difference;
-      distortion[CompositePixelChannel]+=difference*difference;
+      alpha=image_phash[channel].srgb_hu_phash[i];
+      beta=reconstruct_phash[channel].srgb_hu_phash[i];
+      difference+=(beta-alpha)*(beta-alpha);
     }
+    distortion[channel]+=difference;
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+    #pragma omp critical (MagickCore_GetPerceptualHashDistortion)
+#endif
+    distortion[CompositePixelChannel]+=difference;
   }
-  image_moments=(ChannelMoments *) RelinquishMagickMemory(image_moments);
-  reconstruct_moments=(ChannelMoments *) RelinquishMagickMemory(
-    reconstruct_moments);
-  if ((IsImageGray(image,exception) != MagickFalse) ||
-      (IsImageGray(reconstruct_image,exception) != MagickFalse))
-    return(MagickTrue);
   /*
     Compute perceptual hash in the HCLP colorspace.
   */
-  phash_image=PerceptualHashImage(image,HCLpColorspace,exception);
-  if (phash_image == (Image *) NULL)
-    return(MagickFalse);
-  image_moments=GetImageMoments(phash_image,exception);
-  phash_image=DestroyImage(phash_image);
-  if (image_moments == (ChannelMoments *) NULL)
-    return(MagickFalse);
-  phash_image=PerceptualHashImage(reconstruct_image,HCLpColorspace,exception);
-  if (phash_image == (Image *) NULL)
-    {
-      image_moments=(ChannelMoments *) RelinquishMagickMemory(image_moments);
-      return(MagickFalse);
-    }
-  reconstruct_moments=GetImageMoments(phash_image,exception);
-  phash_image=DestroyImage(phash_image);
-  if (reconstruct_moments == (ChannelMoments *) NULL)
-    {
-      image_moments=(ChannelMoments *) RelinquishMagickMemory(image_moments);
-      return(MagickFalse);
-    }
-  for (i=0; i < 7; i++)
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+  #pragma omp parallel for schedule(static,4)
+#endif
+  for (channel=0; channel < MaxPixelChannels; channel++)
   {
-   ssize_t
-      channel;
-      
-    /*
-      Compute sum of moment differences squared.
-    */
-    for (channel=0; channel < MaxPixelChannels; channel++)
+    double
+      difference;
+
+    register ssize_t
+      i;
+
+    difference=0.0;
+    for (i=0; i < MaximumNumberOfImageMoments; i++)
     {
       double
-        difference;
+        alpha,
+        beta;
 
-      if ((fabs(image_moments[channel].ellipse_intensity) < MagickEpsilon) ||
-          (fabs(reconstruct_moments[channel].ellipse_intensity) < MagickEpsilon))
-        continue;
-      difference=log10(fabs(reconstruct_moments[channel].I[i]))-
-         log10(fabs(image_moments[channel].I[i]));
-      distortion[channel]+=difference*difference;
-      distortion[CompositePixelChannel]+=difference*difference;
+      alpha=image_phash[channel].hclp_hu_phash[i];
+      beta=reconstruct_phash[channel].hclp_hu_phash[i];
+      difference+=(beta-alpha)*(beta-alpha);
     }
+    distortion[channel]+=difference;
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+    #pragma omp critical (MagickCore_GetPerceptualHashDistortion)
+#endif
+    distortion[CompositePixelChannel]+=difference;
   }
-  image_moments=(ChannelMoments *) RelinquishMagickMemory(image_moments);
-  reconstruct_moments=(ChannelMoments *) RelinquishMagickMemory(
-    reconstruct_moments);
+  /*
+    Free resources.
+  */
+  reconstruct_phash=(ChannelPerceptualHash *) RelinquishMagickMemory(
+    reconstruct_phash);
+  image_phash=(ChannelPerceptualHash *) RelinquishMagickMemory(image_phash);
   return(MagickTrue);
 }
 
-
 static MagickBooleanType GetRootMeanSquaredDistortion(const Image *image,
   const Image *reconstruct_image,double *distortion,ExceptionInfo *exception)
 {
@@ -1288,9 +1280,8 @@ MagickExport MagickBooleanType GetImageDistortion(Image *image,
   if (image->debug != MagickFalse)
     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
   if (metric != PerceptualHashErrorMetric)
-    if ((reconstruct_image->columns != image->columns) ||
-        (reconstruct_image->rows != image->rows))
-      ThrowBinaryException(ImageError,"ImageSizeDiffers",image->filename);
+    if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
+      ThrowBinaryException(ImageError,"ImageMorphologyDiffers",image->filename);
   /*
     Get image distortion.
   */
@@ -1426,11 +1417,10 @@ MagickExport double *GetImageDistortions(Image *image,
   if (image->debug != MagickFalse)
     (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",image->filename);
   if (metric != PerceptualHashErrorMetric)
-    if ((reconstruct_image->columns != image->columns) ||
-        (reconstruct_image->rows != image->rows))
+    if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
       {
         (void) ThrowMagickException(exception,GetMagickModule(),ImageError,
-          "ImageSizeDiffers","`%s'",image->filename);
+          "ImageMorphologyDiffers","`%s'",image->filename);
         return((double *) NULL);
       }
   /*
@@ -1581,6 +1571,10 @@ MagickExport MagickBooleanType IsImagesEqual(Image *image,
     mean_error,
     mean_error_per_pixel;
 
+  size_t
+    columns,
+    rows;
+
   ssize_t
     y;
 
@@ -1588,16 +1582,17 @@ MagickExport MagickBooleanType IsImagesEqual(Image *image,
   assert(image->signature == MagickSignature);
   assert(reconstruct_image != (const Image *) NULL);
   assert(reconstruct_image->signature == MagickSignature);
-  if ((reconstruct_image->columns != image->columns) ||
-      (reconstruct_image->rows != image->rows))
-    ThrowBinaryException(ImageError,"ImageSizeDiffers",image->filename);
+  if (ValidateImageMorphology(image,reconstruct_image) == MagickFalse)
+    ThrowBinaryException(ImageError,"ImageMorphologyDiffers",image->filename);
   area=0.0;
   maximum_error=0.0;
   mean_error_per_pixel=0.0;
   mean_error=0.0;
+  rows=MagickMax(image->rows,reconstruct_image->rows);
+  columns=MagickMax(image->columns,reconstruct_image->columns);
   image_view=AcquireVirtualCacheView(image,exception);
   reconstruct_view=AcquireVirtualCacheView(reconstruct_image,exception);
-  for (y=0; y < (ssize_t) image->rows; y++)
+  for (y=0; y < (ssize_t) rows; y++)
   {
     register const Quantum
       *restrict p,
@@ -1606,12 +1601,11 @@ MagickExport MagickBooleanType IsImagesEqual(Image *image,
     register ssize_t
       x;
 
-    p=GetCacheViewVirtualPixels(image_view,0,y,image->columns,1,exception);
-    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,reconstruct_image->columns,
-      1,exception);
+    p=GetCacheViewVirtualPixels(image_view,0,y,columns,1,exception);
+    q=GetCacheViewVirtualPixels(reconstruct_view,0,y,columns,1,exception);
     if ((p == (const Quantum *) NULL) || (q == (Quantum *) NULL))
       break;
-    for (x=0; x < (ssize_t) image->columns; x++)
+    for (x=0; x < (ssize_t) columns; x++)
     {
       register ssize_t
         i;
@@ -1637,10 +1631,13 @@ MagickExport MagickBooleanType IsImagesEqual(Image *image,
           continue;
         distance=fabs(p[i]-(double) GetPixelChannel(reconstruct_image,
           channel,q));
-        mean_error_per_pixel+=distance;
-        mean_error+=distance*distance;
-        if (distance > maximum_error)
-          maximum_error=distance;
+        if (distance >= MagickEpsilon)
+          {
+            mean_error_per_pixel+=distance;
+            mean_error+=distance*distance;
+            if (distance > maximum_error)
+              maximum_error=distance;
+          }
         area++;
       }
       p+=GetPixelChannels(image);
@@ -1757,9 +1754,9 @@ MagickExport Image *SimilarityImage(Image *image,const Image *reference,
   assert(exception->signature == MagickSignature);
   assert(offset != (RectangleInfo *) NULL);
   SetGeometry(reference,offset);
-  *similarity_metric=1.0;
-  if ((reference->columns > image->columns) || (reference->rows > image->rows))
-    ThrowImageException(ImageError,"ImageSizeDiffers");
+  *similarity_metric=MagickMaximumValue;
+  if (ValidateImageMorphology(image,reference) == MagickFalse)
+    ThrowImageException(ImageError,"ImageMorphologyDiffers");
   similarity_image=CloneImage(image,image->columns-reference->columns+1,
     image->rows-reference->rows+1,MagickTrue,exception);
   if (similarity_image == (Image *) NULL)
@@ -1822,12 +1819,17 @@ MagickExport Image *SimilarityImage(Image *image,const Image *reference,
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
       #pragma omp critical (MagickCore_SimilarityImage)
 #endif
+      if ((metric == NormalizedCrossCorrelationErrorMetric) ||
+          (metric == UndefinedErrorMetric))
+        similarity=1.0-similarity;
       if (similarity < *similarity_metric)
         {
           offset->x=x;
           offset->y=y;
           *similarity_metric=similarity;
         }
+      if (metric == PerceptualHashErrorMetric)
+        similarity=MagickMin(0.01*similarity,1.0);
       if (GetPixelReadMask(similarity_image,q) == 0)
         {
           SetPixelBackgoundColor(similarity_image,q);