26、卷积码与Turbo码解码技术详解

卷积码与Turbo码解码技术详解

1. 硬判决解码与软判决解码对比

在加性高斯白噪声（AWGN）信道中，对卷积码的硬判决解码和软判决解码进行比较，发现与线性分组码的情况类似，软判决解码比硬判决解码有大约2 - 3 dB的性能优势。

2. Turbo码与迭代解码

Shannon的随机编码定理表明，要达到信道容量的码必须具有随机或接近随机的分布以及大的码块长度。然而，对于随机生成的码，解码复杂度会随着码块长度呈指数增长。因此，缺乏结构和大码块长度使得对达到容量的码进行最大似然（ML）解码变得不切实际。

Turbo编码是一种通过大长度的伪随机交织器连接两个简单码来生成具有接近随机结构和大码块长度的码的方法。由于生成的码基于简单码的组合，因此可以通过基于其组成码解码的迭代方案进行解码。这种解码方法称为迭代解码或Turbo解码，它不是最优的，但对于许多码，经过几次迭代后，其性能接近ML解码。

2.1 Turbo编码器

并行级联卷积码（PCCC），也称为Turbo码，由Berrou等人（1993）和Berrou和Glavieux（1996）提出。基本的Turbo编码器如图10.23所示，它并行使用两个递归系统卷积编码器，第二个编码器前有一个交织器。两个递归系统卷积编码器可以相同也可以不同。Turbo编码器的标称输出速率为Re = 1/3，但通过对二进制卷积编码器输出的奇偶校验位进行删余，可以实现更高的速率，如1/2或2/3。

交织器通常选择为块伪随机交织器，它在将信息序列的比特送入第二个编码器之前对其进行重新排序。使用两个递归卷积编码器和交织器会产生一个低重量码字很少的码。虽然这并不一定意味着级联码的自由距离特别大，但交织器与两个编码器的结合使得码字的最近邻相对较少，即码字相对稀疏。因此，Turbo码的编码增益部分得益于这一特性（即交织导致的最近邻码字数量的减少，称为重数）。

一个标准的Turbo码由组成码（通常相似）和交织模式（通常用TI表示）完全描述。组成码通常是非常简单的递归系统卷积码（RSCC）。

递归系统卷积码是系统卷积码，即输入比特直接作为编码比特的一部分出现，其中奇偶校验位由递归（反馈）线性滤波器生成。例如，一个速率为1/2的RSCC，输入比特直接作为编码比特的一部分，奇偶校验位由线性反馈移位寄存器生成。

组成码主要是速率为1/2的递归系统卷积码，通常用编码器前馈连接的八进制表示与反馈连接的八进制表示的比率来描述。例如，一个27/31 RSC编码器的反馈和前馈连接分别为31（二进制11001）和27（二进制10111）。

2.2 示例：517 RSCC编码

以下是一个MATLAB程序，用于生成当一般二进制序列由图10.24所示的517 RSCC编码时的奇偶校验位。

function [c_sys ,c_pc]=RSCC_57_Encoder(u);
% RSCC_57_Encoder  Encoder for 517 RSCC
% [c_sys,c_pc]=RSCC_57_Encoder(u)
% returns c_sys the systematic bits and 
% c_pc, the parity check bits of the code 
% when input is u and the encoder is 
% initiated at 0-state.
u = [0 1 1 1 0 0 1 0 0 1 1 0 0 1 0 0 1 1 1 1];
L = length(u);
I = 1;
% Initializing the values of the shift register:
r1 = 0;
r2 = 0;
r3 = 0;
while I <= L
    u_t = u(I);
    % Generating the systematic bits:
    c1(I) = u_t;
    % Updating the values of the shift register:
    rLt = mod(mod(r3 + r2,2) + u_t,2);
    r3 = r2;
    r2 = r1;
    r1 = rLt;
    % Generating the parity check bits:
    c2(I) = mod(r1 + r3,2);
    I = I + 1;
end
c_sys = c1;
c_pc = c2;

3. 卷积码的MAP解码 - BCJR算法

BCJR算法以Bahl、Cocke、Jelinek和Raviv（1974）的名字命名，是一种逐符号最大后验概率（MAP）卷积码解码算法。在该算法中，解码器使用MAP算法对每个输入符号进行解码，而不是寻找最可能的输入序列。

假设卷积编码器的状态集用S表示，由于编码器输入可以是0或1，从阶段i - 1到阶段i的状态转移（即从σi - 1∈S到σi∈S）可以是由于ui = 0或ui = 1。令Se表示对应于ui = -B（-B = 0, 1）的所有(σi - 1, σi)对的集合。

逐符号最大后验解码器接收解调器输出y = (y1, Y2, … , YN)，并根据此观测值使用最大后验规则对ui进行解码：

[
u_i = \arg\max_{u_i\in{0,1}}P(u_i|Y) = \arg\max_{u_i\in{0,1}}\frac{P(u_i,Y)}{P(Y)} = \arg\max_{u_i\in{0,1}}P(u_i,Y) = \arg\max_{-B\in{0,1}}\sum_{(\sigma_{i - 1},\sigma_i)\in S_{-B}}p(\sigma_{i - 1},\sigma_i,y)
]

为了进行最大后验解码，需要计算αi - 1(σi - 1)、βi(σi)和γi(σi - 1, σi)的值。

3.1 前向递归计算αi(σi)

