# Need FFT Code for Matlab (not built in)

889 views (last 30 days)

Show older comments

Does anyone have FFT code for Matlab? I don't want to use the built-in Matlab function.

##### 5 Comments

Yasin Arafat Shuva
on 19 Jul 2021

clear variables

close all

%%Implementing the DFT in matrix notation

xn = [0 2 4 6 8]; % sequence of x(n)

N = length(xn); % The fundamental period of the DFT

n = 0:1:N-1; % row vector for n

k = 0:1:N-1; % row vecor for k

WN = exp(-1i*2*pi/N); % Wn factor

nk = n'*k; % creates a N by N matrix of nk values

WNnk = WN .^ nk; % DFT matrix

Xk = xn * WNnk; % row vector for DFT coefficients

disp('The DFT of x(n) is Xk = ');

disp(Xk)

magXk = abs(Xk); % The magnitude of the DFT

%%Implementing the inverse DFT (IDFT) in matrix notation

n = 0:1:N-1;

k = 0:1:N-1;

WN = exp(-1i*2*pi/N);

nk = n'*k;

WNnk = WN .^ (-nk); % IDFS matrix

x_hat = (Xk * WNnk)/N; % row vector for IDFS values

disp('and the IDFT of x(n) is x_hat(n) = ');

disp(real(x_hat))

% The input sequence, x(n)

subplot(3,1,1);

stem(k,xn,'linewidth',2);

xlabel('n');

ylabel('x(n)')

title('plot of the sequence x(n)')

grid

% The DFT

subplot(3,1,2);

stem(k,magXk,'linewidth',2,'color','r');

xlabel('k');

ylabel('Xk(k)')

title('DFT, Xk')

grid

% The IDFT

subplot(3,1,3);

stem(k,real(x_hat),'linewidth',2,'color','m');

xlabel('n');

ylabel('x(n)')

title('Plot of the inverse DFT, x_{hat}(n) = x(n)')

grid

### Answers (4)

Youssef Khmou
on 15 Mar 2013

hi John

Yes there are many versions of Discrete Fourier Transform :

% function z=Fast_Fourier_Transform(x,nfft)

%

% N=length(x);

% z=zeros(1,nfft);

% Sum=0;

% for k=1:nfft

% for jj=1:N

% Sum=Sum+x(jj)*exp(-2*pi*j*(jj-1)*(k-1)/nfft);

% end

% z(k)=Sum;

% Sum=0;% Reset

% end

% return

this is a part of the code , try my submission it contains more details :

##### 4 Comments

Walter Roberson
on 6 Sep 2017

Tobias
on 11 Feb 2014

Edited: Walter Roberson
on 6 Sep 2017

For anyone searching an educating matlab implementation, here is what I scribbled together just now:

function X = myFFT(x)

%only works if N = 2^k

N = numel(x);

xp = x(1:2:end);

xpp = x(2:2:end);

if N>=8

Xp = myFFT(xp);

Xpp = myFFT(xpp);

X = zeros(N,1);

Wn = exp(-1i*2*pi*((0:N/2-1)')/N);

tmp = Wn .* Xpp;

X = [(Xp + tmp);(Xp -tmp)];

else

switch N

case 2

X = [1 1;1 -1]*x;

case 4

X = [1 0 1 0; 0 1 0 -1i; 1 0 -1 0;0 1 0 1i]*[1 0 1 0;1 0 -1 0;0 1 0 1;0 1 0 -1]*x;

otherwise

error('N not correct.');

end

end

end

##### 3 Comments

Jan Afridi
on 6 Sep 2017

Edited: Jan Afridi
on 6 Nov 2017

##### 1 Comment

Jan Suchánek
on 6 Dec 2018

Hello i would like to ask you what had to be different if you would like to use it for harmonical function.

I tried to implement your solution but it was working only for random vectors but not for any type of harmonical.

i was using this for the harmonical function.

______________________

t2=0:1:255;

h1=sin(t2);

x=h1;

X11 = myFFT(x);

______________________

it gave me this error:

Error using *

Inner matrix dimensions must agree.

Error in myFFT (line 18)

X = [1 0 1 0; 0 1 0 -1i; 1 0 -1 0;0 1 0 1i]*[1 0 1 0;1

0 -1 0;0 1 0 1;0 1 0 -1]*x;

Error in myFFT (line 7)

Xp = myFFT(xp);

Error in myFFT (line 7)

Xp = myFFT(xp);

Error in myFFT (line 7)

Xp = myFFT(xp);

Error in myFFT (line 7)

Xp = myFFT(xp);

Error in myFFT (line 7)

Xp = myFFT(xp);

Error in myFFT (line 7)

Xp = myFFT(xp);

Error in ZBS_projekt (line 246)

X11 = myFFT(x);

Ryan Black
on 11 Mar 2019

Edited: Ryan Black
on 26 Aug 2020

Function I wrote that optimizes the DFT or iDFT for BOTH prime and composite signals...

The forward transform is triggered by -1 and takes the time-signal as a row vector:

[Y] = fftmodule(y,-1)

The inverse transform auto-normalizes by N, is triggered by 1 and takes the frequency-signal as a row vector:

[y] = fftmodule(Y,1)

Methodology:

The algorithm decimates to N's prime factorization following the branches and nodes of a factor tree. In formal literature this may be referred to as Mixed Radix FFT, but its really just recursive decimation of additive groups and this method is easily derivable via circular convolutions :)

At the prime tree level, algorithm either performs a naive DFT or if needed performs a single Rader's Algorithm Decomposition to (M-1), zero-pads to power-of-2, then proceeds to Rader's Convolution routine (not easy to derive). Finally it upsamples through the origianlly strucutured nodes and branches incorporating twiddle factors for the solution.

##### 0 Comments

### See Also

### Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!