★ FEATURED ARTICLE 简介围绕16QAM星座图与误码率仿真整理的MATLAB学习素材面向正在学习数字调制、通信系统仿真与BER性能评估的本科生、研究生及工程入门者可用于理解幅度相位联合调制、AWGN信道下端到端链路搭建与解调判决流程。压缩包共8个文件均为.m脚本与函数文件整体约3KB涵盖主流程调用、星座图绘制、四进制与二进制转换、调制解调、随机二进制序列生成及比特成形等模块便于按功能拆解阅读和二次修改。已有1927人学习下载说明其在通信仿真入门场景中具有一定参考价值。读者可据此复现16QAM信号映射、星座点分布观察、接收端最大似然检测以及不同信噪比下的误码率曲线并借助模块化脚本梳理从数据源到性能评估的完整仿真思路适合课程实验、自学练习和代码对照。1. 星座图仿真不是画点16QAM 的误码率要从映射、加噪、判决三处对账16QAM 星座图仿真和误码率仿真在 MATLAB 里十分钟就能跑出一条曲线但曲线和理论差 3 dB 的时候很多人就卡住了——不知道问题出在映射表、噪声方差还是判决门限。这份仿真资源想解决的就是这件事从比特到复数符号的格雷映射、平均符号能量归一化、AWGN 加噪、最近邻判决、蒙特卡洛统计一条链路跑通最后能和 erfc 理论曲线对上。适合正在做数字通信课程设计、刚进链路级仿真门槛的开发者也适合要给自家接收机算法做误码性能基线的人。下面按「映射表怎么建 → 星座图怎么画 → 误码率怎么统计 → 哪里容易翻车 → 怎么扩到衰落信道」的顺序拆开。2. 16QAM 映射表怎么建格雷码、归一化与 MATLAB 查表实现这一章是整个仿真的地基。映射表建错后面画出来的星座图看着也是 16 个点、看着也挺对称但误码率曲线会莫名其妙比理论好或者比理论差而且这种错很难靠肉眼发现。血泪经验是先把映射和归一化写死成一张常量表再往下写别的。2.1 为什么 16QAM 必须配格雷码16QAM 每个符号携带 4 个比特16 个星座点排成 4×4 的方阵。符号在传输中被打错是常态关键在于一个符号错了会连带错几个比特。如果采用自然二进制顺序排列相邻星座点之间可能相差 2 个甚至 3 个比特一个符号判决滑到隔壁就一次性错 23 bit采用格雷码排列任意两个相邻星座点的比特表示只差 1 位符号错一次最多错 1 bit。这个差别在仿真里体现得非常直接。不加格雷映射时同样的 Es/N0 下 BER 会比格雷映射高出一截高信噪比区间尤其明显——因为这时候错误主要来自「滑到最近邻点」而不是「滑到远处」。有人拿自然映射跑出来的曲线去和教科书上的公式比发现怎么也对不上原因就在这里教科书的公式默认格雷映射。我的做法是把「2 bit → 4 个电平」的映射单独抽出来做成常量表。I 路 2 bit 对应 4 个电平Q 路 2 bit 对应 4 个电平两路正交叠加就是 16 个点。2 bit 格雷码的顺序是00 → 01 → 11 → 10对应电平-3, -1, 1, 3这样相邻电平-3 和 -1、-1 和 1、1 和 3对应的标号都只差一位。2.2 归一化把平均符号能量压到 1电平取±1, ±3的时候每个符号的平均能量是(1199)/4 5单维I/Q 两维加起来是 10。如果不做归一化平均符号能量 Es 10你按 Es 1 去算噪声方差实际加进去的噪声就偏大 10 dB曲线整体左移看起来「仿真比理论好 10 dB」这是个很典型的假象。标准做法是整张表除以sqrt(10)让 Es 归一化到 1。这样后面写噪声方差的时候可以直接用N0 1/(Es/N0 线性值)不用每次去乘系数减少一处出错点。提示归一化只对电平表做一次不要在调制函数里每次调用都现算sqrt(10)也不要在解调端用另一套电平值。一旦调制和解调的电平表不一致判决门限就会整体偏BER 曲线会在高信噪比处卡在一个错误平台上不往下掉。2.3 MATLAB 查表实现4 比特拆 I/Q 两路下面是调制与加噪的两个核心函数。整个实现只有十几行但没有循环、没有bi2de靠的是向量化索引跑几百万比特也在秒级。function s qam16_mod(bits) % 16QAM 格雷映射调制每 4 bit 一组前 2 bit 定 I 路后 2 bit 定 Q 路 % 输入 bits 长度必须是 4 的整数倍输出为列向量复符号平均能量 Es 1 bits bits(:); tab [-3 -1 3 1] / sqrt(10); % 格雷码电平表索引 1..4 对应标号 00,01,10,11 idxI bits(1:4:end)*2 bits(2:4:end); % I 路 2bit 标号 0..3 idxQ bits(3:4:end)*2 bits(4:4:end); % Q 路 2bit 标号 0..3 s tab(idxI 1) 1j * tab(idxQ 1); % MATLAB 索引从 1 开始所以要 1 end function r awgn_esn0(s, esn0_dB) % 按 Es/N0 加复高斯白噪声前提是 s 的平均符号能量已经归一化到 1 esn0 10^(esn0_dB/10); n0 1 / esn0; % Es 1 时N0 1/(Es/N0) n sqrt(n0/2) * (randn(size(s)) 1j*randn(size(s))); % 每维方差 N0/2 r s n; endtab的顺序是[-3 -1 3 1]而不是[-3 -1 1 3]这是整段代码里最容易写错的地方。索引 1、2、3、4 分别对应标号00、01、10、11按格雷码规则标号10应该落在3、标号11落在1所以表里第 3、4 项是反的。写成[-3 -1 1 3]程序照样跑星座图照样好看但误码率会比理论差一截——这就是那种「跑得通但对不上」的典型故障。awgn_esn0里n0 1/esn0成立的前提就是上一节的归一化。每维方差取n0/2是因为复基带噪声的总功率是n0实部虚部各分一半用randn生成实数再拼成复数正好符合这个分布。qam16_mod的输入要求是 4 的倍数这个约束要在调用方保证。常见做法是randi([0 1], 4*N, 1)直接生成N 是符号数。如果比特源来自外部文件或者编码器输出长度不整除 4 的时候需要补零补了多少位要在误码统计时相应扣掉否则统计出来的 BER 会偏低。3. 星座图可视化把加噪前后的点画在一张图上星座图这部分的验收标准不是「好看」是「能不能一眼看出判决边界和噪声扩散的关系」。画对了你能从图上直接读出为什么某个信噪比下开始大量误判画错了就是一堆糊在一起的散点。3.1 Eb/N0 到 Es/N0那 6.02 dB 别省论文和教材上的误码率曲线横轴几乎都是 Eb/N0每比特能量噪声比而加噪的时候需要的是 Es/N0每符号能量噪声比。两者的关系是Es/N0 Eb/N0 10*log10(k)16QAM 的 k 4也就是要加 6.02 dB。这一步漏掉的话仿真曲线相对理论会整体偏移 6 dB形状一模一样位置全错。有个取巧但很稳的习惯在脚本里统一用 Eb/N0 作为外层循环变量进加噪函数之前转成 Es/N0函数内部只认 Es/N0永远不在加噪函数里做比特能量换算。这样「换算」这件事在整条链路上只出现一次出错面就收窄了。k 4; % 16QAM 每符号比特数 ebn0_dB 12; esn0_dB ebn0_dB 10*log10(k); % 12 6.02 18.02 dB r awgn_esn0(qam16_mod(randi([0 1], 4000, 1)), esn0_dB);3.2 四档信噪比的散点图取 10、16、22、30 dB 四档2×2 排开画能很直观地看到「噪声把点吹散 → 点云开始越过判决线 → 点云基本收拢」的过程。30 dB 那一档基本就是干净星座图用来确认映射表是对的10 dB 那一档基本糊成一团用来确认仿真没有算错噪声。figure(Position, [100 100 900 900]); esn0List [10 16 22 30]; for ii 1:4 subplot(2,2,ii); b randi([0 1], 4*2000, 1); s qam16_mod(b); r awgn_esn0(s, esn0List(ii)); plot(real(s), imag(s), o, MarkerSize, 5, LineWidth, 1); % 理想星座点 hold on; plot(real(r), imag(r), ., MarkerSize, 6); % 加噪后的接收点 grid on; axis equal; axis([-1.6 1.6 -1.6 1.6]); xlabel(同相分量 I); ylabel(正交分量 Q); title(sprintf(Es/N0 %d dB, esn0List(ii))); hold off; endaxis equal不能省否则 I/Q 两轴比例不一致看起来像矩形星座会误导判断。坐标范围取±1.6是因为归一化后最外圈点的坐标是±3/√10 ≈ ±0.949留出约 0.65 的余量正好能看清噪声尾巴。3.3 判决域与门限线16QAM 的最优判决AWGN 下就是 I、Q 两路各自独立做最近邻门限落在-2/√10、0、2/√10三条竖线Q 路同理是三条横线。把这几条线画到散点图上噪声点和判决边界的相对位置一目了然。th [-2 0 2] / sqrt(10); for v th xline(v, --, Color, [0.4 0.4 0.4]); yline(v, --, Color, [0.4 0.4 0.4]); end门限值取 ±2 而不是 ±1 或 ±3是因为相邻电平和为 -4、0、4中点就是 ±2。这一点在解调函数里必须和调制表严格对应很多「BER 曲线在高信噪比下有错误平台」的问题根子就在门限算成了±1/√10或者±3/√10。如果要在图上标格雷码标号方便对照可以用text在每点旁边写上 4 位二进制串注意标号顺序仍然是00, 01, 10, 11对应电平-3, -1, 3, 1别按大小顺序标。4. 误码率仿真蒙特卡洛循环的帧结构与误差控制星座图只能定性误码率曲线才是定量结论。这一章的重点是循环怎么写帧长、停止条件、统计口径三样东西定错一样跑一晚上出来的数据就白费。4.1 帧长怎么定20000 符号是个折中误码率仿真的基本结构是「生成比特 → 调制 → 加噪 → 解调 → 比对」。有两种循环写法一是外层固定信噪比、内层一直跑二是把整条链路的比特数固定住。我一般用前者内层按帧跑每帧 20000 个符号8 万比特。20000 这个数字不是玄学是内存和统计独立性的折中。一帧 8 万比特双精度复数符号占 20000×16 字节 ≈ 320 KB加上中间变量一帧不到 2 MB跑一万帧也不会撑爆内存同时每帧重新randi生成比特、重新randn生成噪声天然保证了帧与帧之间统计独立。帧长如果取得太小比如 100 符号循环次数暴涨MATLAB 的函数调用开销占比过高取得太大比如 1000 万符号单帧内存几十 MB而且一旦中途出错前面全白跑。20000 这个量级比较稳。4.2 错误数达标比固定比特数更省时间高信噪比区间的 BER 可能低到 1e-5 甚至 1e-6。如果你固定跑 100 万比特那在 1e-5 处只期望错 10 个统计涨落能有 30% 以上曲线会抖得像心电图跑到 1e-6 就更夸张可能一个错都没有直接断线。我的处理是按错误数停累计错误达到 100 个就停同时设一个总比特上限比如 2000 万防止低信噪比时跑不完。100 个错误对应的相对标准差约 10%画到图上曲线比较平滑。targetErr 100; % 目标错误数 maxBit 2e7; % 总比特上限防低 SNR 死循环 nErr 0; nBit 0; while nErr targetErr nBit maxBit b randi([0 1], 4*20000, 1); % 一帧 20000 符号 s qam16_mod(b); r awgn_esn0(s, esn0_dB); bh qam16_demod(r); nErr nErr sum(bh ~ b); nBit nBit numel(b); end ber nErr / nBit;maxBit要设是因为 02 dB 这一段 BER 在 0.1 以上100 个错误几百比特就到了但如果你把曲线往外推到负信噪比做演示没有上限的话循环可能跑很久。这两个参数按需调整目标错误数一般取 100200往下降到 50 曲线就开始抖了。4.3 与 erfc 理论曲线对齐完整链路矩形 16QAM 在 AWGN 下的精确符号错误率是Ps 3Q(x) - 2.25Q²(x)其中x sqrt(0.8·Eb/N0 线性值)Q(x) 0.5·erfc(x/√2)。格雷映射下比特错误率近似为Ps/4。这个公式在中等以上信噪比区间和仿真的吻合度很好可以直接拿来当验收尺子。function bits qam16_demod(r) % 16QAM 格雷映射解调I/Q 各自最近邻判决门限 ±2/sqrt(10) 和 0 tab [-3 -1 3 1] / sqrt(10); th [-2 0 2] / sqrt(10); invGray [0 1 3 2]; % 排序位置 - 格雷标号 idxI invGray(sum(real(r) th, 2) 1); % 0..3 idxQ invGray(sum(imag(r) th, 2) 1); n numel(r); bits zeros(4*n, 1); bits(1:4:end) floor(idxI/2); % I 路高位 bits(2:4:end) mod(idxI, 2); % I 路低位 bits(3:4:end) floor(idxQ/2); bits(4:4:end) mod(idxQ, 2); end % ---- 主循环与理论曲线对比 ---- ebn0 0:2:14; berSim zeros(size(ebn0)); for ii 1:numel(ebn0) esn0_dB ebn0(ii) 10*log10(4); nErr 0; nBit 0; while nErr 100 nBit 2e7 b randi([0 1], 4*20000, 1); r awgn_esn0(qam16_mod(b), esn0_dB); bh qam16_demod(r); nErr nErr sum(bh ~ b); nBit nBit numel(b); end berSim(ii) nErr / nBit; end ebn0Lin 10.^(ebn0/10); x sqrt(0.8 * ebn0Lin); Qf 0.5 * erfc(x/sqrt(2)); berTheory (3*Qf - 2.25*Qf.^2) / 4; semilogy(ebn0, berSim, o-, ebn0, berTheory, k--); grid on; xlabel(Eb/N0 (dB)); ylabel(BER); legend(仿真, 理论, Location, southwest);invGray [0 1 3 2]是tab的逆向查找表排序位置 0、1、2、3 分别对应电平-3, -1, 1, 3而这三个电平在tab里的标号是 0、1、3、2。写错这张表的表现是星座图完全正常但 BER 曲线会比理论差——因为解调出来的比特顺序和调制时对不上相当于人为增加了错误。bits(1:4:end)这一组赋值的顺序必须和qam16_mod里bits(1:4:end)的取用顺序一一对应I 高位、I 低位、Q 高位、Q 低位。这两处只要有一处顺序换了整条曲线就废了而且因为错的是比特位置而不是符号星座图看不出来。5. 避坑仿真的 BER 和理论曲线对不上先查这五处曲线对不上理论九成的问题集中在这五个地方。按「现象 → 原因 → 解决」列出来下次直接照单排查。5.1 映射与归一化环节的三个坑坑一曲线形状对整体平移了 6 dB。现象是仿真曲线和理论曲线平行但始终差大约 6 dB 的位置。 原因是 Eb/N0 和 Es/N0 换算漏了。16QAM 每符号 4 bit10*log10(4) 6.02 dB加噪函数拿到的是 Es/N0如果直接把 Eb/N0 传进去实际噪声就小了 6 dB仿真显得「更好」。 解决是把换算固定在主循环里做一次加噪函数只接受 Es/N0函数名里带上esn0字样写代码的时候就不会混。坑二高信噪比处 BER 跌到某个值就不动了形成错误平台。现象是曲线在 1e-3 或 1e-4 附近走平再加大功率也不降。 原因是调制表和判决门限不匹配。最常见的是tab写成[-3 -1 1 3]少了格雷交换或者门限写成±1/√10、±3/√10。这时候符号本身没有错但比特映射错位等价于恒定数量的固定错误。 解决是单独写一段自检调制一个已知的 4 比特看输出的星座点坐标是不是落在预期位置再把坐标喂回解调看比特能不能原样出来。这个往返测试十行代码比盯曲线快得多。坑三归一化系数用错曲线整体左移约 10 dB。现象是仿真比理论好一大截而且信噪比越高差距越明显。 原因是电平表没除以sqrt(10)Es 实际是 10 而不是 1但加噪函数仍按 Es1 算N0 1/(Es/N0)等于噪声按 1/10 的功率在加。 解决办法是把归一化写进电平表定义那一行注释里标明Es 1全脚本不再出现第二个sqrt(10)的除法。5.2 信噪比与统计环节的两个坑坑四低信噪比区曲线剧烈抖动相邻点能差一个数量级。现象是 04 dB 那几个点上下跳。 原因是错误数太少。低 SNR 下虽然 BER 高但如果总比特数只有几万统计涨落依然很大更常见的是有人在高 SNR 段仍固定跑 1e5 比特一个错误都没有画出来直接断线。 解决办法是用错误数作为停止条件目标 100 个错同时设总比特上限。另外画图前把ber 0的点剔除否则semilogy会因为 log(0) 报警告并且断线。坑五换了随机种子曲线整体跳了一截。现象是两次运行结果差异超出预期。 原因是rng没固定或者固定在了错误的位置——比如在循环内部每次都rng(1)那每一帧的比特和噪声完全一样等于只有一帧的样本统计量严重不足。 解决是在脚本开头rng(2024)固定一次循环内部绝不重置。同时记住固定种子只保证可复现不改善统计精度曲线是否可信还是要看错误数。注意蒙特卡洛仿真的结果是随机变量两条曲线在低信噪比区有 10%20% 的相对差异属于正常范围。判断「仿真错了」之前先把每个点的错误数打印出来看一眼错误数不到 100 的点不值得较真。6. 从 AWGN 推到瑞利一套能反复用的验证套路AWGN 跑通只是及格线真正能说明接收机抗不抗打得看衰落信道。改动其实只有两行乘一个复高斯信道系数接收端做理想信道估计除回去。function r rayleigh_chan(s, esn0_dB) % 平坦瑞利衰落 AWGN接收端假设理想信道估计完美 CSI h (randn(size(s)) 1j*randn(size(s))) / sqrt(2); % CN(0,1) esn0 10^(esn0_dB/10); n0 1 / esn0; n sqrt(n0/2) * (randn(size(s)) 1j*randn(size(s))); r h .* s n; r r ./ h; % 理想均衡噪声被 h 除等效方差变成 n0/|h|^2 end关键在最后一行。除以h之后信号恢复成s但噪声方差变成了n0/|h|²——深衰落时|h|很小噪声被放大得很厉害这就是瑞利信道误码率曲线在高 SNR 处斜率变缓的原因。这个「除回去」的动作必须在加噪之后做顺序反了就成了先把噪声除掉再加噪等于白送增益。验证套路我一般按三步走。第一步先在 AWGN 下把仿真和erfc理论对齐误差控制在 0.5 dB 以内这一步确认映射、归一化、判决都没问题。第二步把信道换成瑞利和瑞利下的理论 BER0.5*(1 - sqrt(Eb/N0/(1Eb/N0)))是 BPSK 的16QAM 要用数值积分或者近似式对比观察高 SNR 处的斜率是否一致。第三步扫不同的targetErr50 / 100 / 200看曲线是否收敛收敛了才敢往报告里写。下面这张表是我每次跑完会对照的检查项检查项合格标准不合格的典型表现星座点往返调制→解调比特零错误固定数量错误BER 有平台Es 归一化平均abs(s).^2≈ 1曲线整体平移约 10 dBEb/Es 换算与理论曲线差 1 dB平行平移 6 dB统计量每点错误数 ≥ 100低 SNR 段剧烈抖动判决门限±2/√10与0高 SNR 处提前走平有个教训印象很深早些年跑一组 16QAM 曲线仿真和理论差 3 dB我盯着噪声方差和 erfc 公式查了大半天最后发现是qam16_demod里invGray写成了[0 1 2 3]。星座图漂亮得挑不出毛病曲线却怎么也贴不上。从那以后我每次改完映射表都强制走一遍「单符号往返测试」——随便挑 4 个比特过一遍调制再解调打印出来对一眼十秒钟的事能省掉一整晚的排查。希望帮到你。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助？