αi(σi)可以通过以下前向递归计算：

[
\alpha_i(\sigma_i) = \sum_{\sigma_{i - 1}\in S}\gamma_i(\sigma_{i - 1},\sigma_i)\alpha_{i - 1}(\sigma_{i - 1}), 1\leq i\leq N
]

如果假设网格图从全零状态开始，前向递归的初始条件为：

[
\alpha_0(\sigma_0) =
\begin{cases}
1, & \sigma_0 = 0 \
0, & \sigma_0 \neq 0
\end{cases}
]

以下是一个MATLAB脚本，用于实现517 RSC编码器的前向递归方程：

function alpha=forward_recursion(gamma);
% FORWARD_RECURSION computing alpha for 517 RSCC
% alpha=forward_recursion( gamma);
% returns alpha in the form of a matrix.
% gamma is a 16XN matrix of gamma_i(sigma_(i-1),sigm_i)
% Assuming gamma is given
N = size(gamma,2);
Ns = 4;
% Initialization:
alpha = zeros(Ns,N);
alpha_0 = 1;
i = 1;
% Time index
simga_i = [1 3];
% Set of states at i=1
alpha(simga_i(1),i) = gamma(1,i);
alpha(simga_i(2),i) = gamma(3,i);
i = 2;
simga_i = [1 2 3 4];
% Set of states at i=2
alpha(simga_i(1),i) = gamma(1,i) * alpha(1,i - 1);
alpha(simga_i(2),i) = gamma(10,i) * alpha(3,i - 1);
alpha(simga_i(3),i) = gamma(3,i) * alpha(1,i - 1);
alpha(simga_i(4),i) = gamma(12,i) * alpha(3,i - 1);
for i = 3:N - 2
    alpha(simga_i(1),i) = gamma(1,i) * alpha(1,i - 1) + gamma(5,i) * alpha(2,i - 1);
    alpha(simga_i(2),i) = gamma(10,i) * alpha(3,i - 1) + gamma(14,i) * alpha(4,i - 1);
    alpha(simga_i(3),i) = gamma(3,i) * alpha(1,i - 1) + gamma(7,i) * alpha(2,i - 1);
    alpha(simga_i(4),i) = gamma(12,i) * alpha(3,i - 1) + gamma(16,i) * alpha(4,i - 1);
