]> granicus.if.org Git - imagemagick/blobdiff - MagickCore/fourier.c
(no commit message)
[imagemagick] / MagickCore / fourier.c
index b53043a5237619aa8a50883d8e5bc407e712ea5a..a649eca2f53c54d77987d1ba80e19eca51fcdeea 100644 (file)
 %                              Software Design                                %
 %                                Sean Burke                                   %
 %                               Fred Weinhaus                                 %
-%                                John Cristy                                  %
+%                                   Cristy                                    %
 %                                 July 2009                                   %
 %                                                                             %
 %                                                                             %
-%  Copyright 1999-2012 ImageMagick Studio LLC, a non-profit organization      %
+%  Copyright 1999-2014 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  %
@@ -43,7 +43,9 @@
   Include declarations.
 */
 #include "MagickCore/studio.h"
+#include "MagickCore/artifact.h"
 #include "MagickCore/attribute.h"
+#include "MagickCore/blob.h"
 #include "MagickCore/cache.h"
 #include "MagickCore/image.h"
 #include "MagickCore/image-private.h"
 #include "MagickCore/log.h"
 #include "MagickCore/memory_.h"
 #include "MagickCore/monitor.h"
+#include "MagickCore/monitor-private.h"
 #include "MagickCore/pixel-accessor.h"
+#include "MagickCore/pixel-private.h"
 #include "MagickCore/property.h"
 #include "MagickCore/quantum-private.h"
+#include "MagickCore/resource_.h"
+#include "MagickCore/string-private.h"
 #include "MagickCore/thread-private.h"
 #if defined(MAGICKCORE_FFTW_DELEGATE)
 #if defined(MAGICKCORE_HAVE_COMPLEX_H)
@@ -99,6 +105,262 @@ typedef struct _FourierInfo
 %                                                                             %
 %                                                                             %
 %                                                                             %
+%     C o m p l e x I m a g e s                                               %
+%                                                                             %
+%                                                                             %
+%                                                                             %
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%
+%  ComplexImages() performs complex mathematics on an image sequence.
+%
+%  The format of the ComplexImages method is:
+%
+%      MagickBooleanType ComplexImages(Image *images,
+%        const ComplexOperator op,ExceptionInfo *exception)
+%
+%  A description of each parameter follows:
+%
+%    o image: the image.
+%
+%    o op: A complex op.
+%
+%    o exception: return any errors or warnings in this structure.
+%
+*/
+MagickExport Image *ComplexImages(const Image *images,const ComplexOperator op,
+  ExceptionInfo *exception)
+{
+#define ComplexImageTag  "Complex/Image"
+
+  CacheView
+    *Ai_view,
+    *Ar_view,
+    *Bi_view,
+    *Br_view,
+    *Ci_view,
+    *Cr_view;
+
+  const char
+    *artifact;
+
+  const Image
+    *Ai_image,
+    *Ar_image,
+    *Bi_image,
+    *Br_image;
+
+  double
+    snr;
+
+  Image
+    *Ci_image,
+    *complex_images,
+    *Cr_image,
+    *image;
+
+  MagickBooleanType
+    status;
+
+  MagickOffsetType
+    progress;
+
+  ssize_t
+    y;
+
+  assert(images != (Image *) NULL);
+  assert(images->signature == MagickSignature);
+  if (images->debug != MagickFalse)
+    (void) LogMagickEvent(TraceEvent,GetMagickModule(),"%s",images->filename);
+  assert(exception != (ExceptionInfo *) NULL);
+  assert(exception->signature == MagickSignature);
+  if (images->next == (Image *) NULL)
+    {
+      (void) ThrowMagickException(exception,GetMagickModule(),ImageError,
+        "ImageSequenceRequired","`%s'",images->filename);
+      return((Image *) NULL);
+    }
+  image=CloneImage(images,images->columns,images->rows,MagickTrue,exception);
+  if (image == (Image *) NULL)
+    return((Image *) NULL);
+  if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
+    {
+      image=DestroyImageList(image);
+      return(image);
+    }
+  complex_images=NewImageList();
+  AppendImageToList(&complex_images,image);
+  image=CloneImage(images,images->columns,images->rows,MagickTrue,exception);
+  if (image == (Image *) NULL)
+    {
+      complex_images=DestroyImageList(complex_images);
+      return(complex_images);
+    }
+  if (SetImageStorageClass(image,DirectClass,exception) == MagickFalse)
+    {
+      complex_images=DestroyImageList(complex_images);
+      return(complex_images);
+    }
+  AppendImageToList(&complex_images,image);
+  /*
+    Apply complex mathematics to image pixels.
+  */
+  artifact=GetImageArtifact(image,"complex:snr");
+  snr=0.0;
+  if (artifact != (const char *) NULL)
+    snr=StringToDouble(artifact,(char **) NULL);
+  Ar_image=images;
+  Ai_image=images->next;
+  Br_image=images;
+  Bi_image=images->next;
+  if ((images->next->next != (Image *) NULL) &&
+      (images->next->next->next != (Image *) NULL))
+    {
+      Br_image=images->next->next;
+      Bi_image=images->next->next->next;
+    }
+  Cr_image=complex_images;
+  Ci_image=complex_images->next;
+  Ar_view=AcquireVirtualCacheView(Ar_image,exception);
+  Ai_view=AcquireVirtualCacheView(Ai_image,exception);
+  Br_view=AcquireVirtualCacheView(Br_image,exception);
+  Bi_view=AcquireVirtualCacheView(Bi_image,exception);
+  Cr_view=AcquireAuthenticCacheView(Cr_image,exception);
+  Ci_view=AcquireAuthenticCacheView(Ci_image,exception);
+  status=MagickTrue;
+  progress=0;
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+    #pragma omp parallel for schedule(static,4) shared(progress,status) \
+      magick_threads(images,complex_images,images->rows,1L)
+#endif
+  for (y=0; y < (ssize_t) images->rows; y++)
+  {
+    register const Quantum
+      *restrict Ai,
+      *restrict Ar,
+      *restrict Bi,
+      *restrict Br;
+
+    register Quantum
+      *restrict Ci,
+      *restrict Cr;
+
+    register ssize_t
+      x;
+
+    if (status == MagickFalse)
+      continue;
+    Ar=GetCacheViewVirtualPixels(Ar_view,0,y,images->columns,1,exception);
+    Ai=GetCacheViewVirtualPixels(Ai_view,0,y,images->columns,1,exception);
+    Br=GetCacheViewVirtualPixels(Br_view,0,y,images->columns,1,exception);
+    Bi=GetCacheViewVirtualPixels(Bi_view,0,y,images->columns,1,exception);
+    Cr=QueueCacheViewAuthenticPixels(Cr_view,0,y,images->columns,1,exception);
+    Ci=QueueCacheViewAuthenticPixels(Ci_view,0,y,images->columns,1,exception);
+    if ((Ar == (const Quantum *) NULL) || (Ai == (const Quantum *) NULL) || 
+        (Br == (const Quantum *) NULL) || (Bi == (const Quantum *) NULL) ||
+        (Cr == (Quantum *) NULL) || (Ci == (Quantum *) NULL))
+      {
+        status=MagickFalse;
+        continue;
+      }
+    for (x=0; x < (ssize_t) images->columns; x++)
+    {
+      register ssize_t
+        i;
+
+      for (i=0; i < (ssize_t) GetPixelChannels(images); i++)
+      {
+        switch (op)
+        {
+          case AddComplexOperator:
+          {
+            Cr[i]=Ar[i]+Br[i];
+            Ci[i]=Ai[i]+Bi[i];
+            break;
+          }
+          case ConjugateComplexOperator:
+          default:
+          {
+            Cr[i]=Ar[i];
+            Ci[i]=(-Bi[i]);
+            break;
+          }
+          case DivideComplexOperator:
+          {
+            double
+              gamma;
+
+            gamma=PerceptibleReciprocal(Br[i]*Br[i]+Bi[i]*Bi[i]+snr);
+            Cr[i]=gamma*(Ar[i]*Br[i]+Ai[i]*Bi[i]);
+            Ci[i]=gamma*(Ai[i]*Br[i]-Ar[i]*Bi[i]);
+            break;
+          }
+          case MagnitudePhaseComplexOperator:
+          {
+            Cr[i]=sqrt(Ar[i]*Ar[i]+Ai[i]*Ai[i]);
+            Ci[i]=atan2(Ai[i],Ar[i])/(2.0*MagickPI)+0.5;
+            break;
+          }
+          case MultiplyComplexOperator:
+          {
+            Cr[i]=QuantumScale*(Ar[i]*Br[i]-Ai[i]*Bi[i]);
+            Ci[i]=QuantumScale*(Ai[i]*Br[i]+Ar[i]*Bi[i]);
+            break;
+          }
+          case RealImaginaryComplexOperator:
+          {
+            Cr[i]=Ar[i]*cos(2.0*MagickPI*(Ai[i]-0.5));
+            Ci[i]=Ar[i]*sin(2.0*MagickPI*(Ai[i]-0.5));
+            break;
+          }
+          case SubtractComplexOperator:
+          {
+            Cr[i]=Ar[i]-Br[i];
+            Ci[i]=Ai[i]-Bi[i];
+            break;
+          }
+        }
+        Ar+=GetPixelChannels(Ar_image);
+        Ai+=GetPixelChannels(Ai_image);
+        Br+=GetPixelChannels(Br_image);
+        Bi+=GetPixelChannels(Bi_image);
+        Cr+=GetPixelChannels(Cr_image);
+        Ci+=GetPixelChannels(Ci_image);
+      }
+    }
+    if (SyncCacheViewAuthenticPixels(Ci_view,exception) == MagickFalse)
+      status=MagickFalse;
+    if (SyncCacheViewAuthenticPixels(Cr_view,exception) == MagickFalse)
+      status=MagickFalse;
+    if (images->progress_monitor != (MagickProgressMonitor) NULL)
+      {
+        MagickBooleanType
+          proceed;
+
+#if defined(MAGICKCORE_OPENMP_SUPPORT)
+        #pragma omp critical (MagickCore_ComplexImages)
+#endif
+        proceed=SetImageProgress(images,ComplexImageTag,progress++,
+          images->rows);
+        if (proceed == MagickFalse)
+          status=MagickFalse;
+      }
+  }
+  Cr_view=DestroyCacheView(Cr_view);
+  Ci_view=DestroyCacheView(Ci_view);
+  Br_view=DestroyCacheView(Br_view);
+  Bi_view=DestroyCacheView(Bi_view);
+  Ar_view=DestroyCacheView(Ar_view);
+  Ai_view=DestroyCacheView(Ai_view);
+  if (status == MagickFalse)
+    complex_images=DestroyImageList(complex_images);
+  return(complex_images);
+}
+\f
+/*
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%                                                                             %
+%                                                                             %
+%                                                                             %
 %     F o r w a r d F o u r i e r T r a n s f o r m I m a g e                 %
 %                                                                             %
 %                                                                             %
@@ -128,10 +390,13 @@ typedef struct _FourierInfo
 #if defined(MAGICKCORE_FFTW_DELEGATE)
 
 static MagickBooleanType RollFourier(const size_t width,const size_t height,
-  const ssize_t x_offset,const ssize_t y_offset,double *fourier)
+  const ssize_t x_offset,const ssize_t y_offset,double *roll_pixels)
 {
   double
-    *roll;
+    *source_pixels;
+
+  MemoryInfo
+    *source_info;
 
   register ssize_t
     i,
@@ -145,9 +410,10 @@ static MagickBooleanType RollFourier(const size_t width,const size_t height,
   /*
     Move zero frequency (DC, average color) from (0,0) to (width/2,height/2).
   */
