基于Matlab实现64QAM调制解调系统仿真（源码）.rar

clear all; close all; %生成信号 [yY,t_ms] = QAM_T5(4096,4000,pi/3); % 加入噪声 SNR = 15; %信噪比 y = awgn(yY,SNR); %加入噪声 y = 0.8*y; %衰减 u = 0.5; %符号偏移 %匹配滤波 Nsym = 6; beta = 1; sampsPerSym = 8; R = 125000; % Data rate fs = R * sampsPerSym; % Sampling frequency fltDelay = Nsym / (2*R); rcrFilt = comm.RaisedCosineReceiveFilter(... 'Shape', 'Square root', ... 'RolloffFactor', beta, ... 'FilterSpanInSymbols', Nsym, ... 'InputSamplesPerSymbol', sampsPerSym, ... 'DecimationFactor', 1); %下变频 Fs = 50000000; %载波NCO输入频率50M P = 2^32; %载波NCO计数上限 Phid = floor(P/50); %载波NCO初始相位步进量 1M Fc = Phid/P*Fs; %初始载波频率1M Phi = 0; %载波NCO初始计数值 ycos = zeros(1,length(y)); %载波NCO输出余弦序列内存预分配 ysin = zeros(1,length(y)); %载波NCO输出正弦序列内存预分配 Ii = zeros(1,length(y)); %I信号变频输出内存预分配 Qi = zeros(1,length(y)); %Q信号变频输出内存预分配 % hr = 0; %载波NCO溢出时刻寄存器 Ii1 = 0; Ii2 = 0; Ii3 = 0; Ii4 = 0; Qi1 = 0; Qi2 = 0; Qi3 = 0; Qi4 = 0; nc = 1; Io0 = 0; Io1 = 0; Io2 = 0; Io3 = 0; Io4 = 0; Qo0 = 0; Qo1 = 0; Qo2 = 0; Qo3 = 0; Qo4 = 0; n = 0; % n = 0; %50M->640k抽样计数器 C = [ 0, 0, 1, 0; %拉格朗日内插滤波器Farrow结构抽头系数表 -1/6, 1,-1/2,-1/3; 0, 1/2, -1, 1/2; 1/6,-1/2, 1/2,-1/6;]; k = 1; %80k符号抽样计数器 etam = 1; %码定时恢复 内插滤波器 NCO计数器η(m) wk = 0.125; %码定时恢复 内插滤波器 NCO相位步进量1/8 uk = 0; %码定时恢复 内插滤波器 小数间隔 ek = 0; %码定时恢复 Gardner算法 误差值 yIkh = 0; %码定时恢复 Gardner算法 I信号采样码中间时刻值 yQkh = 0; %码定时恢复 Gardner算法 Q信号采样码中间时刻值 de = 0.8; %码定时恢复 Gardner算法改进 最小门限 c1 = 0; %码定时恢复 Gardner算法改进 超时计数器 cp = 0.8; %码定时恢复 环路积分器p参数 cd = 0; %码定时恢复 环路积分器d参数 % 星座图 轨道（能级由低到高），点（90°旋转对称点），象限（点旋转经过的4象限），坐标（x，y） S(1,1,:,:) = [ [ 1, 1];[-1, 1];[-1,-1];[ 1,-1] ]; S(2,1,:,:) = [ [ 3, 1];[-1, 3];[-3,-1];[ 1,-3] ]; S(2,2,:,:) = [ [ 1, 3];[-3, 1];[-1,-3];[ 3,-1] ]; S(3,1,:,:) = [ [ 3, 3];[-3, 3];[-3,-3];[ 3,-3] ]; S(4,1,:,:) = [ [ 5, 1];[-1, 5];[-5,-1];[ 1,-5] ]; S(4,2,:,:) = [ [ 1, 5];[-5, 1];[-1,-5];[ 5,-1] ]; S(5,1,:,:) = [ [ 5, 3];[-3, 5];[-5,-3];[ 3,-5] ]; S(5,2,:,:) = [ [ 3, 5];[-5, 3];[-3,-5];[ 5,-3] ]; S(6,1,:,:) = [ [ 7, 1];[-1, 7];[-7,-1];[ 1,-7] ]; S(6,2,:,:) = [ [ 5, 5];[-5, 5];[-5,-5];[ 5,-5] ]; S(6,3,:,:) = [ [ 1, 7];[-7, 1];[-1,-7];[ 7,-1] ]; S(7,1,:,:) = [ [ 7, 3];[-3, 7];[-7,-3];[ 3,-7] ]; S(7,2,:,:) = [ [ 3, 7];[-7, 3];[-3,-7];[ 7,-3] ]; S(8,1,:,:) = [ [ 7, 5];[-5, 7];[-7,-5];[ 5,-7] ]; S(8,2,:,:) = [ [ 5, 7];[-7, 5];[-5,-7];[ 7,-5] ]; S(9,1,:,:) = [ [ 7, 7];[-7, 7];[-7,-7];[ 7,-7] ]; Err = 0; %载波同步 单次误差 (-pi/4,pi/4) yik = 0; %载波同步 I信号判决坐标 yqk = 0; %载波同步 Q信号判决坐标 phir = 0; %载波同步 累计误差 PD算法改进 c0 = 0; %载波同步 PD算法 计数器 Pp = 2; %载波同步 PD算法 环路积分p参数 Pi = 0.1; %载波同步 PD算法 环路积分i参数 Pd = 4; %载波同步 PD算法 环路积分d参数 Pp2 = 2; %载波同步 DD算法 环路积分p参数 Pi2 = 0.005; %载波同步 DD算法 环路积分i参数 Pd2 = 4; %载波同步 DD算法 环路积分d参数 Ap = 1; %自动增益控制 CMA算法 环路积分p参数 Ai = 0.1; %自动增益控制 CMA算法 环路积分i参数 Ad = 8; %自动增益控制 CMA算法 环路积分d参数 Gain = 1; %自动增益控制 CMA算法 增益值 s = 0; %码定时恢复 载波同步 自动增益控制 状态机 状态变量 PD = 0; %状态机 载波同步 PD算法 状态变量 dk = 4; %状态机 码定时恢复 误差阈值 da = 0.4; %状态机 自动增益控制 误差阈值 dp = 0.005; %状态机 载波同步 PD算法 误差阈值 Ek = dk; %状态机 码定时恢复 误差值 Ea = da; %状态机 自动增益控制 误差值 Ep = dp; %状态机 载波同步 误差值 % parpool('local',6); for h = 1 : length(y) if h > 1; Phid(h) = Phid(h-1); % Phid(h) = Phid(1); end Phi = Phi + Phid(h); ycos(h) = cos(2*pi*mod(Phi/P,1)); ysin(h) = sin(2*pi*mod(Phi/P,1)); Ii(h) = y(h).*ycos(h); Qi(h) = y(h).*(-ysin(h)); % 积分 Ii1 = Ii1 + Ii(h); Qi1 = Qi1 + Qi(h); % 25倍抽取，梳状滤波 nr = n; if mod(h,25) == 0 nc = nc + 1; Io0(nc) = Ii1; Io1(nc) = Io0(nc) - Io0(nc-1); % Io2(nc) = Io1(nc) - Io1(nc-1); % Io3(nc) = Io2(nc) - Io2(nc-1); % Io4(nc) = Io3(nc) - Io3(nc-1); Qo0(nc) = Qi1; Qo1(nc) = Qo0(nc) - Qo0(nc-1); % Qo2(nc) = Qo1(nc) - Qo1(nc-1); % Qo3(nc) = Qo2(nc) - Qo2(nc-1); % Qo4(nc) = Qo3(nc) - Qo3(nc-1); if mod(nc,2) == 0 n = n + 1; Io(n) = Io1(nc)/12; Qo(n) = Qo1(nc)/12; end end %匹配滤波 if n > nr && n > 48 yi = step(rcrFilt,Io(n-48:n)'); yq = step(rcrFilt,Qo(n-48:n)'); yI(n-48) = yi(49); yQ(n-48) = yq(49); end %%加入符号偏移%% if n > nr && n > 48+3 xIm(n-49-2) = (((yI(n-48)*C(4,1)+yI(n-48-1)*C(4,2)+yI(n-48-2)*C(4,3)+yI(n-48-3)*C(4,4))*u+(yI(n-48)*C(3,1)+yI(n-48-1)*C(3,2)+yI(n-48-2)*C(3,3)+yI(n-48-3)*C(3,4)))*u+(yI(n-48)*C(2,1)+yI(n-48-1)*C(2,2)+yI(n-48-2)*C(2,3)+yI(n-48-3)*C(2,4)))*u+(yI(n-48)*C(1,1)+yI(n-48-1)*C(1,2)+yI(n-48-2)*C(1,3)+yI(n-48-3)*C(1,4)); xQm(n-49-2) = (((yQ(n-48)*C(4,1)+yQ(n-48-1)*C(4,2)+yQ(n-48-2)*C(4,3)+yQ(n-48-3)*C(4,4))*u+(yQ(n-48)*C(3,1)+yQ(n-48-1)*C(3,2)+yQ(n-48-2)*C(3,3)+yQ(n-48-3)*C(3,4)))*u+(yQ(n-48)*C(2,1)+yQ(n-48-1)*C(2,2)+yQ(n-48-2)*C(2,3)+yQ(n-48-3)*C(2,4)))*u+(yQ(n-48)*C(1,1)+yQ(n-48-1)*C(1,2)+yQ(n-48-2)*C(1,3)+yQ(n-48-3)*C(1,4)); end % 能量检测 if n > nr && n > 48+3 if n > 51 + 32 E(n-49-2) = 0; for j = 0 : 31; E(n-49-2) = E(n-49-2) + abs(xIm(n-49-2-j))/32 + abs(xQm(n-49-2-j))/32; end else E(n-49-2) = 0; for j = 0 : n-52; E(n-49-2) = E(n-49-2) + abs(xIm(n-49-2-j))/32 + abs(xQm(n-49-2-j))/32; end end end % 状态控制 if n > nr && n > 48+3 if s == 0 if E(n-49-2) >= 1 s = 1; end elseif s == 1 if E(n-49-2) < 1 s = 0; PD = 0; else if k > 15 && Ek < dk s = 2; end end elseif s == 2 if E(n-49-2) < 1 s = 0; PD = 0; else if k > 15 && Ek >= dk s = 1; elseif k > 31 %&& Ea < da if PD == 0 s = 4; PD = 1; c0 = 0; elseif Ea < da s = 3; end end end elseif s == 3 if E(n-49-2) < 1 s = 0; PD = 0; else if k > 15 && Ek >= dk s = 1; elseif k > 31 && Ea >= da s = 2; end end elseif s == 4 if E(n-49-2) < 1 s = 0; PD = 0; else if Ep < dp PD = 2; s = 2; end end