% Exercise 6.3.1. Radix-2 FFT.
clc; clear; close all;
%% a. Generate x[n]:
N = 64; % Number of Samples.
n = 0:N-1;
% x = [ones(1,9) zeros(1,55)];
x = randn(1,64) + j*randn(1,64);
% Compute its 64-point DFT:
X = my_DFT(x);
% Plot the real and imaginary parts of the DFT of x[n]:
k = 0:N-1;
figure('Name','Exercise 6.3.1. Radix-2 FFT');
subplot(2,1,1);
stem(k,real(X));
title('\Ree\{X[k]\}');
xlabel('Sample Number k');
axis tight;
grid on;
hold on;
subplot(2,1,2);
stem(k,imag(X),'r');
title(['\Imm\{X[k]\}']);
xlabel('Sample Number k');
axis tight;
grid on;
hold on;
%% b. Use the Radix-2 FFT formula to compute the DFT of the given signal.
X1 = radix2fft(x);
% The following code is the body of radix2fft.m function.
% It is repeatd here for ease of reference.
%
% Generate sequences g[n] and h[n].
% g = zeros(1,32);
% h = zeros(1,32);
% for i=1:32
% g(i) = x(2*i-1); % N/2 sequence made of the even samples of x[n]: x[0], x[2], x[4], ..., x[30].
% h(i) = x(2*i); % N/2 sequence made of the odd samples of x[n]: x[1], x[3], x[5], ..., x[31].
% end
%
% % Compute their 32-point DFT's:
% G1 = my_DFT(g);
% H1 = my_DFT(h);
%
% % Create the periodic extensions of these 2 DFT's:
% H = [H1 H1];
% G = [G1 G1];
%
% % Create the phasor vector:
% W = exp(-j*2*pi*k/Samples);
%
% % Apply the formula.
% X1 = G + W.*H;
% Now compare the results by plotting them altogether.
subplot(2,1,1);
stem(k,real(X1),'b.');
title('\Ree\{X[k]\} using myDFT() (circles) and radix2fft() (dots)');
xlabel('Sample Number k');
axis tight;
grid on;
subplot(2,1,2);
stem(k,imag(X1),'r.');
title(['\Imm\{X[k]\} using myDFT() (circles) and radix2fft() (dots)']);
xlabel('Sample Number k');
axis tight;
grid on;