-  roll=(double *) AcquireQuantumMemory((size_t) height,width*sizeof(*roll));
-  if (roll == (double *) NULL)
+  source_info=AcquireVirtualMemory(height,width*sizeof(*source_pixels));
+  if (source_info == (MemoryInfo *) NULL)
     return(MagickFalse);
+  source_pixels=(double *) GetVirtualMemoryBlob(source_info);
   i=0L;
   for (y=0L; y < (ssize_t) height; y++)
   {
@@ -163,16 +429,17 @@ static MagickBooleanType RollFourier(const size_t width,const size_t height,
       else
         u=((x+x_offset) > ((ssize_t) width-1L)) ? x+x_offset-(ssize_t) width :
           x+x_offset;
-      roll[v*width+u]=fourier[i++];
+      source_pixels[v*width+u]=roll_pixels[i++];
     }
   }
-  (void) CopyMagickMemory(fourier,roll,height*width*sizeof(*roll));
-  roll=(double *) RelinquishMagickMemory(roll);
+  (void) CopyMagickMemory(roll_pixels,source_pixels,height*width*
+    sizeof(*source_pixels));
+  source_info=RelinquishVirtualMemory(source_info);
   return(MagickTrue);
 }
 
 static MagickBooleanType ForwardQuadrantSwap(const size_t width,
-  const size_t height,double *source,double *destination)
+  const size_t height,double *source_pixels,double *forward_pixels)
 {
   MagickBooleanType
     status;
@@ -188,22 +455,24 @@ static MagickBooleanType ForwardQuadrantSwap(const size_t width,
     Swap quadrants.
   */
   center=(ssize_t) floor((double) width/2L)+1L;
-  status=RollFourier((size_t) center,height,0L,(ssize_t) height/2L,source);
+  status=RollFourier((size_t) center,height,0L,(ssize_t) height/2L,
+    source_pixels);
   if (status == MagickFalse)
     return(MagickFalse);
   for (y=0L; y < (ssize_t) height; y++)
     for (x=0L; x < (ssize_t) (width/2L-1L); x++)
-      destination[width*y+x+width/2L]=source[center*y+x];
+      forward_pixels[y*width+x+width/2L]=source_pixels[y*center+x];
   for (y=1; y < (ssize_t) height; y++)
     for (x=0L; x < (ssize_t) (width/2L-1L); x++)
-      destination[width*(height-y)+width/2L-x-1L]=source[center*y+x+1L];
+      forward_pixels[(height-y)*width+width/2L-x-1L]=
+        source_pixels[y*center+x+1L];
   for (x=0L; x < (ssize_t) (width/2L); x++)
-    destination[-x+width/2L-1L]=destination[x+width/2L+1L];
+    forward_pixels[-x+width/2L-1L]=forward_pixels[x+width/2L+1L];
   return(MagickTrue);
 }
 
 static void CorrectPhaseLHS(const size_t width,const size_t height,
-  double *fourier)
+  double *fourier_pixels)
 {
   register ssize_t
     x;
@@ -213,7 +482,7 @@ static void CorrectPhaseLHS(const size_t width,const size_t height,
 
   for (y=0L; y < (ssize_t) height; y++)
     for (x=0L; x < (ssize_t) (width/2L); x++)
-      fourier[y*width+x]*=(-1.0);
+      fourier_pixels[y*width+x]*=(-1.0);
 }
 
 static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
@@ -224,8 +493,8 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
     *phase_view;
 
   double
-    *magnitude_source,
-    *phase_source;
+    *magnitude_pixels,
+    *phase_pixels;
 
   Image
     *magnitude_image,
@@ -234,12 +503,16 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
   MagickBooleanType
     status;
 
-  register ssize_t
-    x;
+  MemoryInfo
+    *magnitude_info,
+    *phase_info;
 
   register Quantum
     *q;
 
+  register ssize_t
+    x;
+
   ssize_t
     i,
     y;
@@ -255,40 +528,45 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
   /*
     Create "Fourier Transform" image from constituent arrays.
   */
-  magnitude_source=(double *) AcquireQuantumMemory((size_t)
-    fourier_info->height,fourier_info->width*sizeof(*magnitude_source));
-  if (magnitude_source == (double *) NULL)
-    return(MagickFalse);
-  (void) ResetMagickMemory(magnitude_source,0,fourier_info->height*
-    fourier_info->width*sizeof(*magnitude_source));
-  phase_source=(double *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->width*sizeof(*phase_source));
-  if (phase_source == (double *) NULL)
+  magnitude_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->width*sizeof(*magnitude_pixels));
+  phase_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->width*sizeof(*phase_pixels));
+  if ((magnitude_info == (MemoryInfo *) NULL) ||
+      (phase_info == (MemoryInfo *) NULL))
     {
+      if (phase_info != (MemoryInfo *) NULL)
+        phase_info=RelinquishVirtualMemory(phase_info);
+      if (magnitude_info != (MemoryInfo *) NULL)
+        magnitude_info=RelinquishVirtualMemory(magnitude_info);
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
-      magnitude_source=(double *) RelinquishMagickMemory(magnitude_source);
       return(MagickFalse);
     }