end
i = N - 1;
simga_i = [1 2];
% Set of states at i=N-1
alpha(simga_i(1),i) = gamma(1,i) * alpha(1,i - 1) + gamma(5,i) * alpha(2,i - 1);
alpha(simga_i(2),i) = gamma(10,i) * alpha(3,i - 1) + gamma(14,i) * alpha(4,i - 1);
i = N;
simga_i = 1;
% Set of states at i=N
alpha(simga_i(1),i) = gamma(1,i) * alpha(1,i - 1) + gamma(5,i) * alpha(2,i - 1);
alpha = [[1 0 0 0]', alpha];

3.2 反向递归计算βi(σi)

βi - 1(σi - 1)可以通过以下反向递归计算：

[
\beta_{i - 1}(\sigma_{i - 1}) = \sum_{\sigma_i\in S}\beta_i(\sigma_i)\gamma_i(\sigma_{i - 1},\sigma_i), 1\leq i\leq N
]

如果假设网格图终止于全零状态，反向递归的边界条件为：

[
\beta_N(\sigma_N) =
\begin{cases}
1, & \sigma_N = 0 \
0, & \sigma_N \neq 0
\end{cases}
]

以下是一个MATLAB脚本，用于实现517 RSC编码器的反向递归方程：

function beta=backward_recursion(gamma);
% BACKWARD_RECURSION computing beta for 517 RSCC
% beta= backward_recursion( gamma);
% beta in the form of a matrix
% gamma is a 16XN matrix of gamma_i(sigma_(i-1),sigm_i)
% Assuming gamma is given
N = size(gamma,2);
Ns = 4;
% Initialization:
beta = zeros(Ns,N);
beta(1,N) = 1;
i = N;
% Time index
simga_i_1 = [1 2];
% Set of states at i=N
beta(simga_i_1(1),i - 1) = gamma(1,i);
beta(simga_i_1(2),i - 1) = gamma(5,i);
i = N - 1;
simga_i_1 = [1 2 3 4];
% Set of states at i=N-1
beta(simga_i_1(1),i - 1) = gamma(1,N) * gamma(1,i);
beta(simga_i_1(2),i - 1) = gamma(1,N) * gamma(5,i);
beta(simga_i_1(3),i - 1) = gamma(5,N) * gamma(10,i);
beta(simga_i_1(4),i - 1) = gamma(5,N) * gamma(14,i);
for i = N - 2:-1:3
    beta(simga_i_1(1),i - 1) = beta(1,i) * gamma(1,i) + beta(3,i) * gamma(3,i);
    beta(simga_i_1(2),i - 1) = beta(1,i) * gamma(5,i) + beta(3,i) * gamma(7,i);
    beta(simga_i_1(3),i - 1) = beta(2,i) * gamma(10,i) + beta(4,i) * gamma(12,i);
    beta(simga_i_1(4),i - 1) = beta(4,i) * gamma(16,i) + beta(2,i) * gamma(14,i);
end
i = 2;
simga_i_1 = [1 3];
% Set of states at i=2
beta(simga_i_1(1),i - 1) = beta(1,i) * gamma(1,i) + beta(3,i) * gamma(3,i);
beta(simga_i_1(2),i - 1) = beta(2,i) * gamma(10,i) + beta(4,i) * gamma(12,i);
i = 1;
simga_i_1 = 1;
% Set of states at i=1
beta_0(simga_i_1(1)) = beta(1,i) * gamma(1,i) + beta(3,i) * gamma(3,i);

3.3 计算γi(σi - 1, σi)

γi(σi - 1, σi)可以表示为：

[
\gamma_i(\sigma_{i - 1},\sigma_i) = P(\sigma_i,Y_i|\sigma_{i - 1}) = p(\sigma_i|\sigma_{i - 1})p(Y_i|\sigma_i,\sigma_{i - 1}) = P(u_i)p(Y_i|u_i) = P(u_i)P(Y_i|c_i)
]

其中，我们利用了状态对(σi - 1, σi)与输入ui之间存在一一对应的事实。上述表达式清楚地表明了γi(σi - 1, σi)对信息序列在时间i的先验概率P(ui)以及p(Yi|ci)（取决于信道特性）的依赖关系。如果信息序列是等概率的（在没有信息可用时通常这样假设），则P(ui = 0) = P(ui = 1) = 1/2。

对于AWGN信道，y = c + n，其中c表示对应于编码序列的调制信号。在n = 2，卷积码是系统码且调制为BPSK的特殊情况下，直接代入可得：

[
\gamma_i(\sigma_{i - 1},\sigma_i) = \frac{1}{\pi N_0}\exp\left(-\frac{(y_s^i)^2 + (y_p^i)^2 + 2E}{N_0}\right)P(u_i)\exp\left(\frac{2y_s^i c_s^i + 2y_p^i c_p^i}{N_0}\right)
]

其中，上标s和p表示y和c的系统分量和奇偶校验分量。

以下是一个MATLAB脚本，用于计算当使用517 RSC码和BPSK调制的AWGN信道时γi的值：

y = [1.2 -0.8 0.3 .02 -0.7 -0.02 2 0 -1 1];
Ns = 4;
% Number of states
E = 1;
% Energy per symbol
EbNO_dB = 1;
% SNR per bit (dB)
% Computing gamma:
gamma = gamma_calc(y,E,EbNO_dB);
% y is the channel output

function gamma = gamma_calc(y,E,EbNO_dB)
%GAMMA_CALC Computes the gamma matrix for a 517 RSCC over AWGN using BPSK.
% GAMMA= GAMMA_CALC(y,E,EbNO_dB)
% y is channle output, E is symbol energy, and EbNO_dB is SNR/bit (in dB).
M = length(y);
N = M/2;
% Length of channel output
% Depth of Trellis
Ns = 4;
% Numner of states
y_s = y(1:2:M - 1);
% Outputs corresponding to systematic bits
y_p = y(2:2:M);
% Outputs corresponding to parity-check bits
% Generation of the 517 RSCC trellis
[c_s c_p] = Rscc_5_7(N);
Eb = 2*E;
NO = Eb*10^(-EbNO_dB/10);
% Noise variance
gamma = zeros(Ns^2,N);
% Initialization of gamma
for i = 1:N
    % Time index
    gamma(1,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
        *exp((2*(y_s(i)*c_s(1,i) + y_p(i)*c_p(1,i)))/NO);
    gamma(3,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
        *exp((2*(y_s(i)*c_s(3,i) + y_p(i)*c_p(3,i)))/NO);
    if i > 1
        gamma(10,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
            *exp((2*(y_s(i)*c_s(10,i) + y_p(i)*c_p(10,i)))/NO);
        gamma(12,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
            *exp((2*(y_s(i)*c_s(12,i) + y_p(i)*c_p(12,i)))/NO);
    end
    if i > 2
        gamma(5,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
            *exp((2*(y_s(i)*c_s(5,i) + y_p(i)*c_p(5,i)))/NO);
        gamma(7,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
            *exp((2*(y_s(i)*c_s(7,i) + y_p(i)*c_p(7,i)))/NO);
        gamma(14,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
            *exp((2*(y_s(i)*c_s(14,i) + y_p(i)*c_p(14,i)))/NO);
        gamma(16,i) = 1/(2*pi*NO)*exp(-(y_s(i)^2 + y_p(i)^2 + 2*E)/NO) ...
            *exp((2*(y_s(i)*c_s(16,i) + y_p(i)*c_p(16,i)))/NO);
    end
end

function [c_s c_p] = Rscc_5_7(N)
%RSCC_5_7 Generates a trellis of depth N for a 517 RSCC.
% [C_S C_P] = RSCC_5_7(N)
% N is the depth of the trellis
% C_S and c_p are systematic and parity check bits.
Ns = 4;
% Number of states
c_s = zeros(Ns^2,N);
% Initiation of the matrix of systematic bits
c_p = zeros(Ns^2,N);
% Initiation of the matrix of parity check bits
for i = 1:N
    % Time index
    c_s(1,i) = 0; c_p(1,i) = 0;
    c_s(3,i) = 1; c_p(3,i) = 1;
    if i > 1
        c_s(10,i) = 0; c_p(10,i) = 0;
        c_s(12,i) = 1; c_p(12,i) = 1;
    end
    if i > 2
        c_s(5,i) = 0;
        c_p(5,i) = 1;
        c_s(7,i) = 1;
        c_p(7,i) = 0;
        c_s(14,i) = 1; c_p(14,i) = 1;
        c_s(16,i) = 0; c_p(16,i) = 0;
    end
end

3.4 后验L值计算

通过BCJR算法可以得到后验L值：

[
L(u_i) = \ln\frac{P(u_i = 1|Y)}{P(u_i = 0|Y)} = \ln\frac{\sum_{(\sigma_{i - 1},\sigma_i)\in S_1}\alpha_{i - 1}(\sigma_{i - 1})\gamma_i(\sigma_{i - 1},\sigma_i)\beta_i(\sigma_i)}{\sum_{(\sigma_{i - 1},\sigma_i)\in S_0}\alpha_{i - 1}(\sigma_{i - 1})\gamma_i(\sigma_{i - 1},\sigma_i)\beta_i(\sigma_i)}
]

这些后验L值也被称为软输出，对于Turbo码的解码至关重要。接受软输入（向量y）并生成软输出的解码器称为软输入软输出（SISO）解码器。基于L(ui)软值的解码规则为：

[
\hat{u}_i =
\begin{cases}
1, & L(u_i) > 0 \
0, & L(u_i) \leq 0
\end{cases}
]

4. 对数域的BCJR算法

上述BCJR算法的一个问题是它在数值上不稳定，特别是当网格图长度较长时。一种替代方法是其对数域版本，即对数APP（对数后验概率）算法或对数MAP算法。

在对数APP算法中，定义α、β和γ的对数：

[
\hat{\alpha} i(\sigma_i) = \ln(\alpha_i(\sigma_i))
]
[
\hat{\beta}_i(\sigma_i) = \ln(\beta_i(\sigma_i))
]
[
\hat{\gamma}_i(\sigma {i - 1},\sigma_i) = \ln(\gamma_i(\sigma_{i - 1},\sigma_i))
]

可以得到以下前向和反向递归关系：

[
\hat{\alpha} i(\sigma_i) = \ln\left(\sum {\sigma_{i - 1}\in S}\exp(\hat{\alpha} {i - 1}(\sigma {i - 1}) + \hat{\gamma} i(\sigma {i - 1},\sigma_i))\right)
]
[
\hat{\beta} {i - 1}(\sigma {i - 1}) = \ln\left(\sum_{\sigma_i\in S}\exp(\hat{\beta} i(\sigma_i) + \hat{\gamma}_i(\sigma {i - 1},\sigma_i))\right)
]

初始条件为：

[
\hat{\alpha}_0(\sigma_0) =
\begin{cases}
0, & \sigma_0 = 0 \
-\infty, & \sigma_0 \neq 0
\end{cases}
]

后验I值计算为：

[
L(u_i) = \ln\left[\sum_{(\sigma_{i - 1},\sigma_i)\in S_1}\exp(\hat{\alpha} {i - 1}(\sigma {i - 1}) + \hat{\gamma} i(\sigma {i - 1},\sigma_i) + \hat{\beta} i(\sigma_i))\right] - \ln\left[\sum {(\sigma_{i - 1},\sigma_i)\in S_0}\exp(\hat{\alpha} {i - 1}(\sigma {i - 1}) + \hat{\gamma} i(\sigma {i - 1},\sigma_i) + \hat{\beta}_i(\sigma_i))\right]
]

为了提高计算效率，引入以下符号：

[
\max^ \left{x,y\right} \approx \ln(e^x + e^y)
]
[
\max^ \left{x,y,z\right} \approx \ln(e^x + e^y + e^z)
]

则递归关系变为：

[
\hat{\alpha} i(\sigma_i) = \max^ \left{\hat{\alpha} {i - 1}(\sigma {i - 1}) + \hat{\gamma} i(\sigma {i - 1},\sigma_i)\right} {\sigma {i - 1}\in S}
]
[
\hat{\beta} {i - 1}(\sigma {i - 1}) = \max^ \left{\hat{\beta}_i(\sigma_i) + \hat{\gamma}_i(\sigma {i - 1},\sigma_i)\right}_{\sigma_i\in S}
]

后验I值为：

[
L(u_i) = \frac{\max^ \left{\hat{\alpha} {i - 1}(\sigma {i - 1}) + \hat{\gamma} i(\sigma {i - 1},\sigma_i) + \hat{\beta} i(\sigma_i)\right} {(\sigma_{i - 1},\sigma_i)\in S_1}}{\max^ \left{\hat{\alpha} {i - 1}(\sigma {i - 1}) + \hat{\gamma} i(\sigma {i - 1},\sigma_i) + \hat{\beta} i(\sigma_i)\right} {(\sigma_{i - 1},\sigma_i)\in S_0}}
]

当x和y（或x, y, z）的值相差较大时，可以使用近似关系：

[
\max^ \left{x,y\right} \approx \max\left{x,y\right}
]
[
\max^ \left{x,y,z\right} \approx \max\left{x,y,z\right}
]

这种近似会导致性能略有下降，得到的算法称为max - Log - APP算法（或max - log - MAP算法）。也可以使用查找表来计算修正项(\ln(1 + e^{-|x - y|}))以加速算法。

以下是一个MATLAB脚本，用于实现max*操作：

function x_s = max_s(x)
%MAX_S Calculates the max* of the elements of vector x.
% X_S = MAX_S(X) calculats the max* of the elements of vetor x by
% recursively performing the max* operation on two components of
% vector x at each recursion.
L = length(x);
x_s = x(1);
if L > 1
    for l = 2:L
        x_s = max_s_2(x_s,x(l));
    end
end

function c = max_s_2(a,b)
%MAX_S_2 Calculates the max* of the two numbers.
% C = MAX_S_2(A,B) calculats the max* of a and b.
c = max(a,b) + log(1 + exp(-abs(a - b)));

以下是一个MATLAB脚本，用于计算当使用517 RSC码和BPSK调制的AWGN信道时γ、α和β的对数域值：

y = [1.2 -0.8 0.3 .02 -0.7 -0.02 2 0 -1 1];
Ns = 4;
% Number of states
E = 1;
% Energy per bit
EbNO_dB = 10;
% SNR/bit in dB
% Calculation of gamma, alpha, and beta in the logarithmic domain:
gamma = gamma_calc(y,E,EbNO_dB);
gamma_t = log(gamma);
[alpha_t alpha_LO] = forward_recursion_log(gamma_t);
beta_t = backward_recursion_log(gamma_t);

function [alpha_t alpha_LO] = forward_recursion_log(gamma_t)
%FORWARD_RECURSION_LOG Calculates alpha in the logarithmic domain for a 517 RSCC.
% [ALPHA_T ALPHA_LO] = FORWARD_RECURSION_LOG(GAMMA_T)
% gamma_t: gamma in the logarithmic domain
% alpha_t: alpha in logarithmic domain.
% Initialization:
Ns = sqrt(size(gamma_t,1)); % Number of states
N = size(gamma_t,2);
alpha_t = ones(Ns,N)*-inf;
alpha_LO = 0;
i = 1;
% Time index
simga_i = [1 3];
% Set of states at i=1
alpha_t(simga_i(1),i) = alpha_LO + gamma_t(1,i);
alpha_t(simga_i(2),i) = alpha_LO + gamma_t(3,i);
i = 2;
simga_i = [1 2 3 4];
% Set of states at i=2
alpha_t(simga_i(1),i) = gamma_t(1,i) + alpha_t(1,i - 1);
alpha_t(simga_i(2),i) = gamma_t(10,i) + alpha_t(3,i - 1);
alpha_t(simga_i(3),i) = gamma_t(3,i) + alpha_t(1,i - 1);
alpha_t(simga_i(4),i) = gamma_t(12,i) + alpha_t(3,i - 1);
for i = 3:N - 2
    alpha_t(simga_i(1),i) = max_s([gamma_t(1,i) + alpha_t(1,i - 1) gamma_t(5,i) + alpha_t(2,i - 1)]);
    alpha_t(simga_i(2),i) = max_s([gamma_t(10,i) + alpha_t(3,i - 1) gamma_t(14,i) + alpha_t(4,i - 1)]);
    alpha_t(simga_i(3),i) = max_s([gamma_t(3,i) + alpha_t(1,i - 1) gamma_t(7,i) + alpha_t(2,i - 1)]);
    alpha_t(simga_i(4),i) = max_s([gamma_t(12,i) + alpha_t(3,i - 1) gamma_t(16,i) + alpha_t(4,i - 1)]);
end
i = N - 1;
simga_i = [1 2];
% Set of states at i=N-1
alpha_t(simga_i(1),i) = max_s([gamma_t(1,i) + alpha_t(1,i - 1) gamma_t(5,i) + alpha_t(2,i - 1)]);
alpha_t(simga_i(2),i) = max_s([gamma_t(10,i) + alpha_t(3,i - 1) gamma_t(14,i) + alpha_t(4,i - 1)]);
i = N;
simga_i = 1;
% Set of states at i=N
alpha_t(simga_i(1),i) = max_s([gamma_t(1,i) + alpha_t(1,i - 1) gamma_t(5,i) + alpha_t(2,i - 1)]);

function beta_t = backward_recursion_log(gamma_t)
%BACKWARD_RECURSION_LOG Calculates beta's in the logarithmic domain for a 517 RSCC.
% BETA_T = BACKWARD_RECURSION_LOG(GAMMA_T)
% gamma_t: gamma in the logarithmic domain
% beta_t: beta in logarithmic domain.
% Initialization:
Ns = sqrt(size(gamma_t,1)); % Number of states
N = size(gamma_t,2);
beta_t = ones(Ns,N)*-inf;
beta_t(1,N) = 0;
i = N;
% Time index
simga_i_1 = [1 2];
% Set of states at i=N
beta_t(simga_i_1(1),i - 1) = beta_t(1,N) + gamma_t(1,i);
beta_t(simga_i_1(2),i - 1) = beta_t(1,N) + gamma_t(5,i);
i = N - 1;
simga_i_1 = [1 2 3 4];
% Set of states at i=N-1
beta_t(simga_i_1(1),i - 1) = gamma_t(1,N) + gamma_t(1,i);
beta_t(simga_i_1(2),i - 1) = gamma_t(1,N) + gamma_t(5,i);
beta_t(simga_i_1(3),i - 1) = gamma_t(5,N) + gamma_t(10,i);
beta_t(simga_i_1(4),i - 1) = gamma_t(5,N) + gamma_t(14,i);
for i = N - 2:-1:3
    beta_t(simga_i_1(1),i - 1) = max_s([beta_t(1,i) + gamma_t(1,i) beta_t(3,i) + gamma_t(3,i)]);
    beta_t(simga_i_1(2),i - 1) = max_s([beta_t(1,i) + gamma_t(5,i) beta_t(3,i) + gamma_t(7,i)]);
    beta_t(simga_i_1(3),i - 1) = max_s([beta_t(2,i) + gamma_t(10,i) beta_t(4,i) + gamma_t(12,i)]);
    beta_t(simga_i_1(4),i - 1) = max_s([beta_t(4,i) + gamma_t(16,i) beta_t(2,i) + gamma_t(14,i)]);
end
i = 2;
simga_i_1 = [1 3];
% Set of states at i=2
beta_t(simga_i_1(1),i - 1) = max_s([beta_t(1,i) + gamma_t(1,i) beta_t(3,i) + gamma_t(3,i)]);
beta_t(simga_i_1(2),i - 1) = max_s([beta_t(2,i) + gamma_t(10,i) beta_t(4,i) + gamma_t(12,i)]);
i = 1;
simga_i_1 = 1;
% Set of states at i=1
beta_LO(simga_i_1(1)) = max_s([beta_t(1,i) + gamma_t(1,i) beta_t(3,i) + gamma_t(3,i)]);

5. Turbo码的迭代解码

由于Turbo码的网格图状态数量巨大，对其进行最优解码是不可能的。Berrou等人（1993）提出了一种次优的迭代解码算法，即Turbo解码算法，该算法能实现非常接近Shannon理论界限的优异性能。

Turbo解码算法基于对数APP或max - log - APP算法的迭代使用。对于速率为1/2的RSCC，对数似然比（LLR）可以写成三个项的和：

[
L(u_i) = L_cy_i^s + L^{(a)}(u_i) + L^{(e)}(u_i)
]

其中，(L_cy_i^s)称为信道I值，表示信道输出中对应于系统比特的影响；(L^{(a)}(u_i))是先验I值，是信息序列先验概率的函数；(L^{(e)}(u_i))表示外在I值或外在信息，是后验I值的一部分，不依赖于先验概率和信道输出的系统信息。

假设二进制信息序列(u = (u_1, U_2, … , U_N))输入到第一个速率为1/2的RSCC，输出的奇偶校验位为(c^P = (c_1^P, c_2^P, … , c_N^P))。信息序列通过交织器得到(u’ = (u_1’, u_2’, … , u_N’))，并输入到第二个编码器生成奇偶校验序列(c’^P = (c_1’^P, c_2’^P, … , c_N’^P))。序列(u)、(c^P)和(c’^P)经过BPSK调制后在高斯信道上传输，对应的输出序列分别为(y^s)、(y^P)和(y’^P)。

第一个组成码的MAP解码器接收对((y^s, y^P))。在第一次迭代中，解码器假设所有比特等概率，因此先验I值设为零。解码器使用对数APP或max - log - APP算法计算后验L值。在第一个组成解码器的输出端，解码器从后验L值中减去信道I值以计算外在I值。这些值用(L_1^{(e)}(u_i))表示，并通过交织器进行置换，然后作为第二个组成解码器的先验I值。第二个解码器还接收(y’^P)和经过交织器置换后的(y^s)。第二个解码器计算外在I值(L_2^{(e)}(u_i))，并通过交织器的逆置换将其提供给第一个解码器，第一个解码器在下一次迭代中使用这些值作为先验I值。这个过程持续进行，直到达到固定的迭代次数或满足某个准则。最后一次迭代后，使用后验I值(L(u_i))进行最终决策。

Turbo解码器的构建块是一个SISO解码器，其输入为(y^s)、(y^P)和(L^{(a)}(u_i))，输出为(L^{(e)}(u_i))和(L(u_i))。在迭代解码中，(L^{(a)}(u_i))由另一个解码器提供的外在I值代替。Turbo解码器的框图如下：

graph LR
    classDef startend fill:#F5EBFF,stroke:#BE8FED,stroke-width:2px
    classDef process fill:#E5F6FF,stroke:#73A6FF,stroke-width:2px
    classDef decision fill:#FFF6CC,stroke:#FFBC52,stroke-width:2px

    A(y^P):::process --> B(SISO 1):::process
    C(y^s):::process --> B
    B --> D(L(u_i)):::process
    B --> E(L_1^{(e)}(u_i)):::process
    E --> F(交织器):::process
    F --> G(SISO 2):::process
    H(y'^P):::process --> G
    I(y^s经交织器):::process --> G
    G --> J(L(u_i)):::process
    G --> K(L_2^{(e)}(u_i)):::process
    K --> L(交织器逆置换):::process
    L --> B

以下是一个MATLAB脚本，用于计算当使用517 RSC码和BPSK调制的AWGN信道时的外在信息(L_e(u_i))和LLR值(L(u_i))：

y = [1.2 -0.8 0.3 .02 -0.7 -0.02 2 0 -1 1];
E = 1;
% Bit energy
EbNO_dB = 10;
% SNR/bit in dB
P = ones(1,length(y)/2)/2;
% A priori probability of each bit initialized to 1/2
[Le L] = extrinsic(y,E,EbNO_dB,P);

function [Le L] = extrinsic(y,E,EbNO_dB,P)
%EXTRINSIC Generates the extrinsic information and output LLR's for a 517 RSCC
% with BPSK modulation over an AWGN channel.
% [L_e L] = extrinsic(y,E,EbNO_dB,P)
% y: channel output sequence
% E: symbol energy
% EbNO_dB: SNR/bit in dB
% P: A priori probabilities of input bits
M = length(y);
% Length of channel output sequence
N = M/2;
% Depth of Trellis
Ns = 4;
% Number of states
y_s = y(1:2:M - 1);
% Systematic components of y
y_p = y(2:2:M);
% Parity-check components of y
[c_s c_p] = Rscc_5_7(N);
% Generation of the 517 RSCC trellis
Eb = 2*E;
NO = Eb*10^(-EbNO_dB/10);
% Noise variance
gamma = gamma_calc(y,E,EbNO_dB);
gamma_t = log(gamma);
[alpha_t alpha_O_t] = forward_recursion_log(gamma_t);
beta_t = backward_recursion_log(gamma_t);
% Initialization for speed:
Le = zeros(1,N);
L = zeros(1,N);
for i = 1:N
    % Time index
    if i == 1
        ul = alpha_O_t + 2*y_p(i)*c_p(3,i)/NO + beta_t(3,i);
        vl = alpha_O_t + 2*y_p(i)*c_p(1,i)/NO + beta_t(1,i);
        Le(i) = max_s(ul) - max_s(vl);
    elseif i == 2
        ul = alpha_t(1,i - 1) + 2*y_p(i)*c_p(3,i)/NO + beta_t(3,i);
        u2 = alpha_t(3,i - 1) + 2*y_p(i)*c_p(12,i)/NO + beta_t(4,i);
        vl = alpha_t(1,i - 1) + 2*y_p(i)*c_p(1,i)/NO + beta_t(1,i);
        v2 = alpha_t(3,i - 1) + 2*y_p(i)*c_p(10,i)/NO + beta_t(2,i);
        Le(i) = max_s([ul u2]) - max_s([vl v2]);
    elseif i == N - 1
        ul = alpha_t(4,i - 1) + 2*y_p(i)*c_p(14,i)/NO + beta_t(2,i);
        vl = alpha_t(1,i - 1) + 2*y_p(i)*c_p(1,i)/NO + beta_t(1,i);
        v2 = alpha_t(2,i - 1) + 2*y_p(i)*c_p(5,i)/NO + beta_t(1,i);
        v3 = alpha_t(3,i - 1) + 2*y_p(i)*c_p(10,i)/NO + beta_t(2,i);
        Le(i) = max_s(ul) - max_s([vl v2 v3]);
    elseif i == N
        vl = alpha_t(1,i - 1) + 2*y_p(i)*c_p(1,i)/NO + beta_t(1,i);
        v2 = alpha_t(2,i - 1) + 2*y_p(i)*c_p(5,i)/NO + beta_t(1,i);
        Le(i) = -max_s([vl v2]);
    else
        ul = alpha_t(1,i - 1) + 2*y_p(i)*c_p(3,i)/NO + beta_t(3,i);
        u2 = alpha_t(2,i - 1) + 2*y_p(i)*c_p(7,i)/NO + beta_t(3,i);
        u3 = alpha_t(3,i - 1) + 2*y_p(i)*c_p(12,i)/NO + beta_t(4,i);
        u4 = alpha_t(4,i - 1) + 2*y_p(i)*c_p(14,i)/NO + beta_t(2,i);
        vl = alpha_t(1,i - 1) + 2*y_p(i)*c_p(1,i)/NO + beta_t(1,i);
        v2 = alpha_t(2,i - 1) + 2*y_p(i)*c_p(5,i)/NO + beta_t(1,i);
        v3 = alpha_t(3,i - 1) + 2*y_p(i)*c_p(10,i)/NO + beta_t(2,i);
        v4 = alpha_t(4,i - 1) + 2*y_p(i)*c_p(16,i)/NO + beta_t(4,i);
        Le(i) = max_s([ul u2 u3 u4]) - max_s([vl v2 v3 v4]);
    end
    % 计算L值
    L(i) = Le(i) + 2*y_s(i)*c_s(3,i)/NO;
end

综上所述，Turbo码及其迭代解码算法在信道编码领域具有重要的地位，通过上述的理论分析和MATLAB代码实现，我们可以更深入地理解其工作原理和实现过程。在实际应用中，可以根据具体需求选择合适的参数和算法，以达到最佳的性能。

卷积码与Turbo码解码技术详解

6. 关键概念总结

为了更好地理解上述内容，下面对一些关键概念和参数进行总结，形成表格如下：
| 概念 | 解释 |
| — | — |
| 硬判决解码 | 对接收信号进行简单的0或1判决，性能相对软判决解码较差 |
| 软判决解码 | 利用接收信号的更多信息进行解码，在AWGN信道中比硬判决解码有2 - 3 dB性能优势 |
| Turbo码 | 通过伪随机交织器连接两个简单码生成，具有接近随机结构和大码块长度 |
| 迭代解码 | 基于组成码解码的迭代方案，性能接近最大似然解码 |
| BCJR算法 | 逐符号最大后验概率解码算法，用于卷积码解码 |
| 对数APP算法 | BCJR算法的对数域版本，提高数值稳定性 |
| max - Log - APP算法 | 对数APP算法的近似版本，计算效率更高 |
| 软输入软输出（SISO）解码器 | 接受软输入并生成软输出，对Turbo码解码至关重要 |
| 信道I值 (L_cy_i^s) | 表示信道输出中对应于系统比特的影响 |
| 先验I值 (L^{(a)}(u_i)) | 信息序列先验概率的函数 |
| 外在I值 (L^{(e)}(u_i)) | 后验I值的一部分，不依赖先验概率和系统信息 |

7. 不同算法性能对比

为了更直观地展示不同解码算法的性能差异，我们可以从复杂度和性能两个维度进行对比，如下表所示：
| 算法 | 复杂度 | 性能 |
| — | — | — |
| 硬判决解码 | 低 | 较差，比软判决解码低2 - 3 dB |
| 软判决解码 | 中 | 较好，在AWGN信道中有优势 |
| BCJR算法 | 高 | 最优解码，但数值不稳定，长码时复杂度高 |
| 对数APP算法 | 中 | 数值稳定，复杂度适中 |
| max - Log - APP算法 | 低 | 近似算法，性能略有下降，但计算效率高 |
| Turbo解码算法 | 中 | 接近Shannon理论界限，性能优异 |

从表中可以看出，不同算法在复杂度和性能上各有优劣。在实际应用中，需要根据具体的场景和需求来选择合适的算法。例如，对于对复杂度要求较高、对性能要求不是特别苛刻的场景，可以选择max - Log - APP算法；而对于对性能要求极高、对复杂度不太敏感的场景，则可以考虑使用BCJR算法。

8. 实际应用中的考虑因素

在实际应用Turbo码和相关解码算法时，还需要考虑以下几个方面：
- 信道特性 ：不同的信道特性（如AWGN信道、衰落信道等）对解码算法的性能有显著影响。在AWGN信道中，上述算法可以较好地发挥作用，但在衰落信道中，可能需要结合其他技术，如分集技术、信道估计技术等，来提高系统的性能。
- 码长和码率 ：码长和码率的选择会影响Turbo码的性能和复杂度。一般来说，码长越长，编码增益越大，但解码复杂度也会相应增加；码率越高，传输效率越高，但纠错能力会降低。因此，需要根据具体的应用场景和需求来选择合适的码长和码率。
- 交织器设计 ：交织器的设计对Turbo码的性能至关重要。一个好的交织器可以有效地减少码字的相关性，提高码的重量分布，从而提高编码增益。在设计交织器时，需要考虑交织器的长度、交织模式等因素。
- 迭代次数 ：在Turbo解码算法中，迭代次数的选择会影响解码性能和复杂度。一般来说，迭代次数越多，解码性能越好，但解码复杂度也会相应增加。在实际应用中，需要根据具体的场景和需求来选择合适的迭代次数，以达到性能和复杂度的平衡。

9. 未来发展趋势

随着通信技术的不断发展，Turbo码及其解码技术也在不断演进。未来，可能会朝着以下几个方向发展：
- 更高性能的码设计 ：研究人员将继续探索新的码结构和编码方法，以进一步提高Turbo码的性能，使其更接近Shannon理论界限。
- 更低复杂度的解码算法 ：为了满足实际应用中对低复杂度的需求，将开发更高效的解码算法，如基于硬件加速的解码算法、近似解码算法等。
- 与其他技术的融合 ：将Turbo码与其他通信技术，如MIMO技术、OFDM技术等相结合，以提高系统的整体性能和频谱效率。
- 适应新的应用场景 ：随着物联网、5G等新兴技术的发展，对通信系统的可靠性、实时性等方面提出了更高的要求。Turbo码及其解码技术将需要适应这些新的应用场景，进行相应的优化和改进。

10. 总结

本文详细介绍了卷积码的硬判决解码和软判决解码的对比，深入探讨了Turbo码的编码原理、迭代解码算法以及相关的解码技术，包括BCJR算法、对数APP算法和max - Log - APP算法等。通过理论分析和MATLAB代码实现，我们可以更深入地理解这些算法的工作原理和实现过程。

在实际应用中，需要根据具体的场景和需求选择合适的解码算法和参数，以达到性能和复杂度的平衡。同时，随着通信技术的不断发展，Turbo码及其解码技术也将不断演进，为未来的通信系统提供更可靠、高效的编码解决方案。

希望本文能够对读者理解卷积码和Turbo码的解码技术有所帮助，也欢迎读者在实际应用中进行更多的探索和实践。

以下是一个简单的mermaid流程图，展示了Turbo码解码的整体流程：

graph LR
    classDef startend fill:#F5EBFF,stroke:#BE8FED,stroke-width:2px
    classDef process fill:#E5F6FF,stroke:#73A6FF,stroke-width:2px
    classDef decision fill:#FFF6CC,stroke:#FFBC52,stroke-width:2px

    A(接收信号):::process --> B(初始化先验I值):::process
    B --> C(第一个SISO解码器):::process
    C --> D(计算外在I值):::process
    D --> E(交织器):::process
    E --> F(第二个SISO解码器):::process
    F --> G(计算外在I值):::process
    G --> H(交织器逆置换):::process
    H --> C(第一个SISO解码器):::process
    I{是否达到迭代次数或准则}:::decision -->|是| J(输出最终判决):::process
    I -->|否| C

这个流程图清晰地展示了Turbo码迭代解码的过程，从接收信号开始，经过多次迭代，直到达到迭代次数或满足某个准则后输出最终判决。通过这个流程图，我们可以更直观地理解Turbo码解码的工作原理。