gmc.m

来自「B3g_phase2_C语言_Matlab程序及说明」· M 代码 · 共 218 行

clear
% Define parameters of MC filter banks
load ptfilter                                                             % load prototype filter coefficients which are saved in h
prototype_filter=h;
Carrier_number=16;                                                        % Number of subcarriers
Sampling_factor=18;                                                       % Factor of upsampling and downsampling
Filter_length=217;                                                        % Length of filters in taps

% Define the slot structure
Spread_factor=16;                                                         % Spread factor
Subslot_number=8;                                                         % Number of subslots
Symbol_number_subslot=6;                                                  % Symbol number in a subslot
Symbol_number_slot=Subslot_number*Symbol_number_subslot;                  % Symbol number in a slot
Pilot_length=24;                                                          % Length of pilot signal in chips
Gaurd_length=8;                                                           % Length of gaurd time in chips
GaurdPilot_length=Gaurd_length+Pilot_length;                              % Total length of gaurd time and pilot signal in chips
SubslotData_length=Symbol_number_subslot*Spread_factor;                   % Length of data in a subslot in chips
Subslot_length=SubslotData_length+GaurdPilot_length;                      % Length of a subslot in chips
SlotData_length=Subslot_number*SubslotData_length;                        % Length of data in a slot in chips
Slot_length=Subslot_number*Subslot_length+GaurdPilot_length;              % Length of a slot in chips 

% Define spreading codes
Walsh=[1 1; 1 -1]/sqrt(2);                                                % Walsh matrix for speading codes 
for k=1:log2(Spread_factor)-1
    Walsh=[Walsh Walsh;Walsh -Walsh]/sqrt(2);
end
Code_channel=Spread_factor;                                               % Number of code channels
Walsh=Walsh(1:Code_channel,:)';                                           % The codes come from Walsh matrix 

% Initialize the scrambler 
PN_RegStateI = [zeros(1,14) 1];  
PN_RegStateQ = [zeros(1,14) 1];  
PN_PolynI = [0 1 0 0 0 1 1 1 0 1 0 0 0 0 1];  
PN_PolynQ = [0 0 1 1 1 0 0 0 1 1 1 1 0 0 1];

%Define the channel
Antenna_number=2;                                                         % Number of receive antennas
Path_number=6;                                                            % Number of paths  
Path_Gain=[0.8084 0.462 0.253 0.259 0.0447 0.01];                         % Profile of channel model
% Path_Gain=[1 1 1 0 0 0]/sqrt(3);                                                 % Profile of channel model

delay=[7 26 38  54 72  90];                                               % Delays of paths
Fc = 2.4e9;                                                               % Carrier frequency 
V =200;                                                                   % moving speed in km/h
Tc = 1/1.28e6/Sampling_factor;                                            % Chip width 
Time_Begin = 0;                                                           % Initializing the time
Phase = 2*pi*rand(1,2);                                                   % Initializing the phase

% Define pilot signal
Pilot_sequence=[1 1 1 1 1 j -1 -j 1 -j 1 -j 1 -1 -1 1 1 -1 1 -1 1 -j -1 j];
for m=1:Path_number
    Pilot_matrix(m,:)=[Pilot_sequence(24-m+2:24) Pilot_sequence(1:24-m+1)];
end
Gaurd_Pilot=[Pilot_sequence(17:24) Pilot_sequence];

% Define Error Code
%%trel=poly2trellis([3 3],[2 7 7; 7 5 2]);                                   % Define trellis
% trel=poly2trellis([4 4],[13 6 13; 6 13 17]);                               % Define trellis
% InputN=2; OutputN=3;
% tblen=9;                                                                   % Traceback length
trel=poly2trellis(9,[561 753]); 
InputN=1; OutputN=2;
tblen=40;
Slot_number=4;
Symbol_bitN=2;                                                                % 2 for QPSK and 4 for 16QAM
Map_16QAM=[3+j*3 3+j 3-j*3 3-j 1+j*3 1+j 1-j*3 1-j -3+j*3 -3+j -3-j*3 -3-j -1+j*3 -1+j -1-j*3 -1-j]/sqrt(10);

msg_L=Symbol_bitN*SlotData_length/OutputN*InputN*Slot_number-InputN*tblen;    % Number of message bits
msg_SlotL=Symbol_bitN*SlotData_length/OutputN*InputN;
code_L=Symbol_bitN*SlotData_length*Slot_number;
code_SlotL=Symbol_bitN*SlotData_length;

% Raised cosine filter for interpolation
rcflt3=rcosfir(0.15,6,3,1);

