%bin22deci.mfunction y=bin22deci(x)%将二进制数转化为十进制数t=size(x,2);y=(t-1:-1:0);y=2.^y;y=x*y';%************************end of file***********************************%comb.m%AWGN加噪声程序function[iout,qout]=comb(idata,qdata,attn)%******************variables*************************%idata:输入I信道数据%qdata:输入Q信道数据%iout输出I信道数据%qout输出Q信道数据%attn:由信噪比导致的衰减系数%******************************************************iout=randn(1,length(idata)).*attn;qout=randn(1,length(qdata)).*attn;iout=iout+idata(1:length(idata));qout=qout+qdata(1:length(qdata));%************************end of file***********************************%crdemapping.m%数据逆映射载波程序function[iout,qout]=crdemapping(idata,qdata,fftlen,nd);%******************variables*************************%idata:输入I信道的数据%qdata:输入Q信道的数据%iout:输出I信道的数据%qout:输出Q信道的数据%fftlen:FFT的长度%nd:OFDM符号数%*****************************************************iout(1:26,:)=idata(2:27,:);qout(1:26,:)=qdata(2:27,:);iout(27:52,:)=idata(39:64,:);qout(27:52,:)=qdata(39:64,:);%********************end of file***************************%crmapping.m%数据映射载波程序function[iout,qout]=crmapping(idata,qdata,fftlen,nd);%******************variables*************************%idata:输入I信道的数据%qdata:输入Q信道的数据%iout:输出I信道的数据%qout:输出Q信道的数据%fftlen:FFT的长度%nd:OFDM符号数%*****************************************************iout=zeros(fftlen,nd);qout=zeros(fftlen,nd);iout(2:27,:)=idata(1:26,:);qout(2:27,:)=qdata(1:26,:);iout(39:64,:)=idata(27:52,:);qout(39:64,:)=qdata(27:52,:);%********************end of file***************************%deci22bin.mfunction y=deci22bin(x,t)%十进制数x转化为二进制数，二进制数至少表示为t位y=zeros(size(x,1),t);for j=1:size(x,1)i=1;while x(j)>=0&i<=ty(j,i)=rem(x(j),2);%x(j)为偶数时,y(j,i)为0;反之为1x(j)=(x(j)-y(j,i))/2;i=i+1;endy(j,:)=y(j,t:-1:1);%倒序排列end%************************end of file***********************************%giins1.m%插入保护间隔程序function[iout,qout]=giins1(idata,qdata,fftlen,gilen,nd);%******************变量*************************%idata:输入I信道数据%qdata:输入Q信道数据%iout:输出I信道数据%qout:输出Q信道数据%fftlen:FFT长度(points)%gilen:保护间隔长度(points)%*****************************************************idata1=reshape(idata,fftlen,nd);qdata1=reshape(qdata,fftlen,nd);idata2=[idata1(fftlen-gilen+1:fftlen,:);idata1];qdata2=[qdata1(fftlen-gilen+1:fftlen,:);qdata1];iout=reshape(idata2,1,(fftlen+gilen)*nd);qout=reshape(qdata2,1,(fftlen+gilen)*nd);%********************end of file***************************%qpskdemod1.m%QPSK解调程序function[demodata]=qpskdemod1(idata,qdata,para,nd,ml)%******************variables*************************%idata:输入I信道数据%qdata:数据Q信道数据%demodata:解调后数据(para-by-nd matrix)%para:并行信道数%nd:符号数%ml:调制数%(QPSK->2 16QAM->4)%*****************************************************demodata=zeros(para,ml*nd);demodata((1:para),(1:ml:ml*nd-1))=idata((1:para),(1:nd))>=0; demodata((1:para),(2:ml:ml*nd))=qdata((1:para),(1:nd))>=0;%************************end of file***********************************%qpskmod1.m%QPSK调制程序function[iout,qout]=qpskmod1(paradata,para,nd,ml)%******************variables*************************%paradata:输入数据%iout:输出数据I%qout:输出数据Q%para:并行信道数%nd:数据数%ml:调制数%(QPSK->2 16QAM->4)%*****************************************************m2=ml./2;paradata2=paradata.*2-1;count2=0;for jj=1:ndisi=zeros(para,1);isq=zeros(para,1);isi=isi+paradata2((1:para),1+count2);isq=isq+paradata2((1:para),2+count2);iout((1:para),jj)=isi;qout((1:para),jj)=isq;count2=count2+ml;end%********************end of file***************************%viterbi.m%viterbi解码程序function[decoder_output,survivor_state,cumulated_metric]=viterbi(G,k,channel_output) %[decoder_output,survivor_state,cumulated_metric]=viterbi(G,k,channel_output)%其中G是一个n行L*k列矩阵，它的每一行决定了从移位寄存器到输入码字的连接方式.%survivor_state是一个矩阵，它显示了通过网格的最优路径，这个矩阵通过一个单独%的函数metric(x,y)给出。

