I'm making in Matlab an 8 channel FFT, according to an online tutorial I found at https://www.youtube.com/watch?v=EsJGuI7e_ZQ&t=957s . It works fine for the same number of inputs as shown in the video (N=8), but as I increase N (say for N=128) it no longer works.
I am not sure why it no longer works as I increase N, but I suspect it has something to do with the twiddle factors. I assumed the twiddle factors would have an exponent increasing by 8 after each stage, but maybe they change in a different way as one increases N?
Here is a picture from the YouTube tutorial, at 15:16
Here is a signal flow diagram of my Matlab code
And below is my actual MatLab code. If one uses N=8 (which you can define at top of code) then output is correct. Using N=128 (or other higher values) produces an incorrect result.
I am very confused about what could be causing this problems with the output. Any thoughts?
Much appreciated,
clear all
% Generate input data sequence and plot
N = 128;
f1 = 10;
num_cycles = 2;
fs = f1*N/num_cycles;
x_time = 0:1/fs:num_cycles/f1-1/fs;
X = sin(x_time*2*pi*f1);
plot(x_time,X);
title('Input Waveform');
% Split inputs into eight channels
X0_0 = X(1:8:N);
X0_1 = X(5:8:N);
X0_2 = X(3:8:N);
X0_3 = X(7:8:N);
X0_4 = X(2:8:N);
X0_5 = X(6:8:N);
X0_6 = X(4:8:N);
X0_7 = X(8:8:N);
% Compute FFT of each channel
X1_0 = fft(X0_0);
X1_1 = fft(X0_1);
X1_2 = fft(X0_2);
X1_3 = fft(X0_3);
X1_4 = fft(X0_4);
X1_5 = fft(X0_5);
X1_6 = fft(X0_6);
X1_7 = fft(X0_7);
% Generate Twiddle factors
Wn=exp(-1i*2*pi/N);
% Produce output of first stage of butterfly
for k=0:(N/8)-1
X2_0(k+1) = X1_0(k+1) + (Wn^(k*8)) * X1_1(k+1);
X2_1(k+1) = X1_0(k+1) + (Wn^(k*8+4)) * X1_1(k+1);
X2_2(k+1) = X1_2(k+1) + (Wn^(k*8)) * X1_3(k+1);
X2_3(k+1) = X1_2(k+1) + (Wn^(k*8+4)) * X1_3(k+1);
X2_4(k+1) = X1_4(k+1) + (Wn^(k*8)) * X1_5(k+1);
X2_5(k+1) = X1_4(k+1) + (Wn^(k*8+4)) * X1_5(k+1);
X2_6(k+1) = X1_6(k+1) + (Wn^(k*8)) * X1_7(k+1);
X2_7(k+1) = X1_6(k+1) + (Wn^(k*8+4)) * X1_7(k+1);
end
% Produce output of second stage of butterfly
for k=0:(N/8)-1
X3_0(k+1) = X2_0(k+1) + (Wn^(k*8)) * X2_2(k+1);
X3_1(k+1) = X2_1(k+1) + (Wn^(k*8+2)) * X2_3(k+1);
X3_2(k+1) = X2_0(k+1) + (Wn^(k*8+4)) * X2_2(k+1);
X3_3(k+1) = X2_1(k+1) + (Wn^(k*8+6)) * X2_3(k+1);
X3_4(k+1) = X2_4(k+1) + (Wn^(k*8)) * X2_6(k+1);
X3_5(k+1) = X2_5(k+1) + (Wn^(k*8+2)) * X2_7(k+1);
X3_6(k+1) = X2_4(k+1) + (Wn^(k*8+4)) * X2_6(k+1);
X3_7(k+1) = X2_5(k+1) + (Wn^(k*8+6)) * X2_7(k+1);
end
% Produce output of third stage of butterfly
for k=0:(N/8)-1
X4_0(k+1) = X3_0(k+1) + (Wn^(k*8)) * X3_4(k+1);
X4_1(k+1) = X3_1(k+1) + (Wn^(k*8+1)) * X3_5(k+1);
X4_2(k+1) = X3_2(k+1) + (Wn^(k*8+2)) * X3_6(k+1);
X4_3(k+1) = X3_3(k+1) + (Wn^(k*8+3)) * X3_7(k+1);
X4_4(k+1) = X3_0(k+1) + (Wn^(k*8+4)) * X3_4(k+1);
X4_5(k+1) = X3_1(k+1) + (Wn^(k*8+5)) * X3_5(k+1);
X4_6(k+1) = X3_2(k+1) + (Wn^(k*8+6)) * X3_6(k+1);
X4_7(k+1) = X3_3(k+1) + (Wn^(k*8+7)) * X3_7(k+1);
end
% Merge X4 into final output
c = 1;
for k = 1:N/8
X5(c:c+7) = [X4_0(k);X4_1(k);X4_2(k);X4_3(k);X4_4(k);X4_5(k);X4_6(k);X4_7(k)];
c = c+8;
end
% Plot expected FFT and Butterfly FFT
figure
matlab_fft=fft(X);
plot(abs(matlab_fft));
title('Matlab FFT');
figure
plot(abs(X5));
title('Butterfly FFT');