极化码MATLAB仿真：编码、SC译码与误码率统计实战

简介：PCode.zip 是一份面向通信工程学习者与科研人员的极化码 MATLAB 仿真源码包，聚焦 Erdal Arikan 提出的极化码编解码原理与性能验证，适合具备一定 MATLAB 与信道编码基础、希望深入理解 5G 极化码实现细节的读者。压缩包共 32 个文件，全部为 .m 脚本，整体约 35KB，涵盖码构造、编码、信道合成与解码等核心模块，如 pdecode、pcencode、initPC、SC_Decoder_ver2 等函数，并配有 MonteCarlo 仿真与误码率绘图脚本，便于观察码率、码长、迭代次数对系统性能的影响。源码注释较为详细，有助于理解比特翻转、LLR 更新与 BEC 信道处理等关键步骤。目前已有 905 人学习下载，可作为极化码理论学习与仿真实验的参考工具。

1. 从一份 PCode.zip 说起：极化码 MATLAB 仿真到底能跑出什么

如果你正在做极化码（Polar Coding）相关的课程设计、论文复现或者算法预研，大概率会遇到一个很现实的问题：理论推导看懂了，SC 译码的递归公式也背下来了，但真要在 MATLAB 里从零搭一套能跑通编码、信道传输、译码、误码率统计的完整链路，光是索引对齐和递归边界就能耗掉两三天。PCode.zip 就是冲着这个场景来的——它把极化码的编码构造、SC 译码（pdecode）以及配套的仿真脚本打包在一起，解压后直接是 MATLAB 工程文件，不需要额外编译，也不依赖通信工具箱里的冷门函数。适合谁？适合已经学过极化码基本原理、需要快速拿到可运行基线代码做对比实验的人；也适合刚接触极化码、想通过读一份能跑通的代码来反向理解 SC 译码递归逻辑的新手。但要注意，这份资源不是“一键出论文图”的黑匣子，参数配置和信道模型需要你自己根据实验目标调整，下面我会把每个环节拆开讲。

2. 极化码编码链路拆解：从信道极化到码字生成

2.1 极化码的核心机制与 MATLAB 实现映射

极化码的理论基础是信道极化：通过递归的信道合并与分裂操作，让一部分比特信道趋于无噪声、另一部分趋于全噪声，然后把信息比特放在“好”信道上。在 MATLAB 里，这个过程的实现通常围绕三个核心变量展开：码长 N、信息比特数 K、以及冻结比特位置集合。PCode.zip 里的编码脚本一般会先根据高斯近似或巴氏参数法计算各子信道的可靠性，再排序选出 K 个最可靠的索引放信息比特，其余置零作为冻结比特。

这里有个容易被忽略的点：极化码的编码本质上是线性分组码，生成矩阵是核矩阵 F 的克罗内克积。MATLAB 里如果直接用kron函数逐级构造，当 N 到 1024 以上时内存占用会明显上升。常见做法是用递归函数代替显式矩阵构造，或者利用极化码编码的蝶形结构做原地运算。PCode.zip 里的实现我翻了一下，用的是递归编码方式，对 N=1024 这种量级比较友好。

% 极化码递归编码核心逻辑（示意） function x = polar_encode(u, N) % u: 输入比特序列，长度 N % N: 码长，必须是 2 的幂 if N == 1 x = u; return; end % 将输入分为前后两半 u1 = u(1:N/2); u2 = u(N/2+1:N); % 递归编码前半部分 x1 = polar_encode(u1, N/2); x2 = polar_encode(u2, N/2); % 蝶形合并：前半为模二加，后半直接传递 x = [mod(x1 + x2, 2), x2]; end

