cycpn.m

ofdm信道估计的matlab程序,很好的程序,可以下载的,大家快点下载啊

clear all;
clc;
Sample_freq=20000000;%采样频率
Df=20000000/64;%载波间隔
Tfft=1/Df;
Tgi=0.8e-6;
ofdmts=4e-6;%ofdm符号周期时间
Ts=1/Sample_freq;
IFFT_length=64;
GI=16;             % guard interval length
carrier_count=52;
bits_per_symbol=2;%QPSK调制
symbols_per_carrier=10;
Tlength=IFFT_length+GI;
N_number=carrier_count*symbols_per_carrier*bits_per_symbol;%总共传输的bit数
carriers=1:carrier_count;

N_snr=40;          % 每比特信噪比
snr=2; 
     
%------------------------------------------------------------
% vector initialization
X=zeros(1,N_number);
X1=[];
X2=[];
X3=[];
X4=[];
X5=[];
X6=[];
X7=[];
Y1=[];
Y2=[];
Y3=[];
Y4=[];
Y5=[];
Y6=[];
Y7=[];
XX=zeros(1,N_number);
dif_bit=zeros(1,N_number);
dif_bit1=zeros(1,N_number);
dif_bit2=zeros(1,N_number);
dif_bit3=zeros(1,N_number);
X=randint(1,N_number);%产生二进制随即序列（非0即1）
%--------------------------------------------------------
%QPSK调制：(1 1)->pi/4;(0 1)->3*pi/4;(0 0)->-3*pi/4;(1,0)->-pi/4;

s=(X.*2-1)/sqrt(2);
sreal=s(1:2:N_number);
simage=s(2:2:N_number);
X1=sreal+j.*simage;


%--------------------------------------------------------
%串并转换
X2=reshape(X1,carrier_count,symbols_per_carrier).';
%---------------------------------------------------------

IFFT_modulation=zeros(symbols_per_carrier,IFFT_length);%由52添0至64
IFFT_modulation(:,carriers)=X2;
%IFFT_modulation(:,conjugate_carriers)=conj(X3);
X3=ifft(IFFT_modulation,IFFT_length,2);
%X5=X4.';
%---------------------------------------------------------
%插入训练序列
%产生PN序列
connection=[1 0 0 1 0 1 ];
n=length(connection);
l=2^n-1;
register=[zeros(1,n-1) 1];
mseq(1)=register(n);
for i = 2:l,
  newregister(1)=mod(sum(connection.*register),2);
    for j =2:n,
       newregister(j)=register(j-1);
    end;
  register = newregister;
  mseq(i)= register(n);
end;
mseq=[mseq,1];
tr_train=ifft(mseq,IFFT_length);
X4=[tr_train;tr_train;X3]
%---------------------------------------------------------
%---------------------------------------------------------
%加保护间隔（循环前缀）
for k=1:(symbols_per_carrier+2);
   for i=1:IFFT_length;
      X6(k,i+GI)=X4(k,i);
   end
   for i=1:GI;
      X6(k,i)=X4(k,i+IFFT_length-GI);    
   end
