ofdm_awgn.m

本源码介绍了OFDM的基本原理的仿真

clear all
fprintf('*********************** Programm is running, please wait ****************************\n');
fprintf('\n*********************** OFDM in AWGN ****************************\n');
fprintf('\n********************* pun.cc R=3/4; random interleaver(1000);16QAM Glay mapping*********\n');
subchannel=256;
K=0;

for SNR=0:2:10
diary ofdm_awgn.txt;

err=0;

if SNR<12 max=1;
elseif SNR<14  max=3;
elseif SNR<16  max=30;
elseif SNR<18  max=300;
else max=3000;
end

while(K<max)

N=3000;
source=round(rand(1,N));

R=3/4;%code rate of punctured convolutional code from (2.1.7)
[outcode_x]=encoder_pc(source,R);

outcode=reshape(outcode_x,4,length(outcode_x)/4); % serial to parallel

code1=outcode(1,:);
code2=outcode(2,:);
code3=outcode(3,:);
code4=outcode(4,:);
% ********************************************** random interleaver *****************************************************
a1=rand(1,length(code1));[b1,I1]=sort(a1);
for i=1:length(code1) inter1(i)=code1(I1(i));end

a2=rand(1,length(code2));[b2,I2]=sort(a2);
for i=1:length(code2) inter2(i)=code2(I2(i));end

a3=rand(1,length(code3));[b3,I3]=sort(a3);
for i=1:length(code3) inter3(i)=code3(I3(i));end

a4=rand(1,length(code4));[b4,I4]=sort(a4);
for i=1:length(code4) inter4(i)=code4(I4(i));end

 %************************************ 16QAM Glay mapping according to IEEE 802.16a *************
d=1;
j=sqrt(-1);
for i=1:length(inter1)   % 4bit/symbol  

    if ( inter1(i)==0 & inter2(i)==1 ) Im(i)=3*d;
        elseif ( inter1(i)==0 & inter2(i)==0 )  Im(i)=d;
        elseif ( inter1(i)==1 & inter2(i)==0 )  Im(i)=-d;
        elseif ( inter1(i)==1 & inter2(i)==1 )  Im(i)=-3*d;
        else fprintf('\n modulation error occurs\n');
    end 

     if ( inter3(i)==0 & inter4(i)==1 ) Qm(i)=3*d;
        elseif ( inter3(i)==0 & inter4(i)==0 )  Qm(i)=d;
        elseif ( inter3(i)==1 & inter4(i)==0 )  Qm(i)=-d;
        elseif ( inter3(i)==1 & inter4(i)==1 )  Qm(i)=-3*d;
        else fprintf('\n modulation error occurs\n');
    end 
mapping(i)=Im(i)+j*Qm(i);
end
   
%**************************************** OFDM and adding CP ***********************************
if (rem(length(mapping),subchannel)==0)
    tt=length(mapping)/subchannel;
else x=[mapping,zeros(1,subchannel-rem(length(mapping), subchannel))];
    tt=length(x)/subchannel;
end
smat=reshape(x,subchannel,tt);
s=ifft(smat);
ts=[s((subchannel-subchannel/4+1):subchannel,:);s];%adding CP

%**************************************** AWGN channel **********************************
transmittedsym=reshape(ts,1,(subchannel*5/4*tt));
% ray=raylrnd(sqrt(1/2),1,(subchannel*5/4*tt));
% chsym=ray.*transmittedsym;
out_ch=awgn(transmittedsym,SNR,'measured');
% out_ch1=out_ch./ray;
out_ch1=reshape(out_ch,subchannel*5/4,tt);
move_cp=out_ch1((subchannel/4+1:end),:);%move cp
out_fft=fft(move_cp);
m_o=reshape(out_fft,1,subchannel*tt);

%************** demodulation of 16QAM ************************************************************
for i=1:length(mapping)
    
if real(m_o(i))>=2*d demap1(i)=0;demap2(i)=1;
elseif ( real(m_o(i))>=0 & real(m_o(i))<2*d ) demap1(i)=0;demap2(i)=0;
elseif  ( real(m_o(i))>=-2*d & real(m_o(i))<0 ) demap1(i)=1;demap2(i)=0;
elseif   real(m_o(i))<-2*d  demap1(i)=1;demap2(i)=1;
else fprintf('\ndemodulation error occur\n');
end 

if imag(m_o(i))>=2*d demap3(i)=0;demap4(i)=1;
elseif ( imag(m_o(i))>=0 & imag(m_o(i))<2*d ) demap3(i)=0;demap4(i)=0;
elseif  ( imag(m_o(i))>=-2*d & imag(m_o(i))<0 ) demap3(i)=1;demap4(i)=0;
elseif   imag(m_o(i))<-2*d  demap3(i)=1;demap4(i)=1;
else fprintf('\ndemodulation error occur\n');
end 

end

%*************************************** deinterleaver*******************************************

for i=1:length(code1)
    deinter1(I1(i))=demap1(i);
end

for i=1:length(code1)
    deinter2(I2(i))=demap2(i);
end
 
for i=1:length(code1)
    deinter3(I3(i))=demap3(i);
end

for i=1:length(code1)
    deinter4(I4(i))=demap4(i);
end

%****************************************** viterbi decoder ****************************************
outdeinter=[deinter1;deinter2;deinter3;deinter4];
outdeinter=reshape(outdeinter,1,4*1000);

[decoder_output]=decoder_pc(outdeinter,R);
decoder_output=decoder_output';

err=err+sum(xor(source,decoder_output));
K=K+1;
end
%************************************** Display the BitError rate***********************************
ber=err/K/N;
fprintf('*********** SNR=%5d dB ****************\n',SNR);
fprintf('************ ber=%8.4d ****************\n',ber);

K=0;
diary off;
end