上面这段代码展示的是极化码编码的递归骨架。u是经过信息位映射后的完整比特序列，长度 N 必须是 2 的幂。递归终止条件是 N=1，此时直接返回。合并步骤里mod(x1 + x2, 2)对应极化码核矩阵 F 的上支路异或操作，x2对应下支路直通。参数 N 决定了递归深度，N=1024 时递归深度为 10 层，每层做一次长度为 N 的向量操作，整体复杂度是 O(N log N)。如果你把 N 设成非 2 的幂，递归会在 N/2 处产生小数索引，MATLAB 会直接报维度错误——这是新手最容易踩的坑之一。

2.2 冻结比特位置选择：三种常用构造法对比

冻结比特集合的选取直接决定极化码在有限码长下的性能。PCode.zip 里默认用的可能是巴氏参数法（Bhattacharyya parameter），但实际做仿真时，高斯近似（GA）和极化权重（PW）也很常见。三者的区别在于：

如果你做的是 AWGN 信道下的仿真，建议把构造方法切到高斯近似或极化权重。PCode.zip 的脚本里通常有一个construct_method变量或者类似的配置项，改成对应字符串即可。注意：切换构造方法后，冻结比特集合会变，译码端的冻结位置必须同步更新，否则译码器会把信息比特当冻结比特处理，误码率直接飙到 0.5 左右——这不是代码 bug，是配置不一致。

% 高斯近似构造冻结比特集合（简化版） function frozen_bits = ga_construct(N, K, snr_dB) % N: 码长, K: 信息比特数, snr_dB: 设计信噪比 sigma = 10^(-snr_dB/20); % 初始化各子信道对数似然比均值 L = 2 / sigma^2; % 递归计算各分裂信道的可靠性（此处省略递归细节） % 最终得到可靠性排序 reliability_rank reliability_rank = compute_reliability(N, L); % 选取可靠性最低的 N-K 个位置作为冻结比特 [~, sorted_idx] = sort(reliability_rank, 'ascend'); frozen_bits = sorted_idx(1:N-K); end

这段代码的关键参数是snr_dB，它决定了构造时假设的信道质量。设计信噪比和仿真信噪比不需要完全一致，但偏差太大会导致构造失配，性能损失在 0.5 dB 以内可以接受，超过 1 dB 就说明构造方法需要换。compute_reliability是递归计算各子信道可靠性的函数，具体实现可以用高斯近似的递归公式，也可以用查表法加速。PCode.zip 里如果自带这个函数，直接调用即可；如果没有，需要自己补上。

3. SC 译码的 pdecode 实现：递归、判决与回溯

3.1 SC 译码的递归结构在 MATLAB 里怎么落地

SC（Successive Cancellation）译码是极化码最基础的译码算法，核心思想是逐比特判决，每个比特的判决依赖前面所有已判决比特。在 MATLAB 里实现 SC 译码，最直观的方式是递归函数，但递归深度随码长对数增长，N=1024 时递归深度为 10，MATLAB 默认递归限制通常够用，但如果 N 继续增大到 4096 以上，可能需要调整递归深度或者改成迭代实现。

PCode.zip 里的pdecode函数大概率是递归版本。SC 译码的递归过程分两个阶段：先递归计算左子节点的 LLR（对数似然比），判决后再递归计算右子节点。左子节点的 LLR 计算用f函数，右子节点用g函数。这两个函数的数值稳定性很关键，尤其是f函数在 LLR 绝对值很大时容易出现数值溢出。

% SC 译码中的 f 和 g 函数（数值稳定版） function out = f_func(a, b) % f(a,b) = 2*atanh(tanh(a/2)*tanh(b/2)) % 数值稳定实现：避免直接计算 tanh 溢出 out = sign(a).*sign(b).*min(abs(a), abs(b)) ... + log(1 + exp(-abs(a+b))) - log(1 + exp(-abs(a-b))); end function out = g_func(a, b, u) % g(a,b,u) = b + (1-2u)*a out = b + (1 - 2*u) * a; end

