目录
1.QAM的调制原理2.QAM的解调原理3.QAM代码4.结果图5.特点6.加星座图的16QAM代码
1.QAM的调制原理
QAM调制原理如下图所示，基带码元波形经过串并转换分成I、Q两路，然后再经过电平转换(00转换成-1,01转换成-3,10转换成1,11转换成3)，再与对应的载波相乘，然后再相加完成QAM的调制。
本次采用的是16QAM，M=16=2^k,因此k=4,
L=2^(k/2)=4
2.QAM的解调原理
QAM的解调原理如下图所示，QAM信号再分为I、Q两路和对应的载波相乘，然后经过低通滤波器后进行抽样判决，判决之后的I、Q路码元进行合并，I路为最终码元序列的奇数位置码元，Q路为最终码元序列的偶数位置码元，恢复出原始的码元序列。
3.QAM代码
clear all; % 清除所有变量close all; % 关闭所有窗口clc; % 清屏%% 基本参数M=40; % 产生码元数 L=100; % 每码元复制L次,每个码元采样次数Ts=0.001; % 每个码元的宽度,即码元的持续时间Rb=1/Ts; % 码元速率1Kdt=Ts/L; % 采样间隔TotalT=M*Ts; % 总时间1t=0:dt:TotalT-dt; % 时间1TotalT2=(M/2)*Ts; % 总时间2t2=0:dt:TotalT2-dt; % 时间2Fs=1/dt; % 采样间隔的倒数即采样频率%% 产生单极性波形wave=randi([0,1],1,M); % 产生二进制随机码,M为码元个数fz=ones(1,L); % 定义复制的次数L,L为每码元的采样点数fz2=ones(1,2*L); % 定义复制的次数2Lx1=wave(fz,:); % 将原来wave的第一行复制L次，称为L*M的矩阵jidai=reshape(x1,1,L*M); % 产生单极性不归零矩形脉冲波形，将刚得到的L*M矩阵，按列重新排列形成1*(L*M)的矩阵%% 产生I、Q路码元a=wave(1:2:end); % I路码元是原始码元奇数位置的码元b=wave(2:2:end); % Q路码元是原始码元偶数位置的码元% 码元电平转换，要将两两码元转换成四种电平，00转换成-1,01转换成-3,10转换成1,11转换成3for i=1:length(a)/2 if a(2*(i-1)+1:2*i)==[0 1] I(i)=-3; elseif a(2*(i-1)+1:2*i)==[0 0] I(i)=-1; elseif a(2*(i-1)+1:2*i)==[1 0] I(i)=1; elseif a(2*(i-1)+1:2*i)==[1 1] I(i)=3; end if b(2*(i-1)+1:2*i)==[0 1] Q(i)=-3; elseif b(2*(i-1)+1:2*i)==[0 0] Q(i)=-1; elseif b(2*(i-1)+1:2*i)==[1 0] Q(i)=1; elseif b(2*(i-1)+1:2*i)==[1 1] Q(i)=3; endendx2=a(fz,:); % 将原来a的第一行复制L次，称为L*(M/2)的矩阵a_mayuan=reshape(x2,1,L*(M/2));% 产生单极性不归零矩形脉冲波形，将刚得到的L*(M/2)矩阵，按列重新排列形成1*(L*(M/2))的矩阵x3=b(fz,:); % 将原来b的第一行复制L次，称为L*(M/2)的矩阵b_mayuan=reshape(x3,1,L*(M/2));% 产生单极性不归零矩形脉冲波形，将刚得到的L*(M/2)矩阵，按列重新排列形成1*(L*(M/2))的矩阵x4=I(fz2,:); % 将原来I的第一行复制2L次，称为2L*(M/4)的矩阵I_mayuan=reshape(x4,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵x5=Q(fz2,:); % 将原来Q的第一行复制2L次，称为2L*(M/4)的矩阵Q_mayuan=reshape(x5,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵%% 码元波形绘制figure(1); % 绘制第1幅图subplot(511); % 窗口分割成5*1的，当前是第1个子图 plot(t,jidai,'LineWidth',2);% 绘制基带码元波形，线宽为2title('基带信号波形'); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT,-0.1,1.1]) % 坐标范围限制subplot(512); % 窗口分割成5*1的，当前是第2个子图 plot(t2,a_mayuan,'LineWidth',2);% I路码元波形，线宽为2title('I路码元波形'); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-0.1,1.1]) % 坐标范围限制subplot(513); % 窗口分割成5*1的，当前是第3个子图 plot(t2,b_mayuan,'LineWidth',2);% 绘制Q路码元波形，线宽为2title('Q路码元波形'); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-0.1,1.1]) % 坐标范围限制subplot(514); % 窗口分割成5*1的，当前是第4个子图 plot(t2,I_mayuan,'LineWidth',2);% 绘制转换后I路电平波形，线宽为2title('电平转换后I路码元波形');% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-4,4]) % 坐标范围限制subplot(515); % 窗口分割成5*1的，当前是第5个子图 plot(t2,Q_mayuan,'LineWidth',2);% 绘制转换后Q路电平波形，线宽为2title('电平转换后Q路码元波形');% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-4,4]) % 坐标范围限制%% QAM调制fc=10000; % 载波频率2kHz zb1=cos(2*pi*fc*t2); % 载波1zb2=-sin(2*pi*fc*t2); % 载波2I_lu=I_mayuan.*zb1; % I路波形 Q_lu=Q_mayuan.*zb2; % Q路波形 qam=I_lu+Q_lu; % QAM调制figure(2); % 绘制第2幅图subplot(411) % 窗口分割成4*1的，当前是第1个子图 plot(t2,I_lu,'LineWidth',2);% 绘制I路的波形 title('I路信号波形') % 标题axis([0,TotalT2,-4,4]); % 坐标范围限制xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(412) % 窗口分割成4*1的，当前是第2个子图 plot(t2,Q_lu,'LineWidth',2);% 绘制Q路的波形 title('Q路信号波形') % 标题axis([0,TotalT2,-4,4]); % 坐标范围限制xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(413) % 窗口分割成4*1的，当前是第3个子图 plot(t2,qam,'LineWidth',2); % 绘制QAM的波形 title('QAM信号波形') % 标题axis([0,TotalT2,-7,7]); % 坐标范围限制xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 信号经过高斯白噪声信道tz=awgn(qam,20); % 信号qam中加入白噪声，信噪比为SNR=20dBsubplot(414); % 窗口分割成4*1的，当前是第4个子图 plot(t2,tz,'LineWidth',2); % 绘制QAM信号加入白噪声的波形axis([0,TotalT2,-7,7]); % 坐标范围设置title('通过高斯白噪声信道后的信号');% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 解调部分jt1=tz.*zb1; % 相干解调，I路乘以相干载波jt2=tz.*zb2; % 相干解调，Q路乘以相干载波figure(3);subplot(511) % 窗口分割成5*1的，当前是第1个子图 plot(t2,jt1,'LineWidth',2) % 绘制I路乘以相干载波后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路乘以相干载波后的信号")% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(512) % 窗口分割成5*1的，当前是第2个子图 plot(t2,jt2,'LineWidth',2) % 绘制Q路乘以相干载波后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("Q路乘以相干载波后的信号")% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 加噪信号经过滤波器% 低通滤波器设计fp=2*Rb; % 低通滤波器截止频率，乘以2是因为下面要将模拟频率转换成数字频率wp=Rb/(Fs/2)b=fir1(30, fp/Fs, boxcar(31));% 生成fir滤波器系统函数中分子多项式的系数% fir1函数三个参数分别是阶数，数字截止频率，滤波器类型% 这里是生成了30阶(31个抽头系数)的矩形窗滤波器[h,w]=freqz(b, 1,512); % 生成fir滤波器的频率响应% freqz函数的三个参数分别是滤波器系统函数的分子多项式的系数，分母多项式的系数(fir滤波器分母系数为1)和采样点数(默认)512lvbo1=2*fftfilt(b,jt1); % 对信号进行滤波，jt1是等待滤波的信号，b是fir滤波器的系统函数的分子多项式系数lvbo2=2*fftfilt(b,jt2); % 对信号进行滤波，jt2是等待滤波的信号，b是fir滤波器的系统函数的分子多项式系数subplot(513); % 窗口分割成5*1的，当前是第3个子图 plot(w/pi*Fs/2,20*log(abs(h)),'LineWidth',2); % 绘制滤波器的幅频响应title('低通滤波器的频谱'); % 标题xlabel('频率/Hz'); % x轴标签ylabel('幅度/dB'); % y轴标签subplot(514) % 窗口分割成5*1的，当前是第4个子图 plot(t2,lvbo1,'LineWidth',2);% 绘制I路经过低通滤波器后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过低通滤波器后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(515) % 窗口分割成5*1的，当前是第5个子图 plot(t2,lvbo2,'LineWidth',2);% 绘制Q路经过低通滤波器后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过低通滤波器后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 抽样判决% 滤波后的波形进行判决，大于2的判决为3,0-2之间的判决为1，-2-0之间判决为-1，小于-2的判决为-3for j=1:length(lvbo1) if lvbo1(j)>=2 I_panjue(j)=3; elseif (lvbo1(j)>0 && lvbo1(j)<2) I_panjue(j)=1; elseif (lvbo1(j)>=-2 && lvbo1(j)<0) I_panjue(j)=-1; elseif lvbo1(j)<-2 I_panjue(j)=-3; endendfor k=1:length(lvbo2) if lvbo2(k)>=2 Q_panjue(k)=3; elseif (lvbo2(k)>0 && lvbo2(k)<2) Q_panjue(k)=1; elseif (lvbo2(k)>=-2 && lvbo2(k)<0) Q_panjue(k)=-1; elseif lvbo2(k)<-2 Q_panjue(k)=-3; endendfigure(4);subplot(611) % 窗口分割成6*1的，当前是第1个子图 plot(t2,I_panjue,'LineWidth',2);% 绘制经过判决的I路信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(612) % 窗口分割成6*1的，当前是第2个子图 plot(t2,Q_panjue,'LineWidth',2);% 绘制经过判决的Q路信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("Q路经过判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 判决之后再处理，取码元中间时刻判决I_zong=[];Q_zong=[];for m=L:2*L:((2*L)*(M/4)) if I_panjue(m)>=2 I_zong=[I_zong,3]; elseif (I_panjue(m)>0 && I_panjue(m)<2) I_zong=[I_zong,1]; elseif (I_panjue(m)>=-2 && I_panjue(m)<0) I_zong=[I_zong,-1]; elseif I_panjue(m)<-2 I_zong=[I_zong,-3]; endendfor n=L:2*L:((2*L)*(M/4)) if Q_panjue(n)>=2 Q_zong=[Q_zong,3]; elseif (Q_panjue(n)>0 && Q_panjue(n)<2) Q_zong=[Q_zong,1]; elseif (Q_panjue(n)>=-2 && Q_panjue(n)<0) Q_zong=[Q_zong,-1]; elseif Q_panjue(n)<-2 Q_zong=[Q_zong,-3]; endendx6=I_zong(fz2,:); % 将原来I_zong的第一行复制2L次，称为2L*(M/4)的矩阵I_wave=reshape(x6,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵x7=Q_zong(fz2,:); % 将原来Q_zong的第一行复制2L次，称为2L*(M/4)的矩阵Q_wave=reshape(x7,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵subplot(613) % 窗口分割成6*1的，当前是第3个子图 plot(t2,I_wave,'LineWidth',2);% 绘制经过中间时刻判决的I路波形axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过中间时刻判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(614) % 窗口分割成6*1的，当前是第4个子图 plot(t2,Q_wave,'LineWidth',2);% 绘制经过中间时刻判决的Q路波形axis([0,TotalT2,-5,5]); % 设置坐标范围title("Q路经过中间时刻判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% I、Q两路恢复% 电平转换成码元，-3转换成01，-1转换成00,1转换成10,3转换成11I_yuan=[];Q_yuan=[];for i=1:length(I_zong) if I_zong(i)==3 I_yuan=[I_yuan,1 1]; elseif I_zong(i)==1 I_yuan=[I_yuan,1 0]; elseif I_zong(i)==-1 I_yuan=[I_yuan,0 0]; elseif I_zong(i)==-3 I_yuan=[I_yuan,0 1]; endendfor i=1:length(Q_zong) if Q_zong(i)==3 Q_yuan=[Q_yuan,1 1]; elseif Q_zong(i)==1 Q_yuan=[Q_yuan,1 0]; elseif Q_zong(i)==-1 Q_yuan=[Q_yuan,0 0]; elseif Q_zong(i)==-3 Q_yuan=[Q_yuan,0 1]; endendcode = [];% 将I路码元为最终输出的奇数位置码元，将Q路码元为最终输出的偶数位置码元for n=1:M if mod(n, 2)~=0 code = [code, I_yuan((n+1)/2)]; else code = [code, Q_yuan(n/2)]; endendx8=wave(fz,:); % 将原来wave的第一行复制L次，称为L*M的矩阵dout=reshape(x8,1,L*M); % 产生单极性不归零矩形脉冲波形，将刚得到的L*M矩阵，按列重新排列形成1*(L*M)的矩阵subplot(615) % 窗口分割成6*1的，当前是第5个子图 plot(t,dout,'LineWidth',2); % 绘制恢复后的基带码元axis([0,TotalT,-0.1,1.1]); % 设置坐标范围title("恢复后的基带码元"); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(616) % 窗口分割成6*1的，当前是第6个子图 plot(t,jidai,'LineWidth',2);% 绘制原始基带码元axis([0,TotalT,-0.1,1.1]); % 设置坐标范围title("原始基带码元"); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 绘制频谱%% 信源频谱T=t(end); % 时间df=1/T; % 频谱分辨率N=length(jidai); % 采样长度f=(-N/2:N/2-1)*df; % 频率范围mf=fftshift(abs(fft(jidai)));%对信源信号采用快速傅里叶变换并移到矩阵中心figure(5);subplot(211); % 窗口分割成2*1的，当前是第1个子图 plot(f,mf,'LineWidth',2); % 绘制信源频谱波形title("基带信号频谱"); % 标题xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签axis([-20000,20000,-inf,inf]);% 坐标范围限制%% QAM信号频谱T2=t2(end); % 时间2df2=1/T2; % 频谱分辨率2N2=length(qam); % 采样长度2f2=(-N2/2:N2/2-1)*df2; % 频率范围2sf=fftshift(abs(fft(qam))); % 对QAM信号采用快速傅里叶变换并将0-fs频谱移动到-fs/2-fs/2subplot(212) % 窗口分割成2*1的，当前是第2个子图 plot(f2,sf,'LineWidth',2) % 绘制QAM调制信号频谱title("QAM信号频谱") % 标题xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签axis([-20000,20000,-inf,inf]);% 坐标范围限制%% 乘以相干载波后的频谱mmf=fftshift(abs(fft(jt1))); % 对相干载波信号采用快速傅里叶变换并移到矩阵中心figure(6) % 绘制第5幅图subplot(211); % 窗口分割成2*1的，当前是第1个子图 plot(f2,mmf,'LineWidth',2) % 画出I路乘以相干载波后的频谱title("I路乘以相干载波后的频谱")xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签mmf2=fftshift(abs(fft(jt2))); % 对相干载波信号采用快速傅里叶变换并移到矩阵中心subplot(212); % 窗口分割成2*1的，当前是第2个子图 plot(f2,mmf2,'LineWidth',2) % 画出Q路乘以相干载波后的频谱title("Q路乘以相干载波后的频谱")xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签%% 经过低通滤波后的频谱dmf=fftshift(abs(fft(lvbo1)));%对低通滤波信号采用快速傅里叶变换并移到矩阵中心figure(7);subplot(211); % 窗口分割成2*1的，当前是第1个子图 plot(f2,dmf,'LineWidth',2) % 画出I路经过低通滤波后的频谱title("I路经过低通滤波后的频谱");xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签dmf2=fftshift(abs(fft(lvbo2)));%对低通滤波信号采用快速傅里叶变换并移到矩阵中心subplot(212); % 窗口分割成2*1的，当前是第2个子图 plot(f2,dmf2,'LineWidth',2) % 画出Q路经过低通滤波后的频谱title("Q路经过低通滤波后的频谱");xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签4.结果图
5.特点
系统的频带利用率高，高速传输。
6.加星座图的16QAM代码
clear all; % 清除所有变量close all; % 关闭所有窗口clc; % 清屏%% 基本参数M=400; % 产生码元数 L=100; % 每码元复制L次,每个码元采样次数Ts=0.001; % 每个码元的宽度,即码元的持续时间Rb=1/Ts; % 码元速率1Kdt=Ts/L; % 采样间隔TotalT=M*Ts; % 总时间1t=0:dt:TotalT-dt; % 时间1TotalT2=(M/2)*Ts; % 总时间2t2=0:dt:TotalT2-dt; % 时间2Fs=1/dt; % 采样间隔的倒数即采样频率%% 产生单极性波形wave=randi([0,1],1,M); % 产生二进制随机码,M为码元个数fz=ones(1,L); % 定义复制的次数L,L为每码元的采样点数fz2=ones(1,2*L); % 定义复制的次数2Lx1=wave(fz,:); % 将原来wave的第一行复制L次，称为L*M的矩阵jidai=reshape(x1,1,L*M); % 产生单极性不归零矩形脉冲波形，将刚得到的L*M矩阵，按列重新排列形成1*(L*M)的矩阵%% 产生I、Q路码元a=wave(1:2:end); % I路码元是原始码元奇数位置的码元b=wave(2:2:end); % Q路码元是原始码元偶数位置的码元% 码元电平转换，要将两两码元转换成四种电平，00转换成-1,01转换成-3,10转换成1,11转换成3for i=1:length(a)/2 if a(2*(i-1)+1:2*i)==[0 1] I(i)=-3; elseif a(2*(i-1)+1:2*i)==[0 0] I(i)=-1; elseif a(2*(i-1)+1:2*i)==[1 0] I(i)=1; elseif a(2*(i-1)+1:2*i)==[1 1] I(i)=3; end if b(2*(i-1)+1:2*i)==[0 1] Q(i)=-3; elseif b(2*(i-1)+1:2*i)==[0 0] Q(i)=-1; elseif b(2*(i-1)+1:2*i)==[1 0] Q(i)=1; elseif b(2*(i-1)+1:2*i)==[1 1] Q(i)=3; endendx2=a(fz,:); % 将原来a的第一行复制L次，称为L*(M/2)的矩阵a_mayuan=reshape(x2,1,L*(M/2));% 产生单极性不归零矩形脉冲波形，将刚得到的L*(M/2)矩阵，按列重新排列形成1*(L*(M/2))的矩阵x3=b(fz,:); % 将原来b的第一行复制L次，称为L*(M/2)的矩阵b_mayuan=reshape(x3,1,L*(M/2));% 产生单极性不归零矩形脉冲波形，将刚得到的L*(M/2)矩阵，按列重新排列形成1*(L*(M/2))的矩阵x4=I(fz2,:); % 将原来I的第一行复制2L次，称为2L*(M/4)的矩阵I_mayuan=reshape(x4,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵x5=Q(fz2,:); % 将原来Q的第一行复制2L次，称为2L*(M/4)的矩阵Q_mayuan=reshape(x5,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵%% 码元波形绘制figure(1); % 绘制第1幅图subplot(511); % 窗口分割成5*1的，当前是第1个子图 plot(t,jidai,'LineWidth',2);% 绘制基带码元波形，线宽为2title('基带信号波形'); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT,-0.1,1.1]) % 坐标范围限制subplot(512); % 窗口分割成5*1的，当前是第2个子图 plot(t2,a_mayuan,'LineWidth',2);% I路码元波形，线宽为2title('I路码元波形'); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-0.1,1.1]) % 坐标范围限制subplot(513); % 窗口分割成5*1的，当前是第3个子图 plot(t2,b_mayuan,'LineWidth',2);% 绘制Q路码元波形，线宽为2title('Q路码元波形'); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-0.1,1.1]) % 坐标范围限制subplot(514); % 窗口分割成5*1的，当前是第4个子图 plot(t2,I_mayuan,'LineWidth',2);% 绘制转换后I路电平波形，线宽为2title('电平转换后I路码元波形');% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-4,4]) % 坐标范围限制subplot(515); % 窗口分割成5*1的，当前是第5个子图 plot(t2,Q_mayuan,'LineWidth',2);% 绘制转换后Q路电平波形，线宽为2title('电平转换后Q路码元波形');% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签axis([0,TotalT2,-4,4]) % 坐标范围限制%% QAM调制fc=10000; % 载波频率2kHz zb1=cos(2*pi*fc*t2); % 载波1zb2=-sin(2*pi*fc*t2); % 载波2I_lu=I_mayuan.*zb1; % I路波形 Q_lu=Q_mayuan.*zb2; % Q路波形 qam=I_lu+Q_lu; % QAM调制figure(2); % 绘制第2幅图subplot(411) % 窗口分割成4*1的，当前是第1个子图 plot(t2,I_lu,'LineWidth',2);% 绘制I路的波形 title('I路信号波形') % 标题axis([0,TotalT2,-4,4]); % 坐标范围限制xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(412) % 窗口分割成4*1的，当前是第2个子图 plot(t2,Q_lu,'LineWidth',2);% 绘制Q路的波形 title('Q路信号波形') % 标题axis([0,TotalT2,-4,4]); % 坐标范围限制xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(413) % 窗口分割成4*1的，当前是第3个子图 plot(t2,qam,'LineWidth',2); % 绘制QAM的波形 title('QAM信号波形') % 标题axis([0,TotalT2,-7,7]); % 坐标范围限制xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 信号经过高斯白噪声信道tz=awgn(qam,20); % 信号qam中加入白噪声，信噪比为SNR=20dBsubplot(414); % 窗口分割成4*1的，当前是第4个子图 plot(t2,tz,'LineWidth',2); % 绘制QAM信号加入白噪声的波形axis([0,TotalT2,-7,7]); % 坐标范围设置title('通过高斯白噪声信道后的信号');% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 解调部分jt1=tz.*zb1; % 相干解调，I路乘以相干载波jt2=tz.*zb2; % 相干解调，Q路乘以相干载波figure(3);subplot(511) % 窗口分割成5*1的，当前是第1个子图 plot(t2,jt1,'LineWidth',2) % 绘制I路乘以相干载波后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路乘以相干载波后的信号")% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(512) % 窗口分割成5*1的，当前是第2个子图 plot(t2,jt2,'LineWidth',2) % 绘制Q路乘以相干载波后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("Q路乘以相干载波后的信号")% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 加噪信号经过滤波器% 低通滤波器设计fp=2*Rb; % 低通滤波器截止频率，乘以2是因为下面要将模拟频率转换成数字频率wp=Rb/(Fs/2)h1=fir1(30, fp/Fs, boxcar(31));% 生成fir滤波器系统函数中分子多项式的系数% fir1函数三个参数分别是阶数，数字截止频率，滤波器类型% 这里是生成了30阶(31个抽头系数)的矩形窗滤波器[h,w]=freqz(h1, 1,512); % 生成fir滤波器的频率响应% freqz函数的三个参数分别是滤波器系统函数的分子多项式的系数，分母多项式的系数(fir滤波器分母系数为1)和采样点数(默认)512lvbo1=2*fftfilt(h1,jt1); % 对信号进行滤波，jt1是等待滤波的信号，h是fir滤波器的系统函数的分子多项式系数lvbo2=2*fftfilt(h1,jt2); % 对信号进行滤波，jt2是等待滤波的信号，h是fir滤波器的系统函数的分子多项式系数subplot(513); % 窗口分割成5*1的，当前是第3个子图 plot(w/pi*Fs/2,20*log(abs(h)),'LineWidth',2); % 绘制滤波器的幅频响应title('低通滤波器的频谱'); % 标题xlabel('频率/Hz'); % x轴标签ylabel('幅度/dB'); % y轴标签subplot(514) % 窗口分割成5*1的，当前是第4个子图 plot(t2,lvbo1,'LineWidth',2);% 绘制I路经过低通滤波器后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过低通滤波器后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(515) % 窗口分割成5*1的，当前是第5个子图 plot(t2,lvbo2,'LineWidth',2);% 绘制Q路经过低通滤波器后的信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过低通滤波器后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 抽样判决% 滤波后的波形进行判决，大于2的判决为3,0-2之间的判决为1，-2-0之间判决为-1，小于-2的判决为-3for j=1:length(lvbo1) if lvbo1(j)>=2 I_panjue(j)=3; elseif (lvbo1(j)>0 && lvbo1(j)<2) I_panjue(j)=1; elseif (lvbo1(j)>=-2 && lvbo1(j)<0) I_panjue(j)=-1; elseif lvbo1(j)<-2 I_panjue(j)=-3; endendfor k=1:length(lvbo2) if lvbo2(k)>=2 Q_panjue(k)=3; elseif (lvbo2(k)>0 && lvbo2(k)<2) Q_panjue(k)=1; elseif (lvbo2(k)>=-2 && lvbo2(k)<0) Q_panjue(k)=-1; elseif lvbo2(k)<-2 Q_panjue(k)=-3; endendfigure(4);subplot(611) % 窗口分割成6*1的，当前是第1个子图 plot(t2,I_panjue,'LineWidth',2);% 绘制经过判决的I路信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(612) % 窗口分割成6*1的，当前是第2个子图 plot(t2,Q_panjue,'LineWidth',2);% 绘制经过判决的Q路信号axis([0,TotalT2,-5,5]); % 设置坐标范围title("Q路经过判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 判决之后再处理，取码元中间时刻判决I_zong=[];Q_zong=[];for m=L:2*L:((2*L)*(M/4)) if I_panjue(m)>=2 I_zong=[I_zong,3]; elseif (I_panjue(m)>0 && I_panjue(m)<2) I_zong=[I_zong,1]; elseif (I_panjue(m)>=-2 && I_panjue(m)<0) I_zong=[I_zong,-1]; elseif I_panjue(m)<-2 I_zong=[I_zong,-3]; endendfor n=L:2*L:((2*L)*(M/4)) if Q_panjue(n)>=2 Q_zong=[Q_zong,3]; elseif (Q_panjue(n)>0 && Q_panjue(n)<2) Q_zong=[Q_zong,1]; elseif (Q_panjue(n)>=-2 && Q_panjue(n)<0) Q_zong=[Q_zong,-1]; elseif Q_panjue(n)<-2 Q_zong=[Q_zong,-3]; endendx6=I_zong(fz2,:); % 将原来I_zong的第一行复制2L次，称为2L*(M/4)的矩阵I_wave=reshape(x6,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵x7=Q_zong(fz2,:); % 将原来Q_zong的第一行复制2L次，称为2L*(M/4)的矩阵Q_wave=reshape(x7,1,(2*L)*(M/4));% 产生单极性不归零矩形脉冲波形，将刚得到的(2*L)*(M/4)矩阵，按列重新排列形成1*((2*L)*(M/4))的矩阵subplot(613) % 窗口分割成6*1的，当前是第3个子图 plot(t2,I_wave,'LineWidth',2);% 绘制经过中间时刻判决的I路波形axis([0,TotalT2,-5,5]); % 设置坐标范围title("I路经过中间时刻判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(614) % 窗口分割成6*1的，当前是第4个子图 plot(t2,Q_wave,'LineWidth',2);% 绘制经过中间时刻判决的Q路波形axis([0,TotalT2,-5,5]); % 设置坐标范围title("Q路经过中间时刻判决后的信号");% 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% I、Q两路恢复% 电平转换成码元，-3转换成01，-1转换成00,1转换成10,3转换成11I_yuan=[];Q_yuan=[];for i=1:length(I_zong) if I_zong(i)==3 I_yuan=[I_yuan,1 1]; elseif I_zong(i)==1 I_yuan=[I_yuan,1 0]; elseif I_zong(i)==-1 I_yuan=[I_yuan,0 0]; elseif I_zong(i)==-3 I_yuan=[I_yuan,0 1]; endendfor i=1:length(Q_zong) if Q_zong(i)==3 Q_yuan=[Q_yuan,1 1]; elseif Q_zong(i)==1 Q_yuan=[Q_yuan,1 0]; elseif Q_zong(i)==-1 Q_yuan=[Q_yuan,0 0]; elseif Q_zong(i)==-3 Q_yuan=[Q_yuan,0 1]; endendcode = [];% 将I路码元为最终输出的奇数位置码元，将Q路码元为最终输出的偶数位置码元for n=1:M if mod(n, 2)~=0 code = [code, I_yuan((n+1)/2)]; else code = [code, Q_yuan(n/2)]; endendx8=code(fz,:); % 将原来code的第一行复制L次，称为L*M的矩阵dout=reshape(x8,1,L*M); % 产生单极性不归零矩形脉冲波形，将刚得到的L*M矩阵，按列重新排列形成1*(L*M)的矩阵subplot(615) % 窗口分割成6*1的，当前是第5个子图 plot(t,dout,'LineWidth',2); % 绘制恢复后的基带码元axis([0,TotalT,-0.1,1.1]); % 设置坐标范围title("恢复后的基带码元"); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签subplot(616) % 窗口分割成6*1的，当前是第6个子图 plot(t,jidai,'LineWidth',2);% 绘制原始基带码元axis([0,TotalT,-0.1,1.1]); % 设置坐标范围title("原始基带码元"); % 标题xlabel('时间/s'); % x轴标签ylabel('幅度'); % y轴标签%% 绘制频谱%% 信源频谱T=t(end); % 时间df=1/T; % 频谱分辨率N=length(jidai); % 采样长度f=(-N/2:N/2-1)*df; % 频率范围mf=fftshift(abs(fft(jidai)));%对信源信号采用快速傅里叶变换并移到矩阵中心figure(5);subplot(211); % 窗口分割成2*1的，当前是第1个子图 plot(f,mf,'LineWidth',2); % 绘制信源频谱波形title("基带信号频谱"); % 标题xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签axis([-20000,20000,-inf,inf]);% 坐标范围限制%% QAM信号频谱T2=t2(end); % 时间2df2=1/T2; % 频谱分辨率2N2=length(qam); % 采样长度2f2=(-N2/2:N2/2-1)*df2; % 频率范围2sf=fftshift(abs(fft(qam))); % 对QAM信号采用快速傅里叶变换并将0-fs频谱移动到-fs/2-fs/2subplot(212) % 窗口分割成2*1的，当前是第2个子图 plot(f2,sf,'LineWidth',2) % 绘制QAM调制信号频谱title("QAM信号频谱") % 标题xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签axis([-20000,20000,-inf,inf]);% 坐标范围限制%% 乘以相干载波后的频谱mmf=fftshift(abs(fft(jt1))); % 对相干载波信号采用快速傅里叶变换并移到矩阵中心figure(6) % 绘制第5幅图subplot(211); % 窗口分割成2*1的，当前是第1个子图 plot(f2,mmf,'LineWidth',2) % 画出I路乘以相干载波后的频谱title("I路乘以相干载波后的频谱")xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签mmf2=fftshift(abs(fft(jt2))); % 对相干载波信号采用快速傅里叶变换并移到矩阵中心subplot(212); % 窗口分割成2*1的，当前是第2个子图 plot(f2,mmf2,'LineWidth',2) % 画出Q路乘以相干载波后的频谱title("Q路乘以相干载波后的频谱")xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签%% 经过低通滤波后的频谱dmf=fftshift(abs(fft(lvbo1)));%对低通滤波信号采用快速傅里叶变换并移到矩阵中心figure(7);subplot(211); % 窗口分割成2*1的，当前是第1个子图 plot(f2,dmf,'LineWidth',2) % 画出I路经过低通滤波后的频谱title("I路经过低通滤波后的频谱");xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签dmf2=fftshift(abs(fft(lvbo2)));%对低通滤波信号采用快速傅里叶变换并移到矩阵中心subplot(212); % 窗口分割成2*1的，当前是第2个子图 plot(f2,dmf2,'LineWidth',2) % 画出Q路经过低通滤波后的频谱title("Q路经过低通滤波后的频谱");xlabel('频率/Hz'); % x轴标签ylabel('幅度'); % y轴标签figure(8); % 绘制第8幅图scatter(I,Q); % 绘制星座图，如果没有16个点，是因为码元数量不够的原因title('QAM星座图'); % 标题% 让坐标轴放在坐标原点ax = gca;ax.XAxisLocation = 'origin';ax.YAxisLocation = 'origin';
