%ofdm1.m
%QPSK調制、AWGN信道下的仿真程序
%***************************variables********************************
para = 52;%并行信道數
fftlen = 64;%FFT長度
noc = 52;%載波數
nd = 6;%每個仿真循環的OFDM數
m1 = 2;%調制水平：QPSK
gilen = 16;%保護間隔長度
ebn0 = 15;%信噪比
sr = 250000;%OFDM symbol rate (250k syombol/s)
br = sr.*m1;%bit rate per carrier
%**************************主循環部分*********************************
nloop = 200;%仿真的循環數
noe1 = 0;%信道解碼前錯誤數據數
nod1 = 0;%信道解碼前傳輸數據數
%*************************發射機**************************************
for iii = 1:nloop
    %****************************數據產生******************************
    xx = ones(1,para*nd*m1);
    %**************************串并轉換*********************************
    paradata = reshape(xx,para,nd*m1);
    %************************QPSK調制**********************************
    [ich,qch] = qpskmod1(paradata,para,nd,m1);
    kmod = 1/sqrt(2);
    ich = ich.*kmod;
    qch = qch.*kmod;
    %數據映射
    [ich1,qch1] = cramapping(ich,qch,fftlen,nd);
    %*******************************IFFT********************************
    x = ich1 + qch1.*i;
    y = ifft(x);
    ich2 = real(y);
    qch2 = imag(y);
    %**************************插入保護間隔******************************
    [ich3,qch3] = giins1(ich2,qch2,fftlen,gilen,nd);
    fftlen2 = fftlen + gilen;
    %**************************attenuation calculation*********************
    spow = sum(ich3.^2+qch3.^2)/nd./para;
    attn = 0.5*spow*sr/br*10.^(-ebn0/10);
    attn = sqrt(attn);
    %**************************AWGN addition******************************
    [ich4,qch4] = comb(ich3,qch3,attn);
    %**********************去除保護間隔***********************************
    [ich5,qch5] = girem1(ich4,qch4,fftlen2,gilen,nd);
    %*************************FFT*****************************************
    rx = ich5 + qch5.*i;
    ry = fft(rx);
    ich6 = real(ry);
    qch6 = imag(ry);
    %載波逆映射
    [ich7,qch7] = crdemapping(ich6,qch6,fftlen,nd);
    %******************************解調**********************************
    ich7 = ich7./kmod;
    qch7 = qch7./kmod;
    [demodata] = qpskdemod1(ich7,qch7,para,nd,m1);
    %*****************************并串轉換*******************************
    demodata1 = reshape(demodata,1,para*nd*m1);
    %****************************bit error rate(BER)*********************
    noe10 = sum(abs(demodata1 - xx));
    nod10 = length(xx);
    
    noe1 = noe10 + noe1;
    nod1 = nod10 + nod1;
    
    fprintf('%d\t%e\t\n',iii,noe10/nod10);
end

ber1 = noe1/nod1;
%********************************輸出結果**********************************
fprintf('%f\t%e\t%d\t\n',ebn0,ber1,nloop);
%********************************繪圖**************************************
%plot(ebn0,ber1);
%*******************************end of file********************************