% Bit error rate
ber=zeros(1,20);
format long
% Main loop

SNR1=5; SNR2=5;
%ERR=[msg_L*100 msg_L*50 msg_L*16 msg_L*4 msg_L  1000 500  200 100];
%ERR=[msg_L*128 msg_L*128 msg_L*64 msg_L*64 msg_L*32 msg_L*32 msg_L*16 msg_L*16  msg_L*8 msg_L*8  msg_L*4 msg_L*4 msg_L*2 msg_L*2];
%ERR=[1000 500  400 300  200 100];
for SNR=SNR1:SNR2
   
    errors=0;
    k=0;                                                                 
    while (errors<200 | k<50)             %ERR(SNR-SNR1+1)  %500*(1-SNR/SNR2/2)^2
        SNR
        k=k+1
      
        msg=randint(Carrier_number,msg_L,2);                                                                               % Random data
        msg=[msg zeros(Carrier_number,InputN*tblen)];
        for nc=1:Carrier_number 
            code(nc,:)=convenc(msg(nc,:),trel);                                                                            % Encode
            code(nc,:)=reshape(vec2mat(code(nc,:),Subslot_number),1,code_L);                                               % Interleaving
        end
        
        for sn=0:Slot_number-1
            % Transmiter
            % Transmiter:Generate single-carrier baseband signals 
            [PN,PN_RegStateI,PN_RegStateQ]=Complex_PNGen(PN_PolynI,PN_PolynQ,PN_RegStateI,PN_RegStateQ,SlotData_length);   % get a complex PN sequence.
            
            for nc=1:Carrier_number 
                s=[]; 
 
                x=((1-2*code(nc,sn*code_SlotL+(1:2:code_SlotL)))+j*(1-2*code(nc,sn*code_SlotL+(2:2:code_SlotL))))/sqrt(2); % QPSK mapping
               
             %   x=Map_16QAM(bi2de(vec2mat(code(nc,sn*code_SlotL+(1:code_SlotL)),4),'left-msb')'+1);                                                          % 16QAM mapping
               
                Tmp_s=reshape((Walsh*reshape(x,Code_channel,Symbol_number_slot))',1,SlotData_length).*PN;                  % spread spectrum modulation and interleving
%                Tmp_s=reshape((Walsh*reshape(x,Code_channel,Symbol_number_slot)),1,SlotData_length).*PN;                  % spread spectrum modulation without interleving

                Index_s=0;                                                                                                 % Insert pilot and gaurd 
                for kk=1:Subslot_number
                    tmp=Tmp_s(Index_s+(1:SubslotData_length));
                    tmp=[Gaurd_Pilot tmp]; 
                    s=[s tmp];
                    Index_s=Index_s+SubslotData_length;
                end
                s=[s Gaurd_Pilot];
                sb(nc,:)=s;              
            end
            
            % Transmiter:Combine subcarrier baseband signals into multicarrier baseband signal
            sb=sb.*exp(j*2*pi*(0:Carrier_number-1)'*Sampling_factor*(0:Slot_length-1)/Carrier_number);
            st=mc_sfb(sb,Slot_length,prototype_filter,Carrier_number,Sampling_factor,Filter_length);
            
            % Channel
            Signal_length=Sampling_factor*(Slot_length-1)+Filter_length+delay(Path_number);                                 % Length of multicarrier baseband signal 
            
            CH_Data = MultiCHannel(Path_Gain,Fc,V,Tc,Signal_length,Time_Begin,Phase);
            Time_Begin = Time_Begin+Signal_length;  

            
            doas=rand(1,Path_number)*pi*2;                                                      
            
            for p=1:Path_number
                ss(p,:)=[zeros(1,delay(p)) st zeros(1,delay(Path_number)-delay(p))];                            
            end
            
            for na=1:Antenna_number
                Sm(na,:)=sum(ss.*(diag(exp(-j*(na-1)*doas))*CH_Data));                                                     % The transmitted signal is passed through multipath channel 
            end
            %    Add Gauss noise
            Rm=awgn(Sm,SNR-10*log10(Spread_factor/Code_channel),'measured');
            %SNR_tmp=10*log10(sum(sum(abs(Rm.^2)))/sum(sum(abs((Rm-Sm).^2))))
         
           
            %    Receiver
            %    Receiver:Separate the received mulicarrier signal into subcarrier signals
            for na=1:Antenna_number
                tmp=mc_afb(Rm(na,:),2*Slot_length,prototype_filter,Carrier_number,Sampling_factor/2,Filter_length);
                rs(na,:,:)=tmp.*exp(-j*2*pi*(0:Carrier_number-1)'*Sampling_factor*(0:2*Slot_length-1)/Carrier_number/2); 
            end 
            %    pause
            %    Receiver: Process each subcarrier separately
            for nc=1:Carrier_number
                tmp=reshape(rs(:,nc,:),Antenna_number,2*Slot_length);
                R=SC_sync(tmp,Subslot_length,Subslot_number,Pilot_sequence,rcflt3);
                %    pause
                %    R=ones(4,1)*sb(nc,:);   
                %    Receiver: Channel estimation 
                Index_R=Gaurd_length;Fades=zeros(Antenna_number,Path_number*(Subslot_number+1));Index_Fades=0;Noise_variance=[];
                for kk=1:Subslot_number+1
                    Tmp_RP=R(:,Index_R+(1:Pilot_length));
                    Tmp_Fades=Tmp_RP*Pilot_matrix'/Pilot_length;
                    Fades(:,Index_Fades+(1:Path_number))=Tmp_Fades(:,1:Path_number);
                    Tmp_Noise = R(:,Index_R+(1:Pilot_length))-Tmp_Fades(:,1:Path_number)*Pilot_matrix;
                    Nv=sum(sum(abs(Tmp_Noise.*Tmp_Noise)))/Pilot_length/Antenna_number;                                    % Estimate of noise variance 
                    Noise_variance=[Noise_variance Nv];
                    Index_R=Index_R+Subslot_length;
                    Index_Fades=Index_Fades+Path_number;
                end
                
                %    Receiver: Equalization in DFT domain
                Index_R=Gaurd_length+Pilot_length; Index_Fades=0;R_EQ=[];
                for kk=1:Subslot_number
                    Tmp_Fades=(Fades(:,Index_Fades+(1:Path_number))+Fades(:,Index_Fades+Path_number+(1:Path_number)))/2;
                    Tmp_Nv=(Noise_variance(kk)+Noise_variance(kk+1))/2;
                    RM=R(:,Index_R+(1:Subslot_length)); 
                    R_DFT=fft(conj(RM)');
                    H_DFT=fft(conj([Tmp_Fades zeros(Antenna_number,Subslot_length-Path_number)])');
                    Tmp_R_EQ=conj(ifft(sum(R_DFT.*conj(H_DFT),2)./(sum(H_DFT.*conj(H_DFT),2)+Tmp_Nv)))';
                    R_EQ=[R_EQ Tmp_R_EQ(1:SubslotData_length)];
                    Index_Fades=Index_Fades+Path_number;
                    Index_R=Index_R+Subslot_length;
                end
                RR=reshape(Walsh'*reshape(R_EQ.*conj(PN),Symbol_number_slot,Spread_factor)',1,SlotData_length);           % with interleaving  
                %RR=reshape(Walsh'*reshape(R_EQ.*conj(PN),Spread_factor,Symbol_number_slot),1,SlotData_length);             % without interleaving
                
                Dem_signal(nc,sn*code_SlotL+(1:code_SlotL))=reshape([real(RR);imag(RR)],1,code_SlotL);      
                %Dem_signal(nc,sn*code_SlotL+(1:code_SlotL))=reshape([real(RR); abs(real(RR))-2/sqrt(10);imag(RR); abs(imag(RR))-2/sqrt(10)],1,code_SlotL);
                
            end   
        end
        for nc=1:Carrier_number
            Dem_signal(nc,:)=reshape(reshape(Dem_signal(nc,:),code_L/Subslot_number,Subslot_number)',1,code_L);  
         %   [x,qcode]=quantiz(sqrt(2)*DeQpsk_signal(nc,:),[-.75 -.5 -.25 0 .25 .5 .75],[7,6,5,4,3,2,1,0]); 
         %   decoded=vitdec(qcode,trel,tblen,'cont','soft',3);
            decoded=vitdec(Dem_signal(nc,:),trel,tblen,'cont','unquant');

            errors=errors+sum(abs(decoded(InputN*tblen+1:end)-msg(nc,1:msg_L)));
        end
        errors
        ber(1,SNR-SNR1+1)=errors/k/msg_L/Carrier_number;
        ber(1,1:SNR-SNR1+1)
    end
%     figure(1)
   
end 
semilogy(SNR1:SNR2,ber(:,1:SNR2-SNR1+1)')
grid
xlabel('SNR of Received Signal(in dB)')
ylabel('Bit Error Rate')
pause(0.2)
save ber_250 ber