-  status=ForwardQuadrantSwap(fourier_info->height,fourier_info->height,
-    magnitude,magnitude_source);
+  magnitude_pixels=(double *) GetVirtualMemoryBlob(magnitude_info);
+  (void) ResetMagickMemory(magnitude_pixels,0,fourier_info->height*
+    fourier_info->width*sizeof(*magnitude_pixels));
+  phase_pixels=(double *) GetVirtualMemoryBlob(phase_info);
+  (void) ResetMagickMemory(phase_pixels,0,fourier_info->height*
+    fourier_info->width*sizeof(*phase_pixels));
+  status=ForwardQuadrantSwap(fourier_info->width,fourier_info->height,
+    magnitude,magnitude_pixels);
   if (status != MagickFalse)
-    status=ForwardQuadrantSwap(fourier_info->height,fourier_info->height,phase,
-      phase_source);
-  CorrectPhaseLHS(fourier_info->height,fourier_info->height,phase_source);
+    status=ForwardQuadrantSwap(fourier_info->width,fourier_info->height,phase,
+      phase_pixels);
+  CorrectPhaseLHS(fourier_info->width,fourier_info->height,phase_pixels);
   if (fourier_info->modulus != MagickFalse)
     {
       i=0L;
       for (y=0L; y < (ssize_t) fourier_info->height; y++)
         for (x=0L; x < (ssize_t) fourier_info->width; x++)
         {
-          phase_source[i]/=(2.0*MagickPI);
-          phase_source[i]+=0.5;
+          phase_pixels[i]/=(2.0*MagickPI);
+          phase_pixels[i]+=0.5;
           i++;
         }
     }
-  magnitude_view=AcquireCacheView(magnitude_image);
-  phase_view=AcquireCacheView(phase_image);
+  magnitude_view=AcquireAuthenticCacheView(magnitude_image,exception);
   i=0L;
   for (y=0L; y < (ssize_t) fourier_info->height; y++)
   {
@@ -304,31 +582,31 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
         default:
         {
           SetPixelRed(magnitude_image,ClampToQuantum(QuantumRange*
-            magnitude_source[i]),q);
+            magnitude_pixels[i]),q);
           break;
         }
         case GreenPixelChannel:
         {
           SetPixelGreen(magnitude_image,ClampToQuantum(QuantumRange*
-            magnitude_source[i]),q);
+            magnitude_pixels[i]),q);
           break;
         }
         case BluePixelChannel:
         {
           SetPixelBlue(magnitude_image,ClampToQuantum(QuantumRange*
-            magnitude_source[i]),q);
+            magnitude_pixels[i]),q);
           break;
         }
         case BlackPixelChannel:
         {
           SetPixelBlack(magnitude_image,ClampToQuantum(QuantumRange*
-            magnitude_source[i]),q);
+            magnitude_pixels[i]),q);
           break;
         }
         case AlphaPixelChannel:
         {
           SetPixelAlpha(magnitude_image,ClampToQuantum(QuantumRange*
-            magnitude_source[i]),q);
+            magnitude_pixels[i]),q);
           break;
         }
       }
@@ -339,7 +617,9 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
     if (status == MagickFalse)
       break;
   }
+  magnitude_view=DestroyCacheView(magnitude_view);
   i=0L;
+  phase_view=AcquireAuthenticCacheView(phase_image,exception);
   for (y=0L; y < (ssize_t) fourier_info->height; y++)
   {
     q=GetCacheViewAuthenticPixels(phase_view,0L,y,fourier_info->height,1UL,
@@ -354,31 +634,31 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
         default:
         {
           SetPixelRed(phase_image,ClampToQuantum(QuantumRange*
-            phase_source[i]),q);
+            phase_pixels[i]),q);
           break;
         }
         case GreenPixelChannel:
         {
           SetPixelGreen(phase_image,ClampToQuantum(QuantumRange*
-            phase_source[i]),q);
+            phase_pixels[i]),q);
           break;
         }
         case BluePixelChannel:
         {
           SetPixelBlue(phase_image,ClampToQuantum(QuantumRange*
-            phase_source[i]),q);
+            phase_pixels[i]),q);
           break;
         }
         case BlackPixelChannel:
         {
           SetPixelBlack(phase_image,ClampToQuantum(QuantumRange*
-            phase_source[i]),q);
+            phase_pixels[i]),q);
           break;
         }
         case AlphaPixelChannel:
         {
           SetPixelAlpha(phase_image,ClampToQuantum(QuantumRange*
-            phase_source[i]),q);
+            phase_pixels[i]),q);
           break;
         }
       }
