MKL FFT - Real input data to comples result

[Marek_C_]

Hello everyone,

I am currently working on doing fft on my wav audio signal. My goal is to recreate Matlab's function called "spectrogram", which does a STFT with moving Hammington window.

My first attempt was to represent my input data as a vector (size of the window, padded with zeros) and output is vector with same size. It was working nice, I got almost the same result as other function in Matlab using fft.

But since Matlab spectrogram is giving me results represented as complex numbers, I want to do the same with MKL using C.

In user guide manual - page 2708, I can perform the fft by using this procedure:

status= DftiComputeForward(desc_handle, xre_in, xim_in, yre_out, yim_out);

So my thinking was to declare xre_in, xim_in, yre_out, yim_out ad arrays of doubles,where:  xre_in is vector of my input data, xim_in is vector with zeros.

FFT part looks like this:

status = DftiCreateDescriptor(&fft_handle, DFTI_DOUBLE, DFTI_REAL, 1, fftpoint);
status = DftiSetValue(fft_handle, DFTI_PLACEMENT, DFTI_NOT_INPLACE);        
status = DftiCommitDescriptor(fft_handle);
status = DftiComputeForward(fft_handle, Fft_input.pVector, input_imag.pVector, Fft_output.pVector, output_imag.pVector);
status = DftiFreeDescriptor(&fft_handle);

Unfortunately my outputs are all zeros after computing, so I guess my process of thinking was wrong. Where have I done mistake?

I am working with Visual Studio 2010, I don't think I have C99 standard, because I cannot declare complex numbers, and dealing with real ones is easier.

Regards,

Marek

[Evgueni_P_Intel]

Hi Marek,

Please (1) change DFTI_REAL to DFTI_COMPLEX in DftiCreateDescriptor and (2) add the following line before DftiCommitDescriptor.

status = DftiSetValue(fft_handle, DFTI_COMPLEX_STORAGE, DFTI_REAL_REAL);

If you care about time taken by calls to MKL FFT routines, please look at examples/dftc/source/basic_dp_real_dft_1d.c (I assume that you use MKL 11) which explains how to compute only the first half of the output. The output will be stored in one array of complex numbers, real and imaginary parts interleaved. The second half of the output can be recovered as follows out=conj(out[fftpoint-k]) (k>=1) and out[0]=conj(out[0]), since FFT of real signal has conjugate-even symmetry. Please see the MKL Reference Manual for more detail.

Thanks,

Evgueni.

[Marek_C_]

Thank you, I was able to somehow solve my problem, by not using my idea I presented here.

However I have one question. Why does MKL FFT do not produce symmetric fft?

Matlab for fftpoint = 4096 creates symmetric fft, but MKL only first half? (2048 samples are result, the rest is mirror image)

One data point in Matlab is represented by two after MKL's fft. Is it possible to change this somewhere?

Regards

Marek

[Evgueni_P_Intel]

Computing only the first half of the output is the documented behavior of MKL real-to-complex FFTs. It is the applications using MKL real-to-complex FFTs that are expected to restore the second half of the output. For example, you may add to your application a small inline function or a macro to this end.

[LUVp]

How can I add a small inline function or macro to restore the second half of the output?