MATLAB: FFT implementation by theselft

dft implementation

Hello. Im studying discrete fourier transform and for that I want to implement it. I mean, this is the code I made but there is a problem
I defined a little discrete signal, x. And I calculate the DFT and then the inverse DFT but I dont get the same signal.
What Im doing wrong?
Thanks in advance.
    x=[0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 -0.1 -0.2 -0.3 -0.4 -0.5 -0.6 -0.7 -0.8 -0.9 -1 -0.9 -0.8 -0.7 -0.6 -0.5 -0.4 -0.3 -0.2 -0.1] 
sum=0
for k=1:40 
     for j=1:40
       sum= sum+x(j)*exp(- 1i * 2 * pi * (1/20)*(j-1)*(k-1))  
     end
     a(k)=sum
      sum=0
  end
n=1:1:40 
sum=0
for k=1:40 
     for j=1:40
       sum= sum+a(j)*exp( 1i * 2 * pi * (1/20)*(j)*(k))  
     end
     b(k)=(1/40)*sum
      sum=0
  end
subplot(1,2,1)
plot(n, x)
subplot(1,2,2)
plot(n,b)

Best Answer

It would be more efficient to use a vector form of the DFT for the coding in matlab. As an example, see below. Note that this method is only best for shorter sequences ,x(n), if you have a larger value for N, then another method would be preferred.
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
Related Question