%===================================================================

%===================================================================
% Alamouti.m
% ==========
% Description: This m-file simulates the SER performance of Alamouti
%              space-time shceme using maximum likelihood decoding
% Usage:  ser_ml=Alamouti(Nr,max_snr,max_err_symbol,symbol_per_frame)
% Inputs:
%         Nr: the number of receiving antenna,either 1 or 2
%         max_snr: simulation from 0 to max_snr
%         max_err_symbol: the simulation will stop when the number of
%                         received error symbols exceeds max_err_symbol
%         symbol_per_frame: the number of symbols to transmit in a data
%                            frame and we assume the channel fading
%                            coefficients will keep the same the data
%                            frame transmitting period.Becasue of 2 Tx
%                            antenna, symbol_per_frame must be even
% Outputs:
%         ser_ml: symbol-error-ratio of each snr
% ======================================================================
function ser_ml=Alamouti(Nr,max_snr,max_err_symbol,symbol_per_frame);
Nr=1;max_snr=12;max_err_symbol=200;symbol_per_frame=24;
%clc;
%發(fā)射天線數(shù)，固定為2
if mod(symbol_per_frame,2) == 0
    Nt = 2;
    tot_frame_tx=zeros(1,max_snr+1);
    tot_err_symbol =zeros(1,max_snr+1);
    tot_frame_tx = zeros(1,max_snr+1);
    for snr=0:max_snr
        while tot_err_symbol(snr+1)<max_err_symbol
            tot_frame_tx(snr+1) = tot_frame_tx(snr+1)+1;
            %生成偽隨機(jī)二進(jìn)制待發(fā)送序列
            psudo_bit_sequence = rand(1,symbol_per_frame * 2)>0.5;
            %二進(jìn)制信道與QPSK信號(hào)的轉(zhuǎn)換
            qpsk_sequence = QpskMapping(psudo_bit_sequence);
            %噪聲能量（已歸一化）
%           sig = sqrt(0.5/(10^(snr/10)));
            sig = sqrt(1/(10^(snr/10)));
            %生成加性高斯噪聲
            addtional_gaussian_noise   =   sig * (randn(Nr,Nt,symbol_per_frame) + j*randn(Nr,Nt,symbol_per_frame));
            %生成具有瑞利分布特性的信道衰落系數(shù)
            H=RayleighCH(Nr,Nt);
            for k=1:2:symbol_per_frame-1
                X=[qpsk_sequence(k) -conj(qpsk_sequence(k+1)); qpsk_sequence(k+1) conj(qpsk_sequence(k))];
                %接收信號(hào)
                R=H*X + addtional_gaussian_noise(:,:,k);
                %接收信號(hào)的合成（利用發(fā)送矩陣的正交性對(duì)信道衰落矩陣進(jìn)行正交化處理），由Alamoui碼原理可知，信號(hào)合成方式因接收天線數(shù)目而異
                switch Nr
                    case 1
                        s0=conj(H(1))*R(1)+H(2)*conj(R(2));
                        s1=conj(H(2))*R(1)-H(1)*conj(R(2));
                    case 2
                        s0=conj(H(1,1))*R(1,1)+H(1,2)*conj(R(1,2))+conj(H(2,1))*R(2,1)+H(2,2)*conj(R(2,2));  
                        s1=conj(H(1,2))*R(1,1)-H(1,1)*conj(R(1,2))+conj(H(2,2))*R(2,1)-H(2,1)*conj(R(2,2));
                end
                S=[s0 s1];
                %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%    
                %          最大似然譯碼          %
                %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                dh = [1+j -1+j -1-j 1-j]/sqrt(2);
                D1=abs(S(1)*[1 1 1 1]-dh).^2;
                [minScale positionmin]=min(D1);
                X_dec1=dh(positionmin);
                D2=abs(S(2)*[1 1 1 1]-dh).^2;
                [minScale positionmin]=min(D2);
                X_dec2=dh(positionmin);
                X_dec=[X_dec1;X_dec2];
                if sum([round(X(1,1))~=round(X_dec1) round(X(2,1))~=round(X_dec2)])
                   tot_err_symbol(snr+1)=tot_err_symbol(snr+1)+1;
                   break
                end
            end
        end
    end        
    %在不同信噪比條件下，求系統(tǒng)SER（誤符號(hào)率）
    ser_ml=tot_err_symbol./(tot_frame_tx * symbol_per_frame);
%    [minScale positionmin] = min(ser_ml);
    semilogy(0:max_snr,ser_ml);
    xlabel('SNR');
    ylabel('Symbol Error Ratio(SER)');
    %axis([0 max_snr 10^(round(log10(minScale)-0.5)) 1]);
else
    disp('symbol_per_frame must be even!');
    ser_ml = 0;
end