FFT bug - lost pixels in frequency space

[brv]

With the latest ImageMagick 6.9.3 it appears that in the forward FFT transform some boundary values are incorrect. Due to this issue, applying -fft, immediatly followed by -ift, does not restore the original picture precisely, even in the HDRI mode.

This can be verified with the following command:

Code: Select all

convert -set colorspace RGB -size 16x16 xc:black -fill white -draw "point 8,8" -write "test.png" -fft \( -clone 0 -auto-level -write testx-fft-0.png +delete \) \( -clone 1 -write testx-fft-1.png +delete \) -ift testx.png

This command creates a 16x16 black image with a single white pixel near the center. DFT of such an image should be a single complex exponent:
(1/N) exp( 2 pi i (u x0+ v y0) ),
where u and v are variable frequencies, and x0,y0 are image coordinates of the white pixel.
This means that the absolute value of this DFT must always be equal to unit, while the phase should drift linearly (modulo 2pi).

In the output of the command above we obtain different pictures however, namely there are bogus columns of pixels at the left and right sides (i=0, i=width-1) that are always filled with zero amplitude and phase (for any source image, actually). And there is an additional bogus pixel at the coordinates i=1,j=0. This causes image artifacts after we apply inverse DFT, so that after applying "-fft -ift" the final image is not the same as the original (this effect, however, should decrease for larger images). Same problems remain for +fft +ift flavors, and do not depend on HDRI or discrete Q16.

I digged the source code and the errors appear in the ForwardQuadrantSwap function in fourier.c.
Here I comment out the buggy code lines replacing them with what I believe is right:

Code: Select all

 /*
    Swap quadrants.
  */
  center=(ssize_t) floor((double) width/2L)+1L;
  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++)
    for (x=0L; x < (ssize_t) (width/2L); 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++)
    for (x=0L; x < (ssize_t) (width/2L); x++)
      forward_pixels[(height-y)*width+width/2L-x-1L]=
        source_pixels[y*center+x+1L];
//  for (x=0L; x < (ssize_t) (width/2L-1L); x++)
//    forward_pixels[-x+width/2L-1L]=forward_pixels[x+width/2L+1L];
  for (x=0L; x < (ssize_t) (width/2L); x++)
    forward_pixels[width/2L-x-1L]=
      source_pixels[x+1L];
  return(MagickTrue);

As can be seen, pixels columns at i=0 and i=width-1 are not copied from the FFT buffer and a similar problem is with the j=0 row.

With this fix, the behaviour of the fft tool becomes as expected, and the original image is reconstructed accurately after "-fft -ift". This can be checked with the same command above.

Developers who are responsible for FFT in imagemagick, please pay attention to this. Hope to find fixed source in git

[magick]

Thanks for the problem report. We can reproduce it and will have a patch to fix it in GIT master branch @ https://github.com/ImageMagick/ImageMagick later today. The patch will be available in the beta releases of ImageMagick @ http://www.imagemagick.org/download/beta/ by sometime tomorrow.

[brv]

So quick, great, thanks!