test.m

锁相环的仿真以及含有文本

clear all;
clc;close all;
fs=10e4; %采样率，Hz
T=0.1; %仿真时间,s
Tsym=1e-3;%符号持续时间
N=fs*T; %仿真数据长度
M=fs*Tsym;%一个符号的长度
L=T/Tsym;
fc=10.1e3;f0=10e3;%信号中心频率为fc;锁相环中心频率为f0.
t=0:1/fs:T-1/fs; %采样时间序列
Delta_Phase=rand(1)*2*pi; %随机初相，Rad
Carrier=exp(j*(2*pi*fc*t+Delta_Phase)); %载波信号

%% 基带信号
PLL_Mode = 1; %定义锁相环的工作模式：单载波为“1”、BPSK调制为“2”、QPSK调制为“3”
if PLL_Mode == 1
I_Data=ones(1,N);
Q_Data=I_Data;
elseif PLL_Mode == 2

I_Bit=randint(1,L)*2-1; %bpsk调制
I_Data=I_Bit(1)*ones(1,M);
for k=2:L
I_Data=[I_Data I_Bit(k)*ones(1,M)]; %生成各采样点数据
end
Q_Data=zeros(1,length(I_Data));
else
    
I_Bit=randint(1,L)*2-1; %bpsk调制
I_Data=I_Bit(1)*ones(1,M);
Q_Bit=randint(1,L)*2-1; %bpsk调制
Q_Data=Q_Bit(1)*ones(1,M);
for k=2:L
I_Data=[I_Data I_Bit(k)*ones(1,M)]; %生成各采样点数据
Q_Data=[Q_Data I_Bit(k)*ones(1,M)];
end
end
Signal_Source=I_Data + j*Q_Data; %数据bit
% figure
% subplot(2,1,1)
% plot(real(Signal_Source))
% subplot(2,1,2)
% plot(imag(Signal_Source))
%% 调制处理
Signal=Signal_Source(1:length(Carrier)).*Carrier; %数据的I、Q路都是信号，都通过
%载波发送，不需与载波I、Q相对应

snr=6;
Signal_Channel=awgn(Signal,snr,'measured');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%以下为锁相环处理过程
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% 参数清零
Signal_PLL=zeros(1,N);
NCO_Freq =ones(1,N)*f0; 
NCO_Phase =zeros(1,N);
NCO_Out=exp(j*NCO_Phase);
Discriminator_Out=zeros(1,N);
Phase_Control=zeros(1,N);
Freq_Control=zeros(1,N);
PLL_Phase_Part=zeros(1,N);
PLL_Freq_Part=zeros(1,N);
Phase_Diff=zeros(1,N);

%环路参数
gain=2000; %锁相环增益k0*k1
ksi=0.707; %阻尼因子，瞬态响应平坦，不过缓，环路带宽较小
Loopwidth=20; %环路带宽,Hz，与pll稳定性有关
digital_sensitivity=20; %数字灵敏度，k1/fs
k1=digital_sensitivity*fs; %nco的增益
k0=gain/k1; %鉴相器后的增益
%以上参数加上采样率即可确定该类锁相环的特性
omega_n=2*Loopwidth/(ksi+1/4/ksi); %自然频率,rad/s
C2=omega_n^2/fs/gain; %将有源比例积分滤波器用双线性变换进行数字化即得C1,C2
C1=2*omega_n*ksi/gain-C2/2;

for n=2:N
%鉴相器
Signal_PLL(n)=Signal_Channel(n)*NCO_Out(n-1);
I_PLL(n)=real(Signal_PLL(n));
Q_PLL(n)=imag(Signal_PLL(n));
if PLL_Mode == 1
%四象限反正切
Discriminator_Out(n)=atan2(Q_PLL(n),I_PLL(n))*k0;
elseif PLL_Mode == 2
Discriminator_Out(n)=I_PLL(n)*Q_PLL(n)/abs(Signal_PLL(n))^2*k0;
else
sgn=cal_sign(I_PLL(n));
Discriminator_Out(n)=(sgn*Q_PLL(n)-sgn*I_PLL(n))/(sqrt(2)*abs(Signal_PLL(n)))*k0;
end

PLL_Phase_Part(n)=Discriminator_Out(n)*C1; %暂态相差
%Discriminator_Out是相位差 
PLL_Freq_Part(n)=Discriminator_Out(n)*C2+PLL_Freq_Part(n-1); %稳态相差，C2可视为暂态相差对稳态相差的加权系数
Phase_Control(n)=PLL_Phase_Part(n)+PLL_Freq_Part(n); %经过滤波的一个采样间隔内的相位差

%数控振荡器，中心频率和初始相位设为0。实际只能调频，通过调频来调相
% NCO_Phase(n)=NCO_Phase(n-1)+Freq_Control(n)/Freq_Sample*k1;
% NCO_Out(n)=exp(-j*NCO_Phase(n));
Phase_Diff(n)=Phase_Control(n)/fs*k1; %加权后一个采样间隔的相位差
NCO_Freq(n)=Phase_Diff(n)*fs; %除以采样间隔，得到pll中心频率与输入信号的频率差
NCO_Phase(n)=NCO_Phase(n-1)+NCO_Freq(n)/fs; %本地载波相位
NCO_Out(n)=exp(-j*NCO_Phase(n)); %本振信号
end
PLL_Freq=PLL_Freq_Part*fs/2/pi*digital_sensitivity; 
%Hz,锁定的频率为归一化角频率，即freq_pll=f0/fs*2*pi
PLL_Phase=Discriminator_Out*180/pi; %度 %PLL_Phase_Part(2:Simulation_Length)*180/pi;
%画图显示结果
figure
subplot(2,2,1),plot(t,PLL_Freq,'*')
grid on
xlabel('时间(s)')
ylabel('输出频率(Hz)')
title('锁相环频率响应曲线')
%axis([1 Simulation_Length -100 100])
subplot(2,2,2),plot(t,PLL_Phase)
grid on
xlabel('时间(s)')
ylabel('相位差(°)')
title('锁相环相位响应曲线')
%axis([1 Simulation_Length -2 2])