# Beat detection and FFT

I am working on a platformer game which includes music with beat detection. I am currently detecting beats by checking for when the current amplitude exceeds a historical sample. This doesn't work well with genres of music, like rock, which have a pretty steady amplitude.

So I looked further and found algorithms splitting the sound into multiple bands using FFT... then I found the Cooley-Tukey FFt algorithm

The only problem I'm having is that I am quite new to audio and I have no idea how to use that to split the signal up into multiple signals.

So my question is :

How do you use a FFT to split a signal into multiple bands ?

Also for the guys interested, this is my algorithm in c# :

// C = threshold, N = size of history buffer / 1024
public void PlaceBeatMarkers(float C, int N)
{
List<float> instantEnergyList = new List<float>();
short[] samples = soundData.Samples;

float timePerSample = 1 / (float)soundData.SampleRate;
int sampleIndex = 0;
int nextSamples = 1024;

// Calculate instant energy for every 1024 samples.
while (sampleIndex + nextSamples < samples.Length)
{

float instantEnergy = 0;

for (int i = 0; i < nextSamples; i++)
{
instantEnergy += Math.Abs((float)samples[sampleIndex + i]);
}

instantEnergy /= nextSamples;

if(sampleIndex + nextSamples >= samples.Length)
nextSamples = samples.Length - sampleIndex - 1;

sampleIndex += nextSamples;
}

int index = N;
int numInBuffer = index;
float historyBuffer = 0;

//Fill the history buffer with n * instant energy
for (int i = 0; i < index; i++)
{
historyBuffer += instantEnergyList[i];
}

// If instantEnergy / samples in buffer < instantEnergy for the next sample then add beatmarker.
while (index + 1 < instantEnergyList.Count)
{
if(instantEnergyList[index + 1] > (historyBuffer / numInBuffer) * C)
beatMarkers.Add((index + 1) * 1024 * timePerSample);
historyBuffer -= instantEnergyList[index - numInBuffer];
historyBuffer += instantEnergyList[index + 1];
index++;
}
}


Well, if your input signal is real (as in, each sample is a real number), the spectrum will be symmetric and complex. Exploiting the symmetry, usually FFT algorithms pack the result by giving you back only the positive half of the spectrum. The real part of each band is in the even samples and the imaginary part in the odd samples. Or sometimes the real parts are packed together in the first half of the response and the imaginary parts in the second half.

In formulas, if X[k] = FFT( x[n] ), you give it a vector i[n] = x[n] , and get an output o[m], then

X[k] = o[2k] + j·o[2k+1]


(although sometimes you get X[k] = o[k] + j·o[k+K/2], where K is the length of your window, 1024 in your example). By the way, j is the imaginary unit, sqrt(-1).

The magnitude of a band is computed as the root of the product of this band with its complex conjugate:

|X[k]| = sqrt( X[k] · X[k]* )


And the energy is defined as the square of the magnitude.

If we call a = o[2k] and b = o[2k+1], we get

X[k] = a + j·b


therefore

E[k] = |X[k]|^2 = (a+j·b)·(a-j·b) = a·a + b·b


Unrolling the whole thing, if you got o[m] as output from the FFT algorithm, the energy in the band k is:

E[k] = o[2k] · o[2k] + o[2k+1] · o[2k+1]


(Note: I used the symbol · to indicate multiplication instead of the usual * in order to avoid confusion with the conjugation operator)

The frequency of the band k, assuming a sampling frequency of 44.1Khz and a window of 1024 samples, is

freq(k) = k / 1024 * 44100 [Hz]


So, for example, your first band k=0 represents 0 Hz, k=1 is 43 Hz, and the last one k=511 is 22KHz (the Nyquist frequency).

I hope this answers your question about how do you get the energy of the signal per band using the FFT.

Addendum: Answering your question in the comment, and assuming you are using the code from the link you posted in the question (The Cooley-Tukey algorithm in C): Let's say you have your input data as a vector of short ints:

// len is 1024 in this example.  It MUST be a power of 2
// centerFreq is given in Hz, for example 43.0
double EnergyForBand( short *input, int len, double centerFreq)
{
int i;
int band;
complex *xin;
complex *xout;
double magnitude;
double samplingFreq = 44100.0;

// 1. Get the input as a vector of complex samples
xin = (complex *)malloc(sizeof(struct complex_t) * len);

for (i=0;i<len;i++) {
xin[i].re = (double)input[i];
xin[i].im = 0;
}

// 2. Transform the signal
xout = FFT_simple(xin, len);

// 3. Find the band ( Note: floor(x+0.5) = round(x) )
band = (int) floor(centerFreq * len / samplingFreq + 0.5);

// 4. Get the magnitude
magnitude = complex_magnitude( xout[band] );

// 5. Don't leak memory
free( xin );
free( xout );

// 6. Return energy
return magnitude * magnitude;
}


My C is a bit rusty (I am mostly coding in C++ nowadays), but I hope I didn´t make any big mistake with this code. Of course if you were interested in the energy of other bands it makes no sense to transform the whole window for each of them, that would be a waste of CPU time. In that case do the transformation once and get all the values you need from xout.

• Oh, I just took a look at the code you linked, it already gives you the results in "complex" form and even provides you with a function to compute the magnitude of a complex number. Then you would only have to compute the square of that magnitude for each element of the output vector, don't need to worry about sorting the results. – CeeJay Mar 15 '11 at 11:39
• As a example if I have all 1024 samples from window 0-1024 and I got them as real values, so no complex part. and I want to calculate the energy in there on the frequency band 43Hz. How would I integrate it then ? (I only need the real part back, the postive part) If you could do it in some pseudocode I'll be in depth of you forever and then I may actually grasp the concept a bit :) – Quincy Mar 15 '11 at 20:34
• The code I wrote is using the C library you linked, which already contains a "complex" structure. This makes the unwrapping I described in my question unnecessary (and the code reflects that) – CeeJay Mar 15 '11 at 21:48