@@ -390,28 +670,34 @@ static MagickBooleanType ForwardFourier(const FourierInfo *fourier_info,
       break;
    }
   phase_view=DestroyCacheView(phase_view);
-  magnitude_view=DestroyCacheView(magnitude_view);
-  phase_source=(double *) RelinquishMagickMemory(phase_source);
-  magnitude_source=(double *) RelinquishMagickMemory(magnitude_source);
+  phase_info=RelinquishVirtualMemory(phase_info);
+  magnitude_info=RelinquishVirtualMemory(magnitude_info);
   return(status);
 }
 
 static MagickBooleanType ForwardFourierTransform(FourierInfo *fourier_info,
-  const Image *image,double *magnitude,double *phase,ExceptionInfo *exception)
+  const Image *image,double *magnitude_pixels,double *phase_pixels,
+  ExceptionInfo *exception)
 {
   CacheView
     *image_view;
 
+  const char
+    *value;
+
   double
-    n,
-    *source;
+    *source_pixels;
 
   fftw_complex
-    *fourier;
+    *forward_pixels;
 
   fftw_plan
     fftw_r2c_plan;
 
+  MemoryInfo
+    *forward_info,
+    *source_info;
+
   register const Quantum
     *p;
 
@@ -425,18 +711,19 @@ static MagickBooleanType ForwardFourierTransform(FourierInfo *fourier_info,
   /*
     Generate the forward Fourier transform.
   */
-  source=(double *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->width*sizeof(*source));
-  if (source == (double *) NULL)
+  source_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->width*sizeof(*source_pixels));
+  if (source_info == (MemoryInfo *) NULL)
     {
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
       return(MagickFalse);
     }
-  ResetMagickMemory(source,0,fourier_info->height*fourier_info->width*
-    sizeof(*source));
+  source_pixels=(double *) GetVirtualMemoryBlob(source_info);
+  ResetMagickMemory(source_pixels,0,fourier_info->height*fourier_info->width*
+    sizeof(*source_pixels));
   i=0L;
-  image_view=AcquireCacheView(image);
+  image_view=AcquireVirtualCacheView(image,exception);
   for (y=0L; y < (ssize_t) fourier_info->height; y++)
   {
     p=GetCacheViewVirtualPixels(image_view,0L,y,fourier_info->width,1UL,
@@ -450,27 +737,27 @@ static MagickBooleanType ForwardFourierTransform(FourierInfo *fourier_info,
         case RedPixelChannel:
         default:
         {
-          source[i]=QuantumScale*GetPixelRed(image,p);
+          source_pixels[i]=QuantumScale*GetPixelRed(image,p);
           break;
         }
         case GreenPixelChannel:
         {
-          source[i]=QuantumScale*GetPixelGreen(image,p);
+          source_pixels[i]=QuantumScale*GetPixelGreen(image,p);
           break;
         }
         case BluePixelChannel:
         {
-          source[i]=QuantumScale*GetPixelBlue(image,p);
+          source_pixels[i]=QuantumScale*GetPixelBlue(image,p);
           break;
         }
         case BlackPixelChannel:
         {
-          source[i]=QuantumScale*GetPixelBlack(image,p);
+          source_pixels[i]=QuantumScale*GetPixelBlack(image,p);
           break;
         }
         case AlphaPixelChannel:
         {
-          source[i]=QuantumScale*GetPixelAlpha(image,p);
+          source_pixels[i]=QuantumScale*GetPixelAlpha(image,p);
           break;
         }
       }
@@ -479,38 +766,47 @@ static MagickBooleanType ForwardFourierTransform(FourierInfo *fourier_info,
     }
   }
   image_view=DestroyCacheView(image_view);
-  fourier=(fftw_complex *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->center*sizeof(*fourier));
-  if (fourier == (fftw_complex *) NULL)
+  forward_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->center*sizeof(*forward_pixels));
+  if (forward_info == (MemoryInfo *) NULL)
     {
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
-      source=(double *) RelinquishMagickMemory(source);
+      source_info=(MemoryInfo *) RelinquishVirtualMemory(source_info);
       return(MagickFalse);
     }
+  forward_pixels=(fftw_complex *) GetVirtualMemoryBlob(forward_info);
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp critical (MagickCore_ForwardFourierTransform)
 #endif
-  fftw_r2c_plan=fftw_plan_dft_r2c_2d(fourier_info->width,fourier_info->width,
-    source,fourier,FFTW_ESTIMATE);
+  fftw_r2c_plan=fftw_plan_dft_r2c_2d(fourier_info->width,fourier_info->height,
+    source_pixels,forward_pixels,FFTW_ESTIMATE);
   fftw_execute(fftw_r2c_plan);
   fftw_destroy_plan(fftw_r2c_plan);
-  source=(double *) RelinquishMagickMemory(source);
-  /*
-    Normalize Fourier transform.
-  */
-  n=(double) fourier_info->width*(double) fourier_info->width;
-  i=0L;
-  for (y=0L; y < (ssize_t) fourier_info->height; y++)
-    for (x=0L; x < (ssize_t) fourier_info->center; x++)
+  source_info=(MemoryInfo *) RelinquishVirtualMemory(source_info);
+  value=GetImageArtifact(image,"fourier:normalize");
+  if ((value == (const char *) NULL) || (LocaleCompare(value,"forward") == 0))
     {
+      double
+        gamma;
+
+      /*
+        Normalize Fourier transform.
+      */
+      i=0L;
+      gamma=PerceptibleReciprocal((double) fourier_info->width*
+        fourier_info->height);
+      for (y=0L; y < (ssize_t) fourier_info->height; y++)
+        for (x=0L; x < (ssize_t) fourier_info->center; x++)
+        {
 #if defined(MAGICKCORE_HAVE_COMPLEX_H)
-      fourier[i]/=n;
+          forward_pixels[i]*=gamma;
 #else
-      fourier[i][0]/=n;
-      fourier[i][1]/=n;
+          forward_pixels[i][0]*=gamma;
+          forward_pixels[i][1]*=gamma;
 #endif
-      i++;
+          i++;
+        }
     }
   /*
     Generate magnitude and phase (or real and imaginary).
@@ -520,19 +816,19 @@ static MagickBooleanType ForwardFourierTransform(FourierInfo *fourier_info,
     for (y=0L; y < (ssize_t) fourier_info->height; y++)
       for (x=0L; x < (ssize_t) fourier_info->center; x++)
       {
-        magnitude[i]=cabs(fourier[i]);
-        phase[i]=carg(fourier[i]);
+        magnitude_pixels[i]=cabs(forward_pixels[i]);
+        phase_pixels[i]=carg(forward_pixels[i]);
         i++;
       }
   else
     for (y=0L; y < (ssize_t) fourier_info->height; y++)
       for (x=0L; x < (ssize_t) fourier_info->center; x++)
       {
-        magnitude[i]=creal(fourier[i]);
-        phase[i]=cimag(fourier[i]);
+        magnitude_pixels[i]=creal(forward_pixels[i]);
+        phase_pixels[i]=cimag(forward_pixels[i]);
         i++;
       }
-  fourier=(fftw_complex *) RelinquishMagickMemory(fourier);
+  forward_info=(MemoryInfo *) RelinquishVirtualMemory(forward_info);
   return(MagickTrue);
 }
 
@@ -541,11 +837,8 @@ static MagickBooleanType ForwardFourierTransformChannel(const Image *image,
   Image *fourier_image,ExceptionInfo *exception)
 {
   double
-    *magnitude,
-    *phase;
-
-  fftw_complex
-    *fourier;
+    *magnitude_pixels,
+    *phase_pixels;
 
   FourierInfo
     fourier_info;
@@ -553,10 +846,15 @@ static MagickBooleanType ForwardFourierTransformChannel(const Image *image,
   MagickBooleanType
     status;
 
+  MemoryInfo
+    *magnitude_info,
+    *phase_info;
+
   size_t
     extent;
 
   fourier_info.width=image->columns;
+  fourier_info.height=image->rows;
   if ((image->columns != image->rows) || ((image->columns % 2) != 0) ||
       ((image->rows % 2) != 0))
     {
@@ -564,43 +862,33 @@ static MagickBooleanType ForwardFourierTransformChannel(const Image *image,
       fourier_info.width=(extent & 0x01) == 1 ? extent+1UL : extent;
     }
   fourier_info.height=fourier_info.width;
-  fourier_info.center=(ssize_t) floor((double) fourier_info.width/2.0)+1L;
+  fourier_info.center=(ssize_t) floor((double) fourier_info.width/2L)+1L;
   fourier_info.channel=channel;
   fourier_info.modulus=modulus;
-  magnitude=(double *) AcquireQuantumMemory((size_t) fourier_info.height,
-    fourier_info.center*sizeof(*magnitude));
-  if (magnitude == (double *) NULL)
+  magnitude_info=AcquireVirtualMemory((size_t) fourier_info.height,
+    fourier_info.center*sizeof(*magnitude_pixels));
+  phase_info=AcquireVirtualMemory((size_t) fourier_info.height,
+    fourier_info.center*sizeof(*phase_pixels));
+  if ((magnitude_info == (MemoryInfo *) NULL) ||
+      (phase_info == (MemoryInfo *) NULL))
     {
+      if (phase_info != (MemoryInfo *) NULL)
+        phase_info=RelinquishVirtualMemory(phase_info);
+      if (magnitude_info == (MemoryInfo *) NULL)
+        magnitude_info=RelinquishVirtualMemory(magnitude_info);
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
       return(MagickFalse);
     }
-  phase=(double *) AcquireQuantumMemory((size_t) fourier_info.height,
-    fourier_info.center*sizeof(*phase));
-  if (phase == (double *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
-      magnitude=(double *) RelinquishMagickMemory(magnitude);
-      return(MagickFalse);
-    }
-  fourier=(fftw_complex *) AcquireQuantumMemory((size_t) fourier_info.height,
-    fourier_info.center*sizeof(*fourier));
-  if (fourier == (fftw_complex *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
-      phase=(double *) RelinquishMagickMemory(phase);
-      magnitude=(double *) RelinquishMagickMemory(magnitude);
-      return(MagickFalse);
-    }
-  status=ForwardFourierTransform(&fourier_info,image,magnitude,phase,exception);
+  magnitude_pixels=(double *) GetVirtualMemoryBlob(magnitude_info);
+  phase_pixels=(double *) GetVirtualMemoryBlob(phase_info);
+  status=ForwardFourierTransform(&fourier_info,image,magnitude_pixels,
+    phase_pixels,exception);
   if (status != MagickFalse)
-    status=ForwardFourier(&fourier_info,fourier_image,magnitude,phase,
-      exception);
-  fourier=(fftw_complex *) RelinquishMagickMemory(fourier);
-  phase=(double *) RelinquishMagickMemory(phase);
-  magnitude=(double *) RelinquishMagickMemory(magnitude);
+    status=ForwardFourier(&fourier_info,fourier_image,magnitude_pixels,
+      phase_pixels,exception);
+  phase_info=RelinquishVirtualMemory(phase_info);
+  magnitude_info=RelinquishVirtualMemory(magnitude_info);
   return(status);
 }
 #endif
@@ -615,7 +903,7 @@ MagickExport Image *ForwardFourierTransformImage(const Image *image,
 #if !defined(MAGICKCORE_FFTW_DELEGATE)
   (void) modulus;
   (void) ThrowMagickException(exception,GetMagickModule(),
-    MissingDelegateWarning,"DelegateLibrarySupportNotBuiltIn","`%s' (FFTW)",
+    MissingDelegateWarning,"DelegateLibrarySupportNotBuiltIn","'%s' (FFTW)",
     image->filename);
 #else
   {
@@ -624,16 +912,19 @@ MagickExport Image *ForwardFourierTransformImage(const Image *image,
 
     size_t
       extent,
+      height,
       width;
 
     width=image->columns;
+    height=image->rows;
     if ((image->columns != image->rows) || ((image->columns % 2) != 0) ||
         ((image->rows % 2) != 0))
       {
         extent=image->columns < image->rows ? image->rows : image->columns;
         width=(extent & 0x01) == 1 ? extent+1UL : extent;
       }
-    magnitude_image=CloneImage(image,width,width,MagickFalse,exception);
+    height=width;
+    magnitude_image=CloneImage(image,width,height,MagickTrue,exception);
     if (magnitude_image != (Image *) NULL)
       {
         Image
@@ -641,7 +932,7 @@ MagickExport Image *ForwardFourierTransformImage(const Image *image,
 
         magnitude_image->storage_class=DirectClass;
         magnitude_image->depth=32UL;
-        phase_image=CloneImage(image,width,width,MagickFalse,exception);
+        phase_image=CloneImage(image,width,height,MagickTrue,exception);
         if (phase_image == (Image *) NULL)
           magnitude_image=DestroyImage(magnitude_image);
         else
@@ -726,7 +1017,7 @@ MagickExport Image *ForwardFourierTransformImage(const Image *image,
                   thread_status;
 
                 thread_status=MagickTrue;
-                if (image->matte != MagickFalse)
+                if (image->alpha_trait == BlendPixelTrait)
                   thread_status=ForwardFourierTransformChannel(image,
                     AlphaPixelChannel,modulus,fourier_image,exception);
                 if (thread_status == MagickFalse)
@@ -791,34 +1082,38 @@ static MagickBooleanType InverseQuadrantSwap(const size_t width,
   /*
     Swap quadrants.
   */
-  center=(ssize_t) floor((double) width/2.0)+1L;
+  center=(ssize_t) floor((double) width/2L)+1L;
   for (y=1L; y < (ssize_t) height; y++)
     for (x=0L; x < (ssize_t) (width/2L+1L); x++)
-      destination[center*(height-y)-x+width/2L]=source[y*width+x];
+      destination[(height-y)*center-x+width/2L]=source[y*width+x];
   for (y=0L; y < (ssize_t) height; y++)
-    destination[center*y]=source[y*width+width/2L];
+    destination[y*center]=source[y*width+width/2L];
   for (x=0L; x < center; x++)
     destination[x]=source[center-x-1L];
   return(RollFourier(center,height,0L,(ssize_t) height/-2L,destination));
 }
 
 static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