n=size(G,1);%检验G的维数if rem(size(G,2),k)~=0error('Size of G and k do not agree')endif rem(size(channel_output,2),n)~=0error('channle output not of the right size')endL=size(G,2)/k;number_of_states=2^((L-1)*k);%产生状态转移矩阵，输出矩阵和输入矩阵for j=0:number_of_states-1for t=0:2^k-1[next_state,memory_contents]=nxt_stat(j,t,L,k);input(j+1,next_state+1)=t;branch_output=rem(memory_contents*G',2);nextstate(j+1,t+1)=next_state;output(j+1,t+1)=bin2deci(branch_output);endendinput;state_metric=zeros(number_of_states,2);depth_of_trellis=length(channel_output)/n;channel_output_matrix=reshape(channel_output,n,depth_of_trellis);survivor_state=zeros(number_of_states,depth_of_trellis+1);[row_survivor col_survivor]=size(survivor_state);%开始非尾信道输出的解码%i为段，j为每一阶段的状态，t为输入for i=1:depth_of_trellis-L+1flag=zeros(1,number_of_states);if i<=Lstep=2^((L-i)*k);elsestep=1;endfor j=0:step:number_of_states-1for t=0:2^k-1branch_metric=0;binary_output=deci2bin(output(j+1,t+1),n);for tt=1:nbranch_metric=branch_metric+metric(channel_output_matrix(tt,i),binary_output(tt));endif((state_metric(nextstate(j+1,t+1)+1,2)>state_metric(j+1,1)+branch_metric)|flag(nextstate(j+1 ,t+1)+1)==0)state_metric(nextstate(j+1,t+1)+1,2)=state_metric(j+1,1)+branch_metric;survivor_state(nextstate(j+1,t+1)+1,i+1)=j;flag(nextstate(j+1,t+1)+1)=1;endendendstate_metric=state_metric(:,2:-1:1);end%开始尾信道输出的解码for i=depth_of_trellis-L+2:depth_of_trellisflag=zeros(1,number_of_states);last_stop=number_of_states/(2^((i-depth_of_trellis+L-2)*k));for j=0:last_stop-1branch_metric=0;binary_output=deci2bin(output(j+1,1),n);for tt=1:nbranch_metric=branch_metric+metric(channel_output_matrix(tt,i),binary_output(tt));endif((state_metric(nextstate(j+1,1)+1,2)>state_metric(j+1,1)+branch_metric)|flag(nextstate(j+1,1 )+1)==0)state_metric(nextstate(j+1,1)+1,2)=state_metric(j+1,1)+branch_metric;survivor_state(nextstate(j+1,1)+1,i+1)=j;flag(nextstate(j+1,1)+1)=1;endendstate_metric=state_metric(:,2:-1:1);end%从最优路径产生解码输出%由段得到状态序列，再由状序列从input矩阵中得到该段的输出state_sequence=zeros(1,depth_of_trellis+1);size(state_sequence);state_sequence(1,depth_of_trellis)=survivor_state(1,depth_of_trellis+1);for i=1:depth_of_trellisstate_sequence(1,depth_of_trellis-i+1)=survivor_state((state_sequence(1,depth_of_trellis+2-i )+1),depth_of_trellis-i+2);endstate_sequence;decoder_output_matrix=zeros(k,depth_of_trellis-L+1);for i=1:depth_of_trellis-L+1dec_output_deci=input(state_sequence(1,i)+1,state_sequence(1,i+1)+1);dec_output_bin=deci2bin(dec_output_deci,k);decoder_output_matrix(:,i)=dec_output_bin(k:-1:1)';enddecoder_output=reshape(decoder_output_matrix,1,k*(depth_of_trellis-L+1));cumulated_metric=state_metric(1,1);%************************end of file***********************************%ofdm1.m%QPSK调制、AWGN信道下仿真程序%**********************初始参数***************************para=52;%并行信道数fftlen=64;%FFT长度noc=52;%载波数nd=6;%每循环中OFDM符号数ml=2;%调制水平:QPSKgilen=16;%保护间隔长度(points)ebn0=2;%信噪比sr=250000;%OFDM symbol rate(250 ksyombol/s)br=sr.*ml;%Bit rate per carrier%**************************主循环部分**************************nloop=200;%仿真的循环数noe1=0;%信道解码前错误数据数nod1=0;%信道解码前传输数据数noe2=0;%信道解码后错误数据数nod2=0;%信道解码后传输数据数%**************************发射机*********************************for iii=1:nloop%**************************数据产生**************************** %信源编码t=[0:pi/25:2*pi];xx=sin(t);init=[-1:.1:1];partition=[-1:.1:.9];predictor=[0 1];encode=dpcmenco(xx,init,partition,predictor);encode2=reshape(encode,51,1);bin=deci22bin(encode2,6);recode=reshape(bin,1,306);%信道编码k0=1;G=[1 0 1 1 0 1 1;1 1 1 1 0 0 1];channelencode=cnv_encd(G,k0,recode);%******************串并转换***********************paradata=reshape(channelencode,para,nd*ml);%reshape:内建功能%**************************QPSK调制***************************** [ich,qch]=qpskmod1(paradata,para,nd,ml);kmod=1/sqrt(2);%sqrt:内建功能ich=ich.*kmod;qch=qch.*kmod;%数据映射[ich1,qch1]=crmapping(ich,qch,fftlen,nd);%*******************IFFT************************x=ich1+qch1.*i;y=ifft(x);%ifft:内建功能ich2=real(y);%real:内建功能(实部)qch2=imag(y);%imag:内建功能（虚部）%*********插入保护间隔**********[ich3,qch3]=giins1(ich2,qch2,fftlen,gilen,nd);fftlen2=fftlen+gilen;%------------------衰减计算-----------------------spow=sum(ich3.^2+qch3.^2)/nd./para;%sum:内建功能attn=0.5*spow*sr/br*10.^(-ebn0/10);attn=sqrt(attn);%---------------AWGN addition---------------[ich4,qch4]=comb(ich3,qch3,attn);%******************去除保护间隔*********[ich5,qch5]=girem1(ich4,qch4,fftlen2,gilen,nd);%******************FFT******************rx=ich5+qch5.*i;ry=fft(rx);ich6=real(ry);qch6=imag(ry);%载波逆映射[ich7,qch7]=crdemapping(ich6,qch6,fftlen,nd);%*****************解调*******************ich7=ich7./kmod;qch7=qch7./kmod;[demodata]=qpskdemod1(ich7,qch7,para,nd,ml);%**************并串转换*****************demodata1=reshape(demodata,1,para*nd*ml);%信道解码[channeldecode]=viterbi(G,k0,demodata1);%信源解码reshapechanneldecode=reshape(channeldecode,51,6);deci=bin22deci(reshapechanneldecode);redeci=reshape(deci,1,51);codebook=[-1:.1:5.3];decode=dpcmdeco(redeci,codebook,predictor);%**************************比特误码率(BER)**************************** noe10=sum(abs(demodata1-channelencode));nod10=length(channelencode);noe20=sum(abs(channeldecode-recode));nod20=length(recode);noe1=noe10+noe1;nod1=nod10+nod1;noe2=noe20+noe2;nod2=nod20+nod2;fprintf('%d\t%e\t%e\n',iii,noe10/nod10,noe20/nod20);endber1=noe1/nod1;ber2=noe2/nod2;%**********************输出结果***************************fprintf('%f\t%e\t%e\t%d\t\n',ebn0,ber1,ber2,nloop);%************************end of file***********************************。