4

I aim to get the DFT of an image in OpenCV.

Using dft function, I'm able to calculate it, and then paint it by calculating its magnitude (then, apply the log and finally normalize it in order to paint values between 0 and 1).

My result is, for the following image, the result I show you (with swap in order to have lower frequencies in the center of the image):

Original picture DFT magnitude

However, if I compare it to the result I obtain using other tools like Halcon, It seems incorrect to my since It seems to have really "high" values (the OpenCV DFT magnitude I mean):

FFT Halcon

I thought it might be for these reasons:

  1. The difference between DFT (at OpenCV) and FFT (Halcon)
  2. The operations I'm performing in order to show the magnitude in OpenCV.

The first one have as problem that it's quite hard for me to analyze, and OpenCV doesn't have a FFT function, as well as Halcon doesn't have a DFT function (if I'm not wrong of course), so I can't compare it directly.

The second one is in which I've been working the most time, but I still don't find the reason if it's there.

There's the code I'm using to paint the magnitude of img (which is my DFT image):

// 1.- To split the image in Re | Im values
Mat planes[] = {Mat_<float>(img), Mat::zeros(img.size(), CV_32F)};

// 2.- To magnitude + phase
split(img, planes);

// Calculate magnitude. I overwrite it, I know, but this is inside a function so it will be never used again, doesn't matter
magnitude(planes[0], planes[1], planes[0]);

// Magnitude Mat
Mat magI = planes[0];

// 3.- We add 1 to all them in order to perform the log
magI += Scalar::all(1);                    // switch to logarithmic scale
log(magI, magI);

// 4.- Swap the quadrants to center frequency
magI = magI(Rect(0, 0, magI.cols & -2, magI.rows & -2));
int cx = magI.cols/2;
int cy = magI.rows/2;

Mat q0(magI, Rect(0, 0, cx, cy));   // Top-Left - Create a ROI per quadrant
Mat q1(magI, Rect(cx, 0, cx, cy));  // Top-Right
Mat q2(magI, Rect(0, cy, cx, cy));  // Bottom-Left
Mat q3(magI, Rect(cx, cy, cx, cy)); // Bottom-Right

// swap quadrants (Top-Left with Bottom-Right)
Mat tmp;
q0.copyTo(tmp);
q3.copyTo(q0);
tmp.copyTo(q3);

// swap quadrant (Top-Right with Bottom-Left)
q1.copyTo(tmp);
q2.copyTo(q1);
tmp.copyTo(q2);

// 5.- Normalize
// Transform the matrix with float values into a
// viewable image form (float between values 0 and 1).
normalize(magI, magI, 0, 1, CV_MINMAX); 

// Paint it
imshow( "Magnitud DFT", magI);

So summarizing: any idea about why do I have this difference between these two magnitudes?

Btc Sources
  • 1,912
  • 2
  • 30
  • 58
  • A guess: The factor is `N` or `N^2`, the number of elements in the _FT. Am I correct? – Avi Ginsburg Jun 01 '15 at 10:55
  • 1
    FFT (Fast Fourier Transform) is just a quick method to compute DFT (Discrete Fourier Transform). The results should be equal up to a small numerical error. As already pointed in the previous comment, one of the algorithms is probably normalising the values and the other one is not. – Hristo Iliev Jun 01 '15 at 10:58
  • I guess Halcon is normalising the result @HristoIliev, but I'm doing it too when painting it in OpenCV (`normalize(magI, magI, 0, 1, CV_MINMAX);`), so shouldn't them be the same? – Btc Sources Jun 01 '15 at 11:00
  • I'm not using a factor when normalising the result in `OpenCV` @AviGinsburg, just the `normalize` function... – Btc Sources Jun 01 '15 at 11:02
  • I have zero experience with OpenCV and Halcon (never even heard of the latter). I'm just noting that many libraries don't normalize the Fourier transform. We're talking about normalizing the transform, not the same as normalizing a vector/matrix. – Avi Ginsburg Jun 01 '15 at 11:05
  • That makes sense, could you tell me how to normalize the complete transform? I mean, should I just sum the `Re` part with the `Im` part of each term (*sqrt* of each one, since It's the norm factor Halcon is using ), and then normalize as I'm doing now? – Btc Sources Jun 01 '15 at 11:16
  • No. That's _not_ the normalization factor I'm talking about. See [this](http://www.mathworks.com/matlabcentral/answers/15770-scaling-the-fft-and-the-ifft) for an example. If `ifft(fft(f(x))) != 1` then there you'll want to divide by that factor to the power of something between 0 and 1. – Avi Ginsburg Jun 01 '15 at 11:31
  • I've readed it @AviGinsburg, and tried with a `1/sqrt(M)` factor. The result seems quite more near to what I was expecting! Thanks! – Btc Sources Jun 01 '15 at 11:49

1 Answers1

5

I'll summarize my comments into an answer.

When one thinks of doing a Fourier transform to work in the inverse domain, the assumption is that doing the inverse transform will return the same function/vector/whatever. In other words, we assume

Inverse transform of a transform yields the original function

This is the case with many programs and libraries (e.g. Mathematica, Matlab/octave, Eigen/unsupported/FFT, etc.). However, with many libraries (FFTW, KissFFT, etc.) this is not the case and there tends to be a scale

Inverse transform of a transform yields a scaled original function

where s is usually the number of elements (m) in the array to the power of something (should be 1 if not scaled in a mismatched fashion in both the transform and the inverse). This is done in order to refrain from iterating over all m elements multiplying by a scale, which is often not important.

That being said, when looking at the scale in the inverse domain, various libraries that do scale the transforms have the liberty to use different scales for the transform and inverse transform. Common scaling pairs for the transform/inverse include {m^-1, m} and {m^-0.5, m^0.5}. Therefore, when comparing results from different libraries, we should be prepared to factors of m (scaled by m^-1 vs. not scaled), m^0.5 (scaled by m^-0.5 vs. not scaled and scaled by m^-1 vs. scaled by m^-0.5) or even other scales if other scaling factors were used.

Note: This scaling factor is not related to normalizing an array, such that all values are [0,1] or that the norm of the array is equal to 1.

Avi Ginsburg
  • 10,323
  • 3
  • 29
  • 56
  • Thanks for the explanation @Avi-Ginsburg , my problem was solved by using a `m^-1` factor instead of normalize the magnitude in OpenCV. – Btc Sources Jun 02 '15 at 09:24