function test_pfb_fft_business

%put freq and phase here
phase=300/180*pi;
%freq = 86119800; %2600*33123
%freq = 110822400; %2600*64*666

freq = 86138000; %2600*33130


%generated signal for ADC
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% amp_diff and DC here is just try to simulate the situation when the I and
% Q stream have different property. Right now, we can leave the amp_diff =1
% and DC = 0, which mean there is not amplitude and phase difference (other than I and Q)
amp_diff=1;
DC = 0;


step = 1/(340787200);
x = [0:step:(2^19-1)*step];

ADC_noise =60;

y = 2048*0.7*sin(x*2*pi*freq+phase);
z = 2048*0.7*cos(x*2*pi*freq+phase);
y_new = 2048*amp_diff*0.7*sin(x*2*pi*freq+phase)+2048*amp_diff*DC;
z_new = 2048*0.7*cos(x*2*pi*freq+phase)+2048*DC;
%z = sin(x*2*pi/4x

y = round(y);
z = round(z);
y_new = round(y_new);
z_new = round(z_new);

%plot(20*log10(abs(fft(complex(y,-z)))))

y_nois = awgn(y,ADC_noise,'measured'); 
z_nois = awgn(z,ADC_noise,'measured'); 
y_new_nois = awgn(y_new,ADC_noise,'measured'); 
z_new_nois = awgn(z_new,ADC_noise,'measured'); 

y_nois = round(y_nois);
z_nois = round(z_nois);
y_new_nois = round(y_new_nois);
z_new_nois = round(z_new_nois);
%hold on
%plot(z,'red')
c_nois=complex(z_nois,-1*y_nois);
c_new_nois=complex(z_new_nois,-1*y_new_nois);

size(c_new_nois)
%size is 1 131072*4

temp_c_new_nois = c_new_nois(1,1:2^17);
temp_c_new_nois_fft = ifft(temp_c_new_nois);
figure
plot(20*log10(abs(temp_c_new_nois_fft)),'.r')
title('matlab 2^17 FFT output')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%ADD PFB_FIR EFFECT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
PFBSize = 11;
TotalTaps = 4;
WindowType = 'hamming';
n_inputs = 0;
nput = 0;
fwidth =1;
a1 = 1;
a2=2;
a3=3;
a4=4;

alltaps = TotalTaps*2^PFBSize;
windowval = transpose(window(WindowType, alltaps));
% plot(windowval)
total_coeffs = windowval .* sinc(fwidth*([0:alltaps-1]/(2^PFBSize)-TotalTaps/2));
%in this case, buf = total_coeff
for i=1:alltaps/2^n_inputs,
    buf(i)=total_coeffs((i-1)*2^n_inputs + nput + 1);
end
coeff1 = buf((a1-1)*2^(PFBSize-n_inputs)+1 : a1*2^(PFBSize-n_inputs));
coeff2 = buf((a2-1)*2^(PFBSize-n_inputs)+1 : a2*2^(PFBSize-n_inputs));
coeff3 = buf((a3-1)*2^(PFBSize-n_inputs)+1 : a3*2^(PFBSize-n_inputs));
coeff4 = buf((a4-1)*2^(PFBSize-n_inputs)+1 : a4*2^(PFBSize-n_inputs));

size(coeff1)
for i=1:1:64
    for j=1:1:2048
        c_new_nois_fir(1,(j+2048*(i-1))) = c_new_nois(1,(j+2048*(i-1)))*coeff1(j)+ c_new_nois(1,(j+2048*(i-1))+2048)*coeff2(j)+c_new_nois(1,(j+2048*(i-1))+4096)*coeff3(j)+c_new_nois(1,(j+2048*(i-1))+6144)*coeff4(j);
       
    end
    i
end

size(c_new_nois_fir)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

c_new_nois_original = c_new_nois;

c_new_nois = c_new_nois_fir;








%c=y;
%plot(20*log10(abs(fft(c_new))))
%figure
%plot(20*log10(abs(fft(c_new_nois))),'.r')


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% do the first FFT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%maybe add lower frequency ( or so call 1/f) type noise later on

fft_matrix_new = zeros(64,2048);

for i=1:1:64
fft_first_new = (ifft(c_new_nois(1,(2048*(i-1)+1):2048*i),2^11));


fft_matrix_new(i,:) = fft_first_new;


end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%transpose and 2nd FFT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
fft_second_new = zeros(1,131072);


for i=1:1:2048
    fft_second_new(1,(64*(i-1)+1):64*i)=ifft(fft_matrix_new(:,i),64);

end

figure
plot(20*log10(abs(fft_matrix_new(1,:))),'.r')
title('first PFB output')
figure

plot(20*log10(abs(fft_second_new)),'.b')
title('final FFT output')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

bin_pos = freq/2600 + 1; 

sprintf('original phase is')
temp_c_new_nois_fft(1,bin_pos)
angle(temp_c_new_nois_fft(1,bin_pos))/pi*180


sprintf('pfb + fft phase is')
fft_second_new(1,bin_pos)
angle(fft_second_new(1,bin_pos))/pi*180

angle(fft_second_new(1,bin_pos+64))/pi*180




% freq =86119800, phase is set to 60, direct 2^17 matlab fft get -60 and PFB+FFT get -93.7598

% freq = 86119800, phase is set to 20, direct 2^17 matlab fft get -20 and PFB+FFT get -53.7598

% freq = 86138000, phase is set to 20, direct 2^17 matlab fft get -20 and PFB+FFT get -132.5052

% freq = 86138000, phase is set to 300, direct 2^17 matlab fft get 60 and PFB+FFT get -52.5052


%%%%% for future simulation %%%%%%%%%%%%%%%%%%


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



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


%select Bin


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



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



%FIR


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



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



%decimation


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



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









% 
% xx = abs(fft(c,2^17))/(2^17);
% for i=1:1:2^17
%     if xx(i)==0
%         xx(i)=1/(2^17);
%     end
% end
% 
% xx_new = abs(fft(c_new,2^17))/(2^17);
% for i=1:1:2^17
%     if xx_new(i)==0
%         xx_new(i)=1/(2^17);
%     end
% end
% 
% 
% yy = real(fft(c,2^17))/(2^17);
% for i=1:1:2^17
%     if yy(i)==0
%         yy(i)=1/(2^17);
%     end
% end
% 
% yy_new = real(fft(c_new,2^17))/(2^17);
% for i=1:1:2^17
%     if yy_new(i)==0
%         yy_new(i)=1/(2^17);
%     end
% end
% 
% 
% 
% 
% zz = imag(fft(c,2^17))/(2^17);
% for i=1:1:2^17
%     if zz(i)==0
%         zz(i)=1/(2^17);
%     end
% end
% 
% zz_new = imag(fft(c_new,2^17))/(2^17);
% for i=1:1:2^17
%     if zz_new(i)==0
%         zz_new(i)=1/(2^17);
%     end
% end
% 
% 
% 
% subplot(6,1,1);
% plot(20*log10((xx)))
% subplot(6,1,2);
% plot(20*log10((xx_new)),'red')
% 
% 
% 
% 
% subplot(6,1,3);
% semilogy(yy)
% 
% subplot(6,1,4);
% semilogy(yy_new,'red')
% 
% subplot(6,1,5);
% semilogy(zz)
% subplot(6,1,6);
% semilogy(zz_new,'red')
% 
% 
% figure
% semilogy(yy_new)
% hold on
% 
% semilogy(zz_new,'red')

end