%%%%%%%%%%%%%%%%%%%%%
%%    Time-Data generation        %%
%%%%%%%%%%%%%%%%%%%%%
clear all
TIME=5;
N=128;
t=linspace(0,TIME,N);
y=2*sin(10*t)-2*sin(25*t+pi/4)+2*cos(40*t+pi/8);
[FileID]=fopen('data.dat','wt+');
fprintf(FileID,'%f\n',y);
fclose('all');

% Ploting the results
kk=0:N-1;
figure
stem(kk,y,'fill')
grid on
axis tight
title('\bf< Gendrated time signal  >','fontsize',12)

%%%%%%%%%%%%%%%%%%%%%
%%         Time-Data read           %%
%%%%%%%%%%%%%%%%%%%%%
clear all;
[FileID]=fopen('data.dat','r');
[x count]= fscanf(FileID,'%f \n');
fclose('all');

%%%%%%%%%%%%%%%%%%%%%
%%      DIT-FFT calculation        %%
%%%%%%%%%%%%%%%%%%%%%
nu=log2(length(x));
N=2^nu;
NM=N-1;

%% Rearrange data by bit reversing
jj=1;
for ii=2:NM
    k=N/2;
     while(k<jj)
           jj=jj-k;
           k=k/2;
     end
      jj=jj+k;
      
      if (ii<jj)
        T=x(jj);
        x(jj)=x(ii);
        x(ii)=T;
    end
end

%% FFT Calculation
for m=1:nu
    LE=2^m;
    LE1=LE/2;
    U=1;
    W=exp(-i*pi/LE1);
    for n=1:LE1
        for p=n:LE:N-1
            pp=p+LE1;
            temp=x(pp)*U;
            x(pp)=x(p)-temp;
            x(p)=x(p)+temp;
        end
     U=U*W;
    end
end

X_real=real(x);
X_imag=imag(x);
X_mag=abs(x);
X_phase=phase(x)*180/pi;

% Ploting the results
kk=0:NM;


figure
subplot(221)
stem(kk,X_real,'fill')
grid on
axis tight
title('\bf< Real >','fontsize',12)

subplot(222)
stem(kk,X_imag,'fill')
grid on
axis tight
title('\bf< Imaginary >','fontsize',12)

subplot(223)
stem(kk,X_mag,'fill')
grid on
axis tight
title('\bf< Magnitude >','fontsize',12)

subplot(224)
stem(kk,X_phase,'fill')
grid on
axis tight
title('\bf< Phase >','fontsize',12)

%%%%%%%%%%%%%%%%%%%%%
%%     Frequency-Data write       %%
%%%%%%%%%%%%%%%%%%%%%
[FileID]=fopen('fft_data.dat','wt+');
for i=1:N
fprintf(FileID,'%f %f\n',X_real(i),X_imag(i));
end
fclose('all');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%
%%     Frequency-Data read       %%
%%%%%%%%%%%%%%%%%%%%%
clear all;
[FileID]=fopen('fft_data.dat','r');
[F count]= fscanf(FileID,'%f\n');
fclose('all');

%%%%%%%%%%%%%%%%%%%
%%    DIT-IFFT calculation     %%
%%%%%%%%%%%%%%%%%%%
N=128;
Freq=reshape(F,2,N)';
x=Freq(:,1)+i*Freq(:,2);
nu=log2(length(x));
NM=N-1;
kk=0:NM;

%% Rearrange data by bit reversing
jj=1;
for ii=2:NM
    k=N/2;
     while(k<jj)
           jj=jj-k;
           k=k/2;
     end
      jj=jj+k;
      
      if (ii<jj)
        T=x(jj);
        x(jj)=x(ii);
        x(ii)=T;
    end
end

%% FFT Calculation
for m=1:nu
    LE=2^m;
    LE1=LE/2;
    U=1;
%    W=exp(-i*pi/LE1);     % For FFT
    W=exp(i*pi/LE1);      % For IFFT
    for n=1:LE1
        for p=n:LE:N-1
            pp=p+LE1;
            temp=x(pp)*U;
            x(pp)=x(p)-temp;
            x(p)=x(p)+temp;
        end
     U=U*W;
    end
end

x=x./N;
X_time=real(x);

% Ploting the results
kk=0:NM;
figure
stem(kk,X_time,'fill')
grid on
axis tight
title('\bf< Time signal by IFFT  >','fontsize',12)