Wednesday, January 8, 2014

Fast Fourier Transform Sample Code in MATLAB

The following matlab code computes and plots the fourier transform of the gaussian function exp(x^2) (the periodic domain is set to -L / 2 to L / 2) in MATLAB:

clear all; close all; clc;

L=20;
n=128;

x2=linspace(-L/2, L/2, n+1);

x=x2(1:n);

k=(2*pi/L)*[0:n/2-1 -n/2:-1];

u=exp(-x.^2);

ut=fft(u);

plot(fftshift(k), abs(fftshift(ut)));

The MATLAB function fft run in O(n * log(n)) complexity and assumes the periodic domain from 0 to 2*pi, therefore the value vector k is scaled for that the frequency spectrum ut corresponds to (-L/2, L/2). Note that n is power of 2 since the reason fft is able to speed up from O(n^2) to O(n*log(n)) using divide and conquest by dividing the computation into two separate computation each time. linspace() generates (n+1) values evenly from -L/2 to L/2 (inclusive, therefore should be n+1 instead of n, but x takes only the first n values from x2). fftshift rearranges the outputs of fft, fft2, and fftn by moving the zero-frequency component to the center of the array.

No comments:

Post a Comment