
简介面向无线通信与信号处理方向的研究者这里提供一套基于STBC的MIMO-OFDM系统误码率分析MATLAB仿真代码。代码实现了M-PSK调制、OFDM调制、STBC编码、AWGN信道传输、STBC译码与OFDM解调全流程通过蒙特卡洛仿真统计不同天线配置下的系统误码率可用于通信课程设计、毕业设计或科研预研。资源包共十个文件其中五个是m脚本分别对应BPSK、QPSK、8PSK、16PSK等调制方式及训练符号生成四个jpg图片为不同天线配置下的误码率曲线对比方便直接查看结果一个txt文件记录关键词说明压缩包整体只有一百六十三KB结构清晰适合快速下载使用。代码运行环境为MATLAB2016b可直接运行出图已有三百五十人下载学习使用者可基于脚本灵活修改调制阶数和天线数目复现多组实验也可进一步扩展信道模型或编码方案作为系统性能研究的基础工具。1. 从误码率曲线反推STBC的MIMO-OFDM链路一份可直接跑的MATLAB仿真包做无线通信课程设计时最烦的就是把基于STBC的MIMO-OFDM系统误码率分析从论文公式变成能跑的仿真。公式里Alamouti矩阵写得很漂亮但一落地到MATLAB就翻车傅里叶变换方向反了、天线维度对不上、误码率曲线卡在一个平台上不去。这份资源把BPSK、QPSK、8PSK、16PSK四种调制以及1x1、2x2、4x4三种天线配置的误码率分析做成了完整代码用MATLAB 2016b可以直接出图。它解决的是“链路能不能闭环”的问题发送端随机比特经M-PSK映射、OFDM调制、频域STBC编码后过AWGN信道接收端做STBC译码和OFDM解调最后统计误码率并画成曲线。正在做MIMO-OFDM仿真、复现论文曲线、或者想抄一份可靠基线的同学可以省下大量调试时间。2. STBC编码与OFDM子载波的接口先弄懂Alamouti矩阵怎么进IFFT2.1 为什么STBC要放在频域而不是时域OFDM系统里STBC不是在时域对整段波形做编码而是对每个子载波上的符号对做空时编码。原因很直接OFDM把频率选择性信道切成了多个平坦子信道所以STBC可以逐子载波地应用接收端译码时也逐子载波处理复杂度低。以经典Alamouti 2x1为例两个连续的频域符号在一对天线上按下面的方式排列% 某个子载波k上的两个已调制符号 s1、s2 % 两个时隙、两根发送天线 txBlock [s1, -conj(s2); ... s2, conj(s1)]; % 行对应时隙1、时隙2列对应天线1、天线2逻辑说明时隙1天线1发 s1天线2发 s2时隙2天线1发 -conj(s2)天线2发 conj(s1)。利用这一对共轭对称关系接收端只需把两时隙收到的信号做简单线性合并就能解出 s1、s2并且信道矩阵被正交化不需要暴力求逆。这就是STBC能拿到分集增益的核心。参数说明txBlock的维度是 2×2如果你在2x2系统里用Alamouti相当于两发两收接收合并的时候多一根接收天线参与最大比合并误码率会比2x1再低一截。而资源中的4x4配置我一般会按两层Alamouti去理解把四根发送天线分成两组每组都走上面的2×2块编码接收端四根天线合并分集。这样扩到4x4不需要重新推导复杂的正交码实现成本最低。2.2 training_symbol.m 到底在同步什么打开压缩包会看到一个training_symbol.m很多人会以为AWGN信道下不需要信道估计训练符号就没什么用。实际上即便信道是纯AWGN、没有频率选择性衰落OFDM接收端也要先定位一个OFDM符号的起点否则FFT窗口错位子载波之间会串扰。training_symbol.m生成的就是用于粗同步的已知序列。% training_symbol.m 的核心逻辑 function trainSym training_symbol(Nfft, M) % 生成频域恒模序列减少峰均比影响 phaseSeq exp(1j * pi/3 * (0:Nfft-1).); trainSym ifft(phaseSeq, Nfft); end逻辑说明频域每个子载波上的幅度都取1只有相位在旋转这样IFFT出来的时域波形能量均匀便于接收端做相关峰检测。接收端把本地生成的trainSym与接收信号做滑动相关峰值位置就是符号边界。原始代码里training_symbol.m被主脚本调用目的就是保证后续OFDM解调时FFT窗口对齐。参数说明Nfft是FFT点数也是子载波总数M在这里不是调制阶数而是为了保持函数签名一致实际可以不用。如果你自己改代码可以把训练符号的相位序列改成伪随机序列效果更稳。3. 拆解 a_Qpsk.m 的完整流程调制、STBC、AWGN、译码每次都在干什么3.1 参数初始化天线数、调制阶数、子载波和CP设置a_Qpsk.m是QPSK对应的主脚本a_Bpsk.m、a_8psk.m、a_16psk.m结构几乎一样区别只在M取值。先看参数区clear; close all; clc; M 4; % 调制阶数: 2BPSK, 4QPSK, 88PSK, 1616PSK numTx 2; % 发送天线数 numRx 2; % 接收天线数 Nsc 64; % OFDM子载波数 Nfft 64; % FFT点数和Nsc一致 Ncp 16; % 循环前缀长度 EbN0_dB 0:2:20; % 比特信噪比扫描点 numMc 1e4; % 每个SNR点蒙特卡洛帧数逻辑说明这里的numTx、numRx可以改成1或4对应1x1、2x2、4x4三种天线配置。Nfft和Nsc相等表示所有子载波都传数据没做导频插空如果你后续要加信道估计得留出导频子载波。Ncp取16只要大于信道最大时延扩展就不用担心子载波间串扰AWGN信道下写8都够但代码里习惯用16。参数说明EbN0_dB是横轴注意它是以每个比特为单位不是每个符号。M-PSK的符号功率相同但每个符号携带的比特数不同所以换算SNR时要用10*log10(log2(M))做偏移。这点很多新手会踩后面避坑章再展开。3.2 数据生成与M-PSK映射% 生成随机比特并映射成M-PSK符号 numBits log2(M) * Nsc * numTx; % 一个“发送帧”包含的总比特数 dataBits randi([0 1], numBits, 1); bitMat reshape(dataBits, log2(M), []); symIdx bi2de(bitMat.); % 每log2(M)个比特转一个十进制数 modSyms pskmod(symIdx, M, 0, gray); % 灰色编码降低相邻符号误判影响逻辑说明randi生成长度为numBits的列向量reshape后每行对应一个PSK符号的二进制下标bi2de把二进制转成0到M-1的整数pskmod再按PSK星座图映射成复数。使用gray参数让相邻星座点只有1个比特不同这样误码率统计时避免“一个符号错误”被计成多个比特错误。参数说明如果改成8PSKM8log2(M)3每个符号3比特改成16PSKM16。这里modSyms是一个列向量长度是Nsc*numTx也就是说把多天线上所有子载波要发的符号全部扁平化生成后面再按天线和子载波重组。3.3 频域STBC编码与IFFT发送% 把调制符号重新排列成 [子载波数, 每个子载波上的符号数] sMatrix reshape(modSyms, Nsc, []); % 每列对应一个天线/时隙的数据 txBlock zeros(2, numTx, Nsc); % 2个时隙, numTx根天线, Nsc个子载波 % 对每个子载波独立做Alamouti编码 for k 1:Nsc s1 sMatrix(k, 1); s2 sMatrix(k, 2); % 两发两收的标准Alamouti块 txBlock(:, 1, k) [s1; -conj(s2)]; txBlock(:, 2, k) [s2; conj(s1)]; end逻辑说明这步是STBC的核心。每个子载波k取出两个已调制符号s1、s2编码成两行对应两个时隙、两列对应两根发送天线。如果你的numTx4则把sMatrix再切分四根天线分两组各做一次Alamouti。1x1时不需要STBC可以直接把txBlock赋成原始符号否则1根天线没法构成正交块。发送端还要把频域块变换到时域加循环前缀才能进信道txTime zeros(2, numTx, Nfft Ncp); % 每根天线每个时隙的时域波形 for slot 1:2 for tx 1:numTx freqData squeeze(txBlock(slot, tx, :)); % 当前天线当前时隙的频域符号 timeData ifft(freqData, Nfft); txTime(slot, tx, :) [timeData(Nfft-Ncp1:end); timeData]; % 加循环前缀 end end逻辑说明IFFT把频域符号转成时域采样点再把最后Ncp个采样复制到头部形成循环前缀。这样即使信号经过多径信道只要时延扩展小于CP长度子载波之间的正交性就保住了。注意这里用numTx做循环变量1x1和4x4都能统一处理。3.4 AWGN信道与接收端FFT% 过AWGN信道按当前SNR加入复数高斯噪声 snr_dB EbN0_dB(idx) 10*log10(log2(M)); % 符号信噪比 snr_lin 10^(snr_dB/10); noiseVar 1 / snr_lin; % 符号功率归一化 rxTime txTime sqrt(noiseVar/2) * (randn(size(txTime)) 1j*randn(size(txTime))); % 接收端去CP、FFT还原到频域 rxFreq zeros(2, numRx, Nfft); for slot 1:2 for rx 1:numRx rBlock squeeze(rxTime(slot, rx, :)); rTime rBlock(Ncp1:end); % 去掉循环前缀 rxFreq(slot, rx, :) fft(rTime, Nfft); end end逻辑说明AWGN信道下信道响应矩阵可以认为是1噪声是复高斯白噪声实部和虚部方差各为noiseVar/2保证叠加后总噪声功率为noiseVar。接收端去CP后做FFT把时域信号还原成频域符号每个子载波上的频域值就是叠加了噪声的编码符号。这里rxFreq的维度是[时隙, 接收天线, 子载波]为STBC译码备好了数据。参数说明如果你的仿真要改成瑞利衰落信道不能直接套这个噪声模型先要把信道矩阵H乘上去然后还要在接收端做信道估计。原资源固定跑AWGN所以噪声模型简单。3.5 STBC译码与误码统计% 对每个子载波做Alamouti译码 estSyms zeros(Nsc, 2); for k 1:Nsc rSlot1 rxFreq(1, :, k).; % 时隙1收到的各接收天线信号 rSlot2 rxFreq(2, :, k).; % 时隙2收到的各接收天线信号 % 在AWGN且H1的前提下等效合成 estSyms(k, 1) rSlot1(1) conj(rSlot2(2)); estSyms(k, 2) rSlot1(2) - conj(rSlot2(1)); % 注意如果numRx2需要把多天线结果求和再合并这里简化展示2x2 end % 解调把估计符号映射回比特 rxSyms pskdemod(estSyms(:), M, 0, gray); rxBits de2bi(rxSyms, log2(M)); rxBits reshape(rxBits, [], 1); ber_current sum(rxBits ~ dataBits) / numBits;逻辑说明Alamouti译码的核心是利用编码时的共轭对称性把信道矩阵对角化。在AWGN且信道系数为1时rSlot1(1)对应的就是s1加噪声rSlot2(2)是conj(s1)加噪声取共轭后两路相加能获得2倍的信噪比增益。对s2则用rSlot1(2)减去conj(rSlot2(1))。pskdemod把复数符号映射回最近的星座点再转换成比特和发送端比较。参数说明这里为了展示原理只写了2x2的合并公式实际代码里numRx4时你要把四根接收天线的路径都累加进去然后归一化。别少累加否则4x4曲线会异常。4. 1x1、2x2、4x4 加 M-PSK 的误码率对比分集增益从哪来、曲线该怎么看4.1 不同天线配置的分集与合并增益STBC在AWGN信道下虽然不能像衰落信道那样提供频率分集但通过多根发送天线重复发送同一组符号接收端做合并后等效信噪比会提升。1x1就是直传没有合并增益2x2相当于两条独立噪声路径经过最大比合并比1x1大约好3dB4x4在理想合并下可以再往低噪声方向推几个dB。不过要注意AWGN下天线数增加带来的主要是阵列增益不是分集增益所以曲线形状不会像瑞利信道那样变陡。资源里给出的图片名也印证了这一点BPSK不同天线、QPSK不同天线、8PSK不同天线、16PSK不同天线每张图都是横轴信噪比、纵轴误码率画了1x1、2x2、4x4三条曲线。运行后能看到相同调制下天线数越多曲线越靠左。4.2 四种M-PSK的曲线走向差异从误码率角度看BPSK和QPSK的符号距离大误码率下降最快8PSK开始星座点变密集同样的信噪比下误码率抬高16PSK最密高信噪比下依然有明显底板。下面这张表是经验趋势具体数值随蒙特卡洛次数波动但相对关系固定调制方式每符号比特星座点间最小距离误码率曲线相对位置BPSK12最低最先掉下去QPSK21.414比BPSK高约1-2dB8PSK30.765比QPSK抬升明显16PSK40.390最高需要更大SNR表格说明如果你把横轴换成每符号信噪比Es/N0BPSK和QPSK曲线几乎重合但资源里很多图直接用Eb/N0所以QPSK比BPSK反而低一点因为QPSK每符号带2比特。这是横轴定义造成的不是系统差异。4.3 蒙特卡洛次数和曲线抖动误码率统计本质上是多次随机实验的频率帧数不够曲线就会抖动。资源代码里外层循环遍历EbN0_dB内层跑numMc次每次独立生成随机比特和噪声。建议至少让每个SNR点统计到50个以上错误比特否则在10^-4以下的点你会看到曲线毛刺。% 外层循环信噪比内层蒙特卡洛的常见框架 for idx 1:length(EbN0_dB) totalErr 0; totalBits 0; for mc 1:numMc [err, bits] run_ber_single(EbN0_dB(idx)); % 单次仿真 totalErr totalErr err; totalBits totalBits bits; end ber(idx) totalErr / totalBits; end逻辑说明这样统计的是累计误码率不会因为单次实验结果为零而变成0。run_ber_single是封装好的单帧处理函数返回错误比特数和总比特数主脚本只负责累加。如果你把numMc从1e4降到1e3低信噪比点可能差异不大但10dB以后曲线会抖成波浪线。参数说明numMc是每个SNR点的帧数不是总帧数。资源里原始代码可能没有单独封装函数而是把所有步骤都写在大循环里但统计思路是一样的。5. 避坑记录跑STBC的MIMO-OFDM仿真时我踩过的五个坑5.1 现象4x4误码率比2x2还高有一次我直接改numTx4、numRx4跑出来的BER曲线竟然比2x2还差完全不符合直觉。原因是STBC编码时我只构造了两根天线的txBlock却把4根天线全用来发送接收端也按4路合并但发送端有两根天线根本没分到编码符号等效信噪比计算全部错乱。 原因STBC编码矩阵维度和实际天线数不匹配或者收发维度的归一化因子没有随numTx变化。 解决先把txBlock的维度改成[时隙, numTx, Nsc]然后用循环给每一对发送天线填Alamouti块。接收合并时再按实际接收天线数累加并除以累加次数对应的归一化系数。我当时调试时在命令行打印size(txBlock)才发现问题。5.2 现象16PSK曲线在高SNR出现地板效应16PSK在SNR到15dB以后误码率停在10^-3上不去。这通常不是代码错而是蒙特卡洛帧数太少导致没有统计到足够多错误比特另一个可能是pskdemod的星座相位偏转没对准或者没有使用gray编码导致一个符号错会连带多个比特错。 原因蒙特卡洛样本量不足或者星座映射顺序和pskmod不匹配。 解决把numMc从1e4提到1e5或者2e5并确认pskmod和pskdemod都用同一个相位偏移和gray选项。我习惯写成pskmod(sym, M, 0, gray)这样相位偏移是0星座轴和文档一致。5.3 现象BPSK和QPSK曲线完全重合但论文里两者差3dB这种情况几乎可以断定横轴用错了。BPSK每符号1比特QPSK每符号2比特如果横轴是Es/N0两曲线确实基本重合如果论文里横轴标的是Eb/N0那QPSK的曲线应当比BPSK左移约3dB。资源里图片可能直接用了SNR作为横轴需要你换算。 原因SNR、Es/N0、Eb/N0三个概念混用。 解决按这个公式换算EsN0_dB EbN0_dB 10*log10(log2(M))。在循环里snr_dB就用EsN0_dB横轴却画成EbN0_dB这样不同调制阶数的曲线才能放在同一张图里比较。5.4 现象OFDM符号边界没对齐误码率一直停在10^-1跑通代码后我试着去掉training_symbol.m结果误码率怎么都下不来。原因很简单接收端不知道每个OFDM符号从哪里开始FFT窗口取错了位置子载波间干扰直接毁掉STBC块。 原因少了训练序列做粗同步或者训练序列长度不够导致相关峰不尖锐。 解决在接收端用training_symbol.m生成本地序列与收到的时域信号做滑动互相关取峰值位置作为FFT窗口起点。注意训练序列要放在发送帧首部不能混在数据里。资源里保留training_symbol.m就是干这个的。5.5 现象仿真跑一个SNR点要十分钟整条曲线要过夜蒙特卡洛循环里嵌套了两层for时隙、天线再加上每个子载波又循环了一遍STBC编码三重循环在MATLAB 2016b里非常慢。我看到有人直接numMc1e6跑完2x2 16PSK的整条曲线等了三个小时。 原因代码写得不面向矩阵子载波循环拖慢了速度。 解决先把STBC编码写成矩阵运算不要用for循环逐子载波处理再把最外层的SNR循环改成parfor并行。如果你只有一台四核电脑parfor能把整条曲线时间缩短到原来的四分之一左右。第六章我会给具体改造方案。6. 把这份代码改造成更实用的BER测试引擎Eb/N0换算与parfor加速蒙特卡洛6.1 让横轴统一成Eb/N0原资源图片的横轴如果标的是SNR不同调制方式的曲线可比性会变弱。我的习惯是在主脚本里保留一个标记变量比如axis_is_EbN0 true然后在信噪比换算处统一处理。EbN0_dB 0:2:20; idx 1:length(EbN0_dB); for k idx EsN0_dB EbN0_dB(k) 10*log10(log2(M)); % 设置AWGN噪声方差时使用EsN0_dB noiseVar 1 / (10^(EsN0_dB/10)); % 跑完单SNR点后ber数组索引仍对应EbN0_dB ber(k) run_ber_single(noiseVar); end逻辑说明把噪声方差的设置从主流程里抽出来run_ber_single只接收noiseVar和调制阶数内部不再管横轴定义。这样你画ber与EbN0_dB的关系时横轴就是统一的比特信噪比。BPSK和QPSK曲线会体现出1比特和2比特的差异16PSK曲线则明显往右移动。参数说明10*log10(log2(M))是把一个符号里的比特数折算成信噪比偏移M越大达到同等Eb/N0所需的Es/N0越高。如果你不这样折算16PSK在低SNR阶段就几乎全错看不出下降趋势。6.2 用parfor替代外层forMATLAB 2016b自带的Parallel Computing Toolbox可以放心用parfor前提是每个迭代之间没有依赖。误码率统计天然满足这个条件每个SNR点独立只需要最后汇总。parfor k 1:length(EbN0_dB) EsN0_dB EbN0_dB(k) 10*log10(log2(M)); noiseVar 1 / (10^(EsN0_dB/10)); ber(k) run_ber_single(noiseVar); end逻辑说明把原先的for改成parforMATLAB会自动把每个SNR点分给不同工作线程。但要注意run_ber_single里不能依赖前一个SNR点的随机数状态否则并行结果和串行结果不一致。解决方法是每次进入函数时用rng(shuffle)重新初始化随机数流或者直接把串行循环里的randi调用原样搬进去不要在循环外共用随机种子。参数说明parfor对内存占用比较敏感如果numMc很大尽量在run_ber_single内部只保存临时变量不要把所有SNR点的调制数据都堆在工作区。我一般会在函数末尾只返回err和bits两个标量。改造完以后我从拿到这份资源开始就养成了一个习惯不管原代码多乱先抽出参数区、再抽出run_ber_single、最后统一横轴和并行化。每次跑通一条MIMO-OFDM的误码率曲线我都强制走一遍这个流程省掉的不只是等待时间还有改一处错一处的心力。希望帮到你。本文还有配套的精品资源点击获取