-  const Image *magnitude_image,const Image *phase_image,fftw_complex *fourier,
-  ExceptionInfo *exception)
+  const Image *magnitude_image,const Image *phase_image,
+  fftw_complex *fourier_pixels,ExceptionInfo *exception)
 {
   CacheView
     *magnitude_view,
     *phase_view;
 
   double
-    *magnitude,
-    *phase,
-    *magnitude_source,
-    *phase_source;
+    *inverse_pixels,
+    *magnitude_pixels,
+    *phase_pixels;
 
   MagickBooleanType
     status;
 
+  MemoryInfo
+    *inverse_info,
+    *magnitude_info,
+    *phase_info;
+
   register const Quantum
     *p;
 
@@ -832,27 +1127,32 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
   /*
     Inverse fourier - read image and break down into a double array.
   */
-  magnitude_source=(double *) AcquireQuantumMemory((size_t)
-    fourier_info->height,fourier_info->width*sizeof(*magnitude_source));
-  if (magnitude_source == (double *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",
-        magnitude_image->filename);
-      return(MagickFalse);
-    }
-  phase_source=(double *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->width*sizeof(*phase_source));
-  if (phase_source == (double *) NULL)
+  magnitude_info=AcquireVirtualMemory((size_t)fourier_info->height,
+    fourier_info->width*sizeof(*magnitude_pixels));
+  phase_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->width*sizeof(*phase_pixels));
+  inverse_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->center*sizeof(*inverse_pixels));
+  if ((magnitude_info == (MemoryInfo *) NULL) ||
+      (phase_info == (MemoryInfo *) NULL) ||
+      (inverse_info == (MemoryInfo *) NULL))
     {
+      if (magnitude_info != (MemoryInfo *) NULL)
+        magnitude_info=RelinquishVirtualMemory(magnitude_info);
+      if (phase_info != (MemoryInfo *) NULL)
+        phase_info=RelinquishVirtualMemory(phase_info);
+      if (inverse_info != (MemoryInfo *) NULL)
+        inverse_info=RelinquishVirtualMemory(inverse_info);
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",
         magnitude_image->filename);
