LDPC Coding with PSK modulation

作者: 程小厮 | 来源:发表于2019-03-18 23:03 被阅读0次

该LDPC软件旨在作为LDPC码基于计算机的仿真的介绍。伪随机不规则低密度奇偶校验矩阵基于Radford M. Neal的程序集，可以在[1]中找到。虽然Neal的收藏品已有详细记载，但在我看来，C源代码仍然势不可挡，特别是如果您不熟悉C语言。我的软件是为MATLAB编写的，它比C语言更具可读性。您可能还想参考[2]中另一个基于MATLAB的LDPC源代码，它具有不同的代码编写风格（实际上Arun在他的日志中有错误） - 似然解码器）。

创建LDPC矩阵

创建LDPC矩阵的步骤（来自Neal的网站）：

选择创建稀疏奇偶校验矩阵的方法之一：（1）Evencol，或（2）Evenboth。

将1添加到没有1或只有1的行中。

如果需要，尽量避免长度为4的循环，其中一对柱在特定行中都有1

Evencol在列之间随机分配相同数量的1，甚至类似于evencol，但也尝试使行具有相同数量的1。

函数 makeLdpc.m 接受5个参数：

M：行数

N：列数

方法：分配1s的方法，0 = Evencol，1 = Evenboth

noCyle：消除长度 - 四个循环

onePerCol：每列1个。

输出是M×N矩阵奇偶校验矩阵H.该函数只能产生R = 1/2 LDPC矩阵。

生成奇偶校验位

使用稀疏LU分解来计算奇偶校验位，利用H的稀疏矩阵属性。设H = [A | B]。我们需要分解A成为LU，其中L是下三角形而U是上三角形。为了减少工作，A必须是非单数的，因此A必须重新排序以在对角线上给出1。重新排序A的步骤（简化来自Neal的网站）：

将L和U设置为全零。

将F设置为H.

对于i = 1到M：

找到第i行和第i列或后一列中的F的非零元素

从i开始重新排列F和H的列，将此元素放入第i列i行

将F的第i列复制到第i行到U的第i列

将F的第i列从第i行复制到L的第i列

将F的第i行添加到第i列中值为1的后续行。

结束

将B设置为重新排列的H的N - M列。

有3种策略可以选择对角线的下一个非零元素：

第一步：从列搜索中选择第i列中第一个找到的非零元素。

Mincol（最小列）：从列i开始选择非零元素，其列中具有最小数量的非零值。

Minprod（最小乘积）：从第i列开始选择非零元素，最小化以下乘积：（1）其行中的非零数减去1，以及（2）其列中的非零数减去1。

令z = Bs，其中s是二进制输入向量。求解c的L（Uc）= z，其中c是奇偶校验向量。假设c是正确的，让u = [c | s]。我们应该有Hu'= 0，其中H是重新排列的原始H.

函数 makeParityChk.m 接受3个参数：

dSource：二进制源

H：奇偶校验矩阵（来自makeLdpc.m）

策略：LU分解策略，0 =第一，1 = Mincol，2 = Minprod，

输出将是：

c：奇偶校验位向量

newH：重新排列H，应该用于编码解码而不是原始的H.

解码

使用迭代置信传播或和积算法（SPA）完成LDPC码解码。提供了四种版本的SPA解码器（在AWGN信道下调制的BPSK）：

硬判决（位翻转）解码器（ decodeBitFlip.m ）

解码0/1消息，如果多数为1，则选择“1”，否则选择“0”。解码器可以用于教程或引入消息传递算法，因为它不使用复杂概率或对数似然函数。与BPSK相比，对于非常低的E b / N 0，预期性能会更差。

概率域SPA解码器（ decodeProbDomain.m ）

基于Gallager的作品[3]。

Log-domain SPA解码器（ decodeLogDomain.m ）

与概率域SPA类似，但使用对数似然而不是概率函数。优点是可以使用加法而不是乘法来完成操作，计算成本更低。

简化的日志域SPA解码器（ decodeLogDomainSimple.m ）

log-domain SPA的修改版本用最小（x）替换Pi（x）。为了进一步简化，对数似然函数可以直接用输入信号波形替换[5]，因此简化的对数域解码器不需要噪声方差信息。

概率域和对数域解码器都需要用于其输入的噪声方差（N 0/2 ）信息。其他解码器的参数包括接收到的噪声信号，H矩阵和迭代次数。

把它们放在一起

下载下面的所有文件并运行 ldpcBer.m ，你会很好！您可能希望首先修改 ldpcBer.m ，运行较小尺寸的LDPC矩阵（即较小的数据帧）或选择您喜欢的解码器。

C-MEX解码器

您需要编译 decodeLogDomainSimpleMex.c ，并用新编译的MEX文件替换任何LDPC解码器。使用Matlab C编译器（在Microsoft Windows XP和Vista下）和Microsoft Visual C ++ 6.0版进行了很好的测试。

ldpcBer.m

makeLdpc.m

makeParityChk.m

decodeBitFlip.m

decodeProbDomain.m

decodeLogDomain.m

decodeLogDomainSimple.m

decodeLogDomainSimpleMex.c

运行LDPC码BER仿真。

