简介这份资源基于MATLAB环境实现16QAM数字通信系统仿真面向通信工程专业学生、科研人员及对调制解调技术感兴趣的开发者重点解决如何用代码完整模拟数字通信链路的问题。资源共包含13个文件全部为M脚本压缩包大小仅5KB但功能覆盖从基带数据生成、16QAM符号映射、上变频到射频、叠加高斯白噪声、下变频恢复基带再到最大似然解调与星座图、误码率曲线绘制等完整流程各功能模块划分清晰便于单独调用和学习。目前已有606人学习使用在同类仿真资源中积累了一定认可度。通过本项目读者既能理解16QAM高阶调制的原理也能掌握噪声对系统性能的影响以及上下变频在无线通信中的作用星座图与误码率曲线脚本可帮助快速评估不同信噪比下的系统性能非常适合课堂教学、自主学习、课程设计与毕业设计参考。此外代码采用模块化结构注释清晰便于二次开发与算法改进。1. 16QAM 仿真不止是 qammod从比特到星座再到中频链路的完整闭环很多人拿到 16QAM 项目会把qammod的输出直接丢进awgn然后画一条 BER 曲线就算交差。这种做法只能证明 MATLAB 自带函数没问题证明不了你对数字通信系统的理解。真正的链路里调制后的符号要插值、成形、上变频到中频接收端再下变频、匹配滤波、做符号同步最容易被忽略的是各级滤波器累计的群延迟它会让解调输出的误码率怎么也收敛不到理论上限。这套仿真的价值在于把整个过程拆成了命名清晰的模块从Qam_modulation.m的比特映射经过rise_cos.m成形、modulate_to_high.m上变频进入generate_noise.m加噪最后用pick_sig.m对齐采样点后交给demodulate_sig.m判决。如果你正在做通信课设或想搞明白 16QAM 各模块接口怎么衔接这套文件是可以直接对着调参的。2. 16QAM 调制映射与格雷编码星座图生成2.1 4 bit 一符号与方形星座的对齐16QAM 的本质是把两个独立的 4 电平信号分别放到 I 路和 Q 路组合出 16 个星座点每个点携带 4 bit 信息。I、Q 两路的幅度通常取-3, -1, 1, 3这样星座图是一个规整的 4×4 方阵。选择格雷编码的目的是让相邻星座点只差 1 个 bit这样在高斯白噪声下发生错误判决时绝大多数误码只落在相邻点比特误码率不会出现突发性跳变。映射前还要做能量归一化。4×4 方阵所有点的平均功率是(4*(3^23^2)8*(3^21^2)4*(1^21^2))/1610所以每个符号要除以sqrt(10)保证平均符号能量为 1。这个归一化直接决定了后面generate_noise.m里噪声方差怎么算如果不除BER 曲线会整体右移几个 dB。2.2 Qam_modulation.m 的格雷映射实现项目里Qam_modulation.m负责把比特流变成复数符号。我一般这样写function sym Qam_modulation(bits) % 输入 bits 是列向量长度必须是 4 的整数倍 bits bits(:); if mod(numel(bits), 4) ~ 0 error(比特数必须是4的整数倍); end % 每 4 bit 一组b1 b2 映射到 I 路b3 b4 映射到 Q 路 b reshape(bits, 4, []).; N size(b, 1); I zeros(N, 1); Q zeros(N, 1); % 格雷码顺序00 01 11 10 amp [-3 -1 1 3]; gray [0 1 3 2]; for n 1:N Ic b(n,1)*2 b(n,2); % 普通二进制序 Qc b(n,3)*2 b(n,4); I(n) amp(gray Ic); % 映射成格雷序对应的幅度 Q(n) amp(gray Qc); end % 归一化16 个星座点平均功率为 10除以 sqrt(10) 后为 1 sym I 1j*Q; sym sym / sqrt(mean(abs(sym).^2)); end代码里gray Ic是逻辑索引因为gray中只有一个元素等于Ic所以amp(gray Ic)取到的是对应格雷序位置的幅度。Ic b1*2b2把两个 bit 变成十进制数比如10是十进制 2而格雷序中第 3 个位置是10对应幅度1这样 16QAM 象限之间的转换不会出现两位同时翻转。2.3 星座图和 I/Q 映射表plot_astrology.m做的事情很简单就是画星座散点。配合归一化后的符号直接用plot打点figure; plot(real(sym), imag(sym), .); grid on; axis equal; xlabel(I); ylabel(Q); title(16QAM 星座图);I/Q 两路的格雷映射关系可以总结成一张表后面调试解调代码时对照用b1 b2 或 b3 b4十进制值格雷序位置对应幅度0001-30112-111331102433. 上下变频与升余弦成形仿真链路里的混频器和滤波器3.1 为什么要让仿真保留中频链路很多仿真省略上下变频直接把基带符号加噪声这样做出的 BER 曲线虽然干净但无法体现无线信道里最重要的频率搬移。modulate_to_high.m存在的意义就是把基带复信号搬到一个数字中频上。数字中频频率fc是相对于采样率的归一化值比如fc0.1表示中频是采样率的 1/10。混频输出的实信号既包含了正频部分也包含了镜像分量接收端下变频后如果不做滤波镜像会直接叠加到基带导致星座图出现一层“雾”。实际仿真中还有一个容易被忽略的点上变频前的基带信号必须做过采样和成形否则频谱副本会混叠。所以完整的发送链路是符号映射 →insert_value插零 →rise_cos成形 →modulate_to_high上变频。3.2 rise_cos.m 生成成形脉冲升余弦滤波器的滚降系数beta决定了带宽和时域旁瓣。beta越大频谱越宽波形衰减越快beta0就是理想低通物理上不可实现。仿真里一般取0.25~0.5。我常用的生成函数function h rise_cos(beta, sps, span) % beta: 滚降系数 0~1 % sps: 每个符号采样点数 % span: 滤波器涵盖的符号数 N span * sps; if mod(N, 2) 0 N N 1; % 保证脉冲长度是奇数中心点明确 end t linspace(-span/2, span/2, N); % 分母加 1e-12 避免在 1/(2*beta) 处出现 0/0 h sinc(t) .* cos(pi * beta * t) ./ (1 - (2 * beta * t).^2 1e-12); h h / sum(h); % 归一化保证直流增益为 1 endt的单位是符号周期不是采样点。sps一般取 8 或 16太低会让频谱混叠太高会拖慢仿真。滤波器的长度是span*sps1群延迟也就基本锁定在这个量级。参数选择可以参考这张表参数取值对链路的影响beta0.25~0.5越大占用带宽越大时域旁瓣衰减越快sps8~16太小时频谱混叠太大时仿真时间成倍增长span6~12越大滤波器越接近理想但群延迟和计算量变大3.3 insert_value 和 modulate_to_high插零与数字上变频insert_value.m做的是零值插值也叫 upsampling。它的作用是给每个符号之间插入sps-1个零让符号序列变成采样率下的波形序列。之后与rise_cos卷积就完成了脉冲成形function y insert_value(x, sps) % x 是符号序列sps 是过采样倍数 y zeros(length(x) * sps, 1); y(1:sps:end) x; endmodulate_to_high.m负责把复基带信号变成实中频信号。常见做法是用复指数exp(j*2*pi*fc*n)与基带信号相乘再取实部function y modulate_to_high(x, fc) % x 是基带复信号fc 是归一化数字中频频率必须小于 0.5 N length(x); n (0:N-1).; lo exp(1j * 2 * pi * fc * n); y real(x .* lo); end取实部后信号的频谱变成关于fc对称的两个边带。接收端下变频用共轭本振相乘可以得到基带和一个2*fc的高频分量后续必须用匹配滤波器把这个高频分量滤掉。3.4 下变频、匹配滤波和 pick_sig 的群延迟对齐接收端最容易出错的一步是信号对齐。发送端和接收端各做一次卷积后总延迟不是 0 而是length(h)个采样点。如果用conv(symUp, h, same)这种截断方式MATLAB 会隐式地把波形中心对齐掩盖掉真实延迟反而把因果性搞乱。我这里用完整卷积并把延迟显式算出来function sym_rx pick_sig(matched, nSym, sps, delay) % matched 是下变频并匹配滤波后的完整卷积输出 % delay 取 length(h)即第一个符号的 1-based 起始位置 sym_rx matched(delay : sps : delay (nSym-1)*sps); end在完整链路中先用conv(symUp, h)得到发送基带再在匹配滤波端用conv(basebandRx, h)。两次卷积后第一个符号的峰值出现在索引length(h)处。下面这段可以做一致性检查test_sym [1; zeros(16,1)]; test_up insert_value(test_sym, 8); test_base conv(test_up, h); test_lo exp(1j*2*pi*fc*(0:length(test_base)-1).); test_if real(test_base .* test_lo); % 直接下变频并匹配滤波观察峰值位置 test_down test_if .* conj(test_lo); test_match conv(test_down, h); [~, delay] max(abs(test_match)); % 这里 delay 应当等于 length(h)4. 高斯白噪声信道下的解调判决与 BER 统计4.1 generate_noise.m 的 SNR 口径generate_noise.m要解决的核心问题不是“怎么生成随机数”而是噪声功率的标定。16QAM 一个符号带 4 bit所以Es/N0和Eb/N0之间差10*log10(4)≈6.02dB很多人直接把EbN0_dB传给噪声函数结果 BER 曲线比理论值偏好几 dB。我在实中频信号上加噪声时按下面的方式换算符号能量function rx generate_noise(tx_if, EsN0_dB, sps) % tx_if 是实中频信号 % EsN0_dB EbN0_dB 10*log10(4) EsN0 10^(EsN0_dB/10); Ps mean(tx_if.^2); % 实信号平均功率 Es 2 * sps * Ps; % 折算到单个符号能量 N0 Es / EsN0; sigma sqrt(N0 / 2); % 实噪声每采样点标准差 rx tx_if sigma * randn(size(tx_if)); end这里的Es折算假设中频信号是平稳的并且不含直流分量。如果你直接在基带符号级加噪声可以省去sps和功率折算直接让噪声方差等于1/(2*EsN0)。两种方式等价的前提是接收端的匹配滤波器对信号和噪声的增益关系一致所以不要把两种口径混着用。4.2 demodulate_sig.m 的最近邻阈值判决16QAM 的硬判决可以按 I/Q 两路分别做阈值再反查格雷码。demodulate_sig.m的典型实现是function bits demodulate_sig(rx_sym) % rx_sym 是复数符号序列 amp2bits [0 0; 0 1; 1 1; 1 0]; % 对应 -3 -1 1 3 I real(rx_sym); Q imag(rx_sym); % 阈值判决 I_amp zeros(size(I)); I_amp(I -2) -3; I_amp(I -2 I 0) -1; I_amp(I 0 I 2) 1; I_amp(I 2) 3; Q_amp zeros(size(Q)); Q_amp(Q -2) -3; Q_amp(Q -2 Q 0) -1; Q_amp(Q 0 Q 2) 1; Q_amp(Q 2) 3; % 幅度反查 bit amps [-3 -1 1 3]; idxI zeros(size(I)); idxQ zeros(size(Q)); for a 1:4 idxI(I_amp amps(a)) a; idxQ(Q_amp amps(a)) a; end bits zeros(length(rx_sym)*4, 1); for n 1:length(rx_sym) bits((n-1)*4 (1:2)) amp2bits(idxI(n), :); bits((n-1)*4 (3:4)) amp2bits(idxQ(n), :); end end阈值取-2, 0, 2而不是星座点中间值-2, 0, 2刚刚好。注意这里用的是最大似然意义上的最近邻判决因为噪声是高斯白噪声时欧氏距离最小等价于概率最大。4.3 project.m 主循环把模块串成完整链路把前面这些模块串起来主仿真脚本大致长这样clear; clc; rng(7); % 基带参数 k 4; % bits/symbol nSym 20000; % 每个信噪比下独立符号数 sps 8; % 每符号采样点 beta 0.35; % 滚降系数 span 6; % 滤波器跨度 fc 0.1; % 归一化中频频率 EbN0_dB 0:2:14; ber zeros(size(EbN0_dB)); % 发射端 bitsTx randi([0 1], k*nSym, 1); symTx Qam_modulation(bitsTx); symUp insert_value(symTx, sps); h rise_cos(beta, sps, span); basebandTx conv(symUp, h); ifTx modulate_to_high(basebandTx, fc); % 循环仿真 for m 1:length(EbN0_dB) EsN0_dB EbN0_dB(m) 10*log10(k); rxIf generate_noise(ifTx, EsN0_dB, sps); n (0:length(rxIf)-1).; basebandRx rxIf .* exp(-1j*2*pi*fc*n); % 下变频 matched conv(basebandRx, h); % 匹配滤波 % 取符号delay length(h) symRx pick_sig(matched, nSym, sps, length(h)); bitsRx demodulate_sig(symRx); ber(m) sum(bitsRx ~ bitsTx) / length(bitsTx); end % 画 BER 曲线 figure; semilogy(EbN0_dB, ber, o-); grid on; xlabel(Eb/N0 (dB)); ylabel(BER);这段代码保留了完整的上下变频和滤波过程。basebandTx用的是完整卷积所以接收端matched也是完整卷积最后用pick_sig从第length(h)个点开始按sps抽样。由于这里用的是单个升余弦滤波器做发送成形和接收匹配滤波频谱是平方关系仿真 BER 会比理论曲线略高这是滤波器选型导致的不是代码 bug。要贴近理论应该把发送和接收都换成根升余弦滤波器。4.4 与理论 BER 对比时注意蒙特卡洛点数画 BER 曲线时理论上限可以直接用berawgnberTheory berawgn(EbN0_DB, qam, 16);仿真点在低信噪比下误码率低需要统计足够的错误比特数。上面的脚本每个信噪比固定跑 20000 个符号也就是 80000 个 bit在Eb/N014dB时可能只出现几个误码曲线会抖动。至少要让累计误码达到 300~500 个再进入下一个信噪比比如改成while 循环误码数不足就增加符号。5. 用星座图和 EVM 定位仿真链路的隐性故障5.1 星座图里的旋转、扩散与偏移plot_astrology.m不只是画点它是最快的排错手段。如果星座图整体顺时针或逆时针旋转说明下变频本振相位不连续或存在残余频偏如果星座点从四个“角”向中心扩散通常是噪声功率标定过大或者匹配滤波器的归一化有问题如果 16 个点错乱地对齐到错误象限则要检查Qam_modulation.m里的gray映射顺序和demodulate_sig.m里的amp2bits是否一一对应。5.2 用 EVM 量化星座质量肉眼判断星座图不够精确可以计算 RMS EVM。EVM 反映的是实际接收符号相对理想星座点的误差矢量大小function evm cal_evm(rx_sym, ref_sym) err rx_sym(:) - ref_sym(:); evm sqrt(mean(abs(err).^2) / mean(abs(ref_sym).^2)) * 100; end参考符号可以用和解调相同的映射表从发送端直接拿也可以用最近邻判决后的星座点。EVM 和系统质量的关系大致如下RMS EVM链路状态 5%信号质量很好BER 接近理论5%~15%轻微噪声或滤波器失配可接受15%~30%明显劣化检查频偏、延迟和噪声口径 30%基本不可用多半是映射表或下变频错误5.3 频偏和延迟不对齐的定位如果星座图旋转但 BER 又不太差说明频偏很小可以用符号间相位差估计残余频偏% 对于重复发送的导频符号按符号周期估计相位增量 phase_diff angle(rx_sym(2:end) .* conj(rx_sym(1:end-1))); freq_offset mean(phase_diff) / (2 * pi * sps);算出的freq_offset对应归一化频偏把它补偿到modulate_to_high的本振频率上即可。延迟不对齐则表现为 BER 不随 SNR 下降或者在某个高信噪比处突然跌不下去。这时不要硬调采样点回到test_sym校准脚本用单脉冲法重新确认pick_sig的起始索引。匹配滤波后的峰值位置如果不在length(h)就把该索引直接作为 delay 传入确保符号采样点落在眼图张开最大处。本文还有配套的精品资源点击获取
