简介：这份QC-LDPC编码MATLAB实现压缩包，面向通信与编码理论初学者，特别是需要动手理解准循环LDPC码构造与编码流程的学习者。包内共3个文件，包含2个m脚本和1个asv自动保存文件，其中QCEncode.m对应编码主程序，QCLDPCBaseH.m用于生成基矩阵或校验矩阵，代码附有详细注释，便于逐步跟踪码字生成、信道模拟与解码迭代过程。资源整体仅4KB，小巧但结构清晰，适合快速上手。已有1254人学习下载。通过学习这份示例，读者能直观掌握QC-LDPC的准循环特性在MATLAB中的实现方式，理解生成矩阵与校验矩阵的配合使用，并可在修改参数后观察不同设置下的性能差异，为后续深入研究LDPC编解码或硬件实现打下基础。
1. QC-LDPC是什么：为什么这份MATLAB资源能省掉你两周的编码时间
做物理层仿真的同行对这类场景应该不陌生：论文里的QC-LDPC校验矩阵写得清清楚楚，落到MATLAB里却四处报错，要么矩阵维度对不上，要么循环移位方向反了，光把基矩阵展开成完整的H矩阵就能耗掉大半天。QC-LDPC，全称准循环低密度奇偶校验码，是5G NR、WiFi 6、卫星高速数传里绕不开的标准信道编码方案，核心思路是用一个很小的基矩阵加循环移位参数，生成超大但极度稀疏的校验矩阵，从而把存储和译码复杂度同时压下来。这份QC-LDPC.zip正是围绕这个需求打包的，里面是完整的MATLAB实现链路，从基矩阵定义、循环移位展开、编译码器构建到误码率曲线绘制，一条龙跑通。适合正在做LDPC相关课程设计或毕业设计的学生，也适合需要快速验证编码增益的物理层算法工程师。
2. 校验矩阵构造：基矩阵展开、循环移位与稀疏化装配
2.1 基矩阵元素、扩展因子与QC结构的对应
QC-LDPC的校验矩阵H不是凭空生成的，它由一个维度为mb×nb的基矩阵B和一个扩展因子Z共同决定。B里每个元素对应一个Z×Z的子矩阵，按照约定，元素为-1时表示全零子矩阵，元素为0时表示单位矩阵，元素为大于0的整数s时表示单位矩阵循环右移s位得到的置换矩阵。这个设计的精妙之处在于，完整的H矩阵维度是(mb×Z)行、(nb×Z)列，实际只需要存储mb×nb个整数。以常见的基矩阵维度12×24、Z=81为例，完整校验矩阵是972×1944，但存储节点上只需要288个整数。
这也解释了QC-LDPC名字里“准循环”三个字的由来：每一行子矩阵块都是单位矩阵的循环移位，整行整列天然具备循环特性。这个结构直接决定了编码可以用移位寄存器实现，译码时的校验节点更新也能按块并行，硬件友好度远高于随机构造的LDPC码。
实践中，基矩阵的行重和列重要控制在合理范围。行重大约在6到10之间比较常见，列重则取决于码率，高码率下列重可能只有2到3。行重太小会导致校验约束过弱，译码收敛慢；行重太大虽然约束强，但校验节点更新的计算量会明显上升。这个权衡在后面的译码参数里还会体现。
2.2 循环移位子矩阵的MATLAB实现
拿到基矩阵后的第一步，是把每个元素展开成对应的Z×Z稀疏矩阵。这里强烈建议直接用MATLAB的sparse类型，而不是用full矩阵。一个Z=512的完整H矩阵如果按稠密存储，光一个子块就是512×512的double数组，几个块加起来内存就吃不消了。
function S = shift_submatrix(shift_val, Z) % 构造 QC-LDPC 的循环移位子矩阵 % shift_val = -1 : 全零方阵 % shift_val = 0 : 单位阵 % shift_val > 0 : 单位阵各列循环右移 shift_val 列 if shift_val < 0 S = sparse(Z, Z); return; end % 列索引做循环平移，行索引保持 1:Z col_idx = mod((0:Z-1) - shift_val, Z) + 1; row_idx = (1:Z)'; S = sparse(row_idx, col_idx, ones(Z,1), Z, Z); end这个实现里有两个容易出错的地方。第一，sparse(i,j,v,m,n)的语法中，i和j必须是向量且维度一致，row_idx用转置列向量、col_idx用行向量会直接报错，要么都取列向量，要么都取行向量。第二，mod运算要先归一化到[0, Z-1]再+1，如果直接用shift_val没用mod包一层，当移位值等于Z时索引会越界。
循环右移的“右”是什么方向并不重要，重要的是编码端和解码端必须用完全相同的约定。常见做法是统一用“列索引减移位值再取模”，也就是上面代码写的这个方向，写进一个公共函数后全工程复用，不要再在别的地方重复实现。
2.3 用三步搭出完整H矩阵并验证行列重
有了移位子矩阵函数，拼装完整校验矩阵就是逐块赋值的问题。最简单直观的写法是双重循环遍历基矩阵的每个元素，按块坐标将子矩阵填入目标位置。
function H = buildH(base_mat, Z) % 基矩阵 base_mat 逐块展开为完整稀疏校验矩阵 [mb, nb] = size(base_mat); H = sparse(mb*Z, nb*Z); for r = 1:mb for c = 1:nb s = base_mat(r,c); if s < 0 continue; % 全零块直接跳过 end idx_r = (r-1)*Z+1 : r*Z; idx_c = (c-1)*Z+1 : c*Z; H(idx_r, idx_c) = shift_submatrix(s, Z); end end end这段代码在基矩阵规模小、Z不超过128时完全够用。逐块赋值虽然中间会产生一些稀疏矩阵的动态填充，但可读性最好，也和手写推导时的思路一一对应，适合先跑通再优化的节奏。
展开完成后，建议立刻做两个验证。第一是行重和列重检查，用sum(H,2)和sum(H,1)统计每行每列非零元素数，核对是否与基矩阵设计一致。第二是H*ones(N,1)的结果，这一步本质是检查每行的非零元素个数，如果出现奇数行重且全1向量乘出来不为0，说明既定矩阵构造时的mod 2约束有问题，后面编码必炸。
% 验证行重列重是否与基矩阵设计一致 row_weights = full(sum(H, 2)); col_weights = full(sum(H, 1)); fprintf('row weight min/max : %d / %d\n', min(row_weights), max(row_weights)); fprintf('col weight min/max : %d / %d\n', min(col_weights), max(col_weights));这里用full包裹是因为sum返回的是稀疏向量，直接min/max有时会得到0值造成错觉。我一般会在这一步额外检查H是否满秩，rank(full(H))只适合小Z场景，大矩阵用spfun配合判断零行零列就够了。这一步验证不通过，就不要往下走编码，否则后面全是白费功夫。
3. 编码实现：生成矩阵、近似下三角快速编码与边界参数
3.1 生成矩阵路线：GF(2)上的高斯消元
拿到校验矩阵H后，最直接的编码方式是求出生成矩阵G，然后信息位乘G得到码字。这里需要在二元域上对H做高斯消元，把它化为系统形式。MATLAB内置的rref是实数域实现，直接拿来用会得到浮点结果，必须自己写GF(2)版本。
function [G, H_sys] = gf2_rref(H) % 在 GF(2) 上对 H 做高斯-若尔当消元，化为系统形式 % 输入 H 需为满秩或可化简为 [P | I] 的结构 [M, N] = size(H); H_aug = [H, speye(M)]; H_aug = full(H_aug); % 小规模时用稠密方便调试 for col = 1:M pivot = find(H_aug(col:end, col), 1) + col - 1; if isempty(pivot) continue; end if pivot ~= col H_aug([col, pivot], :) = H_aug([pivot, col], :); end % 消去其它行的该列 rows_to_elim = find(H_aug(:, col)); rows_to_elim(rows_to_elim == col) = []; for r = rows_to_elim' H_aug(r, :) = mod(H_aug(r, :) + H_aug(col, :), 2); end end H_sys = H_aug(:, 1:M); % 系统部分 P = H_aug(:, M+1:end); % 校验部分 G = mod([P', speye(N-M)], 2); % 生成矩阵 end这个函数有两个代价需要提前知道。第一，高斯消元会把原来极度稀疏的H变成接近稠密的矩阵，G更是几乎全稠密，存储量从O(N)涨到O(N^2)，码长超过2000时内存就很紧张。第二，消元过程如果遇到H不满秩，pivot会提前消失，此时得到的G并不能保证H*G'等于零，需要额外验证。
验证方法是随机生成信息序列，编码后立刻检查mod(H*cw', 2)是否全零，这是编码正确性的金标准。我一般会在生成G之后强制跑一遍这个自检，不通过就重新审视H的构造。
3.2 快速编码路线：Richardson-Urbanke分块三角化
生成矩阵路线写起来简单，但码长一大就跑不动，工程上更多采用Richardson-Urbanke提出的近似下三角编码法。这个方法的思路是对H的行列做重排，把它分成六个子块：
H = [A B T; C D E]
其中T是下三角方阵且对角线全为1。利用分块结构，编码分两步算校验位p1和p2，关键是T的逆用前向代入完成，不显式求逆矩阵，复杂度从O(N^2)降到O(N)。
具体推导可以概括为两个公式。令s为信息位行向量，先计算辅助向量p1 = (ET^(-1)A + C)s'，再计算p2 = T^(-1)(As' + Bp1)。整个过程只需要矩阵乘法和前向代入，不涉及对T求逆，MATLAB里用反斜杠运算符配合稀疏矩阵也能直接处理。
3.3 编码函数骨架与边界参数检查
下面这个编码器骨架是模拟项目X里最常用的封装，把H的分块信息和信息向量接进来，返回完整码字。
function codeword = qc_ldpc_encode(info_bits, H, block_struct) % info_bits : 信息位列向量，长度 K = (nb - mb)*Z % H : 完整校验矩阵 % block_struct : 包含 A,B,C,D,E,T 及 phi_inv 的分块结构体 % 边界检查 K = (size(H,2) - size(H,1)); assert(length(info_bits) == K, '信息位长度与码率不匹配'); % 分块编码 s = info_bits(:); p1 = mod(block_struct.phi_inv * ... mod(block_struct.E * mod(block_struct.T \ ... (block_struct.A * s), 2) + block_struct.C * s, 2), 2); p2 = mod(block_struct.T \ ... (block_struct.A * s + block_struct.B * p1), 2); codeword = [s; p1(:); p2(:)]; codeword = mod(codeword, 2); % 保险起见最终再取一次模 end这里有一个容易忽略的点：反斜杠运算符在稀疏矩阵上默认用的是LU或Cholesky分解，并不保证结果是严格GF(2)的，所以每次运算后都要强制mod取模。另一个边界参数是码率，QC-LDPC的码率由基矩阵维度决定，R=(nb-mb)/nb，在调用编码函数前应该先算好这个值，避免后面做BER仿真时信噪比换算出错。
我的一般习惯是：先用高斯消元版本在小Z下验证算法逻辑，跑通全零自检和单次仿真后，再切到快速编码版本做长码仿真。这样两头都有对照，出了问题能快速定位是编码逻辑还是矩阵构造的问题。
4. 译码实现：对数域BP、最小和近似与迭代控制
4.1 消息传递与LLR的符号约定
QC-LDPC的译码目前工程上几乎全是置信传播（BP）类算法，核心是变量节点和校验节点之间迭代交换软信息。直接的概率域乘法会涉及大量连乘，数值容易下溢，所以实际实现都转到对数域，用LLR（对数似然比）表示消息。
这里最关键的坑是符号约定。LLR定义为log(P(b=0)/P(b=1))还是反过来，由调制映射决定。BPSK调制里，如果0映射为+1、1映射为-1，那么接收信号r对应的LLR是2r/σ²，这个符号方向必须全程一致。如果调制映射或者LLR定义反过来，整条译码路径的符号处理都要跟着变。
校验节点的更新规则是BP的核心。在校验节点m上，传给变量节点n的消息是所有输入到m的其它变量节点消息的“符号乘积乘以最小值”，这就是最小和近似。相比标准BP的tanh域运算，最小和近似损失了一些精度，但运算量大幅下降，实际中通过乘一个归一化因子alpha来补偿，alpha通常在0.7到0.8之间。
4.2 对数域BP译码的MATLAB代码
下面这段代码是译码器的骨架实现，结构清晰，适合先跑通再优化。
function decoded = bp_decode(llr_in, H, max_iter, alpha) % llr_in : 接收符号的 LLR 列向量，长度 N = nb*Z % H : M x N 稀疏校验矩阵 % max_iter : 最大迭代次数 % alpha : 最小和近似归一化因子 [M, N] = size(H); [cn_idx, vn_idx] = find(H); edge_cnt = length(cn_idx); % 消息数组 msg_v2c = llr_in(vn_idx); % 变量节点到校验节点 msg_c2v = zeros(edge_cnt, 1); % 校验节点到变量节点 for it = 1:max_iter % —— 校验节点更新 —— for m = 1:M edges = find(cn_idx == m); if isempty(edges), continue; end mag = min(abs(msg_v2c(edges))); % 最小幅度 sgn_all = prod(sign(msg_v2c(edges))); % 总符号 for e = edges' % 去掉当前边自身符号，等价于其它边的符号乘积 msg_c2v(e) = alpha * mag * sgn_all * sign(msg_v2c(e)); end end % —— 变量节点更新 —— for n = 1:N edges = find(vn_idx == n); if isempty(edges), continue; end total_llr = llr_in(n) + sum(msg_c2v(edges)); for e = edges' msg_v2c(e) = total_llr - msg_c2v(e); end end % —— 硬判决与早停 —— llr_total = llr_in; for e = 1:edge_cnt llr_total(vn_idx(e)) = llr_total(vn_idx(e)) + msg_c2v(e); end hard = double(llr_total >= 0); if mod(H * hard(:), 2) == 0 decoded = hard(:); return; end end % 达到最大迭代次数仍未收敛，按当前LLR硬判决输出 decoded = double(llr_total >= 0); end代码逻辑上要抓住三个要点。第一，校验节点更新用了“总符号乘以当前边符号”来抠掉自身，这个技巧比每次都单独算其它边的乘积快得多，但前提是sign返回±1，数学上正好抵消。第二，变量节点更新是线性叠加，LLR域的消息合并就是求和，不需要任何近似的操作。第三，早停条件用mod(H*hard, 2) == 0，这要求H和hard都是整数类型，double下数值太大时会有浮点误差，建议对H用logical存储。
4.3 迭代次数、归一化因子与早停策略
max_iter的选择直接影响性能与仿真耗时的平衡。BPSK加AWGN场景下，迭代次数从5增加到20，性能提升明显；从20增加到50，增益通常只有0.1到0.2dB，但耗时翻倍不止。做蒙特卡洛仿真时我一般设20，最后一条曲线需要在低误码率区看错误地板时，再临时拉大到50。
alpha的取值相对敏感。alpha偏大（超过0.8）会让消息过度自信，中低信噪比下误码率反而变差；偏小（低于0.7）则收敛太慢，需要更多迭代次数补偿。常见做法是alpha=0.75起步，如果H的行重都比较大（比如8以上），可以试0.8；如果行重偏小（4到6），0.7更稳。
早停策略值得多说一句。如果每条码字都迭代满max_iter，仿真100万帧的时间成本非常高。实际实现里先算mod(H*hard, 2)判断是否满足所有校验方程，满足就提前退出，不满足才继续迭代。高信噪比下大部分码字1到2次迭代就收敛了，早停能省掉80%以上的译码时间。
5. 避坑指南：QC-LDPC仿真中的五个高频翻车点
5.1 坑一：移位方向不一致导致误码率平台
现象：编码端和解码端各写了一个子矩阵生成函数，单测都通过，但整套仿真跑下来误码率在某个值附近出现平台，再怎么增加迭代次数都降不下去。
原因：编码时用的循环右移，解码时按循环左移展开H，两边实际上用的是不同的校验矩阵，相当于编码码字对不上解码约束。
解决：抽一个公共函数，像前面第2章的shift_submatrix那样，全工程只认这一份。如果从网上下载的代码改过方向，务必同时替换编码和解码两处。
5.2 坑二：稀疏矩阵动态赋值引发内存爆炸
现象：Z=128时还能跑，调到Z=512后buildH函数执行时间暴涨，内存直接耗尽，MATLAB卡死报Out of Memory。
原因：逐块稀疏赋值H(idx_r, idx_c) = subM时，稀疏矩阵每插入一批新索引都可能触发内部重组，插入位置不连续时填充量急剧上升，最终变成近似稠密存储。
解决：改成一次性组装。先收集所有非零块的i、j、v向量，再用sparse(i,j,v,m,n)一步生成。
function H = buildH_fast(base_mat, Z) % 一次性组装方案，避免稀疏矩阵动态赋值 [mb, nb] = size(base_mat); I = []; J = []; V = []; for r = 1:mb for c = 1:nb s = base_mat(r,c); if s < 0, continue; end rows = (r-1)*Z + (1:Z)'; cols = mod((0:Z-1) - s, Z) + 1 + (c-1)*Z; I = [I; repmat(rows, 1)]; % 保持向量 J = [J; cols(:)]; V = [V; ones(Z,1)]; end end H = sparse(I, J, V, mb*Z, nb*Z); end注意I、J、V三个向量长度必须一致，且I和J中不能有零下标，sparse函数会直接报错。这个写法的性能瓶颈转移到向量拼接上，Z=512时依然流畅。
5.3 坑三：LLR符号约定与全零自检
现象：仿真结果低信噪比时看起来正常，高信噪比时误码率反而上升，曲线呈“U”型。
原因：LLR定义与调制映射不匹配，解码器内部认为高LLR对应比特1，但实际调制把高LLR对应比特0，高信噪比下符号反转的错误无法被迭代纠正。
解决：先做全零码字自检。发送全零码字、BPSK调制后加噪声解调出LLR，译码结果必须仍为全零。如果译码输出全1，把llr_in取反再试，直到自检通过再跑BER。
5.4 坑四：信噪比定义混乱导致曲线偏移
现象：仿真曲线和理论曲线形状相同，但整体偏移3dB或一个固定倍数，看起来像编码增益不对。
原因：蒙特卡洛仿真里用的信噪比是Es/N0还是Eb/N0没有厘清。BPSK下一个符号携带1个比特，但考虑了码率R后，Es/N0和Eb/N0之间差10*log10(R)，R是编码率。如果按码率折算后仍然用Es/N0直接画图，曲线整体会向右（或向左）平移。
解决：仿真时先明确横轴。公式是esn0_db = ebn0_db + 10*log10(code_rate)。同时加噪用的噪声方差要基于Es/N0计算，不能拿Eb/N0的方差直接套。
5.5 坑五：迭代不收敛与数值饱和
现象：译码函数在部分帧上始终无法通过校验方程的早停判断，最后硬判决输出后误码率比预期高一个数量级。
原因：低信噪比时LLR输入本身噪声很大，迭代过程中消息幅度逐步增大，超过double的表示范围后出现Inf或NaN，再往后全是NaN，译码器完全失效。
解决：对消息幅度做饱和钳制，最简单的做法是把llr_in限制在±20以内，每次校验节点更新后同样对msg_c2v钳制。另一种做法是检测消息中出现NaN时提前终止该帧，按当前硬判决输出，避免NaN在迭代中扩散。
6. 验证与进阶：把蒙特卡洛仿真跑通并调出好性能
6.1 蒙特卡洛仿真框架
跑误码率曲线是验证整套QC-LDPC实现最直接的手段。框架本身不复杂，但统计细节容易漏。下面这段BPSK+AWGN仿真的骨架，噪声方差和LLR映射与前面的约定保持一致。
EbN0_dB = 0:0.5:3.5; code_rate = (nb - mb) / nb; max_iter = 20; nFrames = 2000; BER = zeros(size(EbN0_dB)); for idx = 1:length(EbN0_dB) esn0_dB = EbN0_dB(idx) + 10*log10(code_rate); noise_var = 10^(-esn0_dB/10); for k = 1:nFrames info = randi([0 1], K, 1); cw = qc_ldpc_encode(info, H, block_struct); mod_sym = 1 - 2*cw; % 0 -> +1, 1 -> -1 noise = sqrt(noise_var/2) * randn(size(mod_sym)); rx = mod_sym + noise; llr = 2*rx / noise_var; % 与调制映射对应的LLR dec = bp_decode(llr, H, max_iter, 0.75); BER(idx) = BER(idx) + sum(dec ~= info) / K; end BER(idx) = BER(idx) / nFrames; end semilogy(EbN0_dB, BER, 'o-');注意噪声方差的处理。BPSK信号幅度为±1时，符号能量Es=1，当Es/N0以dB给出后，噪声方差是10^(-(Es/N0)/10)，而由于BPSK每符号仅1比特，这里不需要再除2。如果你用的是QPSK或更高阶调制，Es和Eb的关系就要按调制阶数和码率重新推导。
6.2 参数微调经验与短环检查
MC仿真跑通后，性能还能从两个方向再抠。第一是归一化因子alpha，0.75往往不是最优，逐点尝试0.72、0.75、0.78三组，选BER最低的。第二是检查H里有没有长度为4的短环。短环会让BP译码消息过度相关，形成错误地板，性能曲线会在某个信噪比后不再下降。检查方法不复杂，对H中的每个非零元素，搜索是否存在另一个位置和它构成4环。发现短环时重新选基矩阵或调整扩展因子Z，通常能直接消掉。
还有一个小技巧：把LLR初始化时的噪声方差做一次缩放（比如乘以0.9），能小幅改善中低信噪比性能。这个操作本质上是对输入软信息做加权，工程里常用来弥补最小和近似的精度损失，但对具体矩阵比较敏感，需要实测确认。
从那以后，我每次拿到新的QC-LDPC码本都强制自己走同一套流程：先画基矩阵结构，再跑全零自检，接着单帧译码看迭代收敛曲线，最后才上蒙特卡洛。这套习惯帮我躲过了无数个莫名其妙的3dB偏移和错误地板。希望这段经验能帮到你。
本文还有配套的精品资源，点击获取