-      magnitude_source=(double *) RelinquishMagickMemory(magnitude_source);
       return(MagickFalse);
     }
+  magnitude_pixels=(double *) GetVirtualMemoryBlob(magnitude_info);
+  phase_pixels=(double *) GetVirtualMemoryBlob(phase_info);
+  inverse_pixels=(double *) GetVirtualMemoryBlob(inverse_info);
   i=0L;
-  magnitude_view=AcquireCacheView(magnitude_image);
+  magnitude_view=AcquireVirtualCacheView(magnitude_image,exception);
   for (y=0L; y < (ssize_t) fourier_info->height; y++)
   {
     p=GetCacheViewVirtualPixels(magnitude_view,0L,y,fourier_info->width,1UL,
@@ -866,27 +1166,27 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
         case RedPixelChannel:
         default:
         {
-          magnitude_source[i]=QuantumScale*GetPixelRed(magnitude_image,p);
+          magnitude_pixels[i]=QuantumScale*GetPixelRed(magnitude_image,p);
           break;
         }
         case GreenPixelChannel:
         {
-          magnitude_source[i]=QuantumScale*GetPixelGreen(magnitude_image,p);
+          magnitude_pixels[i]=QuantumScale*GetPixelGreen(magnitude_image,p);
           break;
         }
         case BluePixelChannel:
         {
-          magnitude_source[i]=QuantumScale*GetPixelBlue(magnitude_image,p);
+          magnitude_pixels[i]=QuantumScale*GetPixelBlue(magnitude_image,p);
           break;
         }
         case BlackPixelChannel:
         {
-          magnitude_source[i]=QuantumScale*GetPixelBlack(magnitude_image,p);
+          magnitude_pixels[i]=QuantumScale*GetPixelBlack(magnitude_image,p);
           break;
         }
         case AlphaPixelChannel:
         {
-          magnitude_source[i]=QuantumScale*GetPixelAlpha(magnitude_image,p);
+          magnitude_pixels[i]=QuantumScale*GetPixelAlpha(magnitude_image,p);
           break;
         }
       }
@@ -894,8 +1194,13 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
       p+=GetPixelChannels(magnitude_image);
     }
   }
+  magnitude_view=DestroyCacheView(magnitude_view);
+  status=InverseQuadrantSwap(fourier_info->width,fourier_info->height,
+    magnitude_pixels,inverse_pixels);
+  (void) CopyMagickMemory(magnitude_pixels,inverse_pixels,fourier_info->height*
+    fourier_info->center*sizeof(*magnitude_pixels));
   i=0L;
-  phase_view=AcquireCacheView(phase_image);
+  phase_view=AcquireVirtualCacheView(phase_image,exception);
   for (y=0L; y < (ssize_t) fourier_info->height; y++)
   {
     p=GetCacheViewVirtualPixels(phase_view,0,y,fourier_info->width,1,
@@ -909,27 +1214,27 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
         case RedPixelChannel:
         default:
         {
-          phase_source[i]=QuantumScale*GetPixelRed(phase_image,p);
+          phase_pixels[i]=QuantumScale*GetPixelRed(phase_image,p);
           break;
         }
         case GreenPixelChannel:
         {
-          phase_source[i]=QuantumScale*GetPixelGreen(phase_image,p);
+          phase_pixels[i]=QuantumScale*GetPixelGreen(phase_image,p);
           break;
         }
         case BluePixelChannel:
         {
-          phase_source[i]=QuantumScale*GetPixelBlue(phase_image,p);
+          phase_pixels[i]=QuantumScale*GetPixelBlue(phase_image,p);
           break;
         }
         case BlackPixelChannel:
         {
-          phase_source[i]=QuantumScale*GetPixelBlack(phase_image,p);
+          phase_pixels[i]=QuantumScale*GetPixelBlack(phase_image,p);
           break;
         }
         case AlphaPixelChannel:
         {
-          phase_source[i]=QuantumScale*GetPixelAlpha(phase_image,p);
+          phase_pixels[i]=QuantumScale*GetPixelAlpha(phase_image,p);
           break;
         }
       }
@@ -943,42 +1248,19 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
       for (y=0L; y < (ssize_t) fourier_info->height; y++)
         for (x=0L; x < (ssize_t) fourier_info->width; x++)
         {
-          phase_source[i]-=0.5;
-          phase_source[i]*=(2.0*MagickPI);
+          phase_pixels[i]-=0.5;
+          phase_pixels[i]*=(2.0*MagickPI);
           i++;
         }
     }
