gusucode.com > m4k20功率分配图 mimo功率优化,OFDM典型仿真程序，贪婪算法仿真程序 > code/test_16qam/test.m

    % Program 4-8
% ofdmce.m
%
% Simulation program to realize OFDM transmission system
%
% Programmed by T.Yamamura and H.Harada
%
% GI CE GI data GI data...(data 6symbols)
% (CE: Chanel estimation symbol, GI Guard interval)
%

%********************** preparation part ***************************

%para=52;    % Number of parallel channel to transmit (points)
para=128;
%fftlen=64;  % FFT length
fftlen=128;  % FFT length
%noc=53;     % Number of carriers
nd=6;       % Number of information OFDM symbol for one loop
knd=1;      % Number of known channel estimation (CE) OFDM symbol
ml=2;       % Modulation level : QPSK
sr=250000;  % OFDM symbol rate (250 ksyombol/s)
br=sr.*ml;  % Bit rate per carrier
%gilen=16;   % Length of guard interval (points)
gilen=32;   % Length of guard interval (points)
%ebn0=3;     % Eb/N0
ebn0=30:5:35;
for m=1:2


%%%%%%%%%%%%% fading initialization %%%%%%%%%%%

tstp=1/sr/(fftlen+gilen); % Time resolution
%tstp=1/sr/(fftlen); % Time resolution
%itau=[0];       % Arrival time for each multipath normalized by tstp 
itau=[0,2];       % Arrival time for each multipath normalized by tstp 
%dlvl1=[0];      % Mean power for each multipath normalized by direct wave.
dlvl1=[0,10];      % Mean power for each multipath normalized by direct wave.
%n0=[6];	        % Number of waves to generate fading n0(1),n0(2)
n0=[6,7];	        % Number of waves to generate fading n0(1),n0(2)
%th1=[0.0];      % Initial Phase of delayed wave
th1=[0.0,0.0];      % Initial Phase of delayed wave
%itnd1=[1000];   % set fading counter        
%itnd1=[1000,2000];   % set fading counter    
itnd1=[1000,2000];   % set fading counter  
itnd21=[1000,2000];
itnd22=[100000,102000];
%now1=1;         % Number of directwave + Number of delayed wave
now1=2;         % Number of directwave + Number of delayed wave
fd=150;         % Maximum Doppler frequency
flat=0;         % Flat or not (see ofdm_fading.m)
itnd0=nd*(fftlen+gilen)*20; % Number of fading counter to skip 
%itnd0=nd*(fftlen)*20; % Number of fading counter to skip 

%************************** main loop part **************************

nloop=1;  % Number of simulation loops

noe = 0;    % Number of error data
nod = 0;    % Number of transmitted data
eop=0;      % Number of error packet
nop=0;      % Number of transmitted packet

%************************** transmitter *****************************
for iii=1:nloop
   
seridata=rand(1,para*nd*ml)>0.5;  %  DC=0

paradata=reshape(seridata,para,nd*ml); %size(51  *  nd*ml)

%-------------- ml modulation ---------------- 

[ich,qch]=qpskmod(paradata,para,nd,ml);
kmod=1/sqrt(2);
ich=ich.*kmod;
qch=qch.*kmod;
ch0=ich+sqrt(-1).*qch;
ch1=reshape(ch0,1,para*nd);

%-------------------stbc--------------------
[st_out1,st_out2]=st2tr(ch1,fftlen);
ch11=reshape(st_out1,para,nd);
ch12=reshape(st_out2,para,nd);
ch10=ch0;


% CE data generation
%kndata=zeros(1,fftlen);
%kndata0=2.*(rand(1,52)>0.5)-1;
%kndata(2:27)=kndata0(1:26);
%kndata(39:64)=kndata0(27:52);
%ceich=kndata; % CE:BPSK
%ceqch=zeros(1,64);
kndata=2.*(rand(1,fftlen)>0.5)-1;
cech=kndata;

%------------- data mapping (DC=0) -----------

