gmqam_type2_ofdm_v2.m

4 16 64 QAM plots for simulated and theoretical results for OFDM system

%  Title       : QAM Constellation plots for 4, 16, 64
%  Author      : goutam.ramamurthy.06<AT>student.lth.se
%  History     : 2009-01-30      created
%  References  :
%  Description :
%
%        1         2         3         4         5         6         7         
%23456789012345678901234567890123456789012345678901234567890123456789012345
%

close all;
clear all;
clc;

%%%%%%%%%% GLOBAL SYSTEM PARAMETERS 
%--Case-Specific Simulation Parameters
idx = 1;
SIMU{idx}.name = '4QAM';
SIMU{idx}.M_ary = 4;
SIMU{idx}.legend = ' Coherent 4QAM AWGN Simulation';
SIMU{idx}.linestyle = '-gd'; % linestyle: green Squares
SIMU{idx}.title ='Ps(Sym Err) curve for 4QAM modulation';
idx = idx +1;

SIMU{idx}.name = '16QAM';
SIMU{idx}.M_ary = 16;
SIMU{idx}.legend = ' Coherent 16QAM AWGN Simulation';
SIMU{idx}.linestyle = 'rd-'; % linestyle: red diamonds
SIMU{idx}.title ='Ps(Sym Err) curve for 16QAM modulation';
idx = idx +1;
 
SIMU{idx}.name = '64QAM';
SIMU{idx}.M_ary = 64;
SIMU{idx}.legend = ' Coherent 64QAM AWGN Simulation';
SIMU{idx}.linestyle = '-md'; % linestyle: maroon diamonds
SIMU{idx}.title ='Ps(Sym Err) curve for 64QAM modulation';

%--Common Parameters 
SNRdB = -2:2:24;       % Es/N0 range for simulation
bSNRdB = zeros(length(SIMU), length(SNRdB));
channel = 'PURE_AWGN'; % PURE_AWGN; RAYL;
ON_noise = 1; % 1 = Noise ON; 0=OFF

%--OFDM Parameters 
N_OFDMsyms = 20;      % Choose the total # of symbols in 1RB to simulateTXRX
BW_LTE = 20e6;        % Total Available LTE S/M bandwidth is 20 MHz
BW_guard = BW_LTE/10; % 10% of Total available BW for OOB emissions 
delta_f = 15e3;       % LTE Sub-Carrier Spacing; specified by 3GPP
                      % also delta_f_low = 7.5e3us used for MBMS
                        
%--Derived Parameters
N_c = floor((BW_LTE - BW_guard)/delta_f);   % 1200 for 20MHz BW
m = ceil(log2(N_c));    % radix value for N=2^m in N-IFFT
N_FFT = 2^m;            % N > N_c & for radix-2 IFFT, S.T. N is 2^m;

%~~~~~~~~~~~~~~

