stbc_ofdm_est.asv

是关于空时编码的nrs

%==========================================================================
% "Channel estimation via training in time domain and tracking in frequcy 
% domain for STBC-OFDM Systems"
% Sun Jie-ting, 08.11.2004, East China Uni. Science and Technology,Shanghai
% Results of SER
%==========================================================================

%clc;
function STBC_OFDM_EST(NoOFDMSymbol,NFFT,G,t_a,miu,training_gap,f_dmax,snr_min,snr_max)

Nt = 2;
Nr = 1;
%NFFT = 128;               % FFT length
%G = 10;                   % Guard interval length

%t_a = 0.8*10^(-6);       % Sampling duration of wireless communication
%miu = 0.05;
%training_gap = 10;
%NoOFDMSymbol = 40;        %Number of OFDM symbol
OFDMSymbolbits = NFFT;      %Number of bits per OFDM symbol conveys
length_data = NoOFDMSymbol * NFFT;  
                            % The total data length

%-------------------------------------------------
% Parameters for Monte Carlo channel
%-------------------------------------------------



symbol_duration = NFFT * t_a;   %OFDM symbol duration
number_of_summations = 40;      % Number of summations for Monte-Carlo method



%f_dmax = 50.0;                  % Maximum Doppler frequency

load h_decimation.am -ascii;
h1_initial = h_decimation;

h2_initial = h_decimation;

N_P  = length(h_decimation);
 


u1 = rand(N_P,number_of_summations); % A random variable
u2 = rand(N_P,number_of_summations); % A random variable
%---------------------------
% Data for STBC demodulation
%---------------------------
data_for_destbc = zeros(2,1,16);
data_qpsk = [1+j;-1+j;-1-j;1-j]/sqrt(2);
for k=1:4
    for m=1:4
        data_for_destbc(:,:,4*(k-1)+m) = [data_qpsk(k);data_qpsk(m)];
    end 
end
%--------------------------
% Optimal Training Sequence
%--------------------------
Training_OFDM_T1 = zeros(1,OFDMSymbolbits);
Training_OFDM_T2 = zeros(1,OFDMSymbolbits);
for i=1:OFDMSymbolbits
    if mod(i,4)==0
        Training_OFDM_T1(i) = data_qpsk(4);
    else
        Training_OFDM_T1(i) = data_qpsk(mod(i,4));
    end
    Training_OFDM_T2(i) = Training_OFDM_T1(i)*exp(j*2*pi*(i-1)*(N_P)/NFFT);
end

%----------------------------
% QPSK modulator in base band
%----------------------------------------------
BitFrameC1=(2*round(rand(Nt,OFDMSymbolbits/2,NoOFDMSymbol))-1)/sqrt(Nt);
QpskDataC1=BitFrameC1(1,:,:)+j*BitFrameC1(2,:,:);

BitFrameC2=(2*round(rand(Nt,OFDMSymbolbits/2,NoOFDMSymbol))-1)/sqrt(Nt);
QpskDataC2=BitFrameC2(1,:,:)+j*BitFrameC2(2,:,:);
%--------------------
% Alamouti space-time modulator in base band
%--------------------
Stbc1=zeros(1,OFDMSymbolbits,NoOFDMSymbol);
Stbc2=zeros(1,OFDMSymbolbits,NoOFDMSymbol);
for k=1:1:OFDMSymbolbits/2
    Stbc1(1,2*k-1,:)=QpskDataC1(1,k,:);
    Stbc1(1,2*k,:)=-conj(QpskDataC2(1,k,:));
    Stbc2(1,2*k-1,:)=QpskDataC2(1,k,:);
    Stbc2(1,2*k,:)=conj(QpskDataC1(1,k,:));
end

%------------------------------------------------------------------------
% Transmitted signal of trasmitt antena 1 to receive antenna (channel:h1)
%------------------------------------------------------------------------
%snr_min =12;
%snr_max =12;

step = 1;
symbol_error_rate = zeros(1,(snr_max-snr_min)/step+1);
for snr = snr_min:step:snr_max; 
    disp('snr = ');
    disp(snr);


h1_frame = [];
rs1MP_frame = [];
initial_time=0;                 % Initial time

for i=1:NoOFDMSymbol;

   if mod(i,training_gap) == 0  | i==1
       disp('Insert Training Sequence');
       OFDM_signal_tem = OFDM_Modulator(Training_OFDM_T1,NFFT,G);
   else 
       OFDM_signal_tem = OFDM_Modulator(Stbc1(1,:,i),NFFT,G);
   end
   % OFDM signal from the first transmitt antenna is created
   
  [h1, t] = MCM_channel_model(u1, initial_time, number_of_summations, symbol_duration, ...,
       f_dmax, h1_initial);
  h1_frame = [h1_frame;h1];
  initial_time = t;
   
   rs1MP = conv(OFDM_signal_tem, h1);
%   H1_i = fft([h1 zeros(1,NFFT-N_P)])
%   pause
%   r3 = Stbc1(1,:,i).*H1_i
%   pause
   % The received signal over multhipath channel is created
%   rs1 = awgn(rs1,snr,'measured','dB');
   rs1MP_frame = [rs1MP_frame; rs1MP];
   clear OFDM_signal_tem;
end;


%------------------------------------------------------------------------
% Transmitted signal of antenna 2 to receive antenna 1 (channel: h2)
%-----------------------------------------------------------------------

h2_frame = [];
rs2MP_frame = [];

