简介:这是一份面向通信与信号处理方向的MATLAB误码率仿真资源,完整实现了正交相移键控(QPSK)调制解调、扩频解扩、快速傅里叶变换(FFT)频偏估计与低密度奇偶校验码(LDPC)编译码的级联系统。资源适用于学习无线通信基带算法、开展误码率仿真实验的工程师与高年级学生,也可作为通信课程设计及毕业设计的参考实现。全套资源共30个文件,其中以18个m脚本源码为主,辅以11个mat数据文件和1个txt说明文件,压缩包整体约130KB,基于MATLAB 2024b编写,代码含中文注释,并附有程序操作视频,方便对照运行与调试,目前已有150人学习浏览。仿真链路完整覆盖发送端与接收端:首先生成随机二进制信息序列,经LDPC编码后进行扩频,再通过QPSK调制并加入频偏,随后送入AWGN信道叠加高斯白噪声;接收端依次完成FFT频偏估计、QPSK解调、解扩和LDPC译码,最终统计误比特数并计算误码率。整体流程清晰,能帮助读者深入理解各模块原理及级联系统的性能折中。
1. QPSK调制解调+扩频解扩+FFT频偏估计+LDPC编译码:这套Matlab仿真链路到底在验证什么
做误码率仿真最怕的不是某个模块跑不通,而是每个模块单独调试都正常、拼成整条链路后结果却对不上:信号发出去了、频谱成形也对,误码率曲线偏偏离理论值差好几个dB,根本不知道是频偏没拉回来、扩频码相位没对齐,还是LDPC译码器的LLR(对数似然比)符号反了。这套Matlab仿真工程把QPSK调制解调、扩频解扩、FFT频偏估计、LDPC编译码四段串成一条完整收发链路,逐点统计误码率,并附带中文注释程序与操作视频。对刚接触物理层链路仿真的新手,它是一份可拆解的模板;对写过单模块仿真的熟手,它的价值在于展示四个模块之间的参数怎么对齐、接口怎么衔接——而这恰恰是光看教科书最容易翻车的地方。
2. QPSK调制解调与Eb/N0换算:先把误码率曲线的横轴定义推准确
2.1 QPSK映射与基带波形生成
QPSK位于整个仿真链路的最前端,输入是编码后的比特流,输出是复基带符号。它的正确性直接影响后续所有环节,因此我拿到工程后第一个动作是单独跑通发端映射与收端判决,确认星座点位置。
工程里QPSK映射的基本写法如下:
% QPSK Gray映射:输入比特流,输出归一化复符号 function sym = qpsk_map(bits) % bits 必须是偶数长度,不足时补零保证成对 if mod(length(bits), 2) bits = [bits(:); 0]; end symPts = [1+1j, -1+1j, -1-1j, 1-1j] / sqrt(2); % 归一化到单位能量 pairs = reshape(bits(:), 2, []).'; % 每两比特一组 idx = pairs(:,1) * 2 + pairs(:,2) + 1; % 00->1, 01->2, ... sym = symPts(idx).'; end逻辑说明:这里的关键是Gray映射,相邻星座点间只差一个比特,在误码率统计时不会因为单次符号判决错误引入双比特错误。/ sqrt(2)这一步把四个星座点归一化到单位能量,即每个符号的平均能量为1,方便后续噪声功率标定。索引计算里的pairs(:,1)*2 + pairs(:,2)把二进制比特组映射为十进制下标,顺序与symPts的排列一一对应。
接收端判决是对每个符号的实部和虚部分别做门限判决,门限取0:
% QPSK硬判决解调:real>0判0,imag>0判0 function bits = qpsk_demap(sym) sym = sym(:); bits = zeros(numel(sym) * 2, 1); bits(1:2:end) = real(sym) < 0; bits(2:2:end) = imag(sym) < 0; end这段代码直接体现QPSK的本质:正交载波上的两个BPSK并行传输。需要注意,以上判决语句成立的前提是接收符号已经完成频偏补偿与相位同步;如果残余频偏让星座旋转了30度,门限判决会瞬间失效,误码率直接掉到接近0.5。这也是为什么工程里要在解调前放FFT频偏估计模块。
2.2 Eb/N0与SNR的换算及噪声注入
误码率曲线的横轴是Eb/N0(每比特能量与噪声功率谱密度之比),但Matlab的awgn函数接收的是SNR,两者换算错一步,整条曲线就会整体平移几个dB,看起来像系统性能差,其实只是坐标系没对齐。
复基带仿真中的换算关系如下:
% Eb/N0(dB) 转 SNR(dB),QPSK每符号携带2比特 EbN0dB = 0:2:12; SNRdB = EbN0dB + 10*log10(2); % QPSK 的 k=2 % 发送符号平均能量 Es = 1,单边噪声功率谱密度 N0 = 10^(-SNRdB/10) N0 = 10.^(-SNRdB/10); noiseSigma = sqrt(N0 / 2); % 复噪声每维方差 tx = qpsk_map(infoBits); rx = tx + noiseSigma * (randn(size(tx)) + 1j*randn(size(tx)));参数说明:10*log10(2)就是QPSK相对BPSK的3dB增益来源,每符号两比特等价于每比特能量被摊薄,同样符号能量下Eb更低。noiseSigma的计算里除以2是因为复基带噪声的实部与虚部各占一半方差。很多仿真里误码率曲线偏优或偏劣几个dB,根因就是这里直接把Eb/N0数值传给了awgn,导致噪声功率被少加或多加。
更深一层,这套仿真如果启用了扩频,是否要把扩频增益叠加进SNR?我的经验是不用。扩频增益已经在解扩的相干累积中体现,收端匹配滤波后的信噪比自动提高了10*log10(SF) dB;如果在加噪时就人为叠加增益,等于同一笔收益算了两次,曲线会虚假偏优。
2.3 解调器还要输出软信息给LDPC
硬判决解调适合单独验证QPSK模块,但一旦接入LDPC译码器,硬判决会把信道信息量化成0/1,丢掉幅度信息,LDPC的纠错能力会大幅下降。工程中LDPC译码需要的是LLR(对数似然比),即每个比特为0相对于1的对数概率比。
QPSK软信息计算常见做法是按实部、虚部分别提取两路软值:
% QPSK软解调:输出实部与虚部位的 LLR function llr = qpsk_soft_demap(rx, N0) % rx: 接收符号,N0: 单边噪声功率谱密度 llr = zeros(2 * numel(rx), 1); llr(1:2:end) = 2 * sqrt(2) * real(rx) / N0; % I路 llr(2:2:end) = 2 * sqrt(2) * imag(rx) / N0; % Q路 end参数说明:sqrt(2)来自星座点幅度与单位能量的换算,2/N0是AWGN信道下LLR的标准形式。这里有个容易搞反的细节:如果QPSK映射时把[1+1j]对应比特00,那么real(rx)为正说明I路更可能为0,LLR为正;上式的符号约定必须与映射表一致,否则LDPC译码器会收敛到反方向,误码率比不译码还差。工程代码里中文注释一般会标明“LLR正值对应比特0”,跑之前先确认这一条。
3. 扩频解扩与FFT频偏估计:载波失步时怎么把频偏拉回来
3.1 m序列扩频与相关解扩的数学关系
扩频模块的作用是把每个符号的功率摊到更宽的频带上,收端通过本地序列做相关累积把信噪比抬回来。这里的核心是选取自相关特性尖锐的伪随机序列,m序列就是最常用的一种。
% 生成 m 序列:寄存器初值与反馈抽头决定序列 function seq = mseq(regInit, tapIdx) % regInit: 寄存器初始状态行向量,如 [1 0 0 1] % tapIdx : 反馈抽头位置索引,如 [1 4] n = numel(regInit); len = 2^n - 1; % m序列周期 seq = zeros(1, len); reg = regInit(:).'; for i = 1:len seq(i) = reg(end); % 输出最高位 fb = mod(sum(reg(tapIdx)), 2); % 抽头异或 reg = [fb, reg(1:end-1)]; % 移位 end end参数说明:n级寄存器生成周期为2^n - 1的m序列,自相关函数在零移位处有峰值len,其他移位处接近-1,这个特性决定了解扩的本质是一个匹配滤波过程。tapIdx的选择不是任意的,不同抽头组合决定序列是否为m序列,选取不当周期会缩短,所以仿真里一般直接用本原多项式对应的抽头组合。
解扩操作就是本地序列与接收序列做相关累积:
% 解扩:每个符号内的扩频码片累加 function symOut = despread(rxChips, pnSeq, sf) % rxChips: 接收码片序列,pnSeq: 本地PN序列,sf: 扩频因子 nSym = floor(numel(rxChips) / sf); symOut = zeros(1, nSym); for k = 1:nSym seg = rxChips((k-1)*sf+1 : k*sf); symOut(k) = sum(seg .* pnSeq) / sf; % 相关累积后归一化 end end逻辑说明:seg .* pnSeq在相位对齐时是能量累加,信号分量幅度增长为原来的sf倍,而独立噪声分量累加后幅度只增长sqrt(sf)倍,信噪比净增益正好是sf倍。除以sf是为了把符号幅度还原到发送量级,方便后续的频偏估计和LLR计算使用同一个噪声标定。
3.2 FFT频偏估计:从粗估计到补偿
接收端面对的典型问题是载波残留频偏,它让星座点在每一个符号周期旋转一个固定角度,时间越长旋转角度越大。FFT频偏估计的思路是:把接收信号变换到频域,找到能量最大的谱线位置,那条谱线对应的频率就是频偏估计值。
% FFT频偏估计:解扩后符号进行FFT,找峰值谱线并插值细化 function fHat = fft_freq_est(sym, fs, nfft) % sym : 解扩后的复符号序列 % fs : 符号率 % nfft: FFT点数,一般取1024或2048 seg = sym(:).'; M = min(nfft, numel(seg)); win = hann(M, 'periodic').'; % 加窗抑制频谱泄漏 spec = fft(seg(1:M) .* win, M); specMag = abs(spec(1:floor(M/2))); [~, iPeak] = max(specMag); % 抛物线插值提精:用峰值与左右两根谱线定位 if iPeak > 1 && iPeak < numel(specMag) left = specMag(iPeak-1); right = specMag(iPeak+1); delta = 0.5 * (left - right) / (left - 2*specMag(iPeak) + right); iPeak = iPeak + delta; end fHat = (iPeak - 1) / M * fs; % 归一化到实际频率 end参数说明:FFT分辨率是fs/M,M越大频偏估计越精细,但参与FFT的符号数越多要求观察时间越长。加窗是为了让主瓣展宽从而用抛物线插值得到更准确的峰值位置,虽然是粗略估计,工程里已够用。峰值谱线索引换算频率时,iPeak-1对应的是数字频率(iPeak-1)/M,再乘以符号率fs才是真实的Hz频偏。
补偿时用复指数相乘:rxComp = rx .* exp(-1j*2*pi*fHat*t),其中t从0递增到符号周期倍数。补偿完成后残余频偏小于fs/(2*M),对几十个符号的短帧来说相位旋转量可以忽略,对长帧则需要分段补偿或配合锁相环消除残余值。
3.3 解扩与频偏估计的先后次序怎么定
这是工程实现里一个隐蔽的坑:先估计频偏还是先解扩,结果差别很大。常见做法是先解扩、再频偏估计,原因在于解扩把信噪比抬高了10*log10(sf)dB,FFT峰值在低信噪比下更容易被识别。如果反过来先对未解扩的宽带信号做FFT,信号功率被摊在整个扩频带宽里,峰值往往淹没在噪声中,可能抓到的是噪声尖峰。
如果链路里频偏过大,导致一个符号周期内就旋转了接近90度,解扩相关峰也会被破坏。这时候要先做粗频偏估计把频偏压到符号率的百分之几以内,再做解扩和精估计。我一般会这么判断:先算一下2*pi*fOffset*T_symbol,如果这个相位旋转量超过0.2弧度,就必须先粗补偿再解扩,否则解扩相关峰增益会掉好几个dB,误码率曲线会出现高信噪比段反而抬头的怪异现象。
4. LDPC编译码与误码率统计:如何让曲线真正反映纠错增益
4.1 LDPC编码器配置:校验矩阵与码率选择
LDPC编码环节的输入是信息比特,输出是带校验位的码字。工程里通常直接使用标准校验矩阵构造编码器,避免从零开始构造稀疏矩阵的复杂度。Matlab通信工具箱提供了现成接口,但参数选型直接决定编码增益。
% 使用 DVB-S.2 标准 1/2 码率校验矩阵构造 LDPC 编码器 H = dvbs2ldpc(0.5); % 码率 1/2,码长 64800 encCfg = ldpcEncoderConfig(H); decCfg = ldpcDecoderConfig(H); infoBits = randi([0 1], 32400, 1); % 信息比特长度 32400 codeBits = ldpcEncode(infoBits, encCfg); % 编码后码长 64800参数说明:dvbs2ldpc(0.5)生成的是DVB-S.2标准下码率1/2的稀疏校验矩阵,码长64800比特。码率越低冗余越多,纠错能力越强但有效吞吐下降。仿真中如果只关心算法验证,也可以直接导入自带的ldpcQuasiCyclicMatrix构造QC-LDPC矩阵,码长可以缩减到几百比特,跑得快但瀑布区特性略差。
编码后需要注意码字与QPSK映射的衔接:codeBits是0/1序列,直接送入qpsk_map时要先转换为双极性,因为星座映射函数内部按0/1处理,不会帮你自动转换符号极性。我见过不少人在这个接口上大意,把映射表顺序搞反,导致发送星座整体旋转了45度,误码率从1e-5抬到1e-2。
4.2 BP译码的参数:迭代次数与外信息
LDPC译码用置信传播算法,核心参数是迭代次数与LLR输入格式。迭代次数太少,外信息还没来得及在变量节点与校验节点之间充分传播,纠错增益发挥不出来;太多则增加仿真时长,通常在30到50次迭代后增益趋于饱和。
% QPSK软值接入LDPC译码器 % 注意:LLR正值对应比特0,必须与qpsk_soft_demap的约定一致 llr = qpsk_soft_demap(rxCompensated, N0); decBits = ldpcDecode(llr, decCfg, 50); % 50次迭代参数说明:ldpcDecode的第三个参数是迭代上限,工程里一般先跑一次完整扫描确认25次迭代是否达到性能饱和。如果相同误码率下平均迭代次数已经贴近上限,说明这个工作点在瀑布区边缘,需要增加迭代次数或换更低码率。
LLR输入的符号约定是LDPC接入最容易出错的地方。qpsk_soft_demap按照“LLR为正对应比特0”输出,ldpcDecode内部默认也按这个约定解释,但如果发端映射做了改动,必须同步调整符号。判断方法很简单:先固定信噪比跑一小段,看译码后误码率是否显著低于硬判决;如果反而更高,多半是LLR符号反了。
4.3 误码率统计流程:一次运行要采多少比特
误码率曲线的可信度取决于统计量。低误码率区域如果样本不足,错误比特数可能是0或个位数,曲线在对数坐标上会出现锯齿甚至断线。常见做法是设定“至少收集100个错误比特”作为停止条件,再叠加最小比特数下限防止高误码率时跑太久。
% 误码率统计主循环 BER = zeros(size(EbN0dB)); for snrIdx = 1:numel(EbN0dB) numErr = 0; numBits = 0; while numErr < 100 && numBits < 1e7 % 双停止条件 infoBits = randi([0 1], 32400, 1); codeBits = ldpcEncode(infoBits, encCfg); txSym = qpsk_map(codeBits); % 加信道、频偏、解扩、频偏估计补偿后,得到 rxSym llr = qpsk_soft_demap(rxSym, N0); decBits = ldpcDecode(llr, decCfg, 50); numErr = numErr + sum(decBits ~= infoBits); numBits = numBits + numel(infoBits); end BER(snrIdx) = numErr / numBits; end参数说明:numErr < 100保证了在误码率1e-5量级也能有足够的错误样本支撑统计置信度,numBits < 1e7是防止过低信噪比下无限循环的兜底。统计对象是infoBits而非codeBits,因为用户关心的是有效信息的错误率;如果统计编码后码字,码率1/2时会把误码率拉低一半,曲线看起来偏优。
另外,每一次蒙特卡洛循环都要重新生成信息比特并完整编码,不要复用上一轮的codeBits。复用会导致前后帧相关,统计得到的误码率失去独立采样意义,严重时曲线出现规律性波动,看起来像存在某种周期干扰。
5. 避坑与排查:五个让仿真结果失真并误导设计的常见错误
5.1 频偏估计分辨率不够,残余频偏让星座旋转
现象:中高信噪比段误码率曲线下降到一定程度后不再下降,形成一个平台,星座图上几个点连成弧线而不是聚成清晰四点。 原因:FFT频偏估计最小分辨率是fs/M,残余频偏虽然小,但长帧末尾相位旋转累积到几十度,QPSK判决门限附近的符号被旋转到错误区间。 解决:先增大FFT点数或对多段数据做平均估计;仍不够就把补偿后的信号再送入一个二次精估计,使用抛物线插值后的分谱线定位;工程里常见做法是粗估计加锁相环细跟,前馈FFT只负责把频偏压到捕获范围内。从那以后我每次跑链路都会先看残余频偏对应相位旋转量在帧尾是否超过30度,超了就必加细估计。
5.2 扩频码相位没对齐,整套扩频形同虚设
现象:开启扩频后误码率与不扩频时几乎一样,甚至更差,且怎么调信噪比曲线都不掉下来。 原因:本地PN序列相位与接收端码相位错开,相关累积时信号被当成噪声平均掉了,处理增益没拿到。 解决:解扩前做码相位搜索,对每个整数码片偏移计算一次相关能量,取能量最大者作为同步相位。搜索步进可以先用1个码片粗搜,再在半码片精度上细搜。工程里的操作视频通常对这一节有单独演示,因为很多人误以为只要扩频码一样就能直接解扩,忽略了码同步这个前提。
5.3 Eb/N0与SNR混用,整条曲线平移几dB
现象:仿真误码率曲线与理论QPSK曲线完全平行,但整体向右偏移2到3dB。 原因:把Eb/N0数值直接当SNR注入awgn函数,或少加了10*log10(2)这一项。更隐蔽的错误是在开启扩频后把处理增益重复叠加进信噪比。 解决:统一用N0=10^(-snrDb/10)手工生成噪声,不依赖awgn的自动功率换算;扩频增益不在加噪处叠加,解扩过程会自动体现。写仿真前先把信噪比链路的各项换算列一张表,逐个符号率、比特率、扩频因子代入验证,能省下一整天的调参时间。
5.4 LDPC迭代次数太少,瀑布区消失
现象:误码率曲线有LDPC和没LDPC差别很小,信噪比到2dB后曲线不再快速下降,像卡住了一样。 原因:BP译码迭代次数不足,外信息没有在网络里充分传播,纠错增益没发挥出来。有时候是用了min-sum近似却追求sum-product的理论性能,两者差0.3到0.5dB。 解决:把迭代次数从10提到30、50,观察瀑布区是否持续左移;确认不再明显移动后固定下来。追求精度时用log-domain BP算法替代min-sum,避免浮点下溢。另一点是信噪比偏高时LLR数值很大,译码器内部容易饱和,必要时对输入LLR做截断归一化。
5.5 误码率统计样本不足,曲线低段抖动断裂
现象:低误码率段曲线异常抖动,多个信噪比点误码率是0,log坐标直接断线,看不出瀑布趋势。 原因:统计比特数固定但比特数不够,比如只跑了1e5比特就收工,误码率1e-4量级时平均只能采到10个错误样本,方差极大。 解决:按错误比特数控制停止条件,至少100个,稳妥取200个;同时设最小统计比特数兜底。工程里如果遇到仿真时间过长,可以对高信噪比段开并行循环,把每个信噪比点的蒙特卡洛循环放独立进程里跑,最后汇总曲线。
6. 三个验证手段:用单模块开关确认每个环节都发挥了作用
6.1 模块旁路对比:把某个环节从链路里摘掉再看误码率
整链路最怕“黑匣子”,曲线对上了不知道是哪个模块起作用,对不上也不知道是谁拖后腿。我的做法是给仿真加一个配置开关结构体,每个模块可以独立关闭:
% 旁路开关:置0关闭对应模块 cfg.enableSpread = 1; % 1-扩频解扩使能 cfg.enableFFT = 1; % 1-频偏估计与补偿使能 cfg.enableLDPC = 1; % 1-LDPC编译码使能 if cfg.enableSpread rxSym = despread(rxChips, pnSeq, sf); else rxSym = rxChips; end对比三组曲线,能直观看出每个模块的边际增益。比如开启LDPC后误码率反而上升,基本可以断定LLR符号约定反了或者交织器缺失;开启频偏估计后曲线不变,说明信道里根本没加频偏,模块在空转。这个手段成本极低,但对定位故障模块效率极高,是我拿到任何仿真工程后第一个执行的动作。
6.2 中间探针:星座图、频域峰值与迭代收敛轨迹
旁路对比能定位是哪个模块出了问题,中间探针则能告诉你问题具体出在哪个环节。调试过程中我在关键节点加探针:发射符号星座图、解扩后星座图、频偏估计谱峰位置、LDPC平均迭代次数。其中最有价值的是LDPC平均迭代次数——如果工作点处处贴着迭代上限,说明译码器一直在极限工作,需要更低码率或更高迭代次数;如果平均迭代次数远低于上限,说明还有余量,可以把迭代次数降下来提速。
6.3 参数扫描脚本化:把验证过程固化成可重复的批次任务
手动改参数来回跑仿真,既容易漏记结果,又容易在改回参数时弄错。我习惯把扩频因子、频偏大小、LDPC迭代次数作为扫描变量,统一封装成一个批次脚本:
% 参数扫描:扩频因子×频偏大小×迭代次数 results = struct(); for sf = [16 32 64] for fOffset = [0 0.01 0.05] for iter = [10 30 50] results(sf, fOffset, iter) = runBER(sf, fOffset, iter); end end end扫描完成后把results结构体保存成文件,后续画图、写报告直接调用,不用重新跑仿真。这套批次脚本跑出来的结果除了主误码率曲线,还能附带每个参数组合的迭代余量统计,用来判断仿真参数是否选在合理区间。从那以后,我每次拿到新工程都强制走一遍单模块开关自检、中间探针核查、参数扫描固化这三个动作,再跑整链路误码率;这套验证手段花不了十几分钟,但能避免你在调参上浪费一整个下午。希望这篇拆解和验证习惯能帮到你。
本文还有配套的精品资源,点击获取