0

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

umg

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

img

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:

img

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:

img

user1532868
  • 25
  • 2
  • 9
  • see [What should be the input and output for an FFT image transformation?](https://stackoverflow.com/a/26734979/2521214) what you want to do is most likely just the `wrap` but without input image we can only guess. Also your result looks like with invalid formatting... I assume your imaginary part is half of image and the rest is real part but the width does not match as you image is skewed ... – Spektre Sep 13 '17 at 06:02
  • I try to upload the input image but I cannot, because the my reputation is not enough. the input image is 256x256, I did the FT mainly though the code, fftw_plan p = fftw_plan_dft_r2c_2d(h, w, d, df, FFTW_ESTIMATE); fftw_execute(p); fftw_destroy_plan(p); – user1532868 Sep 13 '17 at 06:10
  • I have updated the description of my question, and add the link to the input image, or you download it here, https://1drv.ms/i/s!AtUbFjG0KnVcjSoCg5pMzYmP3jMN – user1532868 Sep 13 '17 at 06:30
  • When I put it into my app (link is in the answer I linked in my first comment) then I got werry similar to your result when I do 1. `DFFT` 2. `|C|` (power or absolute value) 3. `wrap`. You need to play with the slider to set desired brightness of the result. So You most likely misinterpret the format of the subresult somwhere. I do not use FFTW but check for the organization of the in/out images how exactly the `Re` and `Im` are stored and how are you using them. PS image must be bmp and you need to click on it in the file selector first to make it work ... – Spektre Sep 13 '17 at 06:39
  • As I know, the FT result in FFTW, is a two double data type (fft_complex), one for real, and another one for imaginary. it based on the row-major order, it means I guess first do the FT on row, and then column. because redundancy data, the row length of result is about half of the input image, but the column length is the same. most importantly all the data are stored in a one-dimensional array. therefore, these rows in result are adjacent, and they located on the consecutive memory. but anyway I still cannot figure it out. – user1532868 Sep 13 '17 at 06:51
  • actually, I do not need to do it myself, fftw do the FT on row and column, and give the final complex, just like a blackbox. I will receive the final DFT result, in which only the top half got data(complex), the rest are zero. Then I use the complex data to calculate the magnitude. – user1532868 Sep 13 '17 at 07:11
  • although the image data stored in 1D array, but FFTW need the image size to implement its algorithm. – user1532868 Sep 13 '17 at 07:16
  • Here got a detail information about the process.http://www.fftw.org/fftw3_doc/Multi_002dDimensional-DFTs-of-Real-Data.html#Multi_002dDimensional-DFTs-of-Real-Data – user1532868 Sep 13 '17 at 07:16
  • you should check the raw DFFT 2D result visually if it is not correct than you are feeding the data in wrong way to the FFTW . if it is correct then the problem is in the next step. – Spektre Sep 13 '17 at 07:17
  • my problem is how to do the fft shift, which move the origin (zero-frequency) to the image center. I confirm that the FT and IFT function are working. because I input my image, did the FT, followed by IFT. I can get the image exactly same as the input image. so the problem is how to shift the origin. – user1532868 Sep 13 '17 at 07:21
  • `DFFT+iDFFT` is not a valid test for this as if the intermediate result is wrong the final one will be good again (in case of format problem). You need to visualize the data after **DFFT**. Also `real array requires padding` did you do it in manner the **FFTW** wants it (but I think that is just for in place where `df==d`)? The page is written in a horrible way to understand (put a programmer to write docs and that is what you get :) ). Care about the wrapping latter ... The result should look like mine Re image without wrapping but with doubled the resolution in x axis – Spektre Sep 13 '17 at 07:26
  • :), I will keep check the data, and figure it out where is wrong, thank you for your patience. – user1532868 Sep 13 '17 at 07:43

2 Answers2

0

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.

Kaveh Vahedipour
  • 3,412
  • 1
  • 14
  • 22
  • I do not have your contact, therefore put my question is here. actually my problem is to shift the origin (zero-frequency) to image center. I have read the documentation of FFTW, I thought I was right. I have checked you code, and found your code do the FT from 1D to 3D, and can shift the origin after FT. when you shift the origin, take the 2D for example, function shift2, it mainly switch the the top-left and bottom-right, top-right and bottom-left. I understand that the low frequency locate on the corner, and should do it like this. – user1532868 Sep 13 '17 at 07:41
  • however, the problem is the output result of FFTW, at least on my input image, only have half data (on the top half), the bottom half are zeros, this can be explained on the website, http://www.fftw.org/fftw3_doc/Multi_002dDimensional-DFTs-of-Real-Data.html#Multi_002dDimensional-DFTs-of-Real-Data , I think it it a different condition. certainly there is another explanation that the output was wrong. – user1532868 Sep 13 '17 at 07:41
  • Where is your original data from? Have you tried to use a math package to get a feel for what you are trying in C++? Take a look at `numpy` and `scipy`. Do you have access to MATLAB? generally, what you want is achieved as follows: ```A=ifftshift(fft2(fftshift(A))); A=gassianfilter(A); A=ifftshift(ifft2(fftshift(A)));``` As long as your side lengths are even `fftshift===ifftshift` – Kaveh Vahedipour Sep 13 '17 at 08:52
0

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

user1532868
  • 25
  • 2
  • 9