简介:QC-LDPC.zip是一套面向通信工程和量子信息初学者的准循环低密度奇偶校验码MATLAB仿真程序包。资源共3个文件、大小仅4KB,主体为两个m脚本和一个asv自动备份文件:QCLDPCBaseH.m负责构造具有稀疏性和准循环结构的校验矩阵,QCEncode.m完成从信息序列到码字的编码过程,asv为编辑器备份可忽略。程序覆盖码字生成、编码、信道加噪、迭代解码的完整链路,并对贝叶斯消息传递算法等关键步骤提供逐步注释,便于对照教材理解抽象原理。已有1254人学习下载。学习者可直接运行观察编解码输出,也可调整校验矩阵尺寸、信道信噪比等参数,分析不同设置下纠错性能差异,将理论中的QC-LDPC构造思想落地为可操作仿真实验,适合作为课程设计或科研入门的参考实现。 做通信仿真的人应该都有过这种经历:为了验证一个低密度奇偶校验码(LDPC)的性能,自己搭随机校验矩阵,再写置信传播译码,结果仿真跑了一整夜,误码率曲线还是像过山车一样抖。问题往往不在译码算法,而在校验矩阵的构造。直到我拿到 QC-LDPC 这套 MATLAB 工程,才发现准循环 LDPC 码把矩阵生成、编码和译码仿真拆得清清楚楚,基矩阵加扩展因子的设计让随机矩阵的“黑匣子”变成了可计算的参数。这是一个适合通信专业学生、算法工程师和预研硬件实现的人反复推敲的仿真包,它能帮你搞懂准循环校验矩阵怎么构造、BP 译码怎么收敛、误码率曲线怎么才算正常。
1. 这套 QC-LDPC 的 MATLAB 仿真包,为什么值得跑一遍
做信道编码的仿真,最怕的不是代码跑不通,而是跑通了却不知道结果对不对。QC-LDPC 这个名字里的 “QC” 指准循环结构,它用一块很小的基矩阵加上一个扩展因子 Z 就能描述完整校验矩阵,相比随机 LDPC 的巨型稀疏矩阵,存储开销和构造复杂度都降了一个量级。这套 MATLAB 工程就是围绕这个特性写的:从基矩阵展开校验矩阵开始,到系统编码、BP 译码、误码率仿真,最后画出和理论值可对比的 BER 曲线。如果你是刚啃完 LDPC 原理、想找一份能跑的代码来对照理论,或者在做 5G 类应用的预研,想评估不同码率和码长下的性能,这套代码能直接当试验台用。它的核心价值不是把库里的函数堆给你,而是让你能一步步盯住每个模块的输入输出,知道哪个参数动了会影响什么。
2. 设计思路与目录结构:准循环结构如何落到 MATLAB 代码
2.1 从随机 LDPC 到 QC-LDPC:为什么必须引入基矩阵和扩展因子
随机 LDPC 码的校验矩阵 H 直接用稀疏矩阵生成,行列重量靠概率约束,矩阵越到后面越难保证没有短环,而且存一个几千行乘几千列的稀疏矩阵,在 MATLAB 里操作一次就要吃掉几百 MB 内存。QC-LDPC 的做法完全不同:定义一块 mb 行 nb 列的基矩阵(Base Matrix),里面每个元素只表示循环移位位数,-1 代表全零子矩阵,非负数代表单位阵循环右移的步数。实际校验矩阵 H 的尺寸是 (mb × Z) 行、(nb × Z) 列,其中 Z 是扩展因子。矩阵规模变大,但描述它只需要 mb × nb 个整数,这就是准循环结构最实在的好处。
从构造角度说,QC-LDPC 的每个子矩阵是单位阵的循环位移,天然保证了每一行每一列重量相同,而且只要基矩阵设计时避开特定模式,展开后就不容易出现 4 环。相比随机构造靠随机反复重试,QC-LDPC 是把“试出来”变成“算出来”。在 MATLAB 里,这种结构用矩阵索引加 circshift 就能快速生成,不需要像 C 语言那样逐个比特操作循环移位寄存器。不过用 MATLAB 写生成函数时,最笨的办法反而最好懂:逐个把子矩阵填进 H 的对应区块,再用 sparse() 转成稀疏矩阵。
2.2 工程目录结构与每个文件的作用
拿到 QC-LDPC.zip 解压后,目录结构是这样的:
QC-LDPC/ ├── src/ │ ├── base_matrix_define.m % 定义基矩阵,返回 mb×nb 的整数矩阵 │ ├── expand_matrix.m % 基矩阵 + Z 扩展成完整校验矩阵 H │ ├── encode_qc.m % 编码器,输出系统码码字 │ ├── decode_bp.m % 对数域置信传播译码器 │ ├── decode_min_sum.m % Min-Sum 近似译码器(二次开发用) │ └── run_simulation.m % 主仿真脚本,遍历多个 Eb/N0 画 BER ├── tests/ │ ├── test_h_matrix.m % 检查 H 每行每列重量、是否有环 │ └── test_encode_decode.m % 单帧编解码自检 └── README.md这就是一个典型的研究型工程布局:src 放核心算法,tests 放验证脚本,run_simulation.m 作为总入口。我拆过不少类似项目,这个结构的优点是把“算法实现”和“仿真入口”分开了。你在改译码算法时不需要动仿真脚本,仿真的参数也不用来回翻代码,都在脚本头部集中定义。
2.3 基矩阵设计:码率、码长和循环移位数的关系
基矩阵定义了码率。假设基矩阵是 mb 行 nb 列,扩展因子是 Z,那么完整校验矩阵的行数 M = mb × Z,列数 N = nb × Z,码率 R = 1 - mb/nb。信息位长度 K = N - M = (nb - mb) × Z。比如一块 3 行 6 列的基矩阵,扩展因子取 128,那么 N = 768,M = 384,K = 384,码率就是 0.5。这是 LDPC 里最常见的码率结构。
基矩阵里的数值大小决定了循环移位步长,取值必须在 0 到 Z-1 之间。如果基矩阵里有个元素是 15,而 Z 只有 8,展开时子矩阵就会越界,这是新手最容易翻车的地方。常见做法是基矩阵元素最大不超过 Z-1,或者干脆让 Z 取固定值然后顺着这个值去设计基矩阵。下面的 3×6 基矩阵是个标准示例,前 3 列对应系统位,后 3 列对应校验位:
% 基矩阵示例:3×6,每行重量 4,每列重量 2 base_matrix = [ 0 1 2 -1 0 1; 1 0 -1 2 1 0; 2 -1 1 0 2 1; ];每一行的 -1 表示这个位置不放置子矩阵,是合法的“零块”。在 QC-LDPC 里,校验位对应的子矩阵尽量构造成下三角或近似下三角,这样编码时可以避免算大矩阵求逆。
3. 核心模块拆解:从校验矩阵生成到编码器译码器的完整链路
3.1 校验矩阵展开:用 sparse 函数生成可用的 H
这一节我们看 expand_matrix.m 的核心代码。它的任务非常明确:输入基矩阵和扩展因子,输出一个稀疏校验矩阵 H。我的做法是先用 for 循环逐块填充,最后转成 sparse 格式,代码可读性强,也方便你打印部分矩阵去和理论值对照。
function H = expand_matrix(base_matrix, Z) % 输入: base_matrix - mb×nb 整数矩阵,-1 表示零块 % Z - 扩展因子,正整数 % 输出: H - (mb*Z) 行 × (nb*Z) 列稀疏矩阵 [mb, nb] = size(base_matrix); H = zeros(mb * Z, nb * Z); for i = 1:mb for j = 1:nb s = base_matrix(i, j); if s < 0 continue; % 全零子矩阵,不填充 end % 生成 Z×Z 的单位阵循环移位矩阵 I = eye(Z); sub_matrix = circshift(I, [0, s]); % 向右循环移位 s 位 % 将子矩阵放入 H 的对应区块 row_idx = (i-1)*Z + 1 : i*Z; col_idx = (j-1)*Z + 1 : j*Z; H(row_idx, col_idx) = sub_matrix; end end H = sparse(H); end这段代码的逻辑很直白:每个非 -1 元素都对应一个循环移位单位阵,circshift(I, [0, s]) 把单位阵的每一列向右循环移动 s 位。比如 s = 0 就是单位阵本身,s = 1 时第一列跑到最后一列,形成右移效果。填充完成后再用 sparse() 压缩存储,后续译码时 find(H) 拿到非零元素的行列坐标,能避免双重循环遍历整个大矩阵。
参数上要注意:circshift 的第二个参数第一个元素是行移位,第二个是列移位,写成 [0, s] 表示行不动、列右移。如果你想做下移,就写 [s, 0]。很多新手在这地方搞反方向,导致校验矩阵的行重列重乱了,译码器却照常运行,最后误码率曲线莫名其妙高出一截。
3.2 编码器实现:用近似下三角结构避免大矩阵求逆
拿到 H 矩阵之后,最简单的编码办法是求生成矩阵 G,但 G 的尺寸是 K × N,码长一大就是几万比特,MATLAB 里直接做二进制高斯消元非常慢,而且生成的 G 是稠密矩阵,存储和运算都吃亏。QC-LDPC 的工程实现更常用近似下三角编码:把 H 的列重排成 [A B; C D] 的下三角形式,利用循环移位结构分块求解校验位。
为了让你理解编码流程,我写一个基于下三角结构的编码函数,假设基矩阵已经设计成校验部分为标准下三角(即每个校验列除了对角线外没有别的校验列连接):
function codeword = encode_qc(info, base_matrix, Z) % 输入: info - 1×K 二进制信息向量 % base_matrix - mb×nb 基矩阵 % Z - 扩展因子 % 输出: codeword - 1×N 二进制码字 [mb, nb] = size(base_matrix); K = (nb - mb) * Z; N = nb * Z; % 初始化码字,前 K 位是信息位 codeword = zeros(1, N); codeword(1:K) = info; % 校验部分按基矩阵最后 mb 列逐块求解 % 每个校验块 P_j 长度为 Z,通过线性方程求解 % 流程:从第 1 个校验块到第 mb 个校验块,顺次计算 for p = 1:mb % 累加与当前校验节点连接的信息位 + 已计算校验位 acc = zeros(1, Z); for j = 1:nb s = base_matrix(p, j); if s < 0 continue; end % 取得对应码字段,按循环移位交换 segment = codeword((j-1)*Z + 1 : j*Z); shifted = circshift(segment, [0, -s]); % 向左移 s 位取模 acc = mod(acc + shifted, 2); end % 该块校验位就是累加和(标准下三角时) codeword(K + (p-1)*Z + 1 : K + p*Z) = acc; end end这里的核心思想:校验矩阵的行代表校验方程,每个校验节点要满足 mod(H * codeword', 2) = 0。标准下三角结构让每个新的校验块只依赖之前的信息位和之前的校验块,所以可以按顺序一块一块解出来,不需要求逆。代码中的 acc 循环累加做的是模二和,circshift 的 -s 左移对应展开时 circshift(I, [0, s]) 的逆操作,保证方程成立。
实际工程里,基矩阵往往不是严格下三角,而只是近似下三角,这时需要多几步行列交换和回代。不过对于学习用途,标准下三角已经足够跑通仿真。如果你是做真实系统,我一般建议直接用高斯消元把校验位部分化为下三角,然后再用这个顺序求解,速度比求 G 矩阵快很多。
3.3 BP 译码器:对数似然比消息传递的核心循环
译码器是整套代码里最容易出逻辑错误的地方,但也是最能体现 LDPC 精髓的部分。BP 译码的核心是变量节点和校验节点之间互相传递对数似然比消息。我习惯先写最原始的双重循环版本,保证和公式一一对应,等跑通之后再优化向量化。
function decoded = decode_bp(channel_llr, H, max_iter) % 输入: channel_llr - 1×N 信道初始对数似然比 % H - sparse 校验矩阵 % max_iter - 最大迭代次数 % 输出: decoded - 1×N 判决码字 [M, N] = size(H); [row, col] = find(H); % 初始化消息矩阵,只存非零位置 V2C = zeros(M, N); C2V = zeros(M, N); % 后验LLR初始化为信道LLR LLR_post = channel_llr; for iter = 1:max_iter % ---- 校验节点更新 ---- % 对每个校验节点 i,收集所有相邻变量节点的消息 for i = 1:M edges_i = col(row == i); % 与校验节点 i 相连的变量节点集合 for idx = 1:length(edges_i) j = edges_i(idx); product = 1; for idx2 = 1:length(edges_i) k = edges_i(idx2); if k ~= j % 使用 tanh 规则相乘 product = product * tanh(V2C(i, k) / 2); end end % 防止数值下溢 if abs(product) < 1e-15 product = sign(product) * 1e-15; end C2V(i, j) = 2 * atanh(product); end end % ---- 变量节点更新 ---- for j = 1:N edges_j = row(col == j); % 与变量节点 j 相连的校验节点集合 for idx = 1:length(edges_j) i = edges_j(idx); sum_msg = channel_llr(j); for idx2 = 1:length(edges_j) k = edges_j(idx2); if k ~= i sum_msg = sum_msg + C2V(k, j); end end V2C(i, j) = sum_msg; end end % ---- 后验LLR与判决 ---- LLR_post = channel_llr + sum(C2V, 1); decoded = double(LLR_post < 0); % LLR<0 判定为 1 % ---- 提前停止:校验子全零则退出 ---- syndrome = mod(decoded * H', 2); if sum(syndrome) == 0 break; end end end这段代码的逻辑说明:校验节点更新用的是 tanh 规则,它等价于概率域乘法的对数表示。每个校验节点给某个变量节点传回的消息,是把除了那个变量节点之外的所有变量消息做乘积再映射回 LLR 域。变量节点更新则简单得多:把信道 LLR 加上所有校验节点传回的消息,唯一要注意的是不能包含自己发出的那一条,不然消息会自我加强,导致译码误收敛。
需要特别说明的是提前停止条件:每轮判决后把 decoded 和 H 相乘,如果得到全零向量,说明这个码字满足所有校验方程,译码成功,可以马上跳出循环。对于没有成功收敛的帧,它会一直迭代到 max_iter。这个条件写错的人很多,有人用成 “若迭代次数超过 20 就强制停止”,那会导致曲线在低信噪比时突然掉下去,实际上只是把不收敛的帧当成了错误帧算进了误码率。
4. 运行仿真:从主脚本到误码率曲线
4.1 设计主仿真脚本:参数集中定义、循环遍历 Eb/N0
主脚本 run_simulation.m 是整个工程的发动机。它把基矩阵、扩展因子、信噪比范围、每个信噪比下的仿真帧数全部放在文件开头,方便你改完参数直接跑。下面这段代码我已经按 BPSK+AWGN 信道的标准模型写好,噪声方差的计算和 LLR 映射都用最不容易出错的形式。
% run_simulation.m clear; clc; % ---- 参数区 ---- Z = 128; base_matrix = [ 0 1 2 -1 0 1; 1 0 -1 2 1 0; 2 -1 1 0 2 1; ]; max_iter = 50; EbN0_dB = 0:0.5:3; num_frames = 100; % 每个 SNR 点仿真帧数 % ---- 初始化系统 ---- H = expand_matrix(base_matrix, Z); [mb, nb] = size(base_matrix); M = mb * Z; N = nb * Z; K = N - M; R = K / N; fprintf('码率: %.3f, 码长: %d, 信息位: %d\n', R, N, K); % ---- 遍历信噪比 ---- ber_results = zeros(size(EbN0_dB)); for snr_i = 1:length(EbN0_dB) eb_n0_lin = 10^(EbN0_dB(snr_i) / 10); noise_var = 1 / (2 * R * eb_n0_lin); % BPSK 下 Es=1 total_err = 0; total_bits = 0; for frame = 1:num_frames % 随机信息位 info = randi([0 1], 1, K); % 编码 code = encode_qc(info, base_matrix, Z); % BPSK 映射: 0 -> -1, 1 -> +1 x = 2 * code - 1; % AWGN 噪声 noise = sqrt(noise_var) * randn(1, N); y = x + noise; % 计算对数似然比 LLR = 2*y/noise_var llr = 2 * y / noise_var; % BP 译码 decoded = decode_bp(llr, H, max_iter); % 统计错误比特(只统计信息位或全部码字都可) err = sum(decoded ~= code); total_err = total_err + err; total_bits = total_bits + N; end ber = total_err / total_bits; ber_results(snr_i) = ber; fprintf('Eb/N0 = %.2f dB, BER = %.2e\n', EbN0_dB(snr_i), ber); end % ---- 画图 ---- semilogy(EbN0_dB, ber_results, 'o-'); grid on; xlabel('Eb/N0 (dB)'); ylabel('BER'); title('QC-LDPC 码性能仿真');这个脚本里有几个细节值得强调。第一,noise_var = 1/(2Reb_n0_lin) 是因为 BPSK 符号能量 Es=1,而 Eb=Es/R,噪声单边功率密度 N0=2noise_var,所以 Eb/N0 = 1/(2Rnoise_var)。如果你用的是 4QAM 或高阶调制,这个公式就要改。第二,LLR 计算用 2y/noise_var,在 AWGN 信道下这是最优软信息,符号方向必须和 BPSK 映射一致:code=0 映射到 -1,code=1 映射到 +1,此时 LLR 为正表示更可能是 1,所以后续判决用 LLR<0 判定为 0。如果你的映射反过来,译码器会完全反相,误码率固定在 50% 左右,这是最常见的“能跑但是结果完全错误”的原因。
4.2 参数配置与仿真时间估算
上面的参数设置适合快速验证,实际跑一次大约几十秒。如果你要得到平滑的 BER 曲线,建议每个信噪比点至少跑几百帧,并且覆盖较宽的信噪比范围。下面是我常用的参数表:
| 参数 | 建议值 | 说明 |
|---|---|---|
| Z | 64~512 | 越大码长越长,性能越好,但仿真越慢 |
| max_iter | 50~100 | BP 迭代次数,小于 20 会损失性能 |
| num_frames | 200~1000 | 每个 SNR 点帧数,BER 低时要加大 |
| EbN0_dB | 0~3 dB | 码率 0.5 的 QC-LDPC 在这个范围出现瀑布区 |
| 噪声模型 | AWGN | BPSK 调制,单径信道 |
经验值是:Z=128,码率 0.5,每帧 768 比特,100 帧,信噪比 6 个点,用 MATLAB R2022 跑,大约需要 5~10 分钟。其中 BP 译码是瓶颈,尤其双重循环版本。如果你想快一点,可以把 max_iter 降到 30,误差不大。Z 越大,每帧需要编码和译码的比特数越多,帧数要相应减少。
4.3 结果解读:什么样的 BER 曲线才算正常
跑完仿真,得到的曲线通常是一个“瀑布”形状:低信噪比时 BER 缓慢下降,到了一个拐点后急剧下降,然后进入平台期。对于码率 0.5、码长几百比特的 QC-LDPC,瀑布区大体在 1~2 dB。如果曲线整体右移或出现拐点后下降缓慢,要么是校验矩阵有短环,要么是译码迭代次数不够,要么是噪声方差算错。
我通常会用一组对比来判断系统是否正常:先跑一个 uncoded BPSK 的 BER 曲线做参考(也就是不编码直接把信息比特映射发送),再跑 QC-LDPC 的 BER 曲线。在 Eb/N0=3 dB 时,未编码系统的 BER 大约在 1e-2 量级,而码率 0.5 的 LDPC 应该降到 1e-4 以下。如果两条曲线几乎重合,那肯定是译码器或编码器出了问题,或者是信噪比映射错了一个因子。
5. 常见问题与避坑:校验矩阵构造和译码收敛的 4 个典型坑
5.1 现象:基矩阵展开后 H 的行重列重忽大忽小
原因:我在用 circshift 填子矩阵时,一开始写成 circshift(I, [s, 0]),导致子矩阵变成了列循环上移。虽然单个子矩阵还是每行每列一个 1,但放到整个 H 里,行列的分布规律被打乱了。更典型的错误是基矩阵里的移位数值超出 0~Z-1 范围,比如 Z=64 时写了 s=64,circshift 会自动模 Z,导致子矩阵和基矩阵描述不一致。
解决:展开后立刻验证 H 的每行每列重量。在测试脚本里加两行代码:
assert(all(sum(H, 2) == sum(H, 2, 'all')/M), '行重异常'); assert(all(sum(H, 1) == sum(H, 1, 'all')/N), '列重异常');而且打印出前两个子矩阵,用 imagesc 看是不是单位阵的移位。后来我把这个检查写进了 expand_matrix.m 的尾部,每次生成后自动断言,从那以后矩阵构造的问题少了大半。
5.2 现象:仿真出来的 BER 永远是 0.5
原因:BPSK 映射和 LLR 符号方向不匹配。比如我曾在代码里写 BPSK 映射 x = 2code - 1,那么 code=0 对应 x=-1,code=1 对应 x=1。LLR 应该是 2y/noise_var,这个值是正的说明更可能 code=1,负的说明更可能 code=0。但我在判决时写成了 decoded = double(LLR_post > 0),等于把判决完全反过来。
解决:统一从定义出发查一遍。最简单的验证方法:把 max_iter 设为 1,noise_var 设成极大,也就是信道几乎全是噪声,此时译码输出应该接近随机猜测,BER 约 0.5。然后逐步减小噪声,观察 BER 是否单调下降。如果下降方向相反,肯定是符号反了。我后来习惯在仿真前先跑一个帧:构造一个全零码字,加少量噪声,打印译码后的信息位,一眼就能看出映射对不对。
5.3 现象:BP 译码器迭代不收敛,提前停止条件形同虚设
原因:校验节点更新里的连乘操作,当 V2C 绝对值较大时,tanh(x/2) 会接近 ±1,连乘后 product 可能为 0,C2V 变成 NaN。NaN 在变量节点更新里会像病毒一样传播,整个迭代后 LLR_post 全是 NaN,decoded 全是 0,校验和恰好全零,提前跳出,误码率看起来是 0,实际上根本没译码。
解决:在做连乘前对所有 V2C 做一次限幅,比如把消息限制在 [-10, 10] 之间。在代码里加一行:
V2C = max(min(V2C, 10), -10);同时把 product 加上防下溢判断,我前面给的代码里已经写了 abs(product) < 1e-15 的处理。如果不做限制,交替相乘会出现 0 × ∞ 的组合,NaN 很难躲开。这个问题不跑大数据量仿真基本不会暴露,一旦暴露就是整片曲线崩掉。
5.4 现象:Min-Sum 译码比 BP 性能差一截,瀑布区右移 1 dB
原因:Min-Sum 把校验节点更新时复杂的 tanh 连乘简化成 找最小绝对值 和 符号乘积,这相当于对校验消息做了一次近似,近似误差在校验节点度数小的时候特别明显。我最初直接用 Min-Sum 跑,没做任何修正,比 BP 差 0.8 dB 左右。
解决:在 Min-Sum 的校验节点更新后乘一个归一化因子或加一个偏移量。归一化因子常见取 0.75~0.9,偏移量常见取 0.5。我用的做法是:
new_msg = sign_product * (min_abs - 0.5); if min_abs < 0.5 new_msg = 0; end这样补偿之后,Min-Sum 的性能能拉到和 BP 差别很小。这个话题我会在下一章展开写详细实现。
6. 二次开发:把 BP 换成 Min-Sum、支持多码率的两步改造
6.1 改造一:在 decode_bp 基础上写出修正 Min-Sum
Min-Sum 的核心是高速译码,很多工程场景用它替代 BP。改造不需要推翻整个文件,只需改校验节点更新那部分。我把原来的 tanh 积换成了符号连乘加找最小绝对值,并加上偏移修正。你可以复制下面这个函数,在原有工程里替换 decode_bp 的接口,输入输出完全一致。
function decoded = decode_min_sum(channel_llr, H, max_iter) % 输入输出与 decode_bp 相同,内部用修正 Min-Sum 更新校验节点 [M, N] = size(H); [row, col] = find(H); V2C = zeros(M, N); C2V = zeros(M, N); offset = 0.5; % 偏移修正因子 for iter = 1:max_iter % 校验节点更新 for i = 1:M edges_i = col(row == i); if isempty(edges_i) continue; end % 收集所有进入该校验节点的消息绝对值 abs_msgs = abs(V2C(i, edges_i)); min1 = min(abs_msgs); sign_prod = prod(sign(V2C(i, edges_i))); for idx = 1:length(edges_i) j = edges_i(idx); others = abs_msgs; others(idx) = Inf; % 排除自身 min2 = min(others); if min2 == Inf min2 = min1; % 度数=1 时用自身 end new_msg = sign_prod * sign(V2C(i,j)) * max(min2 - offset, 0); C2V(i, j) = new_msg; end end % 变量节点更新、后验LLR判决、提前停止均与BP相同 % 此处省略重复部分,实际工程可直接复制 decode_bp 后半段 end end这个函数在计算每个输出消息时,排除当前目标变量节点的消息,在剩余的绝对值中找最小的那个,然后减去偏移量并截断到非负。sign_prod 乘以 sign(V2C(i,j)) 的写法要注意:整体符号等于所有输入符号的乘积,也等于目标消息符号乘以其他消息符号乘积,所以最终符号是 sign_prod * sign(V2C(i,j))。我的经验是,offset 取 0.5 在大多数场合效果不错,如果码重量较大,可以尝试 0.7 或 0.9,用几次仿真对比选最优。
改造完后,把 run_simulation.m 里的 decode_bp 换成 decode_min_sum,然后对比两条 BER 曲线。你会发现曲线基本重合,但单帧译码时间明显降低。我经常用两个函数跑同一组信噪比,画成对比图,做汇报时非常有说服力。
6.2 改造二:把基矩阵参数化,一个仿真入口支持多种码率
原工程里 base_matrix 是写死的常量,想换码率必须手动改数组。我建议把基矩阵定义独立成函数,入参直接传目标码率,函数内部返回一个满足约束的基矩阵。下面是一个简单示例,生成码率 0.5 和 0.75 两种基矩阵:
function base_matrix = get_base_matrix(rate) % 根据目标码率返回一个 QC-LDPC 基矩阵 % 实际场景需要从标准文档或设计算法生成,这里只做演示 switch rate case 0.5 base_matrix = [ 0 1 2 -1 0 1; 1 0 -1 2 1 0; 2 -1 1 0 2 1; ]; case 0.75 % 4×16 基矩阵,码率 1-4/16=0.75 base_matrix = [ 0 1 2 -1 0 5 1 -1 2 0 3 -1 1 0 -1 4; 1 -1 0 2 1 0 4 1 -1 2 0 3 -1 1 0 -1; 2 1 -1 0 -1 2 0 -1 1 3 1 0 4 -1 2 1; -1 0 1 2 3 -1 1 0 -1 1 2 -1 0 1 0 2; ]; otherwise error('不支持的码率'); end end然后在主脚本里用 rate 变量控制一切:
rate = 0.75; Z = 128; base_matrix = get_base_matrix(rate); H = expand_matrix(base_matrix, Z); [mb, nb] = size(base_matrix); M = mb * Z; N = nb * Z; K = N - M;这样做的好处是,你想评估不同码率时只需要改 rate 和对应基矩阵,编码、译码、仿真代码完全复用。注意,不同码率下最小和译码的偏移量最优值可能不同,最好分别调一次。 我从这个工程里学到的一条教训是:把基矩阵写成配置,而不是埋在仿真脚本里。有一次我为了验证码率 0.75 的性能,直接在原脚本里替换了 base_matrix 数组,结果后面跑其他实验时误用了旧矩阵,输出结果对不上。从那以后我每次改造都会强制先跑一遍 tests/test_h_matrix.m,确认矩阵尺寸和码率一致,再动仿真。做编码仿真,变量名和矩阵维度的错误往往比算法错误更隐蔽,形成一套固定的自检流程,比事后排查省心得多。希望这篇解析能帮你把 QC-LDPC 这套工程跑起来,也让你在改代码时少踩几个我踩过的坑。
本文还有配套的精品资源,点击获取