-  magnitude_view=DestroyCacheView(magnitude_view);
   phase_view=DestroyCacheView(phase_view);
-  magnitude=(double *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->center*sizeof(*magnitude));
-  if (magnitude == (double *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",
-        magnitude_image->filename);
-      magnitude_source=(double *) RelinquishMagickMemory(magnitude_source);
-      phase_source=(double *) RelinquishMagickMemory(phase_source);
-      return(MagickFalse);
-    }
-  status=InverseQuadrantSwap(fourier_info->width,fourier_info->height,
-    magnitude_source,magnitude);
-  magnitude_source=(double *) RelinquishMagickMemory(magnitude_source);
-  phase=(double *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->width*sizeof(*phase));
-  if (phase == (double *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",
-        magnitude_image->filename);
-      phase_source=(double *) RelinquishMagickMemory(phase_source);
-      return(MagickFalse);
-    }
-  CorrectPhaseLHS(fourier_info->width,fourier_info->width,phase_source);
+  CorrectPhaseLHS(fourier_info->width,fourier_info->height,phase_pixels);
   if (status != MagickFalse)
     status=InverseQuadrantSwap(fourier_info->width,fourier_info->height,
-      phase_source,phase);
-  phase_source=(double *) RelinquishMagickMemory(phase_source);
+      phase_pixels,inverse_pixels);
+  (void) CopyMagickMemory(phase_pixels,inverse_pixels,fourier_info->height*
+    fourier_info->center*sizeof(*phase_pixels));
+  inverse_info=RelinquishVirtualMemory(inverse_info);
   /*
     Merge two sets.
   */
@@ -988,10 +1270,11 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
        for (x=0L; x < (ssize_t) fourier_info->center; x++)
        {
 #if defined(MAGICKCORE_HAVE_COMPLEX_H)
-         fourier[i]=magnitude[i]*cos(phase[i])+I*magnitude[i]*sin(phase[i]);
+         fourier_pixels[i]=magnitude_pixels[i]*cos(phase_pixels[i])+I*
+           magnitude_pixels[i]*sin(phase_pixels[i]);
 #else
-         fourier[i][0]=magnitude[i]*cos(phase[i]);
-         fourier[i][1]=magnitude[i]*sin(phase[i]);
+         fourier_pixels[i][0]=magnitude_pixels[i]*cos(phase_pixels[i]);
+         fourier_pixels[i][1]=magnitude_pixels[i]*sin(phase_pixels[i]);
 #endif
          i++;
       }
