Scaling performed by gpu::dft with OpenCV in C++

Question

I want to use a GPU-accelerated algorithm, to perform a fast and memory saving dft. But, when I perform the gpu::dft, the destination matrix is scaled as it is explained in the documentation. How I can avoid this problem with the scaling of the width to dft_size.width / 2 + 1? Also, why is it scaled like this? My Code for the DFT is this:

cv::gpu::GpuMat d_in, d_out;
d_in = in;
d_out.create(d_in.size(), CV_32FC2 );
cv::gpu::dft( d_in, d_out, d_in.Size );

where in is a CV_32FC1 matrix, which is 512x512.

The best solution would be a destination matrix which has the size d_in.size and the type CV_32FC2.

Answer 1

This is due to complex conjugate symmetry that is present in the output of an FFT. Intel IPP has a good description of this packing (the same packing is used by OpenCV). The OpenCV dft function also describes this packing.

So, from the gpu::dft documentation we have:

If the source matrix is complex and the output is not specified as real, the destination matrix is complex and has the dft_size size and CV_32FC2 type.

So, make sure you pass a complex matrix to the gpu::dft function if you don't want it to be packed. You will need to set the second channel to all zeros:

Mat realData;

// ... get your real data...

Mat cplxData = Mat::zeros(realData.size(), realData.type());

vector<Mat> channels;
channels.push_back(realData);
channels.push_back(cplxData);

Mat fftInput;
merge(channels, fftInput);

GpuMat fftGpu(fftInput.size(), fftInput.type());
fftGpu.upload(fftInput);

// do the gpu::dft here...

There is a caveat though...you get about a 30-40% performance boost when using CCS packed data, so you will lose some performance by using the full-complex output.

Hope that helps!

Answer 2

Scaling is done for obtaining the result within the range of +/- 1.0. This is the most useful form for most applications that need to deal with frequency representation of the data. For retrieving a result which is not scaled just don't enable the DFT_SCALE flag.

Edit

The width of the result is scaled, because it is symmetric. So all you have to do is append the former values in a symmetric fashion.

The spectrum is symmetric, because at half of the width the sampling theorem is fulfilled. For example a 2048 point DFT for a signal source with a samplerate of 48 kHz can only represent values up to 24 kHz and this value is represented at half of the width.

Also for reference take a look at Spectrum Analysis Using the Discrete Fourier Transform.