4

I'm playing with some toy code, to try to verify that I understand how discrete fourier transforms work in OpenCV. I've found a rather perplexing case, and I believe the reason is that the flags I'm calling cv::dft() with, are incorrect.

I start with a 1-dimensional array of real-valued (e.g. audio) samples. (Stored in a cv::Mat as a column.)

I use cv::dft() to get a complex-valued array of fourier buckets.

I use cv::dft(), with cv::DFT_INVERSE, to convert it back.

I do this several times, printing the results. The results seem to be the correct shape but the wrong magnitude.

Code:

cv::Mat samples(1, 2, CV_64F);
samples.at<double>(0, 0) = -1;
samples.at<double>(0, 1) = 1;
std::cout << "samples(" << samples.type() << "):" << samples << std::endl;
for (int i = 0; i < 3; ++i) {
  cv::Mat buckets;
  cv::dft(samples, buckets, cv::DFT_COMPLEX_OUTPUT);
  samples = cv::Mat();
  cv::dft(buckets, samples,
          cv::DFT_INVERSE | cv::DFT_COMPLEX_INPUT | cv::DFT_REAL_OUTPUT);

  std::cout << "buckets(" << buckets.type() << "):" << buckets << std::endl;
  std::cout << "samples(" << samples.type() << "):" << samples << std::endl;
}

Output:

samples(6):[-1, 1]
buckets(14):[0, 0, -2, 0]
samples(6):[-2, 2]
buckets(14):[0, 0, -4, 0]
samples(6):[-4, 4]
buckets(14):[0, 0, -8, 0]
samples(6):[-8, 8]

I would have expected the above output to repeat. E.g. [-1, 1], [0, 0, -1, 0], .... Instead, the magnitude is doubling on each round-trip.

Is my understanding wrong? Or am I using the wrong flags? Etc.

kd8azz
  • 613
  • 5
  • 21
  • DFT_inverse(DFT) often multiplies your result by the length of the array. This comes down to the exact definitions used of the DFT. You schould pay attention to the factor in front which can either be N, sqrt(N) or 1 – Unlikus Feb 13 '22 at 21:35
  • When you look at the formulas on https://de.wikipedia.org/wiki/Diskrete_Fourier-Transformation You see that the inverse transformation has the factor 1/N before. This would make Inverse_DFT(DFT) = id. Your inverse transformation probably does not divide by N, so you get N times of what you expected. – Unlikus Feb 13 '22 at 21:46
  • 1
    You can add the `DFT_SCALE` flag to make the inverse transformation divide by N – Unlikus Feb 13 '22 at 21:50
  • If you copy your comment to an answer, I'll accept it, so you can get internet points. :) – kd8azz Feb 13 '22 at 21:51

1 Answers1

3

The inverse DFT in opencv will not scale the result by default, so you get your input times the length of the array. This is a common optimization, because the scaling is not always needed and the most efficient algorithms for the inverse DFT just use the forward DFT which does not produce the scaling. You can solve this by adding the cv::DFT_SCALE flag to your inverse DFT.

Some libraries scale both forward and backward transformation with 1/sqrt(N), so it is often useful to check the documentation (or write quick test code) when working with Fourier Transformations.

Unlikus
  • 1,419
  • 10
  • 24