initial_time=0;                 % Initial time
for i=1:NoOFDMSymbol;

   if mod(i,training_gap) == 0 | i==1
       disp('Insert Training Sequence');
       OFDM_signal_tem = OFDM_Modulator(Training_OFDM_T2,NFFT,G);
   else 
       OFDM_signal_tem = OFDM_Modulator(Stbc2(1,:,i),NFFT,G);
   end
   

   % OFDM signal from the second antenna is created
   [h2, t] = MCM_channel_model(u2, initial_time, number_of_summations, symbol_duration, ...,
       f_dmax, h2_initial);
   h2_frame = [h2_frame;h2];
   initial_time = t;
%   H2_i = fft([h2 zeros(1,NFFT-N_P)])
%   pause
%   r4 = Stbc2(1,:,i).*H2_i
%   pause   
   rs2MP = conv(OFDM_signal_tem, h2);
   % The received signal over multhipath channel is created
%   rs2 = awgn(rs2,snr,'measured','dB');
   rs2MP_frame = [rs2MP_frame; rs2MP];
   clear OFDM_signal_tem;
end;


%-----------------
% Recever: OFDM demodulator, channel estimation
%-----------------


for i=1:NoOFDMSymbol;
    disp('Number of Frame');
    disp(i);

    rsBeforeAwgn_i = rs1MP_frame(i,:) + rs2MP_frame(i,:);
%    rs_i =   rsBeforeAwgn_i;
   rs_i = awgn(rsBeforeAwgn_i,snr,'measured','dB');
   Demodulated_rs_i = OFDM_Demodulator(rs_i,NFFT,NFFT,G);

    if mod(i,training_gap)==0 | i==1
       disp('Training Estimation');
       p1 = zeros(1,N_P);
       p2 = zeros(1,N_P);       
       for l = 1:N_P
           for k = 1:NFFT
               p1(l) = p1(l) + Demodulated_rs_i(k)*conj(Training_OFDM_T1(k))*exp(j*2*pi*(k-1)*(l-1)/NFFT);
               p2(l) = p2(l) + Demodulated_rs_i(k)*conj(Training_OFDM_T2(k))*exp(j*2*pi*(k-1)*(l-1)/NFFT);               
           end
       end
       h1_training_estimation = 1/NFFT*p1;
       h2_training_estimation = 1/NFFT*p2;
       H1_i_estimation = fft([h1_training_estimation zeros(1,NFFT-N_P)]);
       H2_i_estimation = fft([h2_training_estimation zeros(1,NFFT-N_P)]);
%       H1_i_estimation = zeros(1,NFFT);
%       H2_i_estimation = zeros(2,NFFT);       

%       for l = 1:NFFT
%           for k = 1:N_P
%               H1_i_estimation(l) = H1_i_estimation(l) + h1_training_estimation(k)*exp(-j*2*pi*l*(k-1)/NFFT);
%               H2_i_estimation(l) = H2_i_estimation(l) + h2_training_estimation(k)*exp(-j*2*pi*l*(k-1)/NFFT);               
%           end
%       end
       H1_i = fft([h1_frame(i,:) zeros(1,NFFT-N_P)]);
       H2_i = fft([h2_frame(i,:) zeros(1,NFFT-N_P)]);
    else

%       H1_i = fft([h1_frame(i,:) zeros(1,NFFT-N_P)]);
%       H2_i = fft([h2_frame(i,:) zeros(1,NFFT-N_P)]);
        H1_i = H1_i_estimation;
        H2_i = H2_i_estimation;
 %--------------------------
 % Demodulated signal 
 %--------------------------
        for l=1:OFDMSymbolbits/2;
           R1=Demodulated_rs_i(2*l-1);
           R2=Demodulated_rs_i(2*l);
           R2conj = conj(R2);
           R_n2l=[R1;R2conj];
           H_n2l = [H1_i(2*l-1) H2_i(2*l-1);conj(H2_i(2*l)) -conj(H1_i(2*l))];
           minNorm = 10000;
           num_min = 0;
           for m=1:16
               dec_norm = norm(R_n2l-H_n2l*data_for_destbc(:,:,m))^2;
               if dec_norm < minNorm;
                   minNorm = dec_norm;
                   num_min = m;
               end
            end
            deStbc = data_for_destbc(:,:,num_min);
            estimated_OFDM = [deStbc(1) -conj(deStbc(2));deStbc(2) conj(deStbc(1))];
            H1_i_estimation(2*l-1) = miu*R1/estimated_OFDM(1,1) - miu*estimated_OFDM(2,1)*H2_i(2*l-1)/estimated_OFDM(1,1) + (1-miu)*H1_i(2*l-1);
            H2_i_estimation(2*l-1) = miu*R1/estimated_OFDM(2,1) - miu*estimated_OFDM(1,1)*H1_i(2*l-1)/estimated_OFDM(2,1) + (1-miu)*H2_i(2*l-1);
            H1_i_estimation(2*l) = miu*R2/estimated_OFDM(1,2) - miu*estimated_OFDM(2,2)*H2_i(2*l)/estimated_OFDM(1,2) + (1-miu)*H1_i(2*l);
            H2_i_estimation(2*l) = miu*R2/estimated_OFDM(2,2) - miu*estimated_OFDM(1,2)*H1_i(2*l)/estimated_OFDM(2,2) + (1-miu)*H2_i(2*l);
            TrSTBC = [Stbc1(1,2*l-1,i);Stbc2(1,2*l-1,i)];
            if deStbc==TrSTBC
                disp('ok');
            else
%                disp('No');
                symbol_error_rate(snr-snr_min+1) = symbol_error_rate(snr-snr_min+) + 1;
            end
        end;
    end;
end;
end;

ser = symbol_error_rate/length_data;
semilogy(snr_min:step:snr_max, ser,'b*');

ylabel('SER');
xlabel('SNR');
data = [snr_min:step:snr_max; ser];
save STBC_OFDM_EST.am data -ascii;
clear all;