%%%新方案，接收端先進行合并，將受到jamming影響的子載波去除

clear all
clc
format long
%本次仿真載頻為2GHz，帶寬1MHz，子載波數64個，cp為8
%子載波間隔為15.625kHz
%一個ofdm符號長度為64us，cp長度為8us
%系統調制級別定為QPSK
%最大doppler頻率為50Hz,信道是慢變的
%多徑信道為4徑，功率延遲譜服從負指數分布~exp(-t/trms),trms=(1/5)*tmax＝1.2us時長，各徑延遲取為delay=[0 3e-6 6e-6]
t0=clock;
cp_length=16;%cp長度為16
N_carrier=64;%OFDM子載波個數
bandwidth=1000000;%系統帶寬為1MHz
t_interval=(1/bandwidth)*N_carrier/(cp_length+N_carrier);%采樣間隔64/72us，加上循環前綴后，采樣率增加
delta_f=bandwidth/N_carrier;%Hz
% SNR_dB=[0 4 8 12 16 20 24 28 32 35];%Eb/N0
% SJR_dB=[-10 -5 -3 0 3 5 10 15 20 25];%符號信干比
SNR_dB=[2 4 6 8 10];
%SNR_dB=[10];
% Pe=1e-3;
% groupnumber1=[4 16 32];
% ber_snr_persjr_ofdm=zeros(length(groupnumber1),length(SNR_dB)); 
% ber_snr_reallocatpow=zeros(length(groupnumber1),length(SNR_dB)); 
% for kk=1:length(groupnumber1)
%     groupsize=N_carrier/groupnumber1(kk);
%     groupnumber=groupnumber1(kk);
    for i=1:length(SNR_dB)%每個SNR點上仿真若干次
    snr=10^(SNR_dB(i)/10);
    error_bit_ofdm=0;
    error_bit_ofdm1=0;
    total_bit_num=0;%發送總比特數統計
%     error_bit_ofdm1=0;%for 二次功率分配
%     total_bit_num1=0;%for 二次功率分配
    loop_num=10; %共仿真1000次

     for l=1:loop_num
    ofdm_symbol_num=12;%每次仿真產生10個ofdm符號,則每次仿真共有200×64個星座映射符號；QPSK調制下，1個星座映射符號包含2個bit
    map_flag=2;
    bit_source=rand(map_flag*N_carrier,ofdm_symbol_num)>0.5;
    [m,n]=size(bit_source);
    total_bit_num=total_bit_num+m*n;
    map_out=map_module(bit_source,map_flag);      
            
   ofdm_modulation_out=sqrt(N_carrier)*ifft(map_out,N_carrier);%作64點逆FFT運算，完成ofdm調制,前面乘系數sqtr(64)是為了保持ifft前后的符號能量不變
   
      ofdm_cp_out=insert_cp(ofdm_modulation_out,cp_length);%插入循環前綴 
      %%%%%%%%以下計算噪聲方差%%%%%%%%%%
            [nn,mm]=size(ofdm_cp_out);
            spow2=0;
            for k=1:nn
                for b=1:mm
                    spow2=spow2+real(ofdm_cp_out(k,b))^2+imag(ofdm_cp_out(k,b))^2;
                end
            end
            spow1=spow2/(nn*mm);%信號平均能量
            
%             map_flag=sum(bit_num_ofdmsig)/(N_carrier*ofdm_symbol_num);
           sgma=sqrt(spow1/(2*snr)/map_flag);
        %%%%信道參數%%%%%%
