fftw shift zero-frequency to image center

I am struggling on the result of FFTW applied on image.

I can do the Fourier transform and inverse Fourier transform on image, I am sure the two functions are working. Now I am trying to analyze the result of Fourier transform and try to apply a filter on it, etc.

Currently I am confusing on the FT result. As I know the the zero-frequency should be on the top-left corner. I should shift it to the image center if I want to check frequency friendly. There are at least two way I have tested. one method is described on the FFTW FAQ webpage, all the pixels in the image were multiplied by (-1)^(i + j) before Fourier transform, where i, and j are the pixel's index.

image was multiplied by (-1)^(i+j), then its magnitude is

but actually ,it should be like this (which I did in imagej) Fourier transform in imagej:

I also tried using my code to switch them after calculating its magnitude, but got the similar wrong result.

void ftShift2D(double * d, int w, int h)
{

    int nw = w/2;
    int nh = h/2;
    if ( w%2 != 0) {nw++;}
    if ( h%2 != 0) {nh++;}

        printf("%d, %d\n", nw, nh);
        for (int i = 0; i < nh ; i++)
        {
                for(int j = 0; j < nw; j++)
                {
                        int t1 = i * w + j;
                        int t2 = (i + nh + 1) * w + j + nw;
                        swapXY<double>(d[t1], d[t2]);
                }

                for (int k = nw; k < w; k++)
                {
                        int t1 = i  * w + k;
                        int t2 = (i + nh　＋ 1 ) * w + k - nw;
                        swapXY<double>(d[t1], d[t2]);
                }
        }
}

I understood the layout of FT result in FFTW, as is described in their website. I also checked the output, this layout can be confirmed, only top half have data, the rest is 0. therefore, it seems explain why the two methods used above cannot work. I am trying to design new method to shift the zero-frequency to image center. could you give me a better way to do it?

I insist on this, aim to use the filter, it will be easier to understand. if the zero-frequency did not be shifted to center, could you guys give me another suggestion that I can, for example, apply a high Gaussian filter on it in frequency domain.

thanks a lot.

the input image is:

The FT and IFT code are

    void Imgr2c(double * d, fftw_complex * df, int w, int h)
{
        fftw_plan p = fftw_plan_dft_r2c_2d(h, w, d, df, FFTW_ESTIMATE);
        fftw_execute(p);

        fftw_destroy_plan(p);

}


void Imgc2r(fftw_complex * in, double * ftReverse, int w, int h)
{
        fftw_plan p = fftw_plan_dft_c2r_2d(h, w, in, ftReverse, FFTW_ESTIMATE);
        fftw_execute(p);

        int l = w * h;
        for (int i = 0; i < l; i++)
        {
                ftReverse[i] /= l;
        }

        fftw_destroy_plan(p);
}

[Edit by Spektre]

This is how it should more or less look like:

回答1:

You may want to check out how I did it here: https://github.com/kvahed/codeare/blob/master/src/matrix/ft/DFT.hpp The functions doing the job are at the bottom. If there are any issues, feel free to contact me personally.

回答2:

I found what is my problem! I did not understand the output layout of DFT in FFTW deeply. After some tests, I found the the dimension in width should w/2 +　1, and the dimension in height did not changed. the code to do FT, IFT, and shift the origin have no problem. what I need to do is changing the dimension when try to access it.

magnitude after DFT

magnitude after shifting the origin Cheers,

Yaowang