LDPC Code

Using MATLAB and 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.

 ldpcBer.m   makeLdpc.m    makeParityChk.m    decodeBitFlip.m   decodeProbDomain.m   decodeLogDomain.m   decodeLogDomainSimple.m Run LDPC codes BER simulation.   Create a R = 1/2 LDPC matrix.   Create parity check bits.   Hard-decision SPA decoder   Probability-domain SPA decoder.   Log-domain SPA decoder.   Simplified log-domain SPA decoder.   C-MEX of simplified log-domain SPA decoder.

Simulation results

• 18 February 2008: Eb/N0 value correction.
• 11 March 2009: C-MEX source code for LDPC decoder added.
• 8 January 2013: Result plot added.

Reference (Drop me an email if you need any of these materials)