rf.m

mimo2x2天线选择系统的全系统matlab程序,先前的是dsp程序.

function RF(real_ch,SNR_dB)

close all;

% If argument is omitted, the channel is simulated
if (nargin < 1)
    real_ch = 0;
end
if (nargin <2)
    SNR_dB = 20;
end

disp('-- Start RF --')

% Load the baseband channels
load('TXOutput');

% Load the global parameters
load('tx_param');

save('SimType','real_ch');
% Build the sent stream
S = [channel1_block; channel2_block];

% Output of the baseband channel
[Nr,Nc] = size(S);
Y = [];

% For plotting
h11 = zeros(1,Nc);

% ---------------------------------------------------
% Simulate a channel in baseband
% ---------------------------------------------------

if (real_ch == 0)
    
    % Time-varying channel
    step = 2*10^(-5);
  
    % Channel Matrix
    
    % H = [1 0; 1 0; 0 1; 0 1];
    
    % H = sqrt(1/2)*ones(4,2);
%    H = sqrt(1/2)*randn(4,2) + j*sqrt(1/2)*randn(4,2);
    H =   [   0.5189 + 0.0943i    0.7957 - 0.0453i;
             -1.1536 - 0.2596i   -0.3614 - 0.7444i;
              0.3960 - 0.7802i    0.6947 + 0.4259i;
              0.8576 - 0.0625i   -0.9832 - 0.3027i
          ];
%       
%     % Real channel
% H = 1/2 * [         0.2223 + 0.2129i   0.0452 - 0.2562i
%                     0.1307 + 0.2811i   0.1384 - 0.2228i
%                     0.0035 - 0.0033i   0.0099 + 0.0217i
%                     0.0046 - 0.0019i   0.0013 + 0.0240i ];

    % Keep a recoed of all the channel matrices
    H_var = zeros(4,2*frames_len);
    
    % Proceed frame by frame
    index = 1;
    for fr=1:nr_frames
                 
        % Apply the channel matrix (MUI)
        H = H + step*(sqrt(1/2)*randn(4,2) + j*sqrt(1/2)*randn(4,2));
        H_var(:,2*fr-1:2*fr) = H;
        
        % Number of previous frames
        if mod(fr-1,4)==0
            Y = [Y H*S(:,index:index+Lf_sync-1)];
            index = index+Lf_sync;
        else
            Y = [Y H*S(:,index:index+Lf-1)];
            index = index+Lf;
        end
        
        h11(index) = H(1,1);
  
    end

    save('ChannelMatrix','H_var');
    
    % See page 6 to 11 in the DC Project Assignement
    % Es_bb Energy of a symbol in baseband
    Es_bb_th = 10 * norm(pulse_shape)^2/L;
    
    % Es energy of a symbol in the upconverted signal
    Es_th = Es_bb_th/2;
    
    % SNR = Es/N0
    % Where Sn(f) = N0/2 power spectral density of the real noise
    N0 = Es_th / 10^(SNR_dB/10);
    sigma_sqr_bb = 2*N0*L;

    disp(['Theoritical Es:             ', num2str(Es_th)])
    disp(['Variance of the bb noise:   ', num2str(sigma_sqr_bb)]);
    disp(['Eb/N0 (SNR)                 ', num2str(10*log10(Es_th/4/N0))]);
%     disp(['Es/N0 (SNR)                 ', num2str(10*log10(Es_th/N0))]);
%     disp('');
    disp(['SNR (bb)                    ', num2str(10*log10(Es_bb_th/sigma_sqr_bb))]);

    % Uncorrelated noise
    Noise = sqrt(sigma_sqr_bb/2) * ( randn(4,Nc) + j*randn(4,Nc) );
    
    % Correlated noise
    % Noise = sqrt(sigma_sqr_bb/2) * ( ones(4,1) * randn(1,Nc) + j * ones(4,1) *randn(1,Nc) );
    
    Y_n = Y + Noise;

else
    
    Y_n = S;
    
end

% ---------------------------------------------------
% Upconversion
% ---------------------------------------------------

N_channels = size(Y_n,1);

% Upconvert
for k=1:N_channels
    Y_up(k,:) = upconv(real(Y_n(k,:)), imag(Y_n(k,:)), Fc_tx, Fs, 0);
end

% Add the sinusoid for frequency offset estimation
index = 0:time_end;
cosinus = 3*cos(2*pi*(Fc_tx/Fs)*index);

if (real_ch==0)
    % Noise:
    sigma_sqr = N0/2*L;
    sigma = sqrt(sigma_sqr);    
    disp(['Variance of the real noise: ', num2str(sigma_sqr)]);
else
    % For the real channel, don't add more noise for the tail and queue!
    sigma = 0;
end

Y_R = [sigma * randn(N_channels,init_len), ones(N_channels,1)*cosinus, Y_up, sigma * randn(N_channels,tail_len)];

% ---------------------------------------------------
% Send through the real channel
% ---------------------------------------------------

if(real_ch == 1)
   
	% Use the real channel
    N = size(Y_R,2);
	
    % Scaling to use full power
	scaling = 1/(1.2*max(max(Y_R))) * 32766;
    	
    figure;
    Nc = 2;
    Nr = 2;
    
	subplot(Nr,Nc,1);
    plot(Y_R(1,:))
	title('Channel 1 (Left)');
	
    subplot(Nr,Nc,2);
    plot(Y_R(2,:))
	title('Channel 2 (Right)');
    
    Y_R = Y_R * scaling;

    
    
	[Xin,no_samples_delivered,buffer_size]=DSK6713_rxtx(N,2,2,Y_R,Fs);
	
    % Scale down by the same factor used for scaling up
    Xin = Xin/scaling;
    
    N_out = size(Xin,2);
    
    % Put Left in Y_R(1) and Right in Y_R(3)
    Y_R = [Xin(1,:) ; zeros(1,N_out); Xin(2,:); zeros(1,N_out)];
	
    subplot(Nr,Nc,3);
    plot(Y_R(1,:))
	title('Channel A [Y_R(1)] (Left)');
	
    subplot(Nr,Nc,4);
    plot(Y_R(3,:))
	title('Channel B [Y_R(3)] (Right)');
	
end

save('ChannelOutput','Y_R');
disp('End')

% ---------------------------------------------------
% Plots
% ---------------------------------------------------

if (real_ch == 0)
	
    figure;
    Nc = 2;
    Nr = 2;

    subplot(Nc,Nr,1);
	
	plot(real(S(1,:)));
	title('Signal sent on ant. 1')
    
	subplot(Nc,Nr,2);
	plot(real(Y(1,:)));
	hold on;
	plot(real(Y_n(1,:)),'g');
	legend('I rec. on ant. 1','Noisy I rec. on ant. 1');
	
	subplot(Nc,Nr,3);
	plot(Y_R(1,:))
	title('Signal rec. on ant. 1')
	
    subplot(Nc,Nr,4);
	plot(abs(h11))
	title('Magnitude of H11')
    
end