%         num=1;
%         delay=[0];
%         trms=1.5e-6;
%         var_pow=10*log10(exp(-delay/trms));
%         
    num=6;

    delay=[0 1e-6 2e-6 3e-6 4e-6 6e-6];
    trms=1.5e-6;
    var_pow=10*log10(exp(-delay/trms));%各徑功率衰減,以dB形式給出
    fd=50;%最大doppler頻率為50Hz
    %t_interval=0.8889e-6;%采樣間隔64/72us，加上循環前綴后，采樣率增加
    counter=10000;%200000000;%各徑信道的采樣點間隔，應該大于信道采樣點數。由以上條件現在信道采樣點數
    count_begin=(l-1+100000000)*5*counter;%每次仿真信道采樣的開始位置
    %%%%%%以上為信道參數%%%%%%
    [passchan_ofdm_symbol,Hk]=multipath_chann(ofdm_cp_out,num,var_pow,delay,fd,t_interval,counter,count_begin,cp_length);
            
    sig_temp=cut_cp(passchan_ofdm_symbol,cp_length);
    ofdm_demodulation_out1=fft(sig_temp,N_carrier)/sqrt(N_carrier);
    HHk=ofdm_demodulation_out1./map_out;
    HHk1=HHk./abs(HHk);
    
    passnoise_ofdm_symbol=add_noise(sgma,passchan_ofdm_symbol);%加入隨機高斯白噪聲，receive_ofdm_symbol為最終接收機收到的ofdm符號塊
  
     
    cutcp_ofdm_symbol=cut_cp(passnoise_ofdm_symbol,cp_length);%去除循環前綴
    
    ofdm_demodulation_out=fft(cutcp_ofdm_symbol,N_carrier)/sqrt(N_carrier);%作128點FFT運算，完成ofdm解調
   %receive_ofdm_symbol=ofdm_demodulation_out;
    receive_ofdm_symbol=ofdm_demodulation_out./HHk;
    receive_ofdm_symbol1=ofdm_demodulation_out./HHk1;
 
    receive_bit_sig=de_map_module(receive_ofdm_symbol,map_flag);
    receive_bit_sig1=de_map_module(receive_ofdm_symbol1,map_flag);
    %以下過程統計接收信號中的錯誤比特數
           [mn,nn]=size(bit_source);
           err_num=sum(sum(rem(bit_source+receive_bit_sig,2)));
           error_bit_ofdm=error_bit_ofdm+err_num;
           err_num1=sum(sum(rem(bit_source+receive_bit_sig1,2)));
           error_bit_ofdm1=error_bit_ofdm1+err_num1;
 %%%%%%%%%以上為不做二次功率分配，原參考文獻方法%%%%%%%%%%%%            
           
           
end%for l=1:loop_num

ber_snr_persjr_ofdm(i)=error_bit_ofdm/total_bit_num;
ber_snr_persjr_ofdm1(i)=error_bit_ofdm1/total_bit_num;
%ber_snr_reallocatpow(kk,i)=error_bit_ofdm1/total_bit_num1;
% fre_eff_persnr(kk,i)=total_bit_num/(loop_num*ofdm_symbol_num*(1/bandwidth)*N_carrier)/bandwidth;
% fre_eff_persnr_reallocatpow(kk,i)=total_bit_num1/(loop_num*ofdm_symbol_num*(1/bandwidth)*N_carrier)/bandwidth;

%save bitalloc_data ber_snr_persjr_ofdm 
end%for i=1:length(SNR_dB)

% end %for kk
elapse_time=etime(clock,t0);
% figure(1)
% % semilogy(SNR_dB_theo,ber_theo(1,:),'b')
% % hold on
semilogy(SNR_dB,ber_snr_persjr_ofdm(1,:),'b-*'),hold on
semilogy(SNR_dB,ber_snr_persjr_ofdm1(1,:),'r-.'),hold off
grid on
% semilogy(SNR_dB,ber_snr_persjr_ofdm(2,:)/2,'g-o')
% semilogy(SNR_dB,ber_snr_persjr_ofdm(3,:)/6,'m-+')
% hold off
% xlabel('SNR(dB)');
% ylabel('仿真結果，實際Pe=1e-3');
% 
% figure(2)
% plot(SNR_dB,fre_eff_persnr(1,:),'r-*'),hold on
% % plot(SNR_dB,fre_eff_persnr(1,:),'r-*',SNR_dB,fre_eff_persnr(2,:),'g-o',SNR_dB,fre_eff_persnr(3,:),'m-+')
% plot(SNR_dB,fre_eff_persnr(2,:),'g-o')
% plot(SNR_dB,fre_eff_persnr(3,:),'m-')
% xlabel('SNR(dB)');
% ylabel('Spectral Efficiency(bps/Hz)');
% hold off