f_func里用的是 min-sum 近似加修正项，比直接套公式稳定得多。sign(a).*sign(b).*min(abs(a),abs(b))是 min-sum 的主项，后面两个 log 项是修正。如果你追求速度、对精度要求不高，可以只保留 min-sum 主项，性能损失大约 0.2~0.3 dB。g_func里的u是左子节点的判决结果，取值 0 或 1，(1-2*u)把 0 映射为 +1、1 映射为 -1，这是极化码译码里的标准操作。

3.2 译码流程中的比特判决与冻结位处理

SC 译码的判决逻辑是：如果当前比特索引属于冻结比特集合，直接判为 0；否则根据 LLR 的符号判决，LLR ≥ 0 判 0，LLR < 0 判 1。这个逻辑本身很简单，但在递归回溯时，判决顺序必须严格从左到右，不能并行化——这是 SC 译码的固有串行特性。

PCode.zip 的pdecode函数里，冻结比特集合通常作为输入参数传入。这里有个实操细节：冻结比特集合的索引方式要和编码端一致。有的实现用 0-based 索引，有的用 1-based，MATLAB 默认是 1-based，但如果你从 Python 或 C 代码移植过来，很容易在这里翻车。我一般会在译码函数入口加一行断言，检查冻结比特索引的最大值是否等于 N，如果不等就直接报错，省得后面误码率异常时到处找原因。

% SC 译码主循环（简化示意） function u_hat = sc_decode(llr, frozen_bits, N) % llr: 接收到的对数似然比序列 % frozen_bits: 冻结比特索引集合（1-based） % N: 码长 u_hat = zeros(1, N); % 递归译码 u_hat = recursive_decode(llr, frozen_bits, N, u_hat); end function u = recursive_decode(llr, frozen_bits, N, u) if N == 1 idx = length(u) + 1; % 当前判决位置 if ismember(idx, frozen_bits) u(idx) = 0; else u(idx) = (llr < 0); % LLR<0 判 1 end return; end % 计算左子节点 LLR a = llr(1:N/2); b = llr(N/2+1:N); llr_left = f_func(a, b); % 递归译码左子节点 u = recursive_decode(llr_left, frozen_bits, N/2, u); % 计算右子节点 LLR u_left = u(end-N/2+1:end); llr_right = g_func(a, b, u_left); % 递归译码右子节点 u = recursive_decode(llr_right, frozen_bits, N/2, u); end

这段代码里ismember(idx, frozen_bits)是冻结比特判断，idx的计算依赖length(u)，这要求递归过程中u的更新顺序严格从左到右。llr_left和llr_right的计算分别对应 SC 译码的两个递归分支。注意u_left的提取位置：u(end-N/2+1:end)取的是最近判决的左半部分比特，这个索引如果写错，译码结果会完全乱掉。我见过有人把end-N/2+1写成end-N/2，结果误码率一直卡在 0.3 下不去，排查了半天才发现是索引差了一位。

4. 仿真脚本配置与误码率统计：参数怎么设才不翻车

4.1 信噪比扫描与蒙特卡洛仿真次数

极化码仿真的核心输出是误码率（BER）或误块率（BLER）随信噪比变化的曲线。PCode.zip 里通常有一个主脚本，负责设置信噪比范围、每个信噪比下的仿真帧数、以及调用编码和译码函数。这里有两个参数直接决定仿真结果的平滑度和可信度：信噪比步长和每点帧数。

信噪比步长一般取 0.5 dB 或 1 dB，步长太大曲线拐点看不清楚，步长太小仿真时间成倍增加。每点帧数取决于目标误码率量级：如果目标 BER 在 10^-3 左右，每点至少 10^4 帧才能保证统计置信度；如果目标 BER 到 10^-5，帧数要加到 10^6 量级。PCode.zip 默认的帧数可能偏少，跑出来的曲线毛刺比较多，建议根据实际需求调大。

