www.gusucode.com > MRC系统误码率仿真源码程序 > MRC系统误码率仿真源码程序/comparison_Alamouti_MRC/alamouti_new.m

    function BER2m = alamouti_new(M, frmLen, numPackets, EbNo)
%   alamouti_new:  Orthogonal space-time block coding for 2xM antenna configurations.
%
%   BER2m = alamouti_new(M, frmLen, numPackets, EbNo) computes the bit-error rate 
%   estimates via simulation for an orthogonal space-time block coded 
%	configuration using two transmit antennas and M receive antennas, where the
%	frame length, number of packets simulated and the Eb/No range of values are
%	given by M, frmLen, numPackets, EbNo parameters respectively.
%
%   The simulation uses the full-rate Alamouti encoding scheme for BPSK 
%   modulated symbols with appropriate receiver combining.
%
%   Suggested parameter values:
%       M = 1 or 2; frmLen = 100; numPackets = 1000; EbNo = 0:2:20;

% Simulation parameters
N = 2;              % Number of transmit antennas
rate = 1; inc = N/rate; repFactor = 2;

% Create BPSK mod-demod objects
P = 2; 	% modulation order 
bpskmod = modem.pskmod('M', P, 'SymbolOrder', 'Gray', 'InputType', 'Integer');
bpskdemod = modem.pskdemod(bpskmod);

% Pre-allocate variables for speed
txEnc = zeros(frmLen/rate, N); r  = zeros(frmLen/rate, M);
H  = zeros(frmLen/rate, N, M);
z = zeros(frmLen, M); z1 = zeros(frmLen/N, M); z2 = z1;
error2m = zeros(1, numPackets); BER2m = zeros(1, length(EbNo));

% Loop over EbNo points
for idx = 1:length(EbNo)
    % Loop over the number of packets
    for packetIdx = 1:numPackets
        data = randi([0 P-1], frmLen, 1);         % data vector per user/channel
        tx = modulate(bpskmod, data);        % BPSK modulation

        % Alamouti Space-Time Block Encoder, G2, full rate
        %   G2 = [s1 s2; -s2* s1*]
        s1 = tx(1:N:end); s2 = tx(2:N:end);
        txEnc(1:inc:end, :) = [s1 s2];
        txEnc(2:inc:end, :) = [-conj(s2) conj(s1)];

        % Create the Rayleigh channel response matrix
        H(1:inc:end, :, :) = (randn(frmLen/rate/repFactor, N, M) + ...
                                1i*randn(frmLen/rate/repFactor, N, M))/sqrt(2);
        %   held constant for repFactor symbol periods
        H(2:inc:end, :, :) = H(1:inc:end, :, :);

        % Received signal for each Rx antenna
        for i = 1:M
            % with power normalization
            r(:, i) = awgn(sum(H(:, :, i).*txEnc, 2)/sqrt(N), EbNo(idx)); 
        end

        % Combiner - assume channel response known at Rx
        hidx = 1:inc:length(H);
        for i = 1:M
            z1(:, i) = r(1:inc:end, i).* conj(H(hidx, 1, i)) + ...
                       conj(r(2:inc:end, i)).* H(hidx, 2, i);

            z2(:, i) = r(1:inc:end, i).* conj(H(hidx, 2, i)) - ...
                       conj(r(2:inc:end, i)).* H(hidx, 1, i);
        end
        z(1:N:end, :) = z1; z(2:N:end, :) = z2;
                
        % ML Detector (minimum Euclidean distance)
        demod2m = demodulate(bpskdemod, sum(z, 2));

        % Determine bit errors
        error2m(packetIdx) = biterr(demod2m, data); 
    end % end of FOR loop for numPackets

    % Calculate BER for current idx
    BER2m(idx) = sum(error2m)/(numPackets*frmLen);

end  % end of for loop for EbNo