ARTICLE DETAIL

资讯详情

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

LDPC MinSum译码算法仿真:从校验矩阵到误码率验证的完整链路

LDPC MinSum译码算法仿真:从校验矩阵到误码率验证的完整链路 简介这是一份面向通信工程与信息论方向学习者的LDPC编译码仿真资源基于MinSum最小和译码算法实现适合初学者理解LDPC码构造、迭代译码及误码性能评估流程。压缩包共13个文件以11个MATLAB脚本.m为主另含2个自动保存备份文件.asv整体仅9KB代码结构轻量。仿真工程涵盖校验矩阵生成、信息比特重排、BPSK调制、AWGN信道模拟、MinSum译码、误码率统计等环节可帮助读者将LDPC编码理论转化为可运行的实验代码并可调整码率、码长、迭代次数和信噪比等参数观察不同条件下的纠错表现。工程内函数分工明确从编码、加噪到译码、统计串成完整链路便于按模块逐段学习通过对比不同信噪比下的BER曲线也能直观感受MinSum算法在复杂度与性能之间的折中。目前已有395人学习适合正在学习信道编码或需要快速搭建LDPC仿真的开发者参考。1. 为什么 MinSum 译码算法是 LDPC 解码器落地的首选在信道解码器的工程选型里LDPC 码的译码时常不是那个性能最强的算法而是最贴近物理实现的那个。MinSum最小和译码算法把置信传播里的非线性 tanh 运算替换成取最小值和符号相乘硬件开销小了一个数量级代价是 0.1~0.3dB 的增益损失。这套 LDPC_encode_decode_minsum_sim 仿真工程正好把一个完整的链路拆开校验矩阵构造、GF(2) 编码、BPSK 调制、AWGN 信道、MinSum 迭代译码和误码率统计。适合刚接触 LDPC 译码的工程师、写通信算法验证脚本的研究生以及想在 FPGA 上实现 MinSum 之前先用 MATLAB 摸清边界的人。2. 校验矩阵构造从随机稀疏矩阵到可用的 LDPC 码LDPC 编码的前提是拿到一个性能达标的校验矩阵 H。H 的构造方式直接决定译码瀑布区的位置和错误平层仿真的第一步也在这里。2.1 为什么 H 要“稀疏且均匀”LDPC 的全称是低密度奇偶校验密度指的是校验矩阵中“1”的比例。如果 H 里每个行列都塞满 1译码时校验方程的约束过多消息传递会出现反馈环路迭代很难收敛。工程上常用规则码或准循环码让每行有相同的“1”个数行重每列也有相同的“1”个数列重。行重决定了校验节点的度列重决定了变量节点的度。从 Tanner 图角度看变量节点和校验节点之间的连线就是 H 中 1 的位置。连线越少迭代时外信息相关度越低误码性能越好。但同时列重不能太小否则最小距离不理想纠错能力不足。常规选择是列重 3、行重 6码率约 1/2。genH.m 里通常就是按这个原则生成随机稀疏矩阵的。生成时有一个关键限制不能出现长度为 4 的短环。Tanner 图中如果两个变量节点被同一个校验节点对连接信息会在 4 步之内回到自身导致迭代增益被抵消。所以 genH.m 在随机放完 1 之后会做一次 girth 检查发现 4-cycle 就重新扰动该位置。2.2 genH.m 的典型生成流程在 MATLAB 里我不会用无脑 rand 打散而是按列重和行重约束逐列填充再用贪心方式消除 4-cyclefunction H genH(M, N, dc, dv) % M: 校验方程数行数 N: 码长列数 % dc: 行重 dv: 列重这里按规则码处理 H zeros(M, N); col_w zeros(1, N); % 当前列重 max_attempts 200; % 按列放置 1尽量让每列权重接近 dv for col 1:N placed 0; while placed dv rows randperm(M, dv - placed); for r rows if col_w(col) dv, break; end if sum(H(r, :)) dc, continue; end % 检查候选位置是否会形成 4-cycle candidate H; candidate(r, col) 1; if has4Cycle(candidate), continue; end H(r, col) 1; col_w(col) col_w(col) 1; placed placed 1; end end end end function flag has4Cycle(A) % 只要存在两行在某两列同时为 1就构成 4-cycle rows size(A, 1); for r1 1:rows-1 r1pos find(A(r1, :)); for r2 r11:rows r2pos find(A(r2, :)); if length(intersect(r1pos, r2pos)) 2 flag true; return; end end end flag false; end这个脚本的逻辑是先给每一列随机分配“1”的位置分配时检查行重不能超过 dc同时用 has4Cycle 排除两个校验节点在相同两个变量节点上同时出现的情况。随机尝试很多次之后如果仍然失败会把矩阵重新打散再试直到生成满足 girth 6 的 H。has4Cycle 的实现虽然暴力但作为离线仿真已经够用。真正大规模构造时应该用渐进边增长PEG算法按列逐个扩展 Tanner 图保证最小环长最大化。做工程验证时可以先用随机构造后面换成 QC-LDPC 或 PEG 矩阵译码端完全不用改。2.3 H2P.m 做了什么事LDPC 译码时我们需要把 H 从“二维矩阵”变成“校验节点到变量节点的连接表”。MATLAB 的稀疏矩阵读写效率高但迭代译码时反复用 find 查找效率很低。H2P.m 的作用就是把稀疏的 H 转换成校验节点索引列表例如一个 cell 数组cn_list{row}保存第 row 个校验方程关联的所有变量节点下标。这样译码循环里直接索引不用每次扫描整行。转换逻辑很直接function cn_list H2P(H) M size(H, 1); cn_list cell(M, 1); for r 1:M cn_list{r} find(H(r, :)); end end如果你拿到这套代码可以把这个函数看成解码器的“数据管道入口”。后面 MinSum 更新的所有消息传递都基于这个列表。这里有个容易踩的坑如果原始 H 中某些行重不平衡cn_list 中会出现空行或超长行译码时需要做长度均衡否则 MATLAB 的隐式扩展会报维度不匹配。3. 编码链路与 GF(2) 运算怎么把信息比特变成信道符号校验矩阵构造好之后编码不是简单把信息比特直接发给信道而是要找到满足 H 约束的码字。这一章的实现对应压缩包里的 ldpc_encode.m、mul_GF2.m、bpsk.m。3.1 从 H 到生成矩阵 G对于线性分组码码字 c 满足 H c^T 0。编码的实质是求解 H 的零空间。常见做法是把 H 通过高斯消元变成系统形式 H [P | I]于是生成矩阵 G [I | P^T]。信息比特 u 编码成 c u G。这里的关键是消元过程在 GF(2) 域上完成MATLAB 的常规加减要改成异或。mul_GF2.m 就是封装了 GF(2) 上的乘加运算function x mul_GF2(a, b) % GF(2) 矩阵相乘 % 输入可以是矩阵输出为模 2 乘法结果 x mod(a * b, 2); end注意这个函数自己并不过滤维度错误调用前需要保证前一个矩阵的列数等于后一个矩阵的行数。实际编码时还需要对变量做换列操作因为高斯消元会打乱列顺序。reorder_bits.m 就是记录这个列置换并在译码完成后把信息比特还原回原来的顺序。3.2 编码并映射为 BPSK 符号ldpc_encode.m 中处理的是“信息位 u 长度 K生成矩阵 G 大小 K×N”乘完后得到码字 c。由于 LDPC 码字往往不完整恰好等于信息位加校验位很多实验中会先截断或者打孔。本工程里没有打孔所以 c 的长度直接就是 N。随后 bpsk.m 把 0/1 码字映射成信道符号function x bpsk(bits) % BPSK 调制0 - 1, 1 - -1 x 1 - 2 * bits(:).; end注意映射时符号方向。常用约定是比特 0 映射为 1比特 1 映射为 -1。MinSum 译码里的信道对数似然比LLR初始化用的是2y/σ^2这里如果映射反了译码结果整体取反误码率会接近 1。经过 AWGN 信道后接收向量为y x noise。噪声方差由信噪比 Eb/N0 换算得到这一步通常放在 ldpc_demo.m 里统一处理。调试时最容易犯的错误是把符号能量当 1 用实际 BPSK 每个符号携带 1 个比特Es/N0 和 Eb/N0 在数值上相等但如果是其他调制就需要注意。3.3 编码前的矩阵编号问题很多第一次看 LDPC 仿真代码的人会被 reorder_bits.m 和 extract_mesg.m 弄晕。它们是为了解决同一个问题编码时的 G 是在 H 经过列置换之后才得到的信息位在码字中的物理位置可能不再是前 K 个比特。extract_mesg.m 的典型做法是根据列置换表把译码结果的对应位置取出来。如果你自己写仿真我建议在编码侧就把置换关系保存成全局变量否则解码端还原信息位时很容错位。以下是我常用的保存方式% 高斯消元前记录列置换顺序 [G, col_perm] makeSystematicG(H); % 编码时先对信息位做列置换或直接乘 G c mod(u * G, 2);到这里编码链路基本闭合接下来就是本工程的核心MinSum 译码。4. MinSum 译码算法从 BP 到最小和的代码级拆解LDPC 译码的基准算法是置信传播BP也叫和积算法Sum-Product。它迭代更新两类消息变量节点发给校验节点的消息以及校验节点回传给变量节点的消息。4.1 为什么用最小和替换 tanhBP 的校验节点更新公式在 LLR 域可以写成L(r_ji) 2 * atanh( ∏(i ∈ N(j) \ i) tanh( L(q_ij) / 2 ) )这个公式在硬件里实现要打三张表tanh 表、atanh 表、乘法器。每个校验节点每轮迭代都要把相邻变量节点的消息全部乘一遍资源开销大而且对于高行重的码字tanh 表深度不够还会引入量化误差。MinSum 的近似来自观测到atanh(tanh(a) * tanh(b)) ≈ sign(a) * sign(b) * min(|a|, |b|)。所以校验节点更新简化为L(r_ji) ( Π(i ≠ i) sign(L(q_ij)) ) * min(i ≠ i) |L(q_ij)|这一下把乘法、查表全部换成比较器和异或门FPGA 上触发器数量大幅下降。代价是分子消息的绝对值偏大误码性能略有衰减。工程上可以在最小值乘一个小于 1 的归一化因子来补偿通常取 0.75这就是“归一化 MinSum”。4.2 ldpc_decode.m 的消息更新流程伪代码描述整个译码过程如下初始化变量节点消息为信道 LLRL(q_j) 2*y_j / sigma^2进入迭代循环每次迭代执行校验节点更新对每个校验方程计算关联变量节点消息的符号积和最小绝对值变量节点更新本节点信道值加所有传入校验消息之和计算后验信息做硬判决如果校验方程 H*c^T 0 全部满足提前退出迭代MATLAB 实现可以写成function decoded ldpc_decode(y, H, cn_list, max_iter, beta) [M, N] size(H); v2c zeros(M, N); % 变量节点到校验节点消息 c2v zeros(M, N); % 校验节点到变量节点消息 llr_ch 2 * y(:) / (2 * sigma2); % 具体噪声参数外部传入 llr_ch y; % 占位实际按 SNR 换算 for iter 1:max_iter % 校验节点更新利用 cn_list 填充 c2v for j 1:M idx cn_list{j}; signs sign(v2c(j, idx)); abs_msg abs(v2c(j, idx)); [min1, min1pos] min(abs_msg); % 找次小值 temp abs_msg; temp(min1pos) inf; min2 min(temp); prod_sign prod(signs); for i idx % 排除当前节点自身 if i min1pos c2v(j, i) beta * prod_sign * sign(v2c(j, i)) * min2; else c2v(j, i) beta * prod_sign * sign(v2c(j, i)) * min1; end end end % 变量节点更新外信息等于总校验消息减去本节点回传消息 for i 1:N col_idx find(H(:, i)); v2c(col_idx, i) llr_ch(i) sum(c2v(col_idx, i)) - c2v(col_idx, i); end % 后验信息与硬判决 post llr_ch sum(c2v, 1); hard post 0; if mod(H * hard(:), 2) 0 break; end end decoded hard(:); end上面这段为了演示结构做了一定简化实际用时会发现两个工程细节第一校验节点更新不能对本身回传的消息求外信息必须排除当前变量节点。上面用min1pos记录全局最小值的位置如果当前变量节点就是全局最小值所在则使用次小值 min2否则使用 min1。忽略这一点会导致自反馈迭代次数多之后错误会扩散。第二变量节点更新时当前节点的外部消息等于总校验消息减去从该节点收到的消息。这里c2v(col_idx, i)是第 i 个变量节点从周围校验节点收到的消息集合sum(c2v(col_idx, i)) - c2v(col_idx, i)表示去掉自身消息后的外信息总和。如果直接用llr_ch sum(c2v)而不减去自身后验信息会包含与决策相关的先验信息造成正反馈。4.3 迭代终止与停止条件ldpc_decode.m 里会在硬判决后判断mod(H*hard, 2) 0全部校验方程满足就退出。这个条件在低信噪比时非常苛刻可能几十次迭代也满足不了。实际工程常用的是双条件迭代次数达到上限或校验方程不满足的个数小于某个阈值。阈值的选择与码长有关。我一般用“残留校验和”作为停止条件即sum(mod(H*hard, 2))。如果连续两次迭代的残留校验和不下降说明进入错误平层再迭代也没有意义可以提前停止。这也是 MinSum 调试时最容易观察到的收敛信号。4.4 参数对性能的影响下面表格总结了 ldpc_decode.m 中几个核心参数的作用和调节方向参数作用调节方向max_iter最大迭代次数决定译码时延和纠错上限低 SNR 区增加FPGA 上常限 8~12 次beta归一化因子补偿 MinSum 的近似误差取 0.5~0.9通常 0.75 起步再扫sigma2信道噪声方差决定初始 LLR 可信度必须与仿真 SNR 匹配否则误码率曲线偏移最小值/次小值校验节点消息绝对值的信任程度不准确时容易在小迭代次数下产生突发错误用这套代码做实验可以把 beta 当成一个扫参变量观察 BER 曲线在 0.01 量级附近的变化。beta 太大会让消息过度自信高 SNR 时反而性能差beta 太小则收敛慢同样迭代次数达不到预期纠错能力。5. 把 MinSum 放在完整链路里调参仿真脚本与 BER 验证有了前面三章的函数ldpc_demo.m 像一个总装车间把 genH、ldpc_encode、bpsk、MinSum、BER 统计衔接起来。这个章节同时起到收束工程链路和给出调试技巧的作用。5.1 主仿真脚本的最小骨架一个能跑通并输出 BER 曲线的最小骨架如下clear; clc; N 576; K 288; % 码长信息位长度码率 1/2 dc 6; dv 3; H genH(N - K, N, dc, dv); cn_list H2P(H); n_iter [1 3 5 8 12 20]; % 扫迭代次数 sigma2s [0.5 0.7 0.9 1.1]; for sigma2 sigma2s ber zeros(size(n_iter)); for idx 1:length(n_iter) max_iter n_iter(idx); errors 0; total 0; for trial 1:200 u randi([0 1], 1, K); c ldpc_encode(u, H); % 内部已转为系统码 x bpsk(c); y x sqrt(sigma2) * randn(size(x)); decoded ldpc_decode(y, H, cn_list, max_iter); decoded extract_mesg(decoded, H); % 还原信息比特 errors errors sum(decoded ~ u); total total K; end ber(idx) errors / total; end semilogy(n_iter, ber, o-); hold on; end上面脚本里的 sigma2 设置是简化的严格做法要按 Eb/N0 换算。AWGN 信道下 BPSK 的实际噪声方差为sigma2 1 / (2 * rate * 10^(EbN0/10))其中 rate K/N。很多初学者把 rate 忘掉导致曲线横坐标直接右移 3dB。5.2 校验 MinSum 实现是否正确的三条经验第一用全零码字做信道初始化。由于编码后全零码字仍是有效码字接收端 y 全为正译码输出应全为 0。如果出现大量 1说明 LLR 符号约定或调制映射反了。第二对比固定迭代次数。不要一上来就扫 SNR先固定 Eb/N0 2.5dB把 max_iter 从 1 扫到 30。理论上误码率曲线应单调下降如果出现震荡多半是码本 H 有 4-cycle或者校验节点更新时没有排除自身消息。第三把全零码字的校验和打出来。在译码每轮迭代后打印sum(mod(H * hard(:), 2))。如果这个数以“5 → 3 → 1 → 0”的方式递减说明消息更新正常如果出现“3 → 5 → 2”。5.3 从仿真到硬件实现时可以带走什么这套 MATLAB 工程的编排方式本身就是一个可移植的译码框架。如果你要在 FPGA 上做归一化 MinSum可以把 cn_list 换成按校验节点深度排列的存储结构把校验节点更新写成流水线最小值/次小值比较改为双路寄存器。消息量化位数建议 6 比特左右beta 用 0.75 时量化误差还能接受若追求更高精度就改成偏移 MinSum把beta*min换成max(min - offset, 0)固定点实现更稳定。另外在仿真链路上做一个有用的收敛性检查把每条 SNR 点的平均迭代次数一并统计出来。平均迭代次数低于最高迭代次数一半时说明 max_iter 余量足够如果平均迭代次数几乎顶满就说明码率或 H 构造需要调整否则在硬件实时性要求下会频繁超时。把 BER 和平均迭代次数两张图画在一起比单独看一条曲线更能判断这个码是否适合实际信道的时延预算。本文还有配套的精品资源点击获取
返回列表