ARM DSPlib CFFT vs. RFFT - whats the deal?

[Yansi]

Hello!

I got to toy around a bit with FFT on ARM, using the CMSIS DSP library. However, I am quite confused how to understand it, as the documentation from ARM is beyond horrible.

Supposed I have a sequence of sampled data (i.e. a real numbers only sequence) and want to look at the signal spectrum, so only say the magnitude is relevant.

I see there are two pools of DSP functions: CFFT (complex fft) and RFFT (supposedly real fft, meaning real input data, complex output data).

I understand that a complex FFT is the more general one, probably requiring more computational power.

But when looking at for example arm_rfft_fast_f32 function, what I see inside, is that it uses the CFFT function (arm_cfft_f32) anyway and then uses some other function (stage_rfft_f32) of unknown purpose to somehow bastardize the output data.

Doesn't that mean the rfft functions will be even slower than the cfft? Wtf! I thought that the real fft should be faster due to half the input data.  But no, to use the rfft function, I still need to stupidly input a complex data vector, with zero imaginary values, which I would just do using the complex fft function anyway.

So what is the purpose of the rfft functions? Have I missed something?

Thanks for any hints,
Y.

[chris_leyson]

You can compute a Real FFT or IFFT of length N with a length N/2 complex FFT. From Wikipedia "Alternatively, it is possible to express an even-length real-input DFT as a complex DFT of half the length (whose real and imaginary parts are the even/odd elements of the original real data), followed by O(N) post-processing operations." There are a lot of neat tricks you can pull off because of symetry.
https://en.wikipedia.org/wiki/Fast_Fourier_transform

[Berni]

Its just how FFT works. It takes complex input and makes complex output.

Even if you feed it a purely real input this soon turns into numbers with both complex and imaginary parts that need to both be dragged on trough the rest of the process. So a purely real input doesn't save half the computation time (Tho some extra optimization is possible). The other trick is to pack a second set of real data into the empty imaginary data, allowing you to run both trough a single FFT operation and then use some extra math to separate them back out.

I'm guessing they figured that sort of optimization is not worth it and just provided a function that feeds your data to the complex one.

There are plenty of oddities. For example i found that on my STM32H7 the floating point FFT is always faster than a integer FFT. Still not sure exactly why but nothing i tried could get the integer one to run as fast as float (perhaps this kept the internet math part of the cpu free to do loop iterations, or extra scaling used to fit the data in has taken too much time to compute), let alone faster. Overall works pretty well tho.

[Yansi]

I have found on many places they mention the even symmetry of the spectrum.  But that is not enough to my understanding, unfortunately.

Taking what chris_leyson quoted from the wikipedia, it seems that okay,  one can cram the real series numbers within a complex series of half length. Okay, but then what? What is the O(N) post-processing required?

From the CMSIS DSP help it is absolutely of no telling how data shall be formatted for input or output of many of the functions.

[Yansi]

So after some more painful digging within the source code mess, I have finally decoded that:

1) The RFFT function uses inside half length CFFT call.
2) The RFFT instance structure contains a pointer (pCfft) to a "internal CFFT" structure as well, that gets initialized using the common arm_rfft_init_xxx call.
3) there are errors in the documentation regarding how to initialize these instances manually, as NO property fftLenBy2 exists within that structure:

Code: [Select]

Use of the initialization function is optional.
                   However, if the initialization function is used, then the instance structure
                   cannot be placed into a const data section. To place an instance structure
                   into a const data section, the instance structure should be manually
                   initialized as follows:
  <pre>
      arm_rfft_instance_q31 S = {fftLenReal, fftLenBy2, ifftFlagR, bitReverseFlagR, twidCoefRModifier, pTwiddleAReal, pTwiddleBReal, pCfft};
      arm_rfft_instance_q15 S = {fftLenReal, fftLenBy2, ifftFlagR, bitReverseFlagR, twidCoefRModifier, pTwiddleAReal, pTwiddleBReal, pCfft};
Now I understand it a bit more.

How one extracts the output spectrum from thehalf length CFFT fed with just real number series?

[voltsandjolts]

To be fair, all code of a reasonable size looks messy when you are not familiar with the underlying algorithm(s) it implements.

[chris_leyson]

"Handbook of Real-Time Fast Fourier Transforms" by Smith & Smith explains the double length algorithm in section 2.4.2
https://doc.lagout.org/science/0_Computer%20Science/2_Algorithms/Handbook%20of%20Real-Time%20Fast%20Fourier%20Transforms_%20Algorithms%20to%20Product%20Testing%20%5BSmith%20%26%20Smith%201995-05-22%5D.pdf

[Yansi]

Thank you, that book seems it'll be very handy. I am yet to make myself understand those notations there, will look into that.

Now I should probably try running some code to verify, if I understood correctly what is happening there and where. 

To be fair, all code of a reasonable size looks messy when you are not familiar with the underlying algorithm(s) it implements.

Yes, that is a fair point. However that is what the documentation should be for - to familiarize yourself with what the code does, without trying to reverse engineer all of it!

[chris_leyson]

Handbook of Real-Time Fast Fourier Transforms by Winthrop W. Smith. & Joanne M. Smith was one book that sprang to mind because they run through all of the steps required. Despite being some 25 years old it's still a handy reference book.
Forgot to mention: For Real data there is also the Hartley transform which doesn't require any complex multiplications and the same transform can be used to get the Inverse Hartley transform. Ralph V. L. Hartley of Hartley oscillator fame.