%[ich1,qch1]=crmapping(ich,qch,fftlen,nd);
ch21=[cech.' ch11];
ch22=[cech.' ch12];
ch20=[cech.' ch10];
%ich2=[ceich.' ich1]; % I-channel transmission data
%qch2=[ceqch.' qch1]; % Q-channel transmission data

%------------------- IFFT  -------------------

%x=ich2+qch2.*i;
%y=ifft(x);
%ich3=real(y);
%qch3=imag(y);

ch30=ifft(ch20);
ch31=ifft(ch21);
ch32=ifft(ch22);

%---------- Gurad interval insertion ---------

%fftlen2=fftlen+gilen;
%[ich4,qch4]= giins(ich3,qch3,fftlen,gilen,nd+1);
[ch40]= giins(ch30,fftlen,gilen,nd+1); 
[ch41]= giins(ch31,fftlen,gilen,nd+1); 
[ch42]= giins(ch32,fftlen,gilen,nd+1); 
fftlen2=fftlen+gilen;

%ch40=reshape(ch30,1,para*(nd+1));
%ch41=reshape(ch31,1,para*(nd+1));
%ch42=reshape(ch32,1,para*(nd+1));
ch41=(1/sqrt(2)).*ch41;
ch42=(1/sqrt(2)).*ch42;
ich40=real(ch40);
qch40=imag(ch40);
ich41=real(ch41);
qch41=imag(ch41);
ich42=real(ch42);
qch42=imag(ch42);

%---------- Attenuation Calculation ----------

spow=sum(ich40.^2+qch40.^2)/nd./para;
attn=0.5*spow*sr/br*10.^(-ebn0(m)/10);
attn=sqrt(attn);
%attn=0;


%********************** fading channel ****************************** 
%If you would like to simulate performance under fading, please remove "*"
%from the following four sentenses
%[ifade,qfade,ramp,rcos,rsin]=sefade(ich4,qch4,itau,dlvl1,th1,n0,itnd1,now1,length(ich4),tstp,fd,flat);
%itnd1 = itnd1+itnd0;  % Updata fading counter
%ich4=ifade;
%qch4=qfade;
[ifade40,qfade40,ramp,rcos,rsin]=sefade(ich40,qch40,itau,dlvl1,th1,n0,itnd1,now1,length(ich40),tstp,fd,flat);
[ifade41,qfade41,ramp,rcos,rsin]=sefade(ich41,qch41,itau,dlvl1,th1,n0,itnd21,now1,length(ich41),tstp,fd,flat);
[ifade42,qfade42,ramp,rcos,rsin]=sefade(ich42,qch42,itau,dlvl1,th1,n0,itnd22,now1,length(ich42),tstp,fd,flat);
itnd1 = itnd1+itnd0;  % Updata fading counter
itnd21 = itnd21+itnd0;
itnd22 = itnd22+itnd0; 
ich40=ifade40;
qch40=qfade40;
ich41=ifade41;
qch41=qfade41;
ich42=ifade42;
qch42=qfade41;
ich4=ich41+ich42;
qch4=qch41+qch42;

%***************************  Receiver  *****************************
%--------------- AWGN addition --------------- 
%[ich5,qch5]=comb(ich4,qch4,attn);
%ch5=ich5+sqrt(-1).*qch5;
[ich50,qch50]=comb(ich40,qch40,attn);
ch50=ich50+sqrt(-1).*qch50;
%[ich51,qch51]=comb(ich41,qch41,attn);
[ich51,qch51]=comb(ich41,qch41,0);
ch51=ich51+sqrt(-1).*qch51;
%[ich52,qch52]=comb(ich42,qch42,attn);
[ich52,qch52]=comb(ich42,qch42,0);
ch52=ich52+sqrt(-1).*qch52;
[ich5,qch5]=comb(ich4,qch4,attn);
ch5=ich5+sqrt(-1).*qch5;

%----Perfect fading compensation for one path fading ----
%If you would like to simulate performance under perfect compensation, please remove "*"
%from the following four sentenses
%ifade2=1./ramp.*(rcos(1,:).*ich5+rsin(1,:).*qch5);
%qfade2=1./ramp.*(-rsin(1,:).*ich5+rcos(1,:).*qch5);
%ich5=ifade2;
%qch5=qfade2;

%****************** Guard interval removal *********
%[ich6,qch6]= girem(ich5,qch5,fftlen2,gilen,nd+1);
ch6 = girem(ch5,fftlen2,gilen,nd+1);
ch60 = girem(ch50,fftlen2,gilen,nd+1);
ch61 = girem(ch51,fftlen2,gilen,nd+1);
ch62 = girem(ch52,fftlen2,gilen,nd+1);
    
%ch6=reshape(ch5,para,(nd+1));
%ch60=reshape(ch50,para,(nd+1));
%ch61=reshape(ch51,para,(nd+1));
%ch62=reshape(ch52,para,(nd+1));

%------------------  FFT  --------------------
%rx=ich6+qch6.*i;
%ry=fft(rx);
%ich7=real(ry);
%qch7=imag(ry);
ch70=fft(ch60);
%ich70=real(ch70);
%qch70=imag(ch70);
ch71=fft(ch61);
%ich71=real(ch71);
%qch71=imag(ch71);
ch72=fft(ch62);
%ich72=real(ch72);
%qch72=imag(ch72);
ch7=fft(ch6);
%ich7=real(ch7);
%qch7=imag(ch7);

%-------------- Fading compensation by CE symbol --------------
ce=1;
[ich70,qch70,h0] = cecompensation (ch20,ch70,ce);
[ich71,qch71,h1] = cecompensation (ch21,ch71,ce);
[ich72,qch72,h2] = cecompensation (ch22,ch72,ce);
h1=h1(:,knd+1:nd+1);
h11=reshape(h1,1,para*nd);
h2=h2(:,knd+1:nd+1);
h22=reshape(h2,1,para*nd);

%---------- CE symbol removal ----------------

ch88=ich71+sqrt(-1).*qch71+ich72+sqrt(-1).*qch72;
ch88=ch88(:,knd+1:nd+1);
ch99=reshape(ch88,1,para*nd);
%ich99=real(ch88);
%qch99=imag(ch88);

ch80=ich70+sqrt(-1).*qch70;
ch80=ch80(:,knd+1:nd+1);
ich90=real(ch80);
qch90=imag(ch80);

%ch8=ich7+sqrt(-1).*qch7;
ch8=ch7;
ch8=ch8(:,knd+1:nd+1);
ch9=reshape(ch8,1,para*nd);
%ich9=real(ch8);
%qch9=imag(ch8);
%----------------signal combine-------------

s=signalcombine21(ch99,h11,h22,fftlen,nd);%signalcombine2t1r：
                                                  %2发射天线一接收天线的信号合并
s1=signalcombine2t1r(ch9,para*nd,h11,h22,fftlen,nd);

%for j=1:para*nd
   %h11(j)=1;
  %  h22(j)=1;
  %end
                                                     
%s=signalcombine2t1r(ch99,para*nd,h11,h22,fftlen,nd);%signalcombine2t1r：
                                                  %2发射天线一接收天线的信号合并

%---------------空时解码--------------------
for n=1:para*nd              %de_st(n) 是空时解码的输出
    de_st(n)=de_st_out(s(n),h11(n),h22(n));    %这里的de_st_out是对16qam的空时解码
end

ch10=reshape(de_st,para,nd);
ich10=real(ch10)./kmod;
qch10=imag(ch10)./kmod;



%----------------- demoduration --------------

ich10_0=ich90./kmod;
qch10_0=qch90./kmod;
[demodata0]=qpskdemod(ich10_0,qch10_0,para,nd,ml);   
[demodata]=qpskdemod(ich10,qch10,para,nd,ml);   

%--------------  error calculation  ----------

demodata0=reshape(demodata0,1,para*nd*ml);
demodata=reshape(demodata,1,para*nd*ml);
noe2=sum(abs(demodata-seridata));
nod2=length(seridata);


% calculating PER
if noe2~=0
   eop=eop+1;
else
   eop=eop;
end   
   eop;
   nop=nop+1;
   
% calculating BER
noe=noe+noe2;
nod=nod+nod2;

%fprintf('%d\t%e\t%d\n',iii,noe2/nod2,eop);
   
end

per(m)=eop/nop;
ber(m)=noe/nod;
ber_stand(m)=(1/2)*(1-1/sqrt(1+1/10^(0.1*ebn0(m))));

%********************** Output result ***************************

%fprintf('%f\t%e\t%e\t%d\t%d\n',ebn0,ber,per,nloop,fd);
  
%fid = fopen('BERofdmce.dat','a');
%fprintf(fid,'%f\t%e\t%e\t%d\t\n',ebn0,ber,per,nloop);
%fclose(fid);

%******************** end of file ***************************
end
semilogy(ebn0,ber,'s-r',ebn0,ber_stand,'g');
%semilogy(ebn0,ber,'s-r');
xlabel('SNR  (dB)');
ylabel('BER');
legend('without Space-Time coding','theory');
%legend('without Space-Time coding');
grid;