OFDM-tongbu.zip_OFDM 同步_OFDM同步算法_OFDM系统同步_ofdm 改进算法_ofdm同步改进

%********************schmidl algorithm******************* close all; clear all; clc; %参数定义 N=128; %FFT/IFFT 变换的点数或者子载波个数（Nu<=N） Ng=N/8; %循环前缀的长度 (保护间隔的长度) Ns=Ng+N; %包括循环前缀的符号长度 L1=64; ts=25.6e-9; lam=1.55e-6; c=3e8; OSNR=3; B=N/ts; N_b=8; %Number of samples per bit; f3db=2*200e3; % Linewidth of the loca laser&source laser. 100e3Hz per laser fc=40e9; M=4; Modu='qam'; %************利用查表法生成复随机序列********************** m=log2(M); rand('state', sum(100*clock));%Initialize rand to a different state each time random=round(rand(N,2)); % bin=vec2mat(random,m);%m columns Dec=bi2de(random,[],'left-msb'); switch Modu case 'psk' syncheader=pskmod(Dec,M,[],'gray').';%QPSK case 'qam' syncheader=qammod(Dec,M,[],'gray').';%QAM end %QAMTable=[7+7i,-7+7i,-7-7i,7-7i]/7; %buf=QAMTable(randint(N,1,4)+1); QAMTable1=1/sqrt(2)*[1+i 1-i -1+i -1-i]; buf1=QAMTable1(randint(L1,1,4)+1); %*************在奇数子载波的位置插入零********************* %x=zeros(N,1); %index = 1; %for n=1:2:N % x(n)=buf(index); % index=index+1; %end; %**************利用IFFT变换生成Schmidl训练符号*************** sch = ifft(syncheader); cp=sch(N-L1+1:N); cp=cp.*buf1; %sch1=[cp sch]; %*****************添加一个空符号以及一个后缀符号************* m=log2(M); rand('state', sum(100*clock));%Initialize rand to a different state each time random=round(rand(N,2)); % bin=vec2mat(random,m);%m columns Dec=bi2de(random,[],'left-msb'); switch Modu case 'psk' syncheader1=pskmod(Dec,M,[],'gray').';%QPSK case 'qam' syncheader1=qammod(Dec,M,[],'gray').';%QAM end %src = QAMTable(randint(N,1,4)+1).'; sym = ifft(syncheader1).'; %第二个训练符号 sig =[zeros(N,1) sch.' sym]; %[保护间隔 第一个训练符号 第二个训练符号] 的形式 %**********************添加循环前缀************************* tx =[sig(N - Ng +1:N,:);sig]; %***********************经过信道*************************** recv2= reshape(tx,1,size(tx,1)*size(tx,2)); recv3=[recv2(1:Ns) cp recv2(Ns+1:3*Ns)]; recv = awgn(recv3,10,'measured'); B_point1nm=(0.1e-9/lam)*c/lam; D=0.8*16000e-12/1e-9; L=length(recv); t1=(1:L)*(ts/N*N_b); df=1/(L*(ts/N*N_b)); f1=(0:(L-1))*df; %generate carrier shift fc1=-N/ts/2; A_fc1=exp(j*2*pi*fc1*t1); A_fc=exp(j*2*pi*fc*t1); % carrier frequency %Define modulation index St=recv.*A_fc; St1=St; % shift OFDM central frequency to 0-Hz r_st=real(St1); i_st=imag(St1); r_st=r_st/sqrt(mean(r_st.*r_st)); i_st=i_st/sqrt(mean(i_st.*i_st)); Vpi=7; % set Vpi value zl=[]; m1=[0.05]; for k=1:length(m1) Vrout=m1(k)/(2*sqrt(2))*r_st; Viout=m1(k)/(2*sqrt(2))*i_st; Vrout=(2*Vpi/pi)*(Vrout); %without predistortion Viout=(2*Vpi/pi)*(Viout); %without predistortion end % MZM modulation Str=cos(pi/(2*Vpi)*Vrout+pi/2).*A_fc; Sti=cos(pi/(2*Vpi)*Viout+pi/2).*A_fc*exp(j*(pi/2)); recv=Str+Sti; % compute the highest frequency fmax=N/ts+fc; f1=(0:(L-1))*df; % frequency in GHz % going through the fiber phi_D=pi*lam^2*D/c*(f1-fmax).^2; A_dis=exp(j*phi_D); f_st=fft(recv)/length(recv); %将信号转化成频域信号，因为色散引起的相位A_dis是在频率上的 f_st=f_st.*A_dis; %--------------------------- %f_st=Filter1.*f_st; %-------------------------- St_dis=length(recv)*ifft(f_st); av1=sqrt(mean(St_dis.*conj(St_dis))); St_dis=St_dis*sqrt(N)/(av1); OSNR_L=10^(OSNR/10); Pn=0.5*(B/B_point1nm)*mean(recv.*conj(recv))/OSNR_L; randn('state',0); NI=randn(2,L); nt=(1/sqrt(2))*sqrt(Pn)*(NI(1,:)+j*NI(2,:)); St_dis=St_dis+nt; f_st_dis=fft(St_dis)/length(St_dis); %------------------------------------ %f_st_dis=Filter1.*f_st_dis; %------------------------------------ St_dis=length(recv)*ifft(f_st_dis); % add phase noise T=(ts/(N_b*N))*length(St_dis); N1=length(St_dis); rand('state',0); t=(1:N1)*T/N1; t=t*1e9 ; Phi=rand(1,N1/2)*2*pi; Phi(N1/2+1)=0; Phi(1)=0; I1=1:(N1/2+1); omega(1:N1/2+1)=2*pi*(I1-1)/T;omega(1)=omega(2); % to avoid omega=0; A1=sqrt(2*pi*f3db*2./(T*omega.^2)); AM=A1.*exp(i*Phi); AM(N1/2+1:N1)=fliplr(conj(AM(2:(N1/2)+1))); Phit=N1*ifft(AM); Phit=real(Phit);dphit=Phit-Phit(1); %plot(t,dphit); St_dis=exp(j*Phit).*St_dis; %*****************balanced Receiver***************************** St_dis=St_dis.*conj(A_fc1); I1=(abs(St_dis+A_fc)).^2; I2=(abs(St_dis+A_fc.*exp(j*pi))).^2; I=I1-I2; Q1=(abs(St_dis+A_fc.*exp(j*pi/2))).^2; Q2=(abs(St_dis+A_fc.*exp(j*pi).*exp(j*pi/2))).^2; Q=Q1-Q2; St_dis=I+j*Q; %*****************计算符号定时***************************** P=zeros(1,length(St_dis)-Ns-L1); R=zeros(1,length(St_dis)-Ns-L1); for d = 1:1:length(St_dis)-Ns-L1 %for m=0:1:L-1 %recv(m+1:m+L)=recv(m+1:m+L).*conj(buf1); recvw(d:d+L1-1)=St_dis(d:d+L1-1).*conj(buf1); Pc(d)=(recvw(d:d+L1-1))*conj(St_dis(d+Ns:d+Ns+L1-1)).'; Pz(d)=sum(abs(recvw(d:d+L1-1)).^2,2)+sum(abs(St_dis(d+Ns:d+Ns+L1-1)).^2,2); %P(d)=Pc(d)/Pz(d); %R(d-Ns/2) = R(d-Ns/2) + power(abs(recv(d+N/2+m)),2); %end end M=power(abs(Pc),2); P1=M./((Pz./2).^2); [peak head]=max(P); %**********************绘图****************************** figure; %plot(Ng,0:0.03:1,':k'); %hold on %plot(N+Ng,0:0.01:1,':k'); %hold on %plot(N+2*Ng,0:0.03:1,':k'); %hold on %plot((N+Ng+L1),0:0.03:1,'.k'); %hold on %plot((N+2*Ng+L1),0:0.03:1,':k'); % hold on %plot((2*N+3*Ng+L1),0:0.01:1,':k'); %hold on % plot((2*N+3*Ng+L1),0:0.03:1,':k'); % hold on %plot((2*N+3*Ng+L1),0:0.03:1,':k'); %hold on %plot((3*N+3*Ng+L1),0:0.01:1,':k'); %hold on t=1:1:576; %plot(t,P(t),'r:'); %axis([0,576,0,1]); % going through the fiber D=0; phi_D=pi*lam^2*D/c*(f1-fmax).^2; A_dis=exp(j*phi_D); f_st=fft(recv)/length(recv); %将信号转化成频域信号，因为色散引起的相位A_dis是在频率上的 f_st=f_st.*A_dis; %--------------------------- %f_st=Filter1.*f_st; %-------------------------- St_dis=length(recv)*ifft(f_st); av1=sqrt(mean(St_dis.*conj(St_dis))); St_dis=St_dis*sqrt(N)/(av1); OSNR_L=10^(OSNR/10); Pn=0.5*(B/B_point1nm)*mean(recv.*conj(recv))/OSNR_L; randn('state',0); NI=randn(2,L); nt=(1/sqrt(2))*sqrt(Pn)*(NI(1,:)+j*NI(2,:)); St_dis=St_dis+nt; f_st_dis=fft(St_dis)/length(St_dis); %------------------------------------ %f_st_dis=Filter1.*f_st_dis; %------------------------------------ St_dis=length(recv)*ifft(f_st_dis); % add phase noise T=(ts/(N_b*N))*length(St_dis); N1=length(St_dis); rand('state',0); t=(1:N1)*T/N1; t=t*1e9 ; Phi=rand(1,N1/2)*2*pi; Phi(N1/2+1)=0; Phi(1)=0; I1=1:(N1/2+1); omega(1:N1/2+1)=2*pi*(I1-1)/T;omega(1)=omega(2); % to avoid omega=0; A1=sqrt(2*pi*f3db*2./(T*omega.^2)); AM=A1.*exp(i*Phi); AM(N1/2+1:N1)=fliplr(conj(AM(2:(N1/2)+1))); Phit=N1*ifft(AM); Phit=real(Phit);dphit=Phit-Phit(1); %plot(t,dphit); St_dis=exp(j*Phit).*St_dis; %*****************balanced Receiver***************************** St_dis=St_dis.*conj(A_fc1); I1=(abs(St_dis+A_fc)).^2; I2=(abs(St_dis+A_fc.*exp(j*pi))).^2; I=I1-I2; Q1=(abs(St_dis+A_fc.*exp(j*pi/2))).^2; Q2=(abs(St_dis+A_fc.*exp(j*pi).*exp(j*pi/2))).^2; Q=Q1-Q2; St_dis=I+j*Q; %*****************计算符号定时***************************** P=zeros(1,length(St_dis)-Ns-L1); R=zeros(1,length(St_dis)-Ns-L1); for d = 1:1:length(St_dis)-Ns-L1 %for m=0:1:L-1 %recv(m+1:m+L)=recv(m+1:m+L).*conj(buf1); recvw(d:d+L1-1)=St_dis(d:d+L1-1).*conj(buf1); Pc(d)=(recvw(d:d+L1-1))*conj(St_dis(d+Ns:d+Ns+L1-1)).'; Pz(d)=sum(abs(recvw(d:d+L1-1)).^2,2)+sum(abs(St_dis(d+Ns:d+Ns+L1-1)).^2,2); %P(d)=Pc(d)/Pz(d); %R(d-Ns/2) = R(d-Ns/2) + power(abs(recv(d+N/2+m)),2); %end end M=power(abs(Pc),2); P=M./((Pz./2).^2); [peak head]=max(P)