GNOME Bugzilla – Bug 492098
[GstFFT] Broken scaling
Last modified: 2007-11-06 11:57:49 UTC
currently GstFFT somehow scales the data given to kiss the wrong way (i.e. not at all). To make it possible to call a fft and on the same data an inverse fft and get the same data again some scaling is necessary, I have that fixed for the floating point data types locally but I'm still fighting with the fixed point implementations.
The missing scaling also resulted in wrong magnitude being calculated.
Ok, after some investigation it seems that the following might be a good solution:
a) Document that behaviour of _fft() and _inverse_fft(), i.e. for floating point that iFFT(FFT(x)) = x*nfft and for fixed point iFFT(FFT(X)) = x/nfft and that one has to apply the required scaling to get iFFT(FFT(x)) = x. For fixed point this will be rather hard though... think of overflows, etc :)
b) Fix the magnitude calculation or what I would currently prefer, drop the magnitude and phase calculation functions. Dropping them seems to make sense as it's a rather special use case that requires them and even the spectrum element has it's own implementation for them as it has special needs. When dropping them it might make sense to add an example on how to calculate the magnitude in the docs though...
What do you think?
Ok, now a bit better explaination ;)
First of all, take a look at http://en.wikipedia.org/wiki/Discrete_Fourier_transform
For floating point this can be taken more or less as a reference...
As can be seen there, you'll need a 1/nfft somewhere to get iFFT(FFT(x))==x behaviour. Our current implementation of the FFT does exactly that for floating point, so to a) get the original output back you have to either divide the FFT output by nfft or the iFFT output (as the FFT is linear it doesn't matter). b), to calculate the magnitude and get the gain for different frequency bins from that you also have to divide the output of the FFT by nfft to get proper results.
For fixed point things are unfortunately a bit more complicated as overflows can happen (sure, they can happen to floating point too but that's rather unlikely here). To quote kissfft upstream:
> Consider the case of a radix-2 fft. At every stage of the FFT there is
> an addition operation that causes word growth of 1 bit.
> So at each stage, there is a division by two to prevent the possibility
> of overflow. With log2(Nfft) stages, there is an overall division by
> 2^log2(Nfft) == Nfft
So that's why you get iFFT(FFT(x)) = x/nfft. For getting the original input again you have to take special care with fixed point samples for this reasons. For example you could do the following:
- Analyse the time domain data and get the highest integer C that won't give a overflow if all samples are multiplied with it.
- Multiply them all with C.
- get FFT(x) = X.
- get iFFT(X) = C*x / nfft and multiply by nfft, divide by C.
This will always give you the original data back as accurate as possible.
So my proposal for 0.10.15 is to simply document this behaviour (it can't get much better for the above mentioned reasons) and remove the magnitude and phase calculation functions as they're pretty useless even for the spectrum element.
> So my proposal for 0.10.15 is to simply document this behaviour (it can't get
> much better for the above mentioned reasons) and remove the magnitude and phase
> calculation functions as they're pretty useless even for the spectrum element.
Let's do that then.
2007-11-06 Sebastian Dröge <firstname.lastname@example.org>
* tests/check/libs/fft.c: (GST_START_TEST):
Remove the magnitude and phase calculation functions as these have very special use cases and can't even be used for the spectrum element. Also adjust the docs to mention some properties of the used FFT implemention, i.e. how the values are scaled. Fixes #492098.