LDPC码及其译码实现

LDPC码及其译码实现一、LDPC码简介LDPC码最早在20世纪60年代由Gallager在他的博士论文中提出，但限于当时的技术条件，缺乏可行的译码算法，此后的35年间基本上被人们忽略，其间由Tanner在1981年推广了LDPC码并给出了LDPC码的图表示，即后来所称的Tanner图。1995年前后MacKay和Neal等人对LDPC码重新进行了研究，提出了可行的译码算法，从而进一步发现了LDPC码所具有的良好性能，迅速引起强烈反响和极大关注。LDPC码（低密度奇偶校验码）本质上是一种线形分组码，它通过一个生成矩阵G将信息序列映射成发送序列，也就是码字序列。对于生成矩阵G，完全等效地存在一个奇偶校验矩阵H，所有的码字序列C构成了H的零空间(nullspace)，即HCT=0。LDPC码的奇偶校验矩阵H是一个稀疏矩阵，相对于行与列的长度，校验矩阵每行、列中非零元素的数目(我们习惯称作行重、列重)非常小，这也是LDPC码之所以称为低密度码的原因。由于校验矩阵H的稀疏性以及构造时所使用的不同规则，使得不同LDPC码的编码二分图(Taner图)具有不同的闭合环路分布。而二分图中闭合环路是影响LDPC码性能的重要因素，它使得LDPC码在类似可信度传播(BeliefProPagation)算法的一类迭代译码算法下，表现出完全不同的译码性能。当H的行重和列重保持不变或尽可能的保持均匀时，我们称这样的LDPC码为正则LDPC码，反之如果列、行重变化差异较大时，称为非正则的LDPC码。根据校验矩阵H中的元素是属于GF(2)还是GF(q)(q=2p)，我们还可以将LDPC码分为二元域或多元域的LDPC码。二、LDPC译码算法2.1、Gallager概率译码算法Gallager当初为了介绍LDPC码，同时还提出了一种迭代的概率译码算法，Gallager概率译码算法，后来在此基础上又发展出了置信度传播译码算法（BPA，也称SPA或者MPA）。假设一个二进制序列是一个LDPC码字，那么这n个比特就要满足由该码字的校验矩阵所确定的一系列的校验方程。并且，包含某一比特iC的校验方程可能不止一个，这些校验方程中的其他比特又可能包含在其他更多的校验方程中。为了直观的表示这种关系，Gallager引入了校验集合树的概念。图2.1所示为某一比特的校验集合树。(1,1)(1,2)(1,3)(2,1)(2,2)(2,3)(3,1)(3,2)(3,3)d图2.1校验集合树根节点表示比特d，和d相连的每一条边表示包含d的一个校验方程，在图2.1中，d包含在3个校验方程中；第一层中的每一线段上的节点表示这一校验方程中除d以外的其他比特，因此包含d的3个校验方程分别是：1,11,21,30dCCCC2,12,22,30dCCCC3,13,23,30dCCCC校验方程中的加法是模2加。dC即为比特d的数值，1,1C即为比特（1,1）的数值。与第一层各节点相连接的每条边同样表示包含该比特的一个校验方程，第二层中的每一线段上的节点同样表示该校验方程中除第一层比特以外的其他比特。以比特（2,2）为例，和比特（2,2）相连接的边有3条，其中一条与本层节点（2,1），（2,3），及根节点d相连，另外两条与第二层中节点u,v,w和x,y,z相连。因此包含比特（2,2）的3个校验方程分别为：2,22,12,30dCCCC2,20uvwCCCC2,20xyzCCCC第三层（图中未画出）及以后的各层依此类推，每个比特都有相应的以该比特为根节点的校验集合树。假设信道是无记忆信道，iy仅与ic及信道噪声有关，考虑根节点d和第一层节点组成的集合，这些比特组成了包含比特d的j个校验方程，每个校验方程由k个比特组成（包含比特d），显然发送的这些比特满足这j个校验方程。因此假设当发送的码字是011(,,,)ncccc时，那么在以上情况下接收到的符号即为011(,,,)nyyyy。这样在传送码字c时，码字中的各比特满足包含比特d的j个校验方程。当接收到相应的符号序列y时，比特d为1的条件概率可以表示为(1|,)dPcyc。同理，比特d为0的条件概率表示为(0|,)dPcyc。又令当不考虑发送比特间的相关性时，d为1的概率表示为(1|)dPcy，它与信道模型有关。令(1|)ddPPcy，ilP表示d的校验集合树第一层中包含d的第i个校验方程的第l个比特为1的概率，那么有：111111(12)[0|,]1(2.1)[1|,]1(12)kjilddlkiddillPPcycPPcycPP根据式2.1，只要知道了图2.1的第一层中各比特为1的概率ilP，就可以算出比特d的条件概率。在其他根节点的校验树里比特d又作为一个节点参与到根节点的概率计算中去，即比特d从与其有关的节点中获取信息计算出概率，再将其计算出的概率信息传出用于计算其他的节点ic，由于在计算其他节点ic时同样会用到计算比特d时所用过的运算，所以可以通过共享计算的中间结果而使计算量大为降低，进而发展出了BPA（也称SPA或MPA）。2.2BP算法（也称SP或MP算法）校验集合树虽然比较直观，但对于每一个节点都有不同的校验集合树，因此在描述并行计算整个码字各比特的后验概率时，使用校验集合树并不方便，因此介绍一种新的图形表示法，Tanner图。对应于图2.1的Tanner图如图2.2所示。该图仅画出部分变量节点和校验节点。Tanner图中变量节点对应于校验集合树中的节点，校验节点对应于校验集合树中的边。(1,1)(1,2)(1,3)(2,1)(2,2)(2,3)(3,1)(3,2)(3,3)d1S2S3S4S图2.2对应于图2.1校验集合树的部分Tanner图由Tanner图可看出，信息是在变量节点和校验节点之间来回传递的，变量节点d将自身的固有信息再加上与它有关的1S，2S，3S校验节点传给它的额外信息一起传递给校验节点4S。同理，校验节点jc也是将自身的固有信息再加上与它有关的除某一变量节点iv以外的其余变量节点所传给jc的额外信息一起传递给变量节点iv，如此信息便在变量节点和校验节点之间来回传递，不断更新变量节点和校验节点的值，所有变量节点和校验节点更新完一次称之为一次迭代结束，直到变量节点译码成功或者达到最大的迭代次数译码输出。因此，完整的BPA（也称SPA或MPA）表述如下：1)初始化：111mlllqpf，011mlmlqq，其中1lf表示的是信道的先验概率。2)校验节点更新：01mlmlmlqqq，01mlmlmlrrr，则有''()/(1)(2.2)mzmlmllLmlrq0(1)/2(2.3)mlmlrr1(1)/2(2.4)mlmlrr3)变量节点更新：''000()/(2.5)mlmllmlmMlmqpr''111()/(2.6)mlmllmlmMlmqpr其中ml是一个使得011mlmlqq的值4)似后验概率更新000()(2.7)lllmlmMlqpr111()(2.8)lllmlmMlqpr同样l是一个使得011llqq的值5)比特判决：如果00.5iq，判决0ix；否则判1ix，（1,2,,iN）若0THx，或者迭代次数达到最大迭代次数，则结束迭代，将x作为译码输出；否则返回到步骤2）继续迭代。2.3Log-BP算法由2.2介绍的BP算法可以看出BP算法比较复杂，且该算法需要很多连乘运算，需要大量的运算时间和消耗大量的硬件资源，不便于硬件实现，因此，又发展出了一种更加简便的对数域的BP算法，简称Log-BP算法。考虑一个二值随机变量x，它的对数似然比()Lx定义为(0)()ln(1)PxLxPx根据对数似然比的定义，令001lnlllfvf01lnmlmlmlqvq01lnlllqvq01lnmlmlmlrur再对以上的BP算法在对数域内进行一系列的化简，最终可获得Log-BP算法的完整表述：1)初始化：0(2.9)mllvv2)校验节点更新：'''1()/2tanh[tanh](2.10)2mlmllLmlvu3)变量节点更新：''0()/(2.11)mllmlmMlmvvu4)似后验概率更新：0()(2.12)llmlmMlvvu5)比特判决：如果0lv，则判0lx；否则判1lx，（1,2,,lN）若0THx，或者迭代次数达到最大迭代次数，则结束迭代，将x作为译码输出；否则返回到步骤2）继续迭代。2.4Min_Sum算法（最小和译码算法）由于Log-BP算法仍然具有连乘与反三角函数，因此仍然比较复杂，进而又发展出了Log-BP算法的近似算法，Min_Sum算法。最小和译码算法的优点是对LDPC码并行的迭代软译码算法，故其运算比BP算法以及Log-BP算法都要简单，它仅采用求最小值以及相加运算，而不需要乘法运算。通过一系列的数学化简与近似运算，可在Log-BP算法的基础上得到Min_Sum算法的完整表述，Min_Sum算法仅步骤2）与Log-BP算法不同，其余的步骤均和Log-BP算法的步骤相同，Min_Sum算法的步骤2）为：2)校验节点更新：''''()/()/min||(2.13)mlmlmllLmllLmluv其中'ml为每一个'mlv的符号，其余步骤及循环、迭代运算均与Log-BP算法一致。三、译码算法的实现在此选取C语言来实现译码算法，并且采用C与Matlab联调的方式来验证译码程序的结果，Matlab里完成m文件和script的编写，C语言与Matlab之间的联调通过Matlab的MEX文件来实现。译码算法程序流程图：InitialRevcodeSigmaHUpdatechecknodeJudgedecodeErr=1ReturndecodeUpdatevariblenodeMAX=100MAX100Err=0图3.1译码算法程序流程图译码程序输入的输入参数为信道接收到的码字，信道高斯噪声方差以及校验矩阵H，按照以上程序流程运行后返回译码码字，其中MAX为最大迭代次数。具体的译码程序和Matlab的MEX文件中各函数的定义，调用以及各自的实现方式如下图：#includemex.h#includeMin_Sum4_head.hVoidMin_Sum4(double,double*){Structnode*RDsave(FILL*,int,int)Structnode*decode(structnode*,structnode*,float)}VoidmexFunction(intnlhs,mxArray*plhs[],intnrhs,constmxArray*prhs[]){mexErrMsgTxt()mwSizemxGetM(mxArray*)mwSizemxGetN(mxArray*)mxArray*mxCreateDoubleMatrix(int,int,mxREAL)double*mxGetPr(mxArray*)VoidMin_Sum4(double,double*)}#includeMin_Sum4_head.hstructnode*RDsave(FILE*,int,int)structnode*Fun_check(structnode*,structnode*);structnode*Intial(structnode*,structnode*,float);structnode*Fun_step2(structnode*,structnode*);structnode*Fun_step3(structnode*,structnode*,structnode*);structnode*Fun_step4(structnode*,structnode*,structnode*);structnode*decode(