% 仿真主循环（信噪比扫描 + 蒙特卡洛） snr_range = 0:0.5:4; % 信噪比范围，步长 0.5 dB num_frames = 1e4; % 每个信噪比下的帧数 N = 1024; % 码长 K = 512; % 信息比特数 ber = zeros(size(snr_range)); for snr_idx = 1:length(snr_range) snr = snr_range(snr_idx); error_count = 0; total_bits = 0; for frame = 1:num_frames % 生成随机信息比特 info_bits = randi([0 1], 1, K); % 编码（含冻结比特插入） tx_bits = polar_encode_with_frozen(info_bits, N, K, snr); % BPSK 调制 + AWGN sigma = 10^(-snr/20); rx_signal = (1 - 2*tx_bits) + sigma * randn(1, N); % 计算 LLR llr = 2 * rx_signal / sigma^2; % SC 译码 u_hat = sc_decode(llr, frozen_bits, N); % 提取信息比特并统计误码 info_hat = u_hat(info_indices); error_count = error_count + sum(info_hat ~= info_bits); total_bits = total_bits + K; end ber(snr_idx) = error_count / total_bits; fprintf('SNR = %.1f dB, BER = %.4e\n', snr, ber(snr_idx)); end

这段代码里sigma = 10^(-snr/20)是 BPSK 调制下噪声标准差的换算公式，llr = 2 * rx_signal / sigma^2是 AWGN 信道下 LLR 的标准计算方式。info_indices是信息比特在完整码字中的索引集合，必须和编码端的映射一致。num_frames设成 1e4 时，SNR=4 dB 附近如果 BER 已经降到 10^-5 以下，统计到的错误数可能只有个位数，曲线尾部会抖动，这时候要么加大帧数，要么在图上标注置信区间。

4.2 冻结比特集合与信息位索引的一致性校验

前面提过冻结比特集合不一致会导致译码崩溃，这里展开说具体怎么校验。编码端和信息位索引的映射关系通常有两种实现方式：一种是先构造完整比特序列（信息位 + 冻结位），再整体编码；另一种是先编码信息比特，再在输出端插入冻结位对应的零。PCode.zip 用的是哪种，需要看polar_encode_with_frozen函数的实现。

不管哪种方式，校验方法是一样的：编码前记录信息比特在完整序列中的位置索引info_indices，译码后从u_hat里按同样的索引提取信息比特。如果info_indices和冻结比特集合frozen_bits满足union(info_indices, frozen_bits)恰好覆盖 1 到 N 且无交集，说明配置正确。我一般会在脚本开头加一段校验代码：

% 校验信息位与冻结位配置一致性 assert(length(info_indices) + length(frozen_bits) == N, ... '信息位与冻结位总数不等于码长 N'); assert(isempty(intersect(info_indices, frozen_bits)), ... '信息位与冻结位存在重叠索引');

这两行断言花不了几毫秒，但能省掉后面几小时的排查时间。血泪经验：有一次跑仿真，BER 曲线在低信噪比下正常，高信噪比下反而上升，最后发现是冻结比特集合里混入了一个信息位索引，导致译码端把信息比特强制判零，高信噪比下错误反而更明显。

5. 避坑与排查：极化码 MATLAB 仿真里最容易翻车的五个点

5.1 现象：误码率曲线在高信噪比下不降反升

原因：冻结比特集合与信息位索引不一致，译码端把信息比特当冻结比特处理，信噪比越高，LLR 越可靠，但强制判零导致的错误越明显。解决：在编码和译码之间加一致性断言，检查info_indices和frozen_bits的并集是否等于 1:N 且无交集。

5.2 现象：N=1024 时仿真跑不动，内存爆掉

原因：显式构造生成矩阵kron多次后矩阵维度爆炸，N=1024 时生成矩阵是 1024×1024，再往上做克罗内克积内存占用呈指数增长。解决：改用递归编码函数，或者利用蝶形结构做原地运算，避免显式存储生成矩阵。

5.3 现象：SC 译码结果全为 0 或全为 1

