LDPC碼的編譯碼算法研究

2010 業(yè)論文LDPC： 專業(yè)班級(jí)： 電子生姓： 號(hào)指導(dǎo)教師： 教師職：教授2010 年6 月2 日PAGE\*ROMANPAGE\*ROMANIV要wwyyksCLDPC,,、衛(wèi)星通等領(lǐng)域可得到廣泛應(yīng)用文章介紹了LDPC,綜述了方法和方法方法別描述了矩陣造和基于矩陣算法對(duì)LDPC快速方法進(jìn)行析方法主要論述了消息傳遞傳播譯分LDPCLDPC就行了仿真同時(shí)對(duì)LDPC方法發(fā)展及應(yīng)用景作了析。文重點(diǎn)基本原理和類原理和類別從基于生成矩陣和基于矩陣詳細(xì)討論了LDPC算U解法RURUMP判決比特翻轉(zhuǎn)算法和加權(quán)比特翻轉(zhuǎn)算法文下關(guān)鍵詞：LDPCLDPC仿真主要關(guān)鍵詞：LDPCgdgsfCtC,ywyyk,sadfblockcodesinnature,andthedecodingperformanceofLDPCismorenearertotheShannonlimit.Withitsbestperformanceandsimpledecoderstructure,LDPCcodeswillbewidelyusedingdgsfCdescribedintwosteps:theconstructionofparity-checkmatrixandtheencodingmethodbasedonparity-checkmatrix.therapidlycodingmethodforLDPCcode.Astodecodingalgorithm,method,Min-Sumdecodingmethod,Bit-FlippingmethodandWeightedBit-Flippingmethodarediscussed.EmulatefortheLDPCcodes.Thedevelopmentandandisasarticle focuses on encoding and decoding algorithms of LDPC codesAccordingtothedifferentmethodsofdecodingalgorithm,andmakesthetheoreticalysCsgdg錄133(LDPC) 446() 63.2(LU) 73.3(RU) MP131314下LDPC能仿真 15仿真軟件簡介（MATLAB&SIMULINK） 15仿真與結(jié)果析 15仿真系統(tǒng)框圖及系統(tǒng)總流程圖 16BF及其改進(jìn)仿真 1719202122PAGEPAGE20引言率(BER)R比特/小錯(cuò)誤率的信息傳輸?shù)囊痪褪切畔鬏斔俾式档蚐hannonShannon圖1.1數(shù)字通信系統(tǒng)基本模型Shannon的信息從通信系統(tǒng)的整體最佳化來研究信息的傳輸和處。比11所示的兩個(gè)階段的信息處和信道Shannon1．1中的信道部分只是信息傳輸所通過媒)010且(0110))010、1；信高斯隨機(jī)在過去十年里移動(dòng)通技術(shù)得了迅猛發(fā)展廣泛至今已發(fā)展了三代第一代移動(dòng)通（1G）擬傳方式行語音通話主要GG用數(shù)字時(shí)址（TDMA）或址（CDMA）現(xiàn)動(dòng)態(tài)尋址功能GSM、CDMA系統(tǒng)代表現(xiàn)了從擬數(shù)字系統(tǒng)跨越第三代移動(dòng)通信（3G）著重現(xiàn)傳統(tǒng)移動(dòng)通與開放式因網(wǎng)融合各個(gè)國家網(wǎng)絡(luò)將融合一個(gè)整體而在移動(dòng)通更新?lián)Q代技術(shù)重要一項(xiàng)本文論述LDPC即一。設(shè)計(jì)工作而且利產(chǎn)品開放式結(jié)構(gòu)容易地功能行擴(kuò)充數(shù)析處理功能十強(qiáng)運(yùn)它譯行仿真。LDPC分組碼因?yàn)榈兔芏绕媾夹ｒ?yàn)碼是一種特殊的線性分組碼，所以本章將首先對(duì)線性分組碼做一個(gè)概述，為討論LDPC碼作鋪墊。0，l，2?，q．1，qP加和乘運(yùn)算下構(gòu)成GF(q。1CNGF(q(C0C，1?CN-1)CGF(qCq進(jìn)制線性q=2。3K2kN線性碼還有如下一些有用的性質(zhì)：1必然包含一個(gè)全零碼字。2：線性碼的最小距離等于其中一個(gè)最輕非零碼字的漢明重量。這一性質(zhì)表明確定線性碼的最小距離(決定檢和)一的分組碼的。3：線性碼中檢的的碼字，且由所有的非零碼字組成。(g0

g， g (N，K)--進(jìn)制碼空的一組，對(duì)任意一1 K-1個(gè)碼字c∈C，存在唯一的表達(dá)形式C=a g+a g+ +a g0 0 1 1 K-1 K-1

（2-1）Kaa

CG0 1 列而成。)LDPCCmnHHCl比值遠(yuǎn)小于1；2()最多只1相同位置1；3、任意無關(guān)盡量大。這樣C長n位長度大約m信息位長度k m。規(guī)則指H重重別等于常dvdc，因我們要求H秩所以率r1-

m dn=1-dv

(2-2))n=10d=dv c

c45C。1111110000001 00011100 00 100100110(2-3)0 0100101010 001001011果H重重不常我們就不規(guī)則LDPC我們可以認(rèn)LDPCLDPC例碼v可以用重量布多項(xiàng)式來描假設(shè)最大重最大重別dmax和vdmaxc

H重布多項(xiàng)式(x)可以表示

(x

dmaxvi2

xi

（2-4）i(x1H重i(x可以表示(x)=

dmaxci2

xi

）i(x(1)=r：i1xdxr

01xdx

）0n而不特CH線性分組校驗(yàn)矩陣做變換到矩陣這校驗(yàn)矩陣但LDPC這種特屬稀疏矩陣這校驗(yàn)矩陣因義必須出這稀疏校驗(yàn)矩陣才有義。圖中變節(jié)點(diǎn)ci 0l．列位1．m-1與校驗(yàn)矩陣即各校驗(yàn)方) ci fj應(yīng)應(yīng)位置元素1其圖第24圈(，4)。圖2.1 LDPC因子圖表示這是因?yàn)榉且?guī)則碼的“波浪效應(yīng)C更度的變化很LDPC碼都之為規(guī)則的。碼都之為規(guī)則的。于Th成矩陣線性分組)m×nH，uFm，cFnH cT =0H,它具有非系統(tǒng)碼的形式，因首先需要對(duì)給定的校驗(yàn)矩陣H進(jìn)行行列變換分解成B的形式其中A為m-維的矩陣，Bmm

u AB 即 0因“-”表示向- Au的逆元，在二進(jìn)制編碼逆元是本身。.2于校驗(yàn)矩陣的編碼算法(U 分解法)Bp對(duì)B矩陣進(jìn)行LU分解，即B=LU為上三角矩LUp-Au，然后通過以下步驟計(jì)算校驗(yàn)位：1．zAm。2UpLyyUpyp。3.3于校驗(yàn)矩陣的編碼算法(RU 算法)所表示。但在實(shí)際情況，校驗(yàn)矩陣往往不能化為下三角矩陣。Richardson和UrbankcRU算法。RU算法僅僅通過行列置換將一般的低密度奇偶校驗(yàn)矩陣化為一個(gè)近似o(m2)o(n+92)，g為行的一小部分成比例。J5rBCT（m-g（m、）g（mgg（（m-g）（m-并Tl。A B T H=C D

(3-2)msx（spppp1 2 1 分別定義一個(gè)向量p長p長m-編碼步驟：1 21、計(jì)算信源向量上正子z =AsTA2PA零2

（3-3）PA=T-1Z (3-4)2 A回代算可以線性時(shí)間內(nèi)得出這個(gè)向量即計(jì)算PA第一個(gè)比特，2然第二個(gè)比特然第個(gè)等等。23、計(jì)算向量PA2z=CsT

EpA

（3-5）B 24p1。F=-ET-1B+D并求出它逆(g2g維高密度p為1pT-F-1z (3-6)1 Bog2)P15、拋棄臨檢西并出新上子z z

(3-7)c A 1也可以在線性間內(nèi)完成。6、求出二檢p，使得上子零2PT=-T-1z (3-8)2 C可以用回代算法在線性間內(nèi)求出。計(jì)算PTPT3,13.2所示。1 23．1計(jì)算PTFF1sT（O（n）矩陣求逆運(yùn)算高密矩陣和相乘作AsT度O（n）備 注稀疏矩陣和相乘-T-1AsTO（n）稀疏矩陣和相乘CsTE（ ）O（n）稀疏矩陣和減法運(yùn)算3．2作作度備 注AsTO（n）稀疏矩陣和相乘BPT1O（n）稀疏矩陣和相乘AsT+BPTO（n）稀疏矩陣和相乘1加法運(yùn)算O（n）稀疏矩陣和相乘-T1(AsT+BPT)12.4(12,3,6)LDPC12345671111289g=2E

1 UET-1B+D= 1 U58列123406751289s0000AsT

1 1 0 TT1sT1 1 0 0TT1sT0 0TCsT

0 0Tpp1 2

T1sTsT0 T1U1T1sTsT0 TPT1pT=0 1 0 0T11sTpT1 0 0 T1T1sT

pT1 0 1 0Tp1 2p p1 2

0 0 0 0 0 0 1 1 0 1 02.4HcT=0。有效法利用了矩陣稀疏性適用于任何基于稀疏矩陣的。概述LDPC具有良好性能LDPC具有良好性能重要原LDPC采用了基于置傳播迭譯法這是一種迭概率譯法是LDPC與傳統(tǒng)糾錯(cuò)重要區(qū)別。enMP法是最主要C法它具有嚴(yán)格數(shù)學(xué)結(jié)構(gòu)和良好性能使用它能對(duì)譯性能做定量分析。LDPC譯法中很多種都可以被歸結(jié)傳遞法集中。傳遞法主要思想就是通過在變量節(jié)點(diǎn)和節(jié)點(diǎn)之間來回傳遞概在圖可以直觀表示出來，開始收變量變量對(duì)每個(gè)變量做判決得一個(gè)碼字再通過矩陣證碼字正確性。如果譯碼成功直到達(dá)到預(yù)先設(shè)定最大迭代傳遞法為了保證傳遞獨(dú)立性每個(gè)收都是從除自身之外其他而碼長都是有限使得可能永遠(yuǎn)收到自身無關(guān)即存在環(huán)d列重為dcdv

v1。個(gè)而所又自自以d一1個(gè)變量c在上一迭代周期中值如下圖所示樹狀圖表示們之間。圖3.3樹LDPC碼有很多種譯碼方法圖傳遞譯碼算迭代過程中傳消不同形式譯碼方法分為硬是比特值稱之為硬判決BF法具有較低譯碼復(fù)雜度易于工程實(shí)現(xiàn)。但是軟判決譯碼比硬判決譯碼在性能上要損失約2-3dB；如果在譯碼過程中傳是的左右軟判決增益而在衰落道中5dBl比特量化譯碼，而軟判決譯碼可以看成無窮多比特量化譯碼。主要的硬判決譯碼算法有比特翻轉(zhuǎn)算法(BF)、加權(quán)的比特翻轉(zhuǎn)算法(WBF)、簡化的最小和算法而軟判決譯碼可以看成無窮多比特量化譯碼。主要的硬判決譯碼算法有比特翻轉(zhuǎn)算法(BF)、加權(quán)的比特翻轉(zhuǎn)算法(WBF)、簡化的最小和算法(Min-sum)、歸一化最小和算法(NormalizedMin．Sum)、偏移量最小和算法(OffsetMin．sum)等。4.2.1 比特翻轉(zhuǎn)算法Gallager在其論文中提出了硬判決譯碼算法，該算法是一種比較簡單而且容易理解的譯碼算法，它對(duì)運(yùn)算量和存儲(chǔ)量的要求都很低，但是其性能相對(duì)比較差。比特翻轉(zhuǎn)算法(BitFlippingAlgorithm)可看成是置信傳播算法的簡化形式，而加權(quán)BF算法的基礎(chǔ)上加上硬判決譯碼系數(shù)，其性能較比特翻轉(zhuǎn)譯碼算法有一定程度的提高比特翻轉(zhuǎn)算法(BitFlippingAlgorithm)是GallagerGallager硬判決的譯碼算法。設(shè)碼字 c=，cc 為發(fā)送序列，經(jīng) 調(diào)0 1 制為序列 =c，c，c ，xc，0i1，rr，r，r

為接0 1 i i 0 1 收 的 實(shí) 數(shù) 向 量 序 列 ， 由 實(shí) 數(shù) 序 列 可 以 得 到 硬 判 決 二 元 向 量 序 列(z，z， ，z ：0 1

，當(dāng)

0時(shí)zi

（4-1）i ，當(dāng)ri

0時(shí)由此得到碼字伴隨式

sss， s ）= zHT，若

0，則說明接收向量

0 1 2 J-1 j滿足第 j個(gè)程；若 ，則接收向量滿足有程，接收碼字 z，譯碼成；若伴隨式為“0”向量時(shí)，接收序列 z有錯(cuò)誤，此時(shí)則需計(jì)算出每個(gè)碼元不滿足程的個(gè)數(shù) f=

f， ，f

=sH，搜索f中的最大值，翻轉(zhuǎn)對(duì)應(yīng)位置的碼元

0 1 z 程，到譯碼成到j(luò)最數(shù)。BF譯碼算法步驟如下：根據(jù)硬判決二元向量序列得到碼字的伴隨式 s，判斷s是為“0，為“0，則譯碼成，則轉(zhuǎn))；算，出其最值f=max{f}，翻轉(zhuǎn)對(duì)位置的碼元z ；j j(1“0譯碼成功跳出否則重復(fù)上述步驟直達(dá)最大迭次數(shù)。由于校驗(yàn)矩陣稀疏矩陣而且一般機(jī)構(gòu)成所以參與每個(gè)校驗(yàn)方程的比特很少且這些比特在碼字上分布很分散那么任一校驗(yàn)方程所含比特要么無錯(cuò)要么以很高概率只有一個(gè)比特錯(cuò)誤F算法就可以有效地進(jìn)行糾錯(cuò)。即使某一校驗(yàn)方程發(fā)生多于一個(gè)錯(cuò)誤犧牲就是譯碼性能所以下面對(duì)于硬判決譯碼算法提出了一種加權(quán)硬判決譯碼算法它是BF.2.2加權(quán)比特翻轉(zhuǎn)譯碼算法在譯碼接收端通過添加一些可信信息BF，譯碼(WBF)算法就是在選擇需要翻變節(jié)點(diǎn)時(shí)候每一個(gè)碼字中不校驗(yàn)方程個(gè)數(shù)最多碼元信道輸出信息作該判決權(quán)重信息。息。yy

,y,...,

x

x

,x,...,x

w

0 1,w,...,w

N10 1 N1

0 1 N

mn

校驗(yàn)矩陣Nm

mn

1表示與校驗(yàn)節(jié)點(diǎn)m相連變節(jié)點(diǎn)mnHmn

F1?M－101?－1：根據(jù)硬判決二元zss“0，“0；(2)計(jì)算y(2)計(jì)算min

ny－1n∈；i1?N－1計(jì)算E=n n)

(2sm

y

max{Ej

翻對(duì)應(yīng)位置z；j)(1“0譯碼成功跳出碼成功跳出否則重復(fù)上述步驟直達(dá)最大迭次數(shù)。加權(quán)譯碼算法是通過軟判決譯碼算法和附加信息來計(jì)算加權(quán)校驗(yàn)信息這種算法雖然比單純BF下LDPC碼的性能仿真仿真軟件簡介（matlab&simulink）在數(shù)學(xué)類科技應(yīng)用軟件中在數(shù)值計(jì)算進(jìn)行矩陣運(yùn)算、繪制函數(shù)和數(shù)據(jù)、實(shí)現(xiàn)算法、創(chuàng)建用戶界面、連接其他編程語言的程序等，主要應(yīng)用于工程計(jì)算、控制設(shè)計(jì)、信號(hào)處理與通訊、圖像處理、信來解算問題要比用在數(shù)學(xué)類科技應(yīng)用軟件中在數(shù)值計(jì)算進(jìn)行矩陣運(yùn)算、繪制函數(shù)和數(shù)據(jù)、實(shí)現(xiàn)算法、創(chuàng)建用戶界面、連接其他編程語言的程序等，主要應(yīng)用于工程計(jì)算、控制設(shè)計(jì)、信號(hào)處理與通訊、圖像處理、信來解算問題要比用C，F(xiàn)ORTRAN等語言完成相同的事情簡捷得多，并且mathwork也吸收了像Maple等軟件的優(yōu)點(diǎn),成為一個(gè)強(qiáng)大的數(shù)學(xué)軟件。包括擁有數(shù)百個(gè)內(nèi)部函數(shù)的主包和三十工包。工包可以分為工包和學(xué)科工工包用來。學(xué)科工包，控制工包，信號(hào)處理工包，通信工包等于類。使用戶內(nèi)部函數(shù)主包件和工包是可可的程序的編程序的用工CommunicationSignalProcessingSimulink軟件展和方法學(xué)致力與更然地抽象事物的特征，尋求使模型研究者更然地參與活動(dòng)的方Simulink提供了一個(gè)系統(tǒng)級(jí)的建模與動(dòng)態(tài)在科學(xué)計(jì)算上的天然仿真與結(jié)果分析中采用的信道中采用的信道是二進(jìn)制輸?shù)氖腔鶐PSK調(diào)制。一般在中要獲得低的誤比特率需要大量的數(shù)據(jù)幀，而在碼長長大量的數(shù)據(jù)幀的計(jì)算要花費(fèi)很多只選定了一些碼長相LDPCLDPCLDPCLDPCGallager4LDPCLDPCGallager4LDPC4238、6、4、86C娜l;主要針8161/2LDPC算各種硬判決譯算基BP4862關(guān);623并比;然后812取0000;其。5.3 仿真系統(tǒng)框圖及系統(tǒng)總流程圖LDPCLDPC51示5-1LDPC譯系統(tǒng)框圖針LDPC52。5-2 5.4 BF仿真5-35-36，3，0。5-3 BF3LDPC(816,3,6)5-3B105。從每個(gè)單純硬判決BF最當(dāng)B善接近6B兩距大幅B一直處于102B剛接近104FF相對(duì)要好一些小于6BFF為101處和10B善接近6B兩距大幅縮小了同一。而縮小了同一。而兩加權(quán)IMMBF和前兩相有了明顯善尤是4B范圍內(nèi)但大4BF變得很接近F最佳講與原始BF相有了一定提高但效果有限實(shí)際號(hào)道譯這些雖然復(fù)雜度小然而實(shí)用。結(jié) 論編碼定理的指引下，LDPC碼作為一種新的編碼方式，由于其校驗(yàn)矩陣具有稀疏的特點(diǎn)，使得它的譯碼復(fù)雜度與碼長呈線性的關(guān)系，而性能卻很接近Shannon編碼定理的指引下，LDPC碼作為一種新的編碼方式，由于其校驗(yàn)矩陣具有稀疏的特點(diǎn)，使得它的譯碼復(fù)雜度與碼長呈線性的關(guān)系，而性能卻很接近Shannon完成了以下幾個(gè)方面（1）分解法。并用簡明例子對(duì)RU算法做了清晰的解釋。（MP比。。碼的性能對(duì)特算法進(jìn)行。作碼做進(jìn)一的研究1、LDPC碼的校驗(yàn)矩陣結(jié)構(gòu)及其優(yōu)化已經(jīng)證明了非規(guī)則碼進(jìn)一的研究。2、LDPC碼長碼的快速編碼由于編碼的復(fù)雜度與碼長二次方成正難以接受的碼長碼的快速編碼也一個(gè)很值得研究的。3、LDPC碼的譯碼算法LDPC碼的譯碼算法影響其性能的一個(gè)重要因素。在存在閉合環(huán)路時(shí)LDPC碼的譯碼性能一個(gè)研究課題LDPC碼信道編碼技術(shù)發(fā)展的一個(gè)新的里程碑，它以卓越的性能和線性致 謝參考文獻(xiàn)賀鶴云.LDPC碼基礎(chǔ)與應(yīng)用.人民郵電出版社,2009袁東風(fēng),張海剛.LDPC碼理論與應(yīng)用,人民郵電出版社,2008,,LDPC,2006（）,2005周建興.MATLAB從入門到精通.人民郵電出版社,2008（美）亨塞爾曼美）利特菲爾德.精通 Matlab7. 清華大學(xué)出版社,2006張德豐.MATLAB/Simulink建模與仿真.電子工業(yè)出版社,2009MATLAB/SIMULINK清華大學(xué)出版社,2008張建國.LDPC碼的應(yīng)用研究[J].通信技術(shù).2003年11期 2李水平,劉玉君,王云鶴.串行級(jí)聯(lián) LDPC碼[J].信息工程大學(xué)學(xué) 報(bào).2004年02期蘇杰,王琳,趙春雨.規(guī)則 LDPC碼在瑞利平坦衰落信道下的設(shè)計(jì)和仿 真[J].電訊技術(shù).2004年05期王鵬,王新梅.LDPC碼的快速編碼研究[J].西安電子科技大學(xué)學(xué) 報(bào).2004年06期仲海梅,王華.4G的編碼技術(shù)LDPC碼新[J].東通 信技術(shù).2004年12期J Chen,A Dholakia, E Eleftheriou. “Reduced-complexitydecodingofLDPCcodes”.2005.代碼%BiterrorrateofBPSKmodulatedLDPCcodesunderAWGNchannelclc;clearall;%LDPCmatrixsize,ratemustbe1/2%Warning:encoding-decodingcanbeverylongforlargeLDPCmatrix!M=1000;N=2000;%MethodforcreatingLDPCmatrix(0=Evencol;1=Evenboth)method=1;%Eliminatelength-4noCycle=1;%Numberof1spercolumnforLDPCmatrixonePerCol=3;%LDPCmatrixreorderstrategy(0=First;1=Mincol;2=Minprod)strategy=2;%EbN0indBEbN0=[00.511.5];%Numberofiter=5;%Numberofframe(Nbitsperframe)frame=10;%MaketheLDPCmatrixH=makeLdpc(M,N,1,1,onePerCol);fori=1:length(EbN0)ber1(i)=0;ber2(i)=0;%Makerandomdata(0/1)dSource=round(rand(M,frame));forj=1:framefprintf('Frame:%d\n',j);%Encodingmessage[c,newH]=makeParityChk(dSource(:,j),H,strategy);u=[c;dSource(:,j)];%BPSKbpskMod=2*u-1;%AdditionalwhitegaussiannoiseN0=1/(exp(EbN0(i)*log(10)/10));tx=bpskMod+sqrt(N0/2)*randn(size(bpskMod));%Decoding(selectdecodingmethod)%vhat=decodeProbDomain(tx,H,newN0,iter);vhat1=decodeLogDomain(tx,newH,N0,iter);vhat2=decodeLogDomainSimple(tx,newH,%vhat=decodeBitFlip(tx,newH,ter);%Getbiterrorrate(forbrevity,BERcalculationincludesparity[num1,rat1]=biterr(vhat1',u);ber1(i)=(ber1(i)+rat1);[num2,rat2]=biterr(vhat2',u);ber2(i)=(ber2(i)+rat2);end%forj%GetaverageofBERber1(i)=ber1(i)/frame;ber2(i)=ber2(i)/frame;end%fori%Plottheresultsemilogy(EbN0,ber1,'o-');hold;semilogy(EbN0,ber2,'o--');gridon;holdoff;functionH=makeLdpc(M,N,method,noCycle,onePerCol)%CreateR=1/2lowdensityparitycheckmatrix%% M :Numberofrow% N :Numberof% method :Methodfordistributingnon-zeroelement% {0}Evencol:Foreachcolumn,place1suniformlyatrandom% {1}Evenboth:Foreachcolumnandrow,place1suniformlyat% noCyle :Length-4cycle% {0}Ignore(donothing)% {1}Eliminate% onePerCol:Numberofonespercolumn%% H :Lowdensityparitycheckmatrix%%%CopyrightBagawanS.Nugroho,2007%%Numberofonesperrow(N/Mratiomustbe2)ifN/M~=2fprintf('Coderatemustbe1/2\n');endonePerRow=(N/M)*onePerCol;fprintf('CreatingLDPCmatrix...\n');switchmethod%Evencolcase{0}%Distribute1suniformlyatrandomwithincolumnfori=onesInCol(:,i)=randperm(M)';end%Createnonzeroelements(1s)indexr=reshape(onesInCol(1:onePerCol,:),N*onePerCol,1);tmp=repmat([1:N],onePerCol,1);c=reshape(tmp,N*onePerCol,1);%CreatesparsematrixHH=full(sparse(r,c,1,M,N));%case{1}%Distribute1suniformlyatrandomwithincolumnfori=onesInCol(:,i)=randperm(M)';end%Createnonzeroelements(1s)indexr=reshape(onesInCol(1:onePerCol,:),N*onePerCol,1);tmp=repmat([1:N],onePerCol,1);c=reshape(tmp,N*onePerCol,1);%Makethenumberof1sbetweenrowsasuniformaspossible%Orderrow[r,ix]=sort(r);%Ordercolumnindexbasedonrowindexfori=1:N*onePerColcSort(i,:)=c(ix(i));end%Createnewrowindexwithuniformweighttmp=repmat([1:M],onePerRow,r=reshape(tmp,N*onePerCol,1);%CreatesparsematrixH%RemoveanyduplicatenonzeroelementsindexusinglogicalS=and(sparse(r,cSort,1,M,N),ones(M,N));H=end%switch%Checkrowsthathaveno1oronlyhaveone1fori=1:Mn=randperm(N);%Addtwo1sifrowhasnoiflength(find(r==i))==0H(i,n(1))=1;H(i,n(2))=1;%Addone1ifrowhasonlyone1elseiflength(find(r==i))==1H(i,n(1))=1;end%for%Ifdesired,eliminateanylength-4cycleifnoCycle==1fori=%Lookforpairofrow-columnforj=(i+1):Mw=and(H(i,:),H(j,:));c1=find(w);lc=length(c1);iflc>1%Iffound,flipone1to0intherowwithlessnumberof1siflength(find(H(i,:)))<length(find(H(j,:)))%Repeattheprocessuntilonlyonecolumnleftforcc=1:lc-1H(j,c1(cc))=0;endelseforcc=1:lc-1H(i,c1(cc))=end%end%end%forjend%forend%iffprintf('LDPCmatrixiscreated.\n');function[c,newH]=makeParityChk(dSource,H,strategy)%GenerateparitycheckvectorbasesonLDPCmatrixHusingsparseLUdecomposition%% dSource:Binarysource(0/1)% H :LDPCmatrix% strategy:Strategyforfindingthenextnon-zerodiagonalelements% {0}First :Firstnon-zerofoundbycolumnsearch% {1}Mincol:Minimumnumberofnon-zerosinlatercolumns% {2}Minprod:Minimumproductof:% -Numberofnon-zerositscolumnminus1% -Numberofnon-zerositsrowminus1%% c :Checkbits%%%CopyrightBagawanS.Nugroho,2007%%Getthematric[M,N]=size(H);%SetanewmatrixFforLUdecompositionF=H;%LUmatricesL=zeros(M,N-M);U=zeros(M,N-M);%Re-ordertheMx(N-M)submatrixfori=1:M%strategy{0=First;1=Mincol;2=Minprod}switchstrategy%Creatediagonallystructuredmatrixusing'First'strategycase{0}%Findnon-zeroelements(1s)forthediagonal[r,c]=find(F(:,i:end));%Findnon-zerodiagonalelementcandidatesrowIndex=find(r==i);%Findthefirstnon-zerocolumnchosenCol=c(rowIndex(1))+(i-%Creatediagonallystructuredmatrixusing'Mincol'strategycase{1}%Findnon-zeroelements(1s)forthediagonal[r,c]=find(F(:,i:end));colWeight=sum(F(:,i:end),1);%Findnon-zerodiagonalelementcandidatesrowIndex=find(r==i);%Findtheminimumcolumnweight[x,ix]=min(colWeight(c(rowIndex)));%Addoffsettothechosenrowindextomatchthedimensionofthe...%originalmatrixFchosenCol=c(rowIndex(ix))+(i-1);%Creatediagonallystructuredmatrixusing'Minprod'strategycase{2}%Findnon-zeroelements(1s)forthediagonal[r,c]=find(F(:,i:end));colWeight=sum(F(:,i:end),1)-1;rowWeight=sum(F(i,:),2)-1;%Findnon-zerodiagonalelementcandidatesrowIndex=find(r==i);%Findtheminimumproduct[x,ix]=min(colWeight(c(rowIndex))*rowWeight);%Addoffsettothechosenrowindextomatchthedimensionofthe...%originalmatrixFchosenCol=c(rowIndex(ix))+(i-1);otherwisefprintf('Pleaseselectcolumnsre-orderingstrategy!\n');end%switch%Re-orderingcolumnsofbothHandtmp1=F(:,i);tmp2=H(:,i);F(:,i)=F(:,chosenCol);H(:,i)=H(:,chosenCol);F(:,chosenCol)=tmp1;H(:,chosenCol)=tmp2;%FilltheLUmatricescolumnbyL(i:end,i)=F(i:end,i);U(1:i,i)=F(1:i,i);%Therewillbenorowsoperationatthelastifi<M%Findthelaterrowswithelementsincolumni[r2,c2]=find(F((i+1):end,i));%Addcurrentrowtothelaterrowswhichhavea1incolumniF((i+r2),:)=mod(F((i+r2),:)+repmat(F(i,:),length(r2),1),2);end%end%fori%FindB.dsourcez=mod(H(:,(N-M)+1:end)*dSource,2);%ParitycheckvectorfoundbysolvingsparseLUc=2);%ReturntherearrangeHnewH=H;fprintf('Messageencoded.\n');functionvHat=decodeBitFlipping(rx,H,iteration)%Hard-decision/bitflippingsumproductalgorithmLDPCdecoder%% rx :Receivedsignalvector(columnvector)% H :LDPCmatrix% iteration:Numberofiteration%% vHat :Decodedvector(0/1)%%%CopyrightBagawanS.Nugroho,2007%[MN]=%Priorhard-decisionci=0.5*(sign(rx')+1);%rji=%AsscociatethecimatrixwithelementsofHqij=H.*repmat(ci,M,1);%Iterationforn=1:iterationfprintf('Iteration:%d\n',%stepfori=%Findinthecolumnc1=find(H(i,:));%Getthesummationofqij\c1(k)fork=1:length(c1)rji(i,c1(k))=mod(sum(qij(i,c1))+qij(i,c1(k)),2);end%forkend%fori%Verticalstepforj=1:N%Findintherowr1=find(H(:,j));%Numberof1sinarownumOfOnes=length(find(rji(r1,j)));fork=1:length(r1)%Updateqij,set'1'formajorityof1selse'0',excludingr1(k)ifnumOfOnes+ci(j)>=length(r1)-numOfOnes+rji(r1(k),j)qij(r1(k),j)=elseqij(r1(k),j)=0;endend%fork%BitdecodingifnumOfOnes+ci(j)>=length(r1)-vHat(j)=1;elsevHat(j)=0;endend%forend%fornfunctionvHat=decodeProbDomain(rx,H,N0,iteration)%Probability-domainsumproductalgorithmLDPCdecoder%% rx :Receivedsignalvector(columnvector)% H :LDPCmatrix% N0 :Noisevariance% iteration:Numberofiteration%% vHat :Decodedvector(0/1)%%%CopyrightBagawanS.Nugroho,2007%[MN]=size(H);%PriorprobabilitiesP1=ones(size(rx))./(1+exp(-2*rx./(N0/2)));P0=1-P1;%InitializationK0=zeros(M,N);K1=zeros(M,N);rji0=zeros(M,rji1=zeros(M,qij0=H.*repmat(P0',M,1);qij1=H.*repmat(P1',M,1);%Iterationforn=1:iterationfprintf('Iteration:%d\n',%Horizontalstepfori=%Findnon-zerosinthecolumnc1=find(H(i,:));fork=1:length(c1)%Getcolumnproductsofdrji=1;forl=ifl~=kdrji=drji*(qij0(i,c1(l))-qij1(i,c1(l)));end%forrji0(i,c1(k))=(1+drji)/2;rji1(i,c1(k))=(1-drji)/2;end%forend%fori%Verticalstepforj=1:N%Findnon-zerosintherowr1=find(H(:,j));fork=1:length(r1)%GetrowproductsofprodOfrij\ri(l)prodOfrij0=1;prodOfrij1=1;forl=ifl~=kprodOfrij0=prodOfrij0*rji0(r1(l),j);prodOfrij1=prodOfrij1*rji1(r1(l),j);end%for%UpdateconstantsK0(r1(k),j)=P0(j)*prodOfrij0;K1(r1(k),j)=P1(j)*prodOfrij1;%Updateqij0andqij1qij0(r1(k),j)=K0(r1(k),j)./(K0(r1(k),j)+K1(r1(k),j));qij1(r1(k),j)=K1(r1(k),j)./(K0(r1(k),j)+K1(r1(k),j));end%fork%UpdateconstantsKi0=P0(j)*prod(rji0(r1,j));Ki1=P1(j)*prod(rji1(r1,j));%GetQjQi0=Ki0/(Ki0+Qi1=Ki1/(Ki0+%DecodeifQi1>vHat(j)=elsevHat(j)=0;end%forend%fornfunctionvHat=decodeLogDomainSimple(rx,H,iteration)%Simplifiedlog-domainsumproductalgorithmLDPCdecoder%% rx :Receivedsignalvector(columnvector)% H :LDPCmatrix% iteration:Numberofiteration%% vHat :Decodedvector(0/1)%%%CopyrightBagawanS.Nugroho,2007%[MN]=size(H);%Priorlog-likelihood(simplified).Minussignisusedfor0/1to-1/1mappingLci=-rx';%InitializationLrji=zeros(M,Pibetaij=zeros(M,N);%AsscociatetheL(ci)matrixwithnon-zeroelementsofHLqij=H.*repmat(Lci,M,1);forn=1:iterationfprintf('Iteration:%d\n',n);%Getthesignandmagnitudeofalphaij=sign(Lqij);betaij=abs(Lqij);%Horizontalstepfori=%Findnon-zerosinthecolumnc1=find(H(i,:));%Gettheminimumofbetaijfork=1:length(c1)%MinimumofminOfbetaij=realmax;forl=ifl~=kifbetaij(i,c1(l))<minOfbetaijminOfbetaij=betaij(i,endendend%forl%Multiplicationalphaij\c1(k)(use'*'sincealphaijare-1/1s)prodOfalphaij=prod(alphaij