end
%---------------------------------------------------------
%x=symbols_per_carrier*(IFFT_length+GI)
%---------------------------------------------------------
%并串转换
X7=reshape(X6.',1,(symbols_per_carrier+2)*(IFFT_length+GI))

%---------------------------------------------------------
%%%%%%%参数设置%%%%%%%%%
freq_delta=100000;       %频偏200k
Ts=1/Sample_freq;
%SNR=10;
%%%%%%%%%%%%%%%%%%%%%%%
Signal_in=X7.*exp(-j*2*pi*freq_delta*[1:length(X7)]*Ts);    %加入频偏的影响

Error_ber=[];%误比特率
Error_ber1=[];
Error_ber2=[];%误比特率
Error_ber3=[];
%Error_ser=[];%误符号率
for SNR=0:snr:N_snr %0:8:40

Signal_in=awgn(Signal_in,SNR,'measured',round(rand*100),'dB')
%-------------------------------------------------------
Y1= Signal_in;%未进行频偏的接收数据
 %------------------------------------------------------------- 
 %未进行频偏的数据解调结果
    %串并变换
   Y111=Y1(2*Tlength+1:length(Y1));
   Y11=reshape(Y111,IFFT_length+GI,symbols_per_carrier).';%先变成80*12，再转置成12*80，恢复成行为符号，列为载波
   
  %去保护间隔
    for k=1:symbols_per_carrier;
       for i=1:IFFT_length;
           Y12(k,i)=Y11(k,i+GI);
       end
    end
   %进行64点fft
     Y13=fft(Y12,IFFT_length,2);%每行的符号进行64点fft  12*64
     Y14=Y13(:,carriers);%去掉尾部12列原补零点， 12*52
%同步
    Pn1=zeros(1,Tlength);
    Pn2=zeros(1,Tlength);
    Pn1=Signal_in(1:Tlength)
    Pn2=Signal_in(Tlength+1:Tlength*2)
    Rt=sum(Pn1.*conj(Pn2)); % 式3-64要求的
    f_freqoff=angle(Rt)/((2*pi)*Tlength*Ts);
    freoff=f_freqoff %频偏
    
    %去除频偏
    Receive_signal=Signal_in.*exp(j*2*pi*freoff*[1:length(X7)]*Ts);
     Y7=Receive_signal(2*Tlength+1:length(Receive_signal));
     %串并变换
   Y6=reshape(Y7,IFFT_length+GI,symbols_per_carrier).';%先变成80*12，再转置成12*80，恢复成行为符号，列为载波
   
  %去保护间隔
    for k=1:symbols_per_carrier;
       for i=1:IFFT_length;
           Y5(k,i)=Y6(k,i+GI);
       end
    end
   %进行64点fft
     Y4=fft(Y5,IFFT_length,2);%每行的符号进行64点fft  12*64
     Y3=Y4(:,carriers);%去掉尾部12列原补零点， 12*52

 
  %并串变换
  YY=reshape(Y3,1,N_number/bits_per_symbol);%去频偏
  YY1=reshape(Y14.',1,N_number/bits_per_symbol);%未去频偏
 
%------------------------------------------------------------
%QPSK解调
   y_real=sign(real(YY));
   y_image=sign(imag(YY));
   y_re=y_real./sqrt(2);
   y_im=y_image./sqrt(2);
   
   y_real1=sign(real(YY1));
   y_image1=sign(imag(YY1));
   y_re1=y_real1./sqrt(2);
   y_im1=y_image1./sqrt(2); 
 
   r00=[];
   r01=[];
   r10=[];
   r11=[];
 
  for k=1:length(y_real);
     r00=[r00,[y_real(k),y_image(k)]];
  end;
  for k=1:length(y_real1);
     r10=[r10,[y_real1(k),y_image1(k)]];
  end;
 
 for k=1:length(y_re);
     r01=[r01,[y_re(k),y_im(k)]];
 end;
 for k=1:length(y_re1);
     r11=[r11,[y_re1(k),y_im1(k)]];
 end;

  %  XX(find(r01>0))=1;
%-------------------------------------------------------------
%计算在不同信噪比下的误比特率并作图

 dif_bit=s-r01; %原发射数据-接收后未进行信道估计的数据
 dif_bit1=s-r11; %原发射数据-接收后进行信道估计的数据

 ber_snr=0;    %纪录误比特数
    for k=1:N_number;
       if dif_bit(k)~=0;
         ber_snr=ber_snr+1;
       end
   end;
 ber_snr1=0;    %纪录误比特数
    for k=1:N_number;
       if dif_bit1(k)~=0;
          ber_snr1=ber_snr1+1;
      end
    end
 
 Error_ber=[Error_ber,ber_snr];
 Error_ber1=[Error_ber1,ber_snr1];
end

BER=zeros(1,length(0:snr:N_snr));
BER1=zeros(1,length(0:snr:N_snr));

BER=Error_ber./N_number;% 错误的比特数/总比特数
BER1=Error_ber1./N_number;
%-------------------------------------------------------------
%-------------------------------------------------------------
 i=0:snr:N_snr;
semilogy(i,BER,'-*r');

hold on;
semilogy(i,BER1,'-og');
hold on;

grid on;

legend('frequency estimation','not Estimation');
hold off;