信息迭代处理时,当Aij为单位矩阵时,L_ij(k)对应的信息就是Aij中第k行第k列的“1”的信息;当Aij为单位矩阵的循环移位矩阵时,L_ij(k)对应的信息则是根据右移系数Iij推算出对应行的Aij第k列的信息。
2)在译码每次迭代过程中,如果对2类节点信息rij(b)和qij(b)、接收初始化信息Pi(1)和Pi(0)及译码中间环节变量分别存储,则需要存储器的数量为校验矩阵中非零矩阵数量的2倍以上;如果将存储器复用,即校验节点更新时,存储校验节点信息,当进入变量节点更新时,将此存储空间再用来存取变量节点信息,同时,不再为初始化信息和中间变量信息开辟存储空间,而将其并入2类节点信息的存储器中,往复利用,可以节省一半的存储空间.事实证明,当LDPC码长为1 536的情况下,优化前,所需空间约为282.6 KB,经复用,可以节省到约122.88 KB,极大地提高了资源利用率。
3)在用DSP实现BP译码算法时,在DSP系统外扩大容量SDRAM、Flash芯片来实现大容量数据存储。对于部分不参与迭代的变量数组,如从解调模块接收的信息数组、译码判决信息数组等,可以存在片外SDRAM中,如图2所示,需要时可以从片外存取,从而进一步节省DSP片内存储空间,
图2 DSP外扩存储器结构图
算法涉及的所有浮点型数组,均按照单精度浮点float型(32位)保存,而不采用64位的双精度浮点型;译码判决序列按照char型(8位)保存,而不采用整形数组(16位).通过上述优化,在原有的基础上进一步节省了DSP片内约60 KB的存储空间。经验证,以上阐述的存储器复用、片外存储变量、数据格式的选择等资源优化方式节省了约3/4的存储空间。
第四章程序验证译码效果
程序介绍:
使用码率为0.5的准循环LDPC码,码矩阵(1255,2510)
迭代30次,输入10帧数据。
先对数据通过生成矩阵G进行编码,进行BPSK调制,送入高斯白噪声信道,通过校验矩阵H进行译码操作。
利用对数域LLR-BP算法译码,并且计算出误码率。
并对只采用BPSK调制解调的误码性能和使用LDPC编解码的误码性能通过实验验证做错出比较
源程序如下:
%% 计算LDPC性能
clear all;
clc;
%------------初始化-------------%
EbN0 =0:0.5:8;%信噪比
iter2=30;%迭代次数
frame=10;%测量的帧数
rate=0.5;%误码
ber3=zeros(length(EbN0));
ber=zeros(length(EbN0));
q=251;
summ=0;
sumn=0;
count=0;
E0=eye(q);%μ¥λ¾ØÕó
EZ=zeros(q);%0¾ØÕó
load H_01_dual_diag.mat H_01; %校验矩阵
load P.mat P; %生成矩阵
[row_h col_h]=size(H_01);
%---------------将循环矩阵转换为实际的校验矩阵--------------------%
for i=1:row_h
for j=1:col_h
if H_01(i,j)==-1
H_0_matrix((i-1)*q+1:i*q,(j-1)*q+1:j*q)=circshift(EZ,H_01(i,j));
else
H_0_matrix((i-1)*q+1:i*q,(j-1)*q+1:j*q)=circshift(E0,H_01(i,j));
end
end
end
H=H_0_matrix;
for i = 1:length(EbN0)
for j = 1:frame
fprintf('Frame : %d', j);
%-------------------编码-------------------%
x=(sign(randn(1,size(P,1)))+1)/2;
y=mod(x*P,2);
u=[y x];
%-------------------调制-------------------%
bpskMod = 2*u - 1;
%---------------加入高斯信道---------------%
N0 = 1/(exp(EbN0(i)*log(10)/10));
tx =bpskMod+sqrt(N0/(2*rate))*randn(size(bpskMod));
for iii=1:2510
if tx(iii)<=0
tx1(iii)=0;
else
tx1(iii)=1;
end
end
comp=xor(u, tx1);
summ=sum(comp);
sumn=summ/2510;
ber(i)=ber(i)+sumn;
%---------------译码---------------%
vhat1 = decodeProbDomain(tx', H, N0,iter2);%
%---------------计算误码---------------%
[num3, rat3] = biterr(vhat1, u);
ber3(i) = (ber3(i) + rat3);
end
ber(i)=ber(i)/frame;
ber3(i)=ber3(i)/frame;
fprintf(' i=%f',i);
fprintf('ber3= %12.10f,ber= %12.10f',ber3(i),ber(i));
if count==1
ber3(i)=0;
elseif ber3(i)<1/(q*10*frame)
ber3(i)=1/(2*q*10*frame);
count=1;
end
end
%---------------误码曲线---------------%
semilogy(EbN0, ber3, '-b+',EbN0, ber, '-r<'); %
下面是通过运行程序,得出了对于两种不同的编解码方法误码率的比较,如下图所示:如图上方红线代表BPSK误码率,蓝线代表LDPC误码率
有图可以得出结论:
同的信噪比下,采用LDPC编解码技术,得到的误码率明显好于采用BPSK调制技术,并且误码率提高一个的量化级,如在信噪比为2db的情况下,采用LDPC技术误码率比BPSK技术提高了10倍以上。
参考文献:
[1]GALLAGER R G.Low-density parity-check codes[M].Boston:MIT Press,1963:70-81.
[2]GALLAGER R G.Low density parity check codes[J].IRE
Trans on Information Theory,1962,8(1):2l-28.
[3]MACKAY D J C,NEAL R M.Near Shannon limit performance of low-density parity-check codes[J].IE Electronics Let-ters,1996,32(18):1645-1646.
[4]FOSSORIER M P.Quasi-cyclic low-density parity-check codes from circulant permutation matrices [J].IEEE Trans on
Inform Theory,2004,50:1788-1793.
[5]LI Zongwang,CHEN Lei,ZENG Lingqi,et al.Efficient encoding of quasi-cyclic low-density parity-check codes[J].IEEETrans on Communications,2006,54(1):36-45.
[6]窦高奇,高俊,刘冰.准循环LDPC码快速编译码算法及
DSP实现[J].解放军理工大学学报,2008,9(4):324-327.
- function vHat = decodeProbDomain(rx, H, N0, iteration)
- % Probability-domain sum product algorithm LDPC decoder
- %
- % rx : Received signal vector (column vector)
- % H : LDPC matrix
- % N0 : Noise variance
- % iteration : Number of iteration
- %
- % vHat : Decoded vector (0/1)
-
- [M N] = size(H);
- % Prior probabilities
- P1 = ones(size(rx))./(1 + exp(-2*rx./(N0/2)));
- P0 = 1 - P1;
- % Initialization
- K0 = zeros(M, N);
- K1 = zeros(M, N);
- rji0 = zeros(M, N);
- rji1 = zeros(M, N);
- qij0 = H.*repmat(P0', M, 1);
- qij1 = H.*repmat(P1', M, 1);
- % Iteration
- for n = 1:iteration
-
- % fprintf('Iteration : %d', n);
-
- % ----- Horizontal step -----
- for i = 1:M
-
- % Find non-zeros in the column
- c1 = find(H(i, :));
-
- for k = 1:length(c1)
-
- % Get column products of drjic1(l)
- drji = 1;
- for l = 1:length(c1)
- if l~= k
- drji = drji*(qij0(i, c1(l)) - qij1(i, c1(l)));
- end
- end % for l
-
- rji0(i, c1(k)) = (1 + drji)/2;
- rji1(i, c1(k)) = (1 - drji)/2;
-
- end % for k
-
- end % for i
-
- % ------ Vertical step ------
- for j = 1:N
-
- % Find non-zeros in the row
- r1 = find(H(:, j));
-
- for k = 1:length(r1)
-
- % Get row products of prodOfriji(l)
- prodOfrij0 = 1;
- prodOfrij1 = 1;
- for l = 1:length(r1)
- if l~= k
- prodOfrij0 = prodOfrij0*rji0(r1(l), j);
- prodOfrij1 = prodOfrij1*rji1(r1(l), j);
- end
- end % for l
-
- % Update constants
- K0(r1(k), j) = P0(j)*prodOfrij0;
- K1(r1(k), j) = P1(j)*prodOfrij1;
-
- % Update qij0 and qij1
- qij0(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 % for k
-
- % Update constants
- Ki0 = P0(j)*prod(rji0(r1, j));
- Ki1 = P1(j)*prod(rji1(r1, j));
-
- % Get Qj
- Qi0 = Ki0/(Ki0 + Ki1);
- Qi1 = Ki1/(Ki0 + Ki1);
-
- % Decode Qj
- if Qi1 > Qi0
- vHat(j) = 1;
- else
- vHat(j) = 0;
- end
- % if rem(H*yo.',2) == 0,
- % break;
- % end
- end % for j
-
- end % for n
复制代码