receive_for_channel3.m

这个程序主要是实现空时编码的matlab编解码仿真。

function [nr_of_errors_in_rate, nr_of_errors_in_rate1, nr_of_errors_in_rate2, theta1, theta2, alpha1, alpha2, sample_nr, indata_est,training1_est,training2_est,i_symbols, q_symbols]=receive_for_channel3(received_vector1,received_vector2,training_sequence1,training_sequence2,data_block,fs,T,fc,pulsetype)

%   [output1, output2]=functionname(input1, input2)
%
%	Variable:	Explanation:
%	nr_of_errors_in_rate    -   Percent of bits with error
%   theta1                  -   the angle estimated by used channel-estimator
%   
%	received_vector     -   received vector from channel (1x1 experiment)
%   training_sequence1  -   used training sequence
%   data_block          -   true data-block sent over channel
%   fs                  -   used sample-frequency
%   T                   -   symbol-time   
%   fc                  -   carrier-frequency
%   pulsetype           -   used pulsetype, 2 for root-raised-cosine
%
%   Short Theoretical Background for the Function:
%
%   Receiving in two antennas!
%   Assuming received_vector is normalized and correct size. Two transmitters used!
%   Used when doing experiment over real acoustic channel the 27th.
%   
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     
%%% Function part of simulation for Space-Time
%%% coding project, group Grey-2001.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
%   Author: Stefan Uppg錼d
%   Date: 27/3-2001
%   Version: 1.0
%   Revision (Name & Date & Comment):
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    
    model = 1;
    training_length1 = length(training_sequence1);
    training_length2 = length(training_sequence2);

    %%% Alamouti coding performed on training-sequence???
    [training_sequence1, s_antenna2]=alamouti(training_sequence1,1);
    [s_antenna1, training_sequence2]=alamouti(training_sequence2,1);
    
    
    receive1 = received_vector1;
    receive2 = received_vector2;
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%% Receiving in receive-antenna 1.
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
    
    %%%
    %%% Make the signal unbiased, removing the signals mean.
    %%%
    
    [receive1]=remove_mean(receive1);
    
    %%%
    %%% Down-conversion and lowpassfiltering.
    %%%
    
    [quad, inphase]=down_converter(receive1,fc,fs,T,(fc/(fs/2)),20);  % Function call for low_pass() /T Grey
    
    figure(2), subplot(2,2,1), plot(1:length(quad),quad);
    title('Quad after down-conv.');
    subplot(2,2,2), plot(1:length(inphase),inphase);
    title('Inphase after down-conv.');
    
     
    %%%%%%%%
    %%%%%%%% Do the frequency-offset compensation if initial-state.
    %%%%%%%%
    
    run_offset = 1;
    if (run_offset == 0)
        [freq_diff, delta_theta, alfa]=freq_offset(inphase, quad, fs, fc)
    end
    
    %freq_diff = 0.0237;
    %delta_theta = 2*pi*(((1/(1-freq_diff/fc))-1)*fc*(1/(fs*(1-(freq_diff/fc)))));
    %delta_theta = 9.3490*10^-6;       % estimated angle-shift/sampel for sinus sampled at 16 kHz.
    delta_theta = -2.00000*10^-5;       % estimated angle-shift/sampel for sinus sampled at 8 kHz 2500 carrier freq.
    
    a = 0:(length(inphase)-1);
    rotation_vector = exp(-j*a*delta_theta);                               % Rotation compensating vector.
     
    compensated_received = complex(inphase,quad) .* rotation_vector;
    inphase = real(compensated_received);
    quad    = imag(compensated_received);
    
    
    
    %%%
    %%% Matchedfiltering and lowpassfiltering.
    %%%
    
    [mf_quad_block]=matched_filter(fs, T, pulsetype, quad);    % Quadrature-part
    [mf_inphase_block]=matched_filter(fs, T, pulsetype, inphase);    % Inphase-part
    
    subplot(2,2,3), plot(1:length(mf_quad_block),mf_quad_block);
    title('Quad after MF');
    subplot(2,2,4), plot(1:length(mf_inphase_block),mf_inphase_block);
    title('Inphase after MF');
    
    
    %%%
    %%% Performing synchronization
    %%%
    
    %[corr_maxi,sample_nri]=synchronization_long(mf_inphase_block, fs, T, training_sequence1, model)  %%% CALLING THE SYNC_LONG FUNCTION!
    %[corr_maxq,sample_nrq]=synchronization_longb(mf_quad_block, fs, T, training_sequence1, model)
    
    [corr_maxi,sample_nri,sum_of_r_times_c1i,sum_of_r_times_d2i]=synchronization2(mf_inphase_block, fs, T, training_sequence1, training_sequence2, model);
    [corr_maxq,sample_nrq,sum_of_r_times_c1q,sum_of_r_times_d2q]=synchronization2(mf_quad_block, fs, T, training_sequence1, training_sequence2, model);
    max_i1 = max(sum_of_r_times_c1i)
    max_i2 = max(sum_of_r_times_d2i)
    max_q1 = max(sum_of_r_times_c1q)
    max_q2 = max(sum_of_r_times_d2q)
    corr_maxi
    sample_nri
    corr_maxq
    sample_nrq
    
    
    if corr_maxq > corr_maxi       % choose the correlation with the biggest max value
        sample_nr = sample_nrq;
    else
        sample_nr = sample_nri;
    end
    
    
    %%%
    %%% Down-Sampling
    %%%
    
    block_length=length(data_block) + length(training_sequence1) + length(training_sequence2);
    [inphase_symbols]=down_sampler(mf_inphase_block,sample_nr,fs,T,block_length);
    [quad_symbols]=down_sampler(mf_quad_block,sample_nr,fs,T,block_length);
    
    figure(3), subplot(1,2,2), plot(mf_inphase_block), hold on
    for(i=0:block_length-1)
        plot(i*(fs*T)+sample_nr,inphase_symbols(i+1),'o');
    end
    hold off
    title('received data and samplepoints. (not yet combined!)');
    subplot(1,2,1), stem([training_sequence1 data_block]), title('Sent data.');
    
    figure(4), subplot(1,2,1), plot(inphase_symbols,quad_symbols,'x'), axis([-2 2 -2 2]), title('Signals after MF and sampled')
    
    %%
    %% Calculating influence of freq_offset. IS INSTEAD DONE AFTER DOWNCONVERTER.
    %%
    
    %[quad_symbols, inphase_symbols, theta]=freq_offset_compensate(quad_symbols,inphase_symbols);
    %theta
    %figure(7), plot(inphase_symbols,quad_symbols,'x'), axis([-6 6 -6 6]), title('Signals after freq-offset compensation.'), grid on
	
    %%%
    %%% Estimating Channel 1 and 2.
    %%%
    
   
    %[alpha1,theta1]=chan_estim7(quad_symbols(1:training_length1),inphase_symbols(1:training_length1),bpsk(training_sequence1,1),1,0.95);
    
    [alpha1,theta1] = channel_estimatorb(quad_symbols(1:training_length1), ...
       inphase_symbols(1:training_length1), training_sequence1, model, 1)
    
    if model == 1
       %[alpha2,theta2]=chan_estim5(quad_symbols(training_length1+1:training_length1+ ...
       %training_length2),inphase_symbols(training_length1+1:training_length1+ ...
       %training_length2), bpsk(training_sequence2,model))
       %[alpha2,theta2]=chan_estim7(quad_symbols(training_length1+1:training_length1+training_length2),...
       %    inphase_symbols(training_length1+1:training_length1+training_length2),bpsk(training_sequence1,1),1,0.95); 
       
       [alpha2,theta2] = channel_estimatorb(quad_symbols(training_length1+1:training_length1+ ...
       training_length2), inphase_symbols(training_length1+1:training_length1+ ...
       training_length2), training_sequence2, model, 2)
    elseif model == 0
        alpha2 = 1;
        theta2 = 0;
    end
    
    q_symbols = quad_symbols;
    i_symbols = inphase_symbols;
    
    % set true values instead:
    %alpha1 = 1;
    %theta1 = 0;
    %alpha2 = 1;
    %theta2 = 0;
    
     %%%
     %%% No combiner used since only 1 antenna is used to send from so far.
     %%%
     
     [combined_signal]=combiner(alpha1(end), theta1(end), alpha2(end), theta2(end), ...
        inphase_symbols, quad_symbols, model);

    figure(4), subplot(1,2,2), plot(combined_signal,'x'), axis([-6 6 -6 6]), title('Signals after combiner, just before detector.')
    grid on
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%% Receiving in receive-antenna 2.
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
    
    %%%
    %%% Make the signal unbiased, removing the signals mean.
    %%%
    
    [receive2]=remove_mean(receive2);
    
    %%%
    %%% Down-conversion and lowpassfiltering.
    %%%
    
    [quad2, inphase2]=down_converter(receive2,fc,fs,T,(fc/(fs/2)),20);
    
    figure(2), subplot(2,2,1), plot(1:length(quad2),quad2);
    title('Quad2 after down-conv.');
    subplot(2,2,2), plot(1:length(inphase2),inphase2);
    title('Inphase2 after down-conv.');
    
     
    %%%%%%%%
    %%%%%%%% Do the frequency-offset compensation if initial-state.
    %%%%%%%%
    
    run_offset = 1;
    if (run_offset == 0)
        [freq_diff2, delta_theta2, alfa2]=freq_offset(inphase2, quad2, fs, fc)
    end
    
    %freq_diff2 = 0.0237;
    %delta_theta2 = 2*pi*(((1/(1-freq_diff/fc))-1)*fc*(1/(fs*(1-(freq_diff/fc)))));
    %delta_theta2 = 9.3490*10^-6;       % estimated angle-shift/sampel for sinus sampled at 16 kHz.
    delta_theta2 = -2.00000*10^-5;       % estimated angle-shift/sampel for sinus sampled at 8 kHz 2500 carrier freq.
    
    a = 0:(length(inphase2)-1);
    rotation_vector2 = exp(-j*a*delta_theta2);                               % Rotation compensating vector.
     
    compensated_received2 = complex(inphase2,quad2) .* rotation_vector2;
    inphase = real(compensated_received2);
    quad    = imag(compensated_received2);
    
    %%%
    %%% Matchedfiltering and lowpassfiltering.
    %%%
    
    [mf_quad_block2]=matched_filter(fs, T, pulsetype, quad2);    % Quadrature-part
    [mf_inphase_block2]=matched_filter(fs, T, pulsetype, inphase2);    % Inphase-part
    
    subplot(2,2,3), plot(1:length(mf_quad_block2),mf_quad_block2);
    title('Quad after MF');
    subplot(2,2,4), plot(1:length(mf_inphase_block2),mf_inphase_block2);
    title('Inphase after MF');
    
    
    %%%
    %%% Performing synchronization
    %%%
    
    
    
    [corr_maxi2,sample_nri2,sum_of_r_times_c1i2,sum_of_r_times_d2i2]=synchronization2(mf_inphase_block2, fs, T, training_sequence1, training_sequence2, model);
    [corr_maxq2,sample_nrq2,sum_of_r_times_c1q2,sum_of_r_times_d2q2]=synchronization2(mf_quad_block2, fs, T, training_sequence1, training_sequence2, model);
    max_i12 = max(sum_of_r_times_c1i2)
    max_i22 = max(sum_of_r_times_d2i2)
    max_q12 = max(sum_of_r_times_c1q2)
    max_q22 = max(sum_of_r_times_d2q2)
    corr_maxi2
    sample_nri2
    corr_maxq2
    sample_nrq2
    
    
    if corr_maxq2 > corr_maxi2       % choose the correlation with the biggest max value
        sample_nr2 = sample_nrq2;
    else
        sample_nr2 = sample_nri2;
    end
    
    
    %%%
    %%% Down-Sampling
    %%%
    
    block_length=length(data_block) + length(training_sequence1) + length(training_sequence2);
    [inphase_symbols2]=down_sampler(mf_inphase_block2,sample_nr2,fs,T,block_length);
    [quad_symbols2]=down_sampler(mf_quad_block2,sample_nr2,fs,T,block_length);
    
    figure(3), subplot(1,2,2), plot(mf_inphase_block2), hold on
    for(i=0:block_length-1)
        plot(i*(fs*T)+sample_nr,inphase_symbols2(i+1),'o');
    end
    hold off
    title('received data and samplepoints. (not yet combined!)');
    subplot(1,2,1), stem([training_sequence1 data_block]), title('Sent data.');
    
    figure(4), subplot(1,2,1), plot(inphase_symbols2,quad_symbols2,'x'), axis([-2 2 -2 2]), title('Signals after MF and sampled')
    
    %%
    %% Calculating influence of freq_offset. IS INSTEAD DONE AFTER DOWNCONVERTER.
    %%
    
    %[quad_symbols2, inphase_symbols2, theta2]=freq_offset_compensate(quad_symbols2,inphase_symbols2);
    %theta2
    %figure(7), plot(inphase_symbols2,quad_symbols2,'x'), axis([-6 6 -6 6]), title('Signals after freq-offset compensation.'), grid on
	
    %%%
    %%% Estimating Channel 3 and 4.
    %%%
    
   
    
    %[alpha3,theta3]=chan_estim7(quad_symbols2(1:training_length1),inphase_symbols2(1:training_length1),bpsk(training_sequence1,1),1,0.95);
    
    [alpha3,theta3] = channel_estimatorb(quad_symbols2(1:training_length1), ...
       inphase_symbols2(1:training_length1), training_sequence1, model, 1)
    
    if model == 1
       %[alpha4,theta4]=chan_estim5(quad_symbols2(training_length1+1:training_length1+ ...
       %training_length2),inphase_symbols2(training_length1+1:training_length1+ ...
       %training_length2), bpsk(training_sequence2,model))
       %[alpha4,theta4]=chan_estim7(quad_symbols2(training_length1+1:training_length1+training_length2),...
       %    inphase_symbols2(training_length1+1:training_length1+training_length2),bpsk(training_sequence1,1),1,0.95); 
       
       [alpha4,theta4] = channel_estimatorb(quad_symbols2(training_length1+1:training_length1+ ...
       training_length2), inphase_symbols2(training_length1+1:training_length1+ ...
       training_length2), training_sequence2, model, 2)
    elseif model == 0
        alpha4 = 1;
        theta4 = 0;
    end
    
    q_symbols2 = quad_symbols2;
    i_symbols2 = inphase_symbols2;
    
    % set true values instead:
    %alpha3 = 1;
    %theta3 = 0;
    %alpha4 = 1;
    %theta4 = 0;
    
     %%%
     %%% No combiner used since only 1 antenna is used to send from so far.
     %%%
     
     [combined_signal2]=combiner(alpha3(end), theta3(end), alpha4(end), theta4(end), ...
        inphase_symbols2, quad_symbols2, model);

    figure(4), subplot(1,2,2), plot(combined_signal,'x'), axis([-6 6 -6 6]), title('Signals after combiner, just before detector.')
    grid on
    
    inphase_symbols1 = real(combined_signal);
    inphase_symbols2 = real(combined_signal2);
    inphase_symbols = inphase_symbols1 + inphase_symbols2;
    
    figure(5), subplot(1,3,1), plot(inphase_symbols1(training_length1+training_length2:end),'x'), subplot(1,3,2), plot(inphase_symbols2(training_length1+training_length2:end),'x'), subplot(1,3,3), plot(inphase_symbols(training_length1+training_length2:end),'x')
    
    %%%
    %%% Detector
    %%%
    [indata_est1,training1_est,training2_est]=detector(inphase_symbols1,training_length1,training_length2,model);
    [indata_est2,training1_est,training2_est]=detector(inphase_symbols2,training_length1,training_length2,model);
    
    
    [indata_est,training1_est,training2_est]=detector(inphase_symbols,training_length1,training_length2,model);
    disp('Received data is:')
    indata_est;
    disp('Sent data is:')
    data_block;
    nr_of_errors1 = sum(abs(data_block-indata_est1))
    nr_of_errors_in_rate1 = nr_of_errors1/length(data_block)
    nr_of_errors2 = sum(abs(data_block-indata_est2))
    nr_of_errors_in_rate2 = nr_of_errors2/length(data_block)
    
    
    nr_of_errors = sum(abs(data_block-indata_est))
    nr_of_errors_in_rate = nr_of_errors/length(data_block)