tic;
for iSIMU = 1:length(SIMU)     % iterating through Simulations
    M_ary = SIMU{iSIMU}.M_ary; % number of constellation points
    k = log2(M_ary);
    Es = (M_ary-1)*2/3;        % Avg. Ener per QAM sym, norm = sqrt(Es)
    d_minsq = 3*log2(M_ary)/(M_ary-1); % for M-QAM
    
    N_qsyms = N_c * N_OFDMsyms ;
    N_bits = N_qsyms * k;      % an integral multiple of k channels
    
    p = (1:sqrt(M_ary)/2);        % In rect QAMs, for M=2^n & n ~ {even I}
    a_MQAM = [-(2*p-1) (2*p-1)];  % alphabets
    
    for iSNRdB = 1:length(SNRdB) % iterating through SNRs
        fprintf('\nEs/N0: %d dB\n', SNRdB(iSNRdB));
        SNR = 10^(SNRdB(iSNRdB)/10); % Convert to non-dB value
        bSNRdB(iSIMU,iSNRdB) = SNRdB(iSNRdB) - 10*log10(k);
        
        %%%%%%%%%% RANDOM SOURCE
        %--Generate a column of random binary bits
        x_txbits = (1-sign(randn(N_bits,1)))/2; % ~Bernoulli source(p=1/2)
        %%%RANDOM SOURCE-ED
        
        %%%%%%%%%% QAM MAPPER
        %--using matlab's modem object generate QAM syms and normalize
        qam_mod = modem.qammod(M_ary);    % create a modulator object
        qam_mod.symbolorder = 'gray' ;  % Enable/disable gray coding
        qam_mod.inputtype = 'bit';
        qam_mod.disp;                   % Print the properties of modulator
        q_txsyms = modulate(qam_mod, x_txbits);  % QAM modulate
        q_txsyms_norm = q_txsyms/sqrt(Es);
        %%% QAM MAPPER -ED

        %%%%%%%%%% OFDM MOD
        Ofdm_sz = N_FFT*N_OFDMsyms;
        O_train = zeros(Ofdm_sz,1);          % pre-declaration for speed
        for i = 1:N_OFDMsyms %N_qsyms/N_c    % iterate over N OFDMsyms
            txsyms_ext = zeros(N_FFT, 1);    % pad zeros to unused FFT bins
            txsyms_ext(1:N_c/2) = q_txsyms_norm((i-1)*N_c+1 : ...
                                  (i-1)*N_c + N_c/2).';
            txsyms_ext(N_FFT-(N_c/2-1):N_FFT) = q_txsyms_norm( ...
                                                (i-1)*N_c+(N_c/2+1):...
                                                (i-1)*N_c+N_c).';
            O_sym = sqrt(N_FFT).*ifft(txsyms_ext, N_FFT);
            O_train(((i-1)*N_FFT)+1:((i-1)*N_FFT)+N_FFT, 1) = O_sym(:,1);

        end
        %%%OFDM MOD -ED         
        
        
        %++++++++++
        %%%%%%%%%% THE CHANNEL TO CONVOLVE
        noise_len = size(q_txsyms_norm);
        re = randn(noise_len,1);
        im = randn(noise_len,1);
        awgn = 1/sqrt(2)*(re + j*im); % WGN with 0dB variance
        
        if strcmp(channel, 'RAYL') ~= 0
            h_t = (randn(noise_len,1)+j*randn(noise_len,1))/sqrt(2);
            rxdata = h_t.*O_train + ...
                           ON_noise*(1/sqrt(SNR))*awgn;
            H = fft(h_t,N_FFT);  %%%% Known Channel
        else % if strcmp(channel, 'PURE_AWGN') ~= 0
            rxdata = q_txsyms_norm + ON_noise*(1/sqrt(SNR))*awgn; % AWGN
        end
        
        %%%CHANNEL CONVOLUT-ED 
        %++++++++++
        
        
        
        %%%%%%%%%% OFDM DEMOD
        O_rxtrain = [];          % pre-declaration
        rxqsyms = [];            % pre-declaration
        for i = 1:N_OFDMsyms % iterate over N OFDMsyms
            offset = ((i-1)*N_FFT);
            %%%%%%%%%% REM CP
            
            %%%REM CP -ED
            O_rxsym = (1/sqrt(N_FFT)).*fft(rxdata(offset+1:offset+N_FFT),N_FFT);
            if strcmp(channel, 'PURE_AWGN') == 0
                O_rxsym =  O_rxsym./H;
            end
            O_rxtrain = [O_rxtrain; O_rxsym];  % DEL: for debug
            rxqsyms = [rxqsyms;
                       O_rxsym(1:N_c/2); 
                       O_rxsym((N_FFT-N_c/2+1):N_FFT)];
        end
        %%%OFDM DEMOD -ED 
        
        
        %%%%%%%%%% QAM DEMAPPER
        %--using matlab's modem object
        qam_rxsyms_denorm = rxdata * sqrt(Es);
        
        %%%%%%%%%% RECT QAM MAPPING
        %--Segregate Re-Im
        q_rxsyms_re = real(qam_rxsyms_denorm);
        q_rxsyms_im = imag(qam_rxsyms_denorm);
        %--Mapping to the nearest alphabet
        rx_re = floor(q_rxsyms_re/2)*2+1;
        rx_im = floor(q_rxsyms_im/2)*2+1;
        %--mapping boundary conditions for RE and IM
        rx_re(find(max(a_MQAM)<rx_re)) = max(a_MQAM);
        rx_re(find(min(a_MQAM)>rx_re)) = min(a_MQAM);
        rx_im(find(max(a_MQAM)<rx_im)) = max(a_MQAM);
        rx_im(find(min(a_MQAM)>rx_im)) = min(a_MQAM);
        q_rxsyms = rx_re+j*rx_im;
        %%%RECT MAPPED -ED

        qam_demod = modem.qamdemod(qam_mod);  % copy props from mod object
        qam_demod.disp;                     % Prints properties of DEMOD
        x_rxbits = demodulate(qam_demod, qam_rxsyms_denorm); % QAM modulate
        %%%QAM DEMAPPER -ED 

       
        %%%%%%%%%% PERF. MEASUREMENTS
        %--Simulated SER and BER 
        N_symerrs = nnz(q_rxsyms(:)-q_txsyms(:));
        SER(iSIMU,iSNRdB) = N_symerrs/N_qsyms;
        %--Simulated BER
        N_biterrs = nnz(x_rxbits(:)-x_txbits(:));
        BER(iSIMU,iSNRdB) = N_biterrs/N_bits;
                
        %%%PERF. MEASURE -ED
            %--theoretical SER
        EbN0 = SNR/k;           % EbN0 = EsN0/k
        Qfn = qfunc(sqrt(d_minsq * EbN0));
        SER_theo(iSIMU,iSNRdB) = (4*(1-(1/sqrt(M_ary)))*Qfn) - ...
                                  (4*(1-(1/sqrt(M_ary)))^2*(Qfn^2));
        %-- Scatter plots
        figure(iSIMU+10);
        plot(qam_rxsyms_denorm,'.k');
        title(SIMU{iSIMU}.legend);
        axis square;
        axis equal;
        grid on;
    end

    if strcmp(channel, 'RAYL') ~= 0
        ber = berfading(bSNRdB(iSIMU,:),'qam',M_ary,1);
        BER_theo(iSIMU, :)= ber(:);
    else % if strcmp(channel, 'PURE_AWGN')~= 0
        BER_theo(iSIMU, iSNRdB) = (1/k)*SER_theo(iSIMU, iSNRdB); 
    end
end
toc;

%%%%%%%%%% SER PLOTS
figure(1)
for iSIMU = 1:length(SIMU)     % iterating through Simulations
    semilogy(SNRdB, SER_theo(iSIMU, :),'b:','LineWidth',2);
    hold on
    semilogy(SNRdB, SER(iSIMU, :), SIMU{iSIMU}.linestyle, 'Linewidth',1);
end
hold off
axis([0 28 10^-6 1]);
grid on
legend('theoretical', 'Simulated');
xlabel('Es/No, dB');
ylabel('SER(P_s)');
title('SER - Ps versus E_s/N_0');