FFT4 discrepancy (DSP-139)

[andg1]

Answers checklist.

• I have read the ESP-DSP documentation and the issue is not addressed there.
• I have read the documentation ESP-IDF Programming Guide and the issue is not addressed there.
• I have updated my ESP-DSP branch (master or release) to the latest version and checked that the issue is present there.
• I have searched the issue tracker for a similar issue and not found a similar issue.

General issue report

Testing FFT4 example but with real adc data is generating incorrect output of the fft. Comparing to the FFT example the radix-2 input vector is populated in a completely different way and it's not well documented in esp-dsp documentation, also there are some comments simply copy-pasted that are wrong and "input/output array with size of N2. An elements located: Re[0],Re[1], , … Re[N-1], any data… up to N2 result of DCT will be stored to this array from 0…N-1. Size of data array must be N*2!!! " suggest that the vector in radix4 should be half populated with real and half img, if so why use x1[i] = 10 * log10f((x1[i * 2 + 0] * x1[i * 2 + 0] + x1[i * 2 + 1] * x1[i * 2 + 1] + 0.0000001) / N); copied from FFT example?. To get the fft working properly I had to use FFT example (so using radix2) and populate the array using format Re[0], Im[0], … Re[N-1], Im[N-1]

[dmitry1945]

@andg1 ,
thank you for the issue.

The FFT radix-2 and FFT radix-4 produce the same output in case if the length is the same.
The FFT R2 and FFT R4 have different initialization functions. To use R4 you should use dsps_fft4r_init_fc32() if you use FFT R2 you should use dsps_fft2r_init_fc32().
Do you use the dsps_fft4r_init_fc32() for FFT-R4?
If yes, could you please provide more details about the error, that I will be able to reproduce the error?

Thank you,
Dmitry

[andg1]

Thank you for the reply!
Yes I'm using radix4 with correct initialization. In "FFT" example the array is populated with Re[0],Im0, .... Re[N-1], Im[N-1] while in "FFT4" even the radix2 array seems to be populated as Re[0],...,Re[N-1],Im[N],...Im[2N-1] which seems wrong watching the documentation and, in fact, i think it is because the output of the "FFT4" is wrong. Even further if the array is populated as Re[0],...,Re[N-1],Im[N],...Im[2N-1] how can x1[i] = 10 * log10f((x1[i * 2 + 0] * x1[i * 2 + 0] + x1[i * 2 + 1] * x1[i * 2 + 1] + 0.0000001) / N); give you correct results? I'm populating the array with data from adc sampling mains 50Hz sine and using the model from "FFT4" gives mostly correct results but not really as using the format from "FFT" example.

Also from
-dsps_fft4r_fc32_ansi_ document description "[inout] data: input/output complex array. An elements located: Re[0], Im[0], … Re[N-1], Im[N-1] result of FFT will be stored to this array. "
-dsps_bit_rev4r_fc32 says "[inout] data: input/ complex array. An elements located: Re[0], Im[0], … Re[N-1], Im[N-1] result of FFT will be stored to this array. "
but after after
-dsps_cplx2real_fc32_ansi_ says "inout] data: Input complex array and result of FFT2R/FFT4R. input has size of 2N, because contains real and imaginary part. result will be stored to the same array. Input1: input[0..N-1], Input2: input[N..2N-1] " that to me seems suggesting the Re[0],...,Re[N-1],Im[N],...Im[2*N-1] above

So if the first 2 input/output array definition is the correct how the array is populated in "FFT4" example is wrong if I understood the documentation correctly

[dmitry1945]

@andg1
now I got the point.
The description for input of dsps_cplx2real_fc32_ansi() has mistake, as you told, it's copy-paste from other function.
Sure, it should be :

 * @brief      Convert FFT result to complex array for real input data
 *
 * Convert FFT result of complex FFT for real input to complex output.
 * This function have to be used if FFT used to process real data.
 * This function use tabels inside and can be used only it dsps_fft4r_init_fc32(...) was
 * called and FFT4 was initialized.
 * The implementation use ANSI C and could be compiled and run on any platform

 * @param[inout] data: Input complex array and result of FFT2R/FFT4R.
 *               input has size of 2*N, because contains real and imaginary part.
 *               result will be stored to the same array.
 *               Input1: input[0..N-1] if the result is complex Re[0], Im[0]....Re[N-1], Im[N-1], 
 *              and input[0...2*n-1] if result is real Re[0], Re[1],...,Re[2*N-1].

Or something like this.

We will fix it soon in next release.

Thank you very much!
Dmitry

[andg1]

So if I'm correct even the way the array is populated in FFT4 is wrong. Instead of
for (int i = 0 ; i < N ; i++) {
x1[i] = x1[i] * wind[i];
x2[i] = x1[i];
}
does not respect the format Re[0], Im[0], .... should be something similar to "FFT" example
for (int i = 0 ; i < N ; i++) {
y_cf[i * 2 + 0] = x1[i] * wind[i];
y_cf[i * 2 + 1] = 0;
}
So the Re part is data from tone generation and Im=0
Is it correct?

Thanks in advance

[dmitry1945]

@andg1
In the dsps_fft4real_main.c example the input is real signal.

real[0...N-1] -> FFT-> dsps_bit_rev4r_fc32() -> Cplx[0...N/2-1] -> dsps_cplx2real_fc32() -> Cplx[0...N/2-1] = complex spectrum of real input signal.

The FFT (R2 and R4) is always working with complex input. But, if you have real input signal that is located in one real array, like float[512], then you can imagine, that you have a complex[256] array, and you will have odd real samples as Re[...] input data, and even input samples as Im[...] input data.
We pass this complex array to FFT, and then we have to make additional operation dsps_cplx2real_fc32_ansi() that will convert result complex array length N/2 Cplx[0...N/2-1] to half of complex array Cplx[0...N-1], that represent input signal Re[0]=In[0], Re[1] = In[1]...Re[N-1] = In[N-1 and Im[0...N-1] = 0];

In simple words,
if you want to make FFT for input real array with N samples, you will need to make complex array with N complex samples where Re[0...N-1] = In[0...N-1], Im[0..N-1] = 0.
The result will be complex spectrum F[0...N-1] where F[0...N/2-1] = F*[N/2...N-1]
For this method, for N input real samples you will use FFT for N samples.

Or...
You can make FFT N/2 points for input real signal of N samples, and make additional operation Cplx2real. This will give you first half of complex spectrum F[0...N/2], and second half you can make as F*[0..N/2]
For this method, for N input real samples you will use FFT for N/2 samples.

Regards,
Dmitry