创建R = 1/2 LDPC矩阵。

创建奇偶校验位。

硬判决SPA解码器

概率域SPA解码器。

Log-domain SPA解码器。

简化的日志域SPA解码器。

s -implified log-domain SPA解码器的 C-MEX 。

仿真结果

This LDPC software is intended as an introduction to LDPC codes computer based simulation. The pseudo-random irregular low density parity check matrix is based on Radford M. Neal’s programs collection, which can be found in [1]. While Neal’s collection is well documented, in my opinion, C source codes are still overwhelming, especially if you are not knowledgeable in C language. My software is written for MATLAB, which is more readable than C. You may also want to refer to another MATLAB based LDPC source codes in [2], which has different flavor of code-writing style (in fact Arun has error in his log-likelihood decoder).

Creating LDPC matrix

Steps for creating LDPC matrix (from Neal's website):

Select one of the method for creating sparse parity check matrix: (1) Evencol, or (2) Evenboth.

Add 1s to rows that have no 1s or only have one 1 in them.

If desired, try to avoid having cycle of lenght-four, where pair of column, both have 1s in particular rows

Evencol randomly distribute the same amount of 1s between columns, evenboth is similar to evencol, but also try to make rows having the same number of 1s.

Function makeLdpc.m accepts 5 parameters:

M: Number of row

N: Number of column

method: Method for distributing 1s, 0 = Evencol, 1 = Evenboth

noCyle: Eliminate length-four cycle

onePerCol: Number of 1s per column.

The output is M x N matrix parity check matrix H. This function can only create R = 1/2 LDPC matrix.

Generating parity check bits

Parity check bits is computed using sparse LU decomposition, utilizing sparse matrix properties of H. Let H = [A|B]. We need to decompose A become LU, where L is lower triangular and U is upper triangular. In order the reduction to work, A must be non-singular, hence A must be reordered to give 1s on the diagonal. Steps for reordering A are (simplification from Neal's website):

Set L and U to all zeros.

Set F to H

for i = 1 to M:

Find non-zero element of F that is in row i and column i, or in the latter column

Rearrange columns of F and H from i onward to put this element in row i column i

Copy column i of F up to row i to column i of U

Copy column i of F from row i to colunm i of L

Add row i of F to later rows that have value 1 in column i.

end

Set B to the N - M column of the rearranged H.

There are 3 strategies to choose the next non-zero element for the diagonal:

First: choose the first found non-zero element from column i onward from column search.

Mincol (minimal column): choose non-zero element from column i onward that has minimum number of non-zeros in its column.

Minprod (minimal product): choose non-zero element from column i onward that minimized the product of: (1) the number of non-zeros in its row minus 1, and (2) the number of non-zeros in its column minus 1.

Let z = Bs, where s is binary input vector. Solve L(Uc) = z for c, where c is parity check vector. Provided c is correct, let u = [c|s]. We should have Hu' = 0, where H is rearranged original H.

Function makeParityChk.m accepts 3 parameters:

dSource: Binary source

H: Parity check matrix (from makeLdpc.m)

strategy: Strategy for LU decomposition, 0 = First, 1 = Mincol, and 2 = Minprod,

and the output will be:

c: Parity check bits vector

newH: Rearrange H, which should be used for encoding-decoding instead of original H.

Decoding

LDPC code decoding is done using iterative belief propagation or sum-product algorithm (SPA). Four versions of SPA decoder (BPSK modulated under AWGN channel) are presented:

Hard-decision (bit-flip) decoder (decodeBitFlip.m)

Decode 0/1 message, choose '1' if the majority is 1, else '0'. The decoder could be used for tutorial or introduction to message passing algorithms since it does not employ complicated probability or log-likelihood function. Expect worse performance compared to BPSK for very low Eb/N0.

Probability-domain SPA decoder (decodeProbDomain.m)

Based on Gallager's works [3].

Log-domain SPA decoder (decodeLogDomain.m)

Similar to probability-domain SPA, but using log-likelihood instead of probability function. The advantage is operations can be done using additions instead of multiplications, which computationaly less expensive.

Simplified log-domain SPA decoder (decodeLogDomainSimple.m)

A modified version of log-domain SPA, replaces Pi(x) with minimum(x). For further simplification, log-likelihood function can be replaced with incoming signal waveform directly [5], hence simplified log-domain decoder does not need noise variance information.

Both probability-domain and log-domain decoder need noise variance (N0/2) information for their input. Other decoder's paramaters include received noisy signal, H matrix and number of iteration.

Putting it all together

Download all the files below and run ldpcBer.m, you'll be good! You might want to modify ldpcBer.m first, running smaller size of LDPC matrix (i.e. smaller data frame) or chose your prefered decoder.

C-MEX decoder

You need to compile decodeLogDomainSimpleMex.c, and replace any of the LDPC decoder with the new compiled MEX file. It was tested well using Matlab C compiler (under Microsoft Windows XP and Vista) and Microsoft Visual C++ version 6.0.

decodeLogDomainSimple.m

decodeLogDomainSimpleMex.c

Run LDPC codes BER simulation.

Create a R = 1/2 LDPC matrix.

Create parity check bits.

Hard-decision SPA decoder