ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

Turbo码MATLAB仿真从零实现:RSC、交织器与迭代译码要点

Turbo码MATLAB仿真从零实现:RSC、交织器与迭代译码要点 Turbo码这名字懂的人知道它靠的是“迭代”而不是“分集”。1993年Berrou等人提出后把逼近香农限这件事从理论变成了工程可触达的现实。用MATLAB实现Turbo码编码译码我这些年反复上手过不下十次每次体会都不一样第一次觉得交织器和迭代译码就是个状态机嵌套后面才意识到真正决定误码率曲线能不能压下去的是软信息的计算精度、交织器对齐关系、还有尾部比特那一点细节。这篇文章我把一套可以完整跑通、能在AWGN信道下直观看到增益的Turbo码编码译码链路写出来给正在做MATLAB仿真、通信系统课设、或者对现代信道编码感兴趣的读者一个可以直接抄作业的参考。代码不是唯一方案但思路和排查方法是可以复用的。1. 内容整体设计与思路拆解1.1 核心需求解析一条Turbo码仿真链路里到底有什么实现Turbo码第一件事不是急着写代码而是先把链路角色分清楚。一条完整的Turbo码仿真链路至少包含四个部分编码器、交织器、信道软信息计算、迭代译码器。这四个部分相互依赖交织器和译码器之间的关系尤其容易出错。编码器端用的是递归系统卷积码简写RSC它不是普通卷积码而是把反馈结构接进移位寄存器让每个信息比特和它的历史纠缠在一起。把两个RSC编码器并行接起来中间用一个交织器把信息比特重新排列就组成Turbo码的并行级联结构。这个结构是Turbo码名字的由来turbo在法语里是“涡轮机”的意思两台分量码像涡轮叶片一样交替做功。译码器端则是两套软输入软输出的MAP译码器彼此交换“外信息”。第一轮译码时第二路还没有任何先验所以外信息初始化为零下一轮第二路把从校验位里“榨”出来的额外信息反交织后交给第一路第一路又把自己的外信息交织后交给第二路。这样迭代几轮两路译码器互相纠偏误码率就一层层降下去。整个过程说起来简单但真正落到代码上最麻烦的不是单个算法而是编号、排列、尺度这些细节任何一个对不齐仿真出来就是一台噪声机。1.2 为什么选择手写而不是直接调用工具箱MATLAB的Communications Toolbox里其实有现成的comm.TurboEncoder和comm.TurboDecoder如果只想快速得到一条BER曲线调用它们可以在十分钟内搞定。但我不建议第一次接触Turbo码的人直接用黑盒模块。原因很简单Turbo码的性能高度依赖软信息的表达、交织器的对齐、迭代过程中外信息的更新方式这些都藏在黑盒内部。你要是不知道里面发生了什么遇到曲线不对的时候根本无从排查。手写实现看起来费时间实际上是在做“还债”的事。你写一次RSC编码循环就会明白状态转移是怎么回事你写一次LLR计算就会明白BPSK映射和噪声方差为什么能决定译码性能你写一次迭代循环就会明白为什么外信息不能直接喂回自己必须在交织域和自然域之间来回切换。而且这些代码以后要移植到C、Verilog或者做实时验证时都是现成的参考。所以这篇文章的方案是用最基础的MATLAB脚本从零实现核心链路不依赖专用通信工具箱。1.3 实验参数怎么定不贪大先求能跑通参数选择直接影响调试难度和曲线置信度。我在代码里常用的配置是信息比特长度K取1024分量码用约束长度3的RSC生成多项式写成(7,5)八进制码率取1/3不打孔调制方式用BPSK信道用AWGN迭代次数固定6次。这套配置的运算量不大K1024时跑几百帧也就几分钟的事而且交织长度足够长能明显看到迭代增益。为什么不选更大的KK越大交织增益越好这也是Turbo码长码性能逼近香农限的原因但K大之后每帧译码运算量线性增长调试时非常痛苦。K1024的错误统计已经够用在2dB附近跑几十帧就能看到BER降到千分之一量级。码率取1/3是因为每个信息比特附带两路校验输出结构最简单不用处理打孔表先跑通再谈码率适配。参数取值说明信息长度 K1024调试速度和交织增益的折中分量码RSC(2,1,2)生成多项式(7,5)状态数4编码和译码都简单码率 R1/3输出系统位校验1校验2调制/信道BPSK / AWGN软信息推导最简洁迭代次数6该参数下4-6次已足够收敛2. Turbo编码器RSC、交织器与整体组合2.1 RSC编码器的手写结构RSC的全称是递归系统卷积码关键有三个字递归、系统、卷积。“系统”是指输出中直接包含信息比特“递归”是指编码器内部状态会把输出反馈回输入端。光说概念比较虚我写一个可以直接跑的最小版本。约束长度3的RSC(7,5)编码器内部只有两级移位寄存器反馈多项式是7八进制二进制111前向多项式是5二进制101。它的输入输出关系是反馈量等于两个寄存器值的异或送入寄存器的新值等于信息比特与反馈值的异或校验位等于当前反馈输入与第二级寄存器的异或。写成MATLAB函数核心循环如下function [sys, par] rsc_encode_75(u) % u: 单个分量编码器的二进制信息序列0/1 N length(u); sys zeros(1, N); par zeros(1, N); reg1 0; % 前一级寄存器 reg2 0; % 后一级寄存器 for k 1:N fb xor(reg1, reg2); % 反馈多项式 1DD^2 d xor(u(k), fb); % 递归输入 sys(k) u(k); % 系统输出 par(k) xor(d, reg2); % 校验输出对应 1D^2 reg2 reg1; % 移位 reg1 d; end end这个函数输出的par就是一路校验比特。如果你用poly2trellis(3,[7 5])在MATLAB里建立trellis概念上等价。自己手写的好处是可以加断点观察寄存器状态理解“递归”到底是怎么把过去的信息卷进当前输出里的。2.2 交织器设计为什么不能用简单分组交织交织器在Turbo码里的作用是把两个分量编码器看到的错误图样尽量打散。如果两路编码器面对同样的错误位置那么两路校验同时无能为力迭代译码也救不回来如果交织后两个分量码的低权重码字错位整体码字的最小距离就会被抬高误码率曲线不会出现地板效应。最省事的交织器是randperm(K)产生的随机交织索引每帧重新生成也行固定下来也行。我在仿真中喜欢先固定一组索引保证每次结果可复现。对于K1024简单随机交织已经能带来明显增益。想再进一步可以用S随机交织器它要求交织前后相距S以内的位置交织后仍然相距S以上。S越大越难生成一般取floor(sqrt(K/2))左右的量级。注意编码端用的交织索引idx和译码端使用的必须是完全同一组编码端把u(idx)交给第二路编码器译码端就必须用llr_sys(idx)作为第二路译码器的系统位输入。反交织索引可以通过一句代码生成inv_idx(idx) 1:K;这行代码的意思是反过来记录“交织索引里的第几位是从自然序的第几号来的”。调试时最容易错的就是这里交织后的外信息必须反交织回自然序才能作为第一路译码器的先验信息。2.3 Turbo编码器整体组装用上面的RSC函数组合成完整编码器输出三路序列原始系统位、第一路校验位、第二路校验位。注意第二个分量编码器的输入是交织后的序列但它的系统位输出我们不直接传输我们只传输原始信息比特和两路校验比特。function c turbo_encode(u, idx) [sys1, par1] rsc_encode_75(u); [~, par2] rsc_encode_75(u(idx)); % 丢弃交织后的系统位 c [sys1; par1; par2]; % 每列对应一个信息比特的编码输出 end这里有一个我必须说明的简化这个编码器没有做归零处理编码完信息序列后寄存器停在某个未知状态。在标准Turbo码里每一路分量编码器都要发送额外的尾比特把寄存器拉回零态以便译码时精确知道起始和结束状态。但仿真初期把结束状态视为“均匀分布”并不影响整体趋势。真到了要求精确性能的时候再加尾比特核心的RSC递推部分不用改。3. 信道软信息计算译码的一切都从这里开始3.1 BPSK映射与LLR推导Turbo译码器需要的是软信息不能先把接收符号硬判成0/1再交给译码器那样会把信道置信度全扔了。AWGN信道下BPSK调制最常用的映射是编码比特0对应发送符号1编码比特1对应发送符号-1。接收端得到r s n噪声n服从均值为0、方差为σ²的高斯分布。对数似然比LLR的定义是L(x) ln( P(x0|r) / P(x1|r) )在这个映射和对称信道假设下化简结果非常干净LLR 2 * r / σ²这个式子很重要重要到值得抄下来贴在显示器边上。它说明软信息不只是一个“带符号的判决值”还隐含着噪声方差的归一化。如果发射符号不是能量为1的±1公式里的系数要做相应修改如果映射极性反了所有LLR符号同步取反编码位输出就会变成全是错误。3.2 从Eb/N0正确换算噪声方差信噪比参数用Eb/N0而不是Es/N0是因为Turbo码有编码冗余每个信息比特的能量被摊到多个信道符号上。对码率R1/3的链路能量关系是Es R * Eb。如果我们把调制符号能量归一化Es1那么给定目标Eb/N0时EbN0 10^(EbN0dB / 10) N0 1 / (R * EbN0) σ² N0 / 2噪声方差为什么是N0/2因为AWGN的双边功率谱密度是N0/2匹配滤波之后实部噪声方差就是这个值。很多仿真曲线偏左或偏右根本原因不是编码器写错而是这里把码率因子漏了。我经常看到有人直接用sqrt(1/(2*EbN0))当噪声标准差这默认了码率R1对1/3码率的Turbo码来说整条曲线会偏移约10*lg(3)4.77dB完全没法看。信道模拟和LLR计算的MATLAB函数可以写成这样function [llr_sys, llr_p1, llr_p2] awgn_llr(c, EbN0dB, R) EbN0 10^(EbN0dB / 10); N0 1 / (R * EbN0); % Es能量归一化为1 sigma sqrt(N0 / 2); s 1 - 2 * c; % 0 - 1, 1 - -1 r s sigma * randn(size(s)); llr_all 2 * r / (sigma^2); llr_sys llr_all(1, :); % 系统位软信息 llr_p1 llr_all(2, :); % 第一路校验软信息 llr_p2 llr_all(3, :); % 第二路校验软信息 end3.3 系统位软信息如何分给两路译码器Turbo译码器有两个分量两个分量都需要系统位的软信息。第一路直接使用自然序的llr_sys第二路必须使用交织后的llr_sys(idx)。原因很直白第二路RSC编码器编码的就是交织后的信息序列所以第二路译码器“看到”的系统位顺序也应该是交织后的。这个细节几乎每一版实现都会踩一次。常见错误是在第二路译码器里忘了交织llr_sys结果两路译码器对不上号迭代再多轮也没有增益。4. Log-MAP译码器三张表递推出软信息4.1 前向、后向和分支度量到底在算什么MAP译码器的目标是给定完整接收序列计算每个信息比特的后验概率。直接枚举所有码字是不可行的但卷积码的网格结构允许我们用前向-后向递推把计算量压下来。网格图上的每个节点代表一个编码器状态每条边代表一个状态转移分支分支上写着输入比特和输出比特。对数域里定义三个量分支度量γ表示从状态s转移到状态s的概率对数包含信道软信息、码字输出和先验信息三部分前向度量α表示从网格起始状态走到当前状态的所有路径概率之和的对数后向度量β表示从当前状态走到网格终点状态的所有路径概率之和的对数。递推时α从前往后扫β从后往前扫。每一时刻的比特LLR等于所有输入为1的分支对应的αγβ做log-sum-exp减去所有输入为0的分支对应的αγβ。这个操作等价于把通过该比特的所有“合法码头路径”的概率都加起来然后比较两类分支谁更有优势。一个容易忽略的点用普通max代替log-sum-exp就是Max-Log-MAP性能损失大约0.3到0.5dB但计算量大幅下降。实际工程中Max-Log-MAP配合外信息缩放系数0.75是很常见的做法。我的建议是先实现普通max版本链路跑通后再改成精确Log-MAP看看性能提升是否和理论对得上。4.2 分支度量的具体写法分支度量不是拍脑袋写的它来自高斯信道下的条件概率展开。对一个分量译码器假设系统位软信息是Ls校验位软信息是Lp当前先验信息是La那么状态转移s→s对应的分支度量可以写成γ(s,s) u * La / 2 (Ls * cu Lp * cp) / 2其中u是该分支对应的输入比特cu和cp是该分支对应的系统位和校验位经过BPSK映射后的±1值。这个式子看起来简单但实现时要用trellis表查出“从状态s输入比特u后下一个状态是谁输出比特是什么”。手写RSC的网格表也可以自己列但用poly2trellis生成会更省事。译码器里预先提取分支信息trellis poly2trellis(3, [7 5]); numStates trellis.numStates; % nextState(s, u1) 表示当前状态s输入u后的下一状态 % outBits(s, u1) 表示当前状态s输入u后的输出比特码在十进制里 for s 1:numStates for u 0:1 idx u 1; nextState(s, idx) trellis.nextStates(s, idx); outBits(s, idx) trellis.outputs(s, idx); end end实际用的时候需要把outBits拆成系统位和校验位两个0/1值再映射为±1。如果trellis.outputs里两个输出比特的编码顺序和你的预期不一致拆位顺序调一下即可这正是自写链路时“自由度”所在也是一开始容易糊的地方。4.3 外信息是怎么从总LLR里剥出来的分量译码器最终算出的总LLR由三部分组成信道给出的系统位信息、上一轮另一路译码器提供的先验信息、本路新榨出来的外信息。用公式写就是L_tot Ls La Le所以新的外信息Le L_tot - Ls - La这里Ls是系统位的信道LLRLa是输入这个分量译码器的先验LLR。为什么必须把Ls和La减掉因为Turbo码迭代的核心理念是“只交换新信息”。如果直接把包含旧信息的总LLR传给另一路相当于拿同一份证据反复投票容易陷入自激误码率反而恶化。我的调试经验是每次迭代后看一眼Le的方差正常的Le应该随着迭代小幅度增大如果第一轮Le就大得离谱多半是先验和信道项没减干净。5. 迭代译码完整流程与仿真脚本5.1 双译码器怎么交替运行整个迭代译码器是一个循环每一轮迭代里两个分量译码器各运行一次。伪码逻辑如下初始化外信息 Le1 全零 for it 1:maxIter % 第一路译码器自然域 Le1 comp_decode(Ls, Lp1, Le1) % 第二路译码器交织域 % 输入的系统位、校验位、先验都要交织 Le2_inter comp_decode(Ls(idx), Lp2, Le1(idx)) % 第二路输出的外信息在交织域反交织回自然域 Le1 Le2_inter(inv_idx) end注意第二路的先验输入是Le1(idx)不是Le1。第一次迭代时Le1全零所以第二路等于在无先验的条件下运行但从第二次开始先验就是另一路刚刚产生的新信息。反交织这步如果写反或者把Le1(idx)误写成Le2_inter(idx)整个循环就会进入一种“驴唇不对马嘴”的震荡BER该降不降。分量译码器comp_decode的输入输出是输入系统LLR、校验LLR、先验LLR输出新的外信息。内部按前向、后向、总LLR的顺序计算最后执行Le L_tot - Ls - La。5.2 迭代次数和停止准则仿真阶段固定迭代次数最省事6次已经能展示Turbo码的主要增益。工程系统里不可能每帧都固定跑6次因为SNR高的时候可能2次迭代就收敛了继续迭代只是浪费功耗和时延。常用办法是加入CRC校验每轮迭代后对硬判决结果做CRC通过就提前退出。这样能显著降低平均迭代次数。不过要注意CRC本身是额外开销会给链路引入很小比例的漏检实际系统要权衡。我们在MATLAB里做教学仿真时固定迭代更简单画图对比迭代次数对性能的影响也更直观。5.3 主仿真脚本框架把上面所有片段串起来主循环的大致结构是这样K 1024; R 1/3; maxIter 6; numFrames 100; snrVals 0:0.5:2; rng(2026); idx randperm(K); inv_idx(idx) 1:K; BER zeros(size(snrVals)); for s 1:length(snrVals) errBits 0; totalBits 0; for f 1:numFrames u randi([0 1], 1, K); c turbo_encode(u, idx); [Ls, Lp1, Lp2] awgn_llr(c, snrVals(s), R); uhat turbo_decode(Ls, Lp1, Lp2, idx, inv_idx, maxIter); errBits errBits sum(uhat ~ u); totalBits totalBits K; end BER(s) errBits / totalBits; end semilogy(snrVals, BER, o-); grid on;这段代码足够在小规模下跑通。真正跑论文级曲线时每个SNR点至少收集几十个错误比特再停否则BER曲线末尾会抖动得非常厉害。可以改成while errBits 50的控制结构帧数上限再设一个值防止死循环。5.4 运算量评估每帧译码复杂度和K * maxIter * numStates * 2^m有关其中m是输入比特数这里为1状态数4所以复杂度约K*6*4*2很低。K1024时普通电脑跑100帧大概几秒到十几秒。向量化可以进一步提速但第一版不建议写复杂向量化代码循环逻辑可读性更高。等K上到4096或8192再去考虑把α和β递推改成矩阵操作、预计算分支度量表。6. 常见问题与排查技巧实录6.1 曲线出现平台迭代不收敛这是Turbo码仿真最常见的问题现象是BER下降到一定程度后不再下降或者干脆在0.5附近游荡。排查顺序我一般是这样先关掉先验把迭代次数设成1看能不能退化成普通卷积码的单次MAP译码结果。如果退化的结果都不对说明分量译码器本身就有问题别急着找迭代的毛病。接下来检查交织索引和反交织索引是否严格互逆。再检查第二路译码器是否用了Ls(idx)和Lp2以及外信息反交织后是不是送给了第一路。最后检查LLR符号极性如果发送映射是0→1而译码分支度量里拿反了所有软信息会反向错误率直接崩溃。6.2 第三路校验位顺序与打孔错位我写过一版带打孔的Turbo码码率从1/3提到1/2也就是两个信息比特里只用三个传输符号。打孔表本来很简单但调试时发现高SNR区域出现地板。事后发现是译码器里把第二路校验位接错了位置第一路和第二路的校验顺序在打孔后对不齐。遇到这类问题建议先去掉打孔用1/3码率跑通全部链路再叠加打孔。每加一个特性就要重新验证一次基本曲线。6.3 外信息尺度太大或太小在Max-Log-MAP实现里计算出的外信息往往偏乐观直接传给下一轮会让误码率曲线变差。解决办法是给外信息乘一个0.75左右的缩放因子这是LTE等实用系统里常用的经验值。如果你的曲线在高SNR下反而变差可以先试试在Le更新处乘0.75。如果是精确Log-MAP缩放问题没那么严重可以不乘。另一个尺度问题是系统位LLR和校验位LLR之间的比例关系。信道噪声方差要同时作用于Ls、Lp1、Lp2如果某一路在校验位里忘了除以σ²这一路相当于用了错误置信度迭代时会对另一路输出误导性外信息。6.4 状态度量溢出或初始化错误我用对数域实现基本不存在指数上溢问题但如果你用的是线性MAPα和β必须做归一化否则迭代几步后数值就爆掉。在对数域里α和β可以任意加一个常数而不影响LLR结果所以不需要每步归一化这也是Log-MAP的优势。初始化方面如果编码端没有归零处理β的结束状态不能只给零态赋值而应该所有状态等概率也就是β初始化为全0的对数均匀值。如果加了尾比特归零α起始只在零态有效β终止只在零态有效。两种模式混用最后的几个比特就会莫名其妙地错。6.5 常见问题速查表现象最可能原因检查方式BER一直约0.5交织/反交织方向反或LLR极性反打印每次迭代外信息方差是否增大曲线低SNR正常高SNR有平台打孔表错位、外信息过估计去掉打孔、加0.75缩放再试最后一个比特总是错编码器未归零但译码结束状态假设了零态检查β初始化的状态假设曲线偏移几dB噪声方差没乘码率R核对Eb/N0换算公式高SNR帧数太少导致BER抖动错误比特数不够设置累计错误数达到50再停6.6 一个小而有效的调试招数调试Turbo码时我最常用的一招是把中间变量可视化。比如第一轮迭代后把第一路和第二路的外信息画成直方图你会发现它们都近似高斯分布均值附近的点很多尾部有少量大值。如果外信息全是零或者全是同一个值说明链路根本没有传递有效信息如果外信息的符号和原始系统位完全一致说明两路译码器已经“串通”提前收敛到了错误路径。这种可视化检查比盯着一堆数值变量直观得多。另外先把信息位设成全1或者单个1跑通极简场景再上随机比特。极简场景下你可以手动算出编码输出和期望LLR快速定位是编码器的问题还是信道模型的问题。如果整个实现只留一个检查点我建议反复检查外信息交换顺序。Turbo码表面上是两个编码器、两个译码器的对称结构实际上因为交织器的存在两路译码器始终在“自然域”和“交织域”之间穿梭。我调试时习惯在第二轮迭代后打印第一路译码器输入先验的均值正常情况下它应该和上一次迭代输出的外信息均值保持一致如果出现跳变基本可以断定交织或反交织没对齐。另一个个人经验是先用Max-Log-MAP把链路跑通再切精确Log-MAP。很多人一上来就追求最优算法结果公式和Bug混在一起调一天都找不到问题。先跑出下降趋势再优化0.3dB才是效率最高的路径。Turbo码这套东西编码器和译码器的代码规模并不大真正的门槛在于对“软信息如何流动”的理解。一旦你在MATLAB里把它跑通后面再看5G的极化码、LDPC的迭代译码思路都会顺畅很多。
返回列表