基于matlab QPSK的锁相，用来进行QPSK通信系统的仿真和实际信号载波同步的提取

%%% QPSK PLL clear all; close all; % loop parameters % digital bandwidth, theta_0 = w_0*T/2 % damping factor, eta = sqrt(2)/2 theta_0=2*pi/1000; eta = 1; k_i=(4*theta_0*theta_0)/(1+2*eta*theta_0+theta_0*theta_0); k_p=(4*eta*theta_0)/(1+2*eta*theta_0+theta_0*theta_0); % at the transimitter n_dat=5000; x1=(floor(2*rand(1,n_dat))-0.5)/0.5+j*((floor(2*rand(1,n_dat))-0.5)/0.5); x2=reshape([x1;zeros(3,n_dat)],1,4*n_dat); h=rcosine(1,4,'sqrt',0.5,6); h=h/max(h); x3=filter(h,1,x2); % Channel x3_ch = x3.*exp(j*0.3); plot(x3_ch(1:4:4*n_dat),'rx') %% rotate the constellation hold off axis('square') grid on % Receiver phi_2_save=zeros(1,4*n_dat); int_1_save=zeros(1,4*n_dat); int_2_save=zeros(1,4*n_dat); dphi_sv=zeros(1,4*n_dat); int_1=0; int_2=0; SCL=h*h'; gg = h/SCL; reg1 = zeros(1,length(gg)); % MF reg reg2 = zeros(1,4); for n = 1:4*n_dat prod(n) = x3_ch(n)*exp(-j*int_2); % compute one input, high rate reg1 = [prod(n) reg1(1:length(gg)-1)]; MF_out(n) = reg1*gg'; if (mod(n,4)==1) det_out(ceil(n/4)) = conj(MF_out(n))*(sign(real(MF_out(n)))+j*sign(imag(MF_out(n)))); % slicer output (Angle difference) ATN(ceil(n/4)) = -atan(imag(det_out(ceil(n/4)))/real(det_out(ceil(n/4)))); % loop filter int_1 = int_1 + k_i*ATN(ceil(n/4)); int_1_save(ceil(n/4)) = int_1; dphi = k_p*ATN(ceil(n/4))+int_1; dphi_sv(ceil(n/4)) = dphi; % hold for 4 samples reg2 = dphi.*ones(1,4); % dss phi_2_save(n) = int_2; int_2 = int_2 + reg2(1); end if(mod(n,4)==2) phi_2_save(n) = int_2; int_2 = int_2 + reg2(2); end if(mod(n,4)==3) phi_2_save(n) = int_2; int_2 = int_2 + reg2(3); end if(mod(n,4)==0) phi_2_save(n) = int_2; int_2 = int_2 + reg2(4); end end x4 = filter(h/SCL,1,x3_ch); % if no PLL output = MF_out(1:4:4*n_dat); % final output figure subplot(221) plot(x3_ch(1:4:4*n_dat),'rx') title('Received Signal Constellation') grid on axis('square') subplot(222) plot(x4(1:4:4*n_dat),'rx') grid on axis('square') % constellation if without PLL title('Matched Filter Ouput without PLL') subplot(223) plot(output(1:1000),'rx') % first 100 samples of the output grid on axis('square') title('First 500 Samples from MF with PLL') subplot(224) plot(output(2000:5000),'rx') % 2000th - 5000th sample of the output grid on axis('square') title('The Rest Ouput from MF') figure subplot(2,1,2) plot(phi_2_save(1:3000)) grid on axis([-20 2500 -0.1 0.45]) title('Output of The Phase Accumulator') subplot(2,1,1) plot(dphi_sv(1:200)) grid on axis([-5 200 -0.025 0.025]) title('Output of The Loop Filter')