Prev: multiple clock design
Next: THD without using an FFT
From: AverySimonsen on 11 Jul 2010 17:30
>bump . . .anyone have some new light to shed?
>
>Thanks
>Carmen

I am aware that this is a really old thread, but seeing as Carmen came back
years after the original posting, I though I would just shed what little
light I can on the matter in case somebody else is still asking the same
question.

I was stuck the same way Carmen was using the exocortex library and wishing
I had a recognizable fft output to look at. Andor's reply put me on the
right track. I did learn about FFT at uni, but that's years ago and I have
barely thought about it since much less used it in practice, so my approach
is highly results-oriented (as in I didn't strive to understand it fully,
just to get the results I'm after).

Keeping this in mind, the following C# method may provide the sought for
output:

//See www.dspguide.com (chapter 8) for a nice explanation of the
concepts that come into play when one needs to do spectral analysis
public static void ComputeFFTPolarMag(float[] samples, float[] fft,
out float dcComponent)
{
int nextPowerOf2 = RoundToNextPowerOf2(samples.Length);
if (fft.Length < nextPowerOf2)
throw new Exception("length of fft result array must be big
enough for next power of two");
Array.Copy(fft, samples, fft.Length);//Copy FFT input into the
array that will hold the output so we don't destroy the input samples
Fourier.RFFT(fft, nextPowerOf2, FourierDirection.Forward);
for (int i = 2; i < nextPowerOf2; i++) fft[i] /= nextPowerOf2 /
2;//Divide all bins by N/2 to obtain amplitude of respective sinusoids
(http://www.dspguide.com/ch8/5.htm)
fft[0] /= nextPowerOf2;//Special treatment of DC component
//Convert from rectangular coordinates to just the magnitude of
polar coordinates - disregarding the phase 'coordinate'
//(http://www.dsprelated.com/showmessage/58198/1.php)
dcComponent = fft[0];
for (int i = 1; i < nextPowerOf2 / 2; i++)//Skip the first pair
as first bin is at index 2 and 3
fft[i - 1] = (float)Math.Sqrt(fft[i * 2] * fft[i * 2] +
fft[i * 2 + 1] * fft[i * 2 + 1]);
for (int i = nextPowerOf2 / 2 - 1; i < fft.Length; i++)//Zero
the rest
fft[i] = 0.0f;
//The fft array now holds plottable fft results from index 0 to
N/2
}

I thought I would include this utility method as it is quite clever. I
can't take credit, though. I don't remember where I dug it up.

public static int RoundToNextPowerOf2(int x)
{
if (x <= 1) return 1;
x--;
x |= x >> 1;
x |= x >> 2;
x |= x >> 4;
x |= x >> 8;
x |= x >> 16;
x++;
return x;
}

Please note that the above source is a quick and untested adaptation of a
larger method specific to my application intended to provide only what is
requested here and nothing more.

Hope to help.
Lars Ole

>
>>Hello
>>
>>Is there anyone that can help me.
>>
>>I don't care what the FFT does in the background, all I'm asking is how
>to
>>impliment it using a 1D array of audio samples.
>>
>>I am using Exocortex.DSP in C#.
>>
>>If someone could respond, refering to elements of C# (array number 1
>etc,
>>instead of A_k etc) that would be most helpful.
>>
>>Thanks in advance
>>
>>Carmen Branje
>>
>>>cbranje(a)gmail.com wrote:
>>>> Hi There
>>>>
>>>> I'm fairly new to FFT and the like and I'm trying to implement the
>>>> Exocortex.DSP for C# and I have it working but I'm not really sure
>how
>>>> to interpret to results?
>>>>
>>>> So I pass a 1024 sized array of floats (they are 16 bit audio
>samples,
>>>> only the first 700 or 800 are actual values, I've zeroed the rest) to
>>>>
>>>> Fourier.FFT(myarray, 512, FourierDirection.Forward)
>>>>
>>>> So I'm wondering a few things . . . .
>>>>
>>>> the length variable (here I've hard coded 512), what is this length?
>>>
>>>At the risk of being obvious, It's the length of the Fourier transform.
>>>Take a step backwards and look at it like this:
>>>
>>>You have N data points: x[0], x[1], ..., x[N-1]. When you take the
>>>Fourier transform of those N points, you are answering the following
>>>question:
>>>
>>>What are the coefficients of the function f(t) which has the form
>>>
>>>f(t) = a_0 + 2 sum_{k=1}^{N/2-1} ( a_k cos(2 pi k / N t) - b_k sin(2 pi
>>>k / N t) ) + a_{N/2} cos(pi t),
>>>
>>>such that f(n) = x[n] ? That is, you are fitting a periodic function
>>>(f(n) = f(n+N)) through the data points x[n]. The DFT computes the
>>>coefficients a_k and b_k (sometimes called real and imaginary part of
>>>the DFT). Another, but equivalent form of f(t) is
>>>
>>>f(t) = sum_{k=0}^{N/2} A_k cos(2 pi k / N t + phi_k).
>>>
>>>The A_k's are called the magnitudes, and the phi_k are called the
>>>phases. Together, they form an estimate of the spectrum of the data
>>>x[n].
>>>
>>>Other names for these coefficients: the a_k's and b_k's are called the
>>>rectangular coordinates of the DFT, while the A_k's and the phi_k's are
>>>called the polar coordinates of the DFT.
>>>
>>>>
>>>> Then how do I interpret the results of the new values in myarray?
>I've
>>>> heard something's about where complex numbers, where two slot
>>>> represent one complex number? When I graph the resulting array, it
>>>> doesn't look like most of the FFT graphs I've seen?
>>>
>>>Usually, FFT routines return the rectangular coordinates of the DFT
>>>(look up the documentation of your routine for the exact format -
>>>remember that real part = a_k and imaginary part = b_k). For graphing
>>>and human interpretation, we prefer the polar coordinates - in
>>>addition, we neglect the phi_k and only concentrate on the A_k's. The
>>>A_k's can be computed through
>>>
>>>A_k = sqrt(a_k^2 + b_k^2).
>>>
>>>Graphs usually plot 20 log(A_k) = 10 log(a_k^2 + b_k^2). This kind of
>>>graph should be more familiar to you.
>>>
>>>> I'd appreciate a layman's response without to much reference to math,
>>>> as it not really my strong point.
>>>>
>>>> Thanks for your help
>>>>
>>>>
>>>> Carmen Branje
>>>> Toronto
>>>
>>>Regards,
>>>Andor
>>>Zurich
>>>
>>>
>>
>>
>>
>


 | 
Pages: 1
Prev: multiple clock design
Next: THD without using an FFT