From: Jim Peters (jim_at_uazu.net)
Date: 2002-06-15 15:50:03
Michal Wallace wrote:
> What I tried to do was create 10 seconds of fake data,
> simulating "pure" brainwaves (such as a steady 6Hz for
> theta), as well as some contrived mixtures. Then, I
> ran an FFT on them to get the spectrogram.
>
> The output bears some resemblance to what I wanted, but
> it's grossly distorted, and I'm not sure what's wrong.
>
> I think the code *MIGHT* be working fine except that it's
> using the wrong scale for the "bins", because it looks like
> it might be right, except it's all squished into the top
> bars and there seems to be a lot of noise.
Here we go ...
When you do an FFT of a block of sound, you can't just FFT it straight
off. Well, you can, but you end up with lots of distortion because
you have effectively cut a rectangle out of the sound, and there are
big clicks at the start and end of it. (And when you give a chunk of
sound to an FFT, it acts as if that sound chunk is looped forever, so
you can imagine what kind of distortion that creates if you don't tidy
up the edges).
So, you need to apply a 'window' to the data before FFTing it. If you
look at BWView, if you click somewhere on the display, you can see the
'window function' being used in the display above. This looks
something like a bell-shape.
So, if you want nice results from your FFT, you need to window the
data first -- which means multiplying it by a 'window function',
effectively flattening off the edges. Don't think you can choose just
any old function, because each different shape has its own particular
frequency response. There are several types of window, but I used the
'Blackman' window in my code.
You mention that you've already looked at the code, so you should have
no trouble finding this (in analysis.c around line 410 onwards). The
multiplier is given by this:
double mag= 0.42 + 0.5 * cos(ang) + 0.08 * cos(2*ang); // Blackman window
Where 'ang' goes from -pi at one end of the chunk to +pi at the other.
It's a while since I've coded in Python, so I'll prototype it in C for
you:
for (a= 0; a<64; a++) {
double ang= (a + 0.5 - 32.0) * M_PI / 32.0;
double mag= 0.42 + 0.5 * cos(ang) + 0.08 * cos(2*ang);
slice[a] *= mag;
}
I think that should be correct.
Note that a 6Hz wave is not going to be very easy to see if you are
only using a 0.25 second window (that's only 1.5 of a wavelength
you're expecting it to recognise!). Bear in mind that the lowest band
in your FFT will be centred on 0Hz, then 4Hz, then 8Hz, 16Hz, 24Hz,
and so on up to 128Hz.
Jim
-- Jim Peters (_)/=\~/_(_) jim_at_uazu.net (_) /=\ ~/_ (_) Uazú (_) /=\ ~/_ (_) http:// B'ham, UK (_) ____ /=\ ____ ~/_ ____ (_) uazu.net
This archive was generated by hypermail 2.1.4 : 2002-07-27 12:28:43 BST