原因：LLR 符号判决逻辑写反，或者f_func的数值稳定实现里符号处理有误。解决：检查u(idx) = (llr < 0)这一行，确认 LLR ≥ 0 判 0、LLR < 0 判 1；同时用简单测试用例（N=2, K=1）手动验证f_func输出。

5.4 现象：仿真曲线毛刺严重，无法判断性能趋势

原因：每点帧数不足，高信噪比下统计到的错误数太少，BER 估计方差大。解决：将num_frames提高到 1e5 或 1e6，或者在图上标注每个点的错误数，错误数少于 100 的点视为不可信。

5.5 现象：切换构造方法后性能反而变差

原因：构造时假设的信道模型与仿真信道不匹配，比如用巴氏参数法构造却在 AWGN 信道下仿真。解决：AWGN 信道下优先用高斯近似或极化权重构造，设计信噪比与仿真信噪比偏差控制在 1 dB 以内。

6. 进阶技巧：把 SC 译码改成 SCL 并验证增益

SC 译码的局限在于逐比特硬判决，没有利用路径度量做候选保留。SCL（Successive Cancellation List）译码保留 L 条候选路径，在译码结束时选路径度量最小的作为输出，在相同码长下能带来 0.5~1 dB 的增益。PCode.zip 里的pdecode是 SC 版本，但改成 SCL 并不需要重写整个递归结构，核心改动是在每层递归返回时维护一个路径列表，而不是单条路径。

具体做法：把u从一维向量改成 L×N 的矩阵，每行是一条候选路径；在冻结比特位置，所有路径都判 0，路径度量不变；在信息比特位置，每条路径分裂成两条（判 0 和判 1），计算路径度量增量，然后按度量排序保留前 L 条。路径度量增量用 LLR 的绝对值加权：判 0 时增量为log(1+exp(-llr))，判 1 时增量为log(1+exp(llr))，实际实现里常用近似max(0, -llr)和max(0, llr)来降低计算量。

% SCL 译码路径度量更新（简化示意） function [paths, metrics] = scl_update(paths, metrics, llr, L, is_frozen) if is_frozen % 冻结比特：所有路径判 0，度量不变 return; end % 信息比特：每条路径分裂为两条 num_paths = size(paths, 1); new_paths = zeros(2*num_paths, size(paths, 2)); new_metrics = zeros(2*num_paths, 1); for i = 1:num_paths % 判 0 分支 new_paths(2*i-1, :) = [paths(i, :), 0]; new_metrics(2*i-1) = metrics(i) + max(0, -llr); % 判 1 分支 new_paths(2*i, :) = [paths(i, :), 1]; new_metrics(2*i) = metrics(i) + max(0, llr); end % 按度量排序，保留前 L 条 [sorted_metrics, sort_idx] = sort(new_metrics, 'ascend'); keep = sort_idx(1:min(L, length(sort_idx))); paths = new_paths(keep, :); metrics = sorted_metrics(1:min(L, length(sort_idx))); end

这段代码里L是列表大小，典型取值 4、8、16，L 越大性能越好但计算量线性增长。max(0, -llr)和max(0, llr)是路径度量增量的近似计算，比精确的 log-sum-exp 快很多，性能损失可以忽略。注意paths的列数在每次更新后增加 1，最终长度等于 N。SCL 译码结束后，选度量最小的路径，按info_indices提取信息比特。

验证 SCL 增益的方法：固定 N、K、构造方法和信噪比范围，分别跑 SC 和 SCL（L=8），对比 BER 曲线。如果 SCL 在 BER=10^-3 处比 SC 有 0.5 dB 以上的增益，说明实现正确；如果增益不明显甚至更差，检查路径度量更新是否在冻结比特位置错误地分裂了路径，或者排序保留逻辑写反了。从那以后我每次改译码算法，都会先用 N=8、K=4 的小规模用例手动验证每一步的路径和度量，确认无误后再上大码长跑仿真。希望帮到你。

本文还有配套的精品资源，点击获取