@@ -1000,30 +1283,36 @@ static MagickBooleanType InverseFourier(FourierInfo *fourier_info,
       for (x=0L; x < (ssize_t) fourier_info->center; x++)
       {
 #if defined(MAGICKCORE_HAVE_COMPLEX_H)
-        fourier[i]=magnitude[i]+I*phase[i];
+        fourier_pixels[i]=magnitude_pixels[i]+I*phase_pixels[i];
 #else
-        fourier[i][0]=magnitude[i];
-        fourier[i][1]=phase[i];
+        fourier_pixels[i][0]=magnitude_pixels[i];
+        fourier_pixels[i][1]=phase_pixels[i];
 #endif
         i++;
       }
-  phase=(double *) RelinquishMagickMemory(phase);
-  magnitude=(double *) RelinquishMagickMemory(magnitude);
+  magnitude_info=RelinquishVirtualMemory(magnitude_info);
+  phase_info=RelinquishVirtualMemory(phase_info);
   return(status);
 }
 
 static MagickBooleanType InverseFourierTransform(FourierInfo *fourier_info,
-  fftw_complex *fourier,Image *image,ExceptionInfo *exception)
+  fftw_complex *fourier_pixels,Image *image,ExceptionInfo *exception)
 {
   CacheView
     *image_view;
 
+  const char
+    *value;
+
   double
-    *source;
+    *source_pixels;
 
   fftw_plan
     fftw_c2r_plan;
 
+  MemoryInfo
+    *source_info;
+
   register Quantum
     *q;
 
@@ -1034,25 +1323,50 @@ static MagickBooleanType InverseFourierTransform(FourierInfo *fourier_info,
   ssize_t
     y;
 
-  source=(double *) AcquireQuantumMemory((size_t) fourier_info->height,
-    fourier_info->width*sizeof(*source));
-  if (source == (double *) NULL)
+  source_info=AcquireVirtualMemory((size_t) fourier_info->height,
+    fourier_info->width*sizeof(*source_pixels));
+  if (source_info == (MemoryInfo *) NULL)
     {
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",image->filename);
       return(MagickFalse);
     }
+  source_pixels=(double *) GetVirtualMemoryBlob(source_info);
+  value=GetImageArtifact(image,"fourier:normalize");
+  if (LocaleCompare(value,"inverse") == 0)
+    {
+      double
+        gamma;
+
+      /*
+        Normalize inverse transform.
+      */
+      i=0L;
+      gamma=PerceptibleReciprocal((double) fourier_info->width*
+        fourier_info->height);
+      for (y=0L; y < (ssize_t) fourier_info->height; y++)
+        for (x=0L; x < (ssize_t) fourier_info->center; x++)
+        {
+#if defined(MAGICKCORE_HAVE_COMPLEX_H)
+          fourier_pixels[i]*=gamma;
+#else
+          fourier_pixels[i][0]*=gamma;
+          fourier_pixels[i][1]*=gamma;
+#endif
+          i++;
+        }
+    }
 #if defined(MAGICKCORE_OPENMP_SUPPORT)
   #pragma omp critical (MagickCore_InverseFourierTransform)
 #endif
   {
     fftw_c2r_plan=fftw_plan_dft_c2r_2d(fourier_info->width,fourier_info->height,
-      fourier,source,FFTW_ESTIMATE);
+      fourier_pixels,source_pixels,FFTW_ESTIMATE);
     fftw_execute(fftw_c2r_plan);
     fftw_destroy_plan(fftw_c2r_plan);
   }
   i=0L;
-  image_view=AcquireCacheView(image);
+  image_view=AcquireAuthenticCacheView(image,exception);
   for (y=0L; y < (ssize_t) fourier_info->height; y++)
   {
     if (y >= (ssize_t) image->rows)
@@ -1063,35 +1377,40 @@ static MagickBooleanType InverseFourierTransform(FourierInfo *fourier_info,
       break;
     for (x=0L; x < (ssize_t) fourier_info->width; x++)
     {
-      switch (fourier_info->channel)
-      {
-        case RedPixelChannel:
-        default:
-        {
-          SetPixelRed(image,ClampToQuantum(QuantumRange*source[i]),q);
-          break;
-        }
-        case GreenPixelChannel:
+      if (x < (ssize_t) image->columns)
+        switch (fourier_info->channel)
         {
-          SetPixelGreen(image,ClampToQuantum(QuantumRange*source[i]),q);
-          break;
-        }
-        case BluePixelChannel:
-        {
-          SetPixelBlue(image,ClampToQuantum(QuantumRange*source[i]),q);
-          break;
-        }
-        case BlackPixelChannel:
-        {
-          SetPixelBlack(image,ClampToQuantum(QuantumRange*source[i]),q);
-          break;
-        }
-        case AlphaPixelChannel:
-        {
-          SetPixelAlpha(image,ClampToQuantum(QuantumRange*source[i]),q);
-          break;
+          case RedPixelChannel:
+          default:
+          {
+            SetPixelRed(image,ClampToQuantum(QuantumRange*source_pixels[i]),q);
+            break;
+          }
+          case GreenPixelChannel:
+          {
+            SetPixelGreen(image,ClampToQuantum(QuantumRange*source_pixels[i]),
+              q);
+            break;
+          }
+          case BluePixelChannel:
+          {
+            SetPixelBlue(image,ClampToQuantum(QuantumRange*source_pixels[i]),
+              q);
+            break;
+          }
+          case BlackPixelChannel:
+          {
+            SetPixelBlack(image,ClampToQuantum(QuantumRange*source_pixels[i]),
+              q);
+            break;
+          }
+          case AlphaPixelChannel:
+          {
+            SetPixelAlpha(image,ClampToQuantum(QuantumRange*source_pixels[i]),
+              q);
+            break;
+          }
         }
-      }
       i++;
       q+=GetPixelChannels(image);
     }
@@ -1099,7 +1418,7 @@ static MagickBooleanType InverseFourierTransform(FourierInfo *fourier_info,
       break;
   }
   image_view=DestroyCacheView(image_view);
-  source=(double *) RelinquishMagickMemory(source);
+  source_info=RelinquishVirtualMemory(source_info);
   return(MagickTrue);
 }
 
@@ -1108,12 +1427,8 @@ static MagickBooleanType InverseFourierTransformChannel(
   const PixelChannel channel,const MagickBooleanType modulus,
   Image *fourier_image,ExceptionInfo *exception)
 {
-  double
-    *magnitude,
-    *phase;
-
   fftw_complex
-    *fourier;
+    *inverse_pixels;
 
   FourierInfo
     fourier_info;
@@ -1121,10 +1436,14 @@ static MagickBooleanType InverseFourierTransformChannel(
   MagickBooleanType
     status;
 
+  MemoryInfo
+    *inverse_info;
+
   size_t
     extent;
 
   fourier_info.width=magnitude_image->columns;
+  fourier_info.height=magnitude_image->rows;
   if ((magnitude_image->columns != magnitude_image->rows) ||
       ((magnitude_image->columns % 2) != 0) ||
       ((magnitude_image->rows % 2) != 0))
@@ -1134,47 +1453,25 @@ static MagickBooleanType InverseFourierTransformChannel(
       fourier_info.width=(extent & 0x01) == 1 ? extent+1UL : extent;
     }
   fourier_info.height=fourier_info.width;
-  fourier_info.center=(ssize_t) floor((double) fourier_info.width/2.0)+1L;
+  fourier_info.center=(ssize_t) floor((double) fourier_info.width/2L)+1L;
   fourier_info.channel=channel;
   fourier_info.modulus=modulus;
-  magnitude=(double *) AcquireQuantumMemory((size_t) fourier_info.height,
-    fourier_info.center*sizeof(*magnitude));
-  if (magnitude == (double *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",
-        magnitude_image->filename);
-      return(MagickFalse);
-    }
-  phase=(double *) AcquireQuantumMemory((size_t) fourier_info.height,
-    fourier_info.center*sizeof(*phase));
-  if (phase == (double *) NULL)
-    {
-      (void) ThrowMagickException(exception,GetMagickModule(),
-        ResourceLimitError,"MemoryAllocationFailed","`%s'",
-        magnitude_image->filename);
-      magnitude=(double *) RelinquishMagickMemory(magnitude);
-      return(MagickFalse);
-    }
-  fourier=(fftw_complex *) AcquireQuantumMemory((size_t) fourier_info.height,
-    fourier_info.center*sizeof(*fourier));
-  if (fourier == (fftw_complex *) NULL)
+  inverse_info=AcquireVirtualMemory((size_t) fourier_info.height,
+    fourier_info.center*sizeof(*inverse_pixels));
+  if (inverse_info == (MemoryInfo *) NULL)
     {
       (void) ThrowMagickException(exception,GetMagickModule(),
         ResourceLimitError,"MemoryAllocationFailed","`%s'",
         magnitude_image->filename);
-      phase=(double *) RelinquishMagickMemory(phase);
-      magnitude=(double *) RelinquishMagickMemory(magnitude);
       return(MagickFalse);
     }
-  status=InverseFourier(&fourier_info,magnitude_image,phase_image,fourier,
-   exception);
+  inverse_pixels=(fftw_complex *) GetVirtualMemoryBlob(inverse_info);
+  status=InverseFourier(&fourier_info,magnitude_image,phase_image,
+    inverse_pixels,exception);
   if (status != MagickFalse)
-    status=InverseFourierTransform(&fourier_info,fourier,fourier_image,
+    status=InverseFourierTransform(&fourier_info,inverse_pixels,fourier_image,
       exception);
-  fourier=(fftw_complex *) RelinquishMagickMemory(fourier);
-  phase=(double *) RelinquishMagickMemory(phase);
-  magnitude=(double *) RelinquishMagickMemory(magnitude);
+  inverse_info=RelinquishVirtualMemory(inverse_info);
   return(status);
 }
 #endif
@@ -1201,12 +1498,12 @@ MagickExport Image *InverseFourierTransformImage(const Image *magnitude_image,
   fourier_image=(Image *) NULL;
   (void) modulus;
   (void) ThrowMagickException(exception,GetMagickModule(),
-    MissingDelegateWarning,"DelegateLibrarySupportNotBuiltIn","`%s' (FFTW)",
+    MissingDelegateWarning,"DelegateLibrarySupportNotBuiltIn","'%s' (FFTW)",
     magnitude_image->filename);
 #else
   {
     fourier_image=CloneImage(magnitude_image,magnitude_image->columns,
-      magnitude_image->rows,MagickFalse,exception);
+      magnitude_image->rows,MagickTrue,exception);
     if (fourier_image != (Image *) NULL)
       {
         MagickBooleanType
@@ -1287,7 +1584,7 @@ MagickExport Image *InverseFourierTransformImage(const Image *magnitude_image,
               thread_status;
 
             thread_status=MagickTrue;
-            if (magnitude_image->matte != MagickFalse)
+            if (magnitude_image->alpha_trait == BlendPixelTrait)
               thread_status=InverseFourierTransformChannel(magnitude_image,
                 phase_image,AlphaPixelChannel,modulus,fourier_image,exception);
             if (thread_status == MagickFalse)