Mario.Tapilouw

Tuesday, June 01, 2010

Discrete Fourier Transform using FFTW library

I was trying to understand how to use fftw library in my software because it was recommended by people and some of my colleague also use this library. This library is known to be the fastest FFT library. Well, that is one of the reason why I choose this and the other reason is because I don't want to build my own FFT function which might be full of bug and unstable. Therefore I prefer using this library instead.

It's quite straight away. Depends on your platform, the suitable package is available in the fftw website. I'm using Windows so I have to follow this guide for installing and using the library. I don't need to explain it here, it's very clear in the website, if you have any problem just send email to me.

There are four basic steps that we have to do when using fftw. I got several good examples for using fftw, one of them is this one. Here I will explain how it works.

1. Variables declaration.

//fft data window
int n=512;

//the complex value for storing the input values and output values
fftw_complex* b1 = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*n);
fftw_complex* b2 = (fftw_complex*)fftw_malloc(sizeof(fftw_complex)*n);

2. Declaring FFT plan

//create a forward fft plan from b1 into b2
fftw_plan p1 = fftw_plan_dft_1d(n, b1, b2, FFTW_FORWARD, FFTW_ESTIMATE);


3. Fill the data into the input array. The arrangement of the complex data is [real][imaginary] So for example for the first data of b1, b1[0][0] is the real part and b[0][1] is the imaginary part and so on.

4. Perform FFT by executing the plan

//executing fftw plan
fftw_execute(p1);

5. Draw the spectrum (if needed). The spectrum is can be obtained from the first n/2 of the data. So we can draw the spectrum by calculating the magnitude.

for (i=0; i < mag =" sqrt(pow(b2[i][0],2)">
6. After finished using the fftw, the library needs to be released from the memory.
//destroy the plan
fftw_destroy_plan(p1);

//destroy the variables
fftw_free(b1);
fftw_free(b2);

Example of the signal that I tried, it's a superposition of two sine waves with different frequencies.



and the result that I obtained after performing FFT is:

That's it, it works with me. Hope it will be useful for you.

Labels: , ,

3 Comments:

Post a Comment

<< Home