CN106788459A - 一种快速收敛的估计ldpc码错误概率的方法 - Google Patents

Abstract

本发明公开了一种快速收敛的估计LDPC码错误概率的方法，解决了传统估计方法收敛速度慢，所需样本数多的问题。本发明为原始偏置分布乘以一个辅助函数，使新的偏置分布产生的噪声仅落入固定区域内，当该区域对应的记录噪声数目的直方图达到平坦时，对应的概率质量函数收敛，错误概率达到精确值。本发明通过较小的译码运算量可获得码字的较低错误概率；本发明减少了仿真所需耗费的译码运算量，缩短了仿真时间。

Description

一种快速收敛的估计LDPC码错误概率的方法

技术领域

本发明属于数字通信差错控制编码领域，具体涉及一种采用LDPC码编码的通信系统的错误概率估计方法。

背景技术

在数字通信系统中，衡量纠错码性能的指标是比特或码字的错误概率，比特或码字的错误概率反映了系统的可靠性。目前，获取编码通信系统错误概率的方法有三大类：解析法，计算机仿真方法及基于硬件平台的仿真方法。针对码长为有限长的LDPC码的性能估计，现有的解析法应用起来非常困难；并且基于硬件平台的仿真方法的开发周期长、成本高、可移植性差，因此仿真方法的应用获得了极大的推广。

目前，蒙特卡罗(Monte Carlo，MC)法是纠错码性能估计中最常用的仿真方法。它被广泛应用于估计通信系统的各个参量，其中最常见的，是估计通信系统的错误概率。但是，LDPC码在高信噪比下的错误概率很低，运用MC法估计该类码的性能时，估计器方差的收敛速度缓慢，需要耗费大量的仿真时间。

针对该问题，为减小仿真时间，研究者们提出了重点采样(Importance Sampling，IS)法，该方法是一种能降低估计器方差的仿真方法，一定程度上解决了MC法在估计低错误概率时遇到的问题，已在LDPC码的性能仿真中获得了广泛应用，缩短了仿真时间。这种方法的关键在于选取一个合适的偏置分布，使导致译码错误的样本点数在仿真总样本数中所占的比例升高，从而降低错误概率估计器的方差，缩短仿真时间。最早应用于LDPC码性能仿真中的IS法是经典IS法，在该方法中，偏置分布的获取依赖于分析码字结构及错误区域特性，当码字结构是非规则结构，且仿真中采用迭代译码时，上述分析变得十分困难，使得经典的IS法不能估计非规则结构的LDPC码的性能。

自适应重点采样(Adaptive Importance Sampling，AIS)方法是上述经典IS方法的一种改进方法，它无需分析LDPC码的码字结构或搜索陷阱集，能够在有限次的迭代过程中自适应地获得趋近于最优的偏置分布，因而能灵活应用于各种类型的LDPC码中。双自适应重点采样(Dual Adaptive Importance Sampling，DAIS)法和快速平坦直方图(FastFlat Histogram，FFH)法是AIS法在纠错码性能估计领域中的发展和应用。该类方法利用一个预先定义的变量来控制偏置分布的更新，使偏置分布逐渐趋近于最优。

发明人在实现本发明的过程中，发现现有技术中至少存在以下缺点和不足：

AIS方法和DAIS方法的收敛速度较慢，导致所需的样本数较多。

为解决上述问题，本发明提供了一种快速收敛的错误概率估计方法。该方法引入一个辅助函数，在DAIS法偏置分布的基础上乘以该辅助函数，人为地提高了仿真预设区域Ω内的噪声被抽取的概率，使得每次迭代过程中从原偏置分布中抽取的噪声仅落入Ω内，最后引入权重以消除两次偏置的影响。当仿真收敛时，得到该区域内概率质量函数的估计值，进一步获得LDPC码的错误概率。该方法与DAIS方法相比，能够加快仿真的收敛速度，减少仿真所需的样本数。

发明内容

本发明提供了一种快速收敛的估计LDPC码错误概率的方法，能加快传统的双自适应重点采样仿真的收敛速度，减少所需的仿真样本数，并且能够精确地估计采用LDPC码编码的通信系统在高信噪比下的错误概率，详见下文描述。

一种快速收敛的估计LDPC码错误概率的方法。所述方法包括以下步骤：

(1)利用蒙特卡罗仿真方法确定参考信噪比SNR_ref；

(2)输入待仿真的信噪比SNR_obj，判断输入信噪比SNR_obj是否小于SNR_ref，当SNR_obj≤SNR_ref时，采用蒙特卡罗法估计系统的错误概率，流程结束；当SNR_obj>SNR_ref时，执行步骤(3)；

(3)定义控制变量V，将V的范围划分为若干个子区间，确定可信区间；定义辅助函数，确定仿真预设区域Ω；

(4)估计在SNR_obj下噪声落入V的第k个子区间且导致译码错误的概率P_err,k,obj；

(5)估计在SNR_obj下导致译码错误的所有噪声中，落入V的第k个子区间的噪声占总数的条件概率P_k|err,obj；

(6)计算在SNR_obj下的错误概率P_err,obj。

所述利用蒙特卡罗仿真方法确定参考信噪比SNR_ref的步骤具体为：

(1.1)初始化SNR_ref为一合适值；

(1.2)采用蒙特卡罗仿真方法估计SNR_ref下系统的错误概率；

(1.3)判断所得错误概率是否在区间[10^-4,10^-3]之内，当错误概率落在区间[10^-4,10^-3]之内时，得到SNR_ref；当错误概率小于10^-4时，SNR_ref减小一个步长ΔSNR，执行步骤(1.2)～(1.3)直至满足判断条件，得到SNR_ref；当错误概率大于10^-3时，SNR_ref增加一个步长ΔSNR，执行步骤(1.2)～(1.3)直至满足判断条件，得到SNR_ref。

所述定义控制变量V，将V的范围划分为若干个子区间，确定可信区间；定义辅助函数，确定仿真预设区间Ω具体为：

(3.1)控制变量其中n为LDPC码的码长，z_l是噪声z的第l个分量，当z_l与调制后的信号符号相反时，H(q_lz_l)＝1，当z_l与调制后的信号符号相同时，H(q_lz_l)＝0；

(3.2)将V的取值范围划分为若干个子区间，SNR_ref下的划分步骤具体为，

(3.2.1)将SNR_ref下的控制变量范围[V_min,ref,V_max,ref]和子区间数目L_ref暂设定为[0,1]和1000，子区间的长度ΔV＝(V_max,ref-V_min,ref)/L_ref；

(3.2.2)运行不包含译码的无约束仿真，获得仿真样本；

(3.2.3)当仿真样本数超过10⁷时，选取第一个噪声数不为零的子区间k对应的V值为V_min,ref，选取最后一个噪声数不为零的子区间k对应的V值为V_max,ref；

(3.2.4)将[V_min,ref,V_max,ref]划分为L_ref＝L₀个子区间，L₀可根据不同的LDPC码的参数设定，并计算ΔV₀＝(V_max,ref-V_min,ref)/L₀；

(3.3)确定可信区间具体为，

(3.3.1)采用包含译码的无约束仿真估计在SNR_ref下噪声落入V的第k个子区间，且导致译码错误的联合概率P_err,k,ref；采用约束仿真估计SNR_ref下在导致译码错误的所有噪声中，落入V的第k个子区间的噪声占总数的条件概率P_k|err,ref；

(3.3.2)选取联合概率P_err,k,ref不为零的一段连续区域；

(3.3.3)在P_err,k,ref不为零的连续区域内，从最左端开始，选取子区间数为δ的区域并暂定为[V^*,V^**]；

(3.3.4)在子区间数为δ的区域内计算错误概率P_err,ref＝P_err,k,ref/P_k|err,ref，与MC的仿真结果进行对比，当错误概率P_err,ref与MC的仿真结果之间的相对误差小于等于10％时，得到可信区间[V^*,V^**]，当相对误差大于10％时，保持区域内的子区间数为δ，增大区间左右端点的值，重复步骤(3.3.4)以修正可信区间[V^*,V^**]；

(3.4)辅助函数其中β>1，γ＝2；Ω为[k_min,L]。

所述估计在SNR_obj下噪声落入V的第k个子区间且导致译码错误的概率P_err,k,obj具体为：

(4.1)初始化不包含译码的无约束仿真中的控制变量范围[V_min,obj,V_max,obj]，子区间长度ΔV及子区间的数目L_obj，具体为，

(4.1.1)将SNR_obj下的控制变量范围[V_min,obj,V_max,obj]和子区间数目L_obj暂设定为[0,V^**]和1000，子区间的长度ΔV＝(V_max,obj-V_min,obj)/L_obj；

(4.1.2)运行不包含译码的无约束仿真，获得仿真样本；

(4.1.3)选取噪声数不为零且序号k最小的子区间对应的V值为V_min,obj，V_max,obj＝V^**；

(4.1.4)将[V_min,obj,V_max,obj]划分为L_obj个子区间，L_obj＝(V_max,obj-V_min,obj)/ΔV₀；

(4.2)运行包含译码的无约束仿真，估计P_err,k,obj。

所述估计在SNR_obj下导致译码错误的所有噪声中，落入V的第k个子区间的噪声占总数的条件概率P_k|err,obj具体为：

(5.1)初始化约束仿真中的控制变量范围[V′_min,obj,V′_max,obj]，子区间长度ΔV及子区间的数目L'_obj；具体为，

(5.1.1)将[V′_min,obj,V′_max,obj]暂设定为[0,V_max,obj]，即[0,V^**]，子区间长度ΔV＝ΔV₀，运行包含译码的约束仿真，获得仿真样本；

(5.1.2)选取第一个噪声数目不等于零的子区间对应的V值为V′_min,obj，则控制变量范围为[V′_min,obj,V^**]；

(5.1.3)计算L'_obj＝(V′_max,obj-V′_min,obj)/ΔV＝(V^**-V′_min,obj)/ΔV₀；

(5.2)运行约束仿真，估计P_k|err,obj。

所述计算在SNR_obj下的错误概率P_err,obj具体为：

在可信区间内计算平均值

所述步骤(3.2.2)、(4.1.2)中的不包含译码的无约束仿真具体为：

1)设定噪声个数i＝0，迭代次数j＝1，第一次迭代过程中需要的总噪声数目为N₁，噪声落入控制变量第k个子区间的概率其中1≤k≤L；

2)选取一个导致译码器出现判决错误的噪声z₀，i增1；

3)采用Metropolis算法生成一个新的噪声，i增1；

4)计算新的噪声落入的控制变量子区间k，k内的噪声数增1；

5)若i<N₁，重复步骤3)～5)；否则，执行步骤6)；

6)采用Berg递归公式计算

7)若输出否则，i＝0，j增1，下一次迭代过程中的总样本数等于本次迭代中总样本数的1.3倍，重复步骤2)至7)，直至输出

所述步骤(4.2)中包含译码的无约束仿真具体为：

1)设定噪声个数i＝0，迭代次数j＝1，迭代过程中需要的总噪声数目为N₁，N₁'，第一次迭代中

2)选取一个导致译码器出现判决错误的噪声z₀，i增1；

3)采用Metropolis算法生成一个新的噪声，i增1；

4)计算新的噪声落入控制变量的子区间号k，k内的噪声数增1；

5)将新的噪声加在全零码字上，送入LDPC译码器；若译码器输出端的错误码字数为1，第k个子区间内的错误帧数增1；若译码器输出端的错误码字数为0，第k个子区间内的错误帧数不变；

6)若i<N₁，重复步骤3)～6)直至不满足判断条件；否则，执行步骤7)；

7)采用Berg递归公式计算

8)若计算Ω区域内的概率和j增1，执行步骤9)；否则，i＝0，j增1，重复步骤2)～8)，直至执行步骤9)；

9)采用Metropolis算法从引入了辅助函数之后的偏置分布中抽取一个新的噪声，i增1；计算新的噪声落入的控制变量的区间号k，k内的噪声数增1；

10)将新的噪声加在全零的LDPC码字上，送入译码器进行判决；若新的噪声导致译码器出现判决错误，k内的错误帧数增1；若未导致错误，k内的错误帧数不变；

11)若i<N₁'，重复步骤9)～11)，直至判断条件成立，执行步骤12)；否则，执行步骤12)；

12)采用Berg递归公式计算进一步计算

13)若输出并计算P_err|k；否则，i＝0，j增1，重复步骤9)～13)，直到Ω区域内的P_k估计值收敛，输出并计算P_err|k；

14)计算P_err,k＝P_kP_err|k。

所述步骤(3.3.1)、(5.1)、(5.2)中约束仿真具体为：

1)设定噪声个数i＝0，迭代次数j＝1，迭代过程中需要的总噪声数目为N₂，N'₂，第一次迭代中

2)选取一个导致译码器出现判决错误的噪声z₀，i增1；

3)采用修改后的Metropolis算法生成一个新的噪声，i增1；

4)计算新的噪声落入控制变量的子区间号k，k内的噪声数增1；

5)若i<N₂，重复步骤3)～5)直至判断条件不成立；否则，执行步骤6)；若判断条件成立；

6)采用Berg递归公式计算

7)若计算Ω区域内的概率和j增1，执行步骤8)；否则，i＝0，j增1，重复步骤2)～7)，直至满足执行步骤8)；

8)采用修改后的Metropolis算法从引入了辅助函数之后的偏置分布中抽取一个新的噪声，i增1；计算新的噪声落入的控制变量的区间号k，k内的噪声数增1；

9)若i<N'₂，重复步骤8)～9)，直至i＝N'₂，执行步骤10)；否则，执行步骤10)；

10)采用Berg递归公式计算进一步计算

11)若输出并计算P_k|err；否则，i＝0，j增1，重复步骤8)～11)，直到Ω区域内的P_k|err估计值收敛，输出并计算P_k|err。

所述的无约束仿真的步骤3)修改后的Metropolis算法具体为：

1)提取噪声个数i，迭代次数j，噪声z^*,(i)，表示噪声维数的l＝1；

2)令噪声的第l维分量的下一个可能状态其中，ε₂为最大步长因子，调整ε₂使得每一次迭代过程中被接受的状态数目除以总的状态数目接近0.234，Δz～N(0,σ²)；

3)产生两个服从U(0,1)均匀分布的随机数α₁与α₂，判断

4)若l等于码长n，得到下一个可能的状态执行步骤5)；否则，l增1，重复步骤2)～4)，直至满足条件l＝n，执行步骤5)；

5)计算与

6)判断E表示能够导致译码器出现判决错误的噪声的集合；

7)输出z^*,(i+1)。

本发明引入一个辅助函数，使得噪声几乎全部落入V值较大的区域，加快了该区域内记录噪声数目的直方图趋于平坦的速度，进而加快了仿真的收敛速度，获得了精确的LDPC码的错误概率。

附图说明

图1是本发明的快速收敛的估计LDPC码错误概率方法的流程图；

图2是估计P_err,k,obj的流程图；

图3是图2中的Metropolis算法的流程图；

图4是估计P_k|err,obj的流程图；

图5是图4中的修改后的Metropolis算法的流程图；

图6是本发明的方法和DAIS方法的直方图在Ω区域内的收敛速度的对比示意图；

图7是MacKay(96,48)LDPC码的仿真结果图。

具体实施方式

下面结合附图对本发明的快速收敛的估计LDPC码错误概率的方法做详细说明。

如图1所示，本发明包括以下六个步骤，

(1)利用蒙特卡罗仿真方法确定参考信噪比SNR_ref；

(2)输入待仿真的信噪比SNR_obj，判断输入信噪比SNR_obj是否小于SNR_ref，当SNR_obj≤SNR_ref时，采用蒙特卡罗法估计系统的错误概率，流程结束；当SNR_obj>SNR_ref时，执行步骤(3)；

(3)定义控制变量V，将V的范围划分为若干个子区间；定义辅助函数，确定仿真预设区域Ω；

(4)估计在SNR_obj下噪声落入V的第k个子区间且导致译码错误的概率P_err,k,obj；

(5)估计在SNR_obj下导致译码错误的所有噪声中，落入V的第k个子区间的噪声占总数的条件概率P_k|err,obj；

(6)计算在SNR_obj下的错误概率P_err,obj。

以下分别介绍上述六个步骤的具体实施步骤。

所述的步骤(1)具体为

(1.1)初始化SNR_ref为一合适值；

(1.2)采用蒙特卡罗仿真方法估计SNR_ref下系统的错误概率；

(1.3)判断所得错误概率是否在区间[10^-4,10^-3]之内，当错误概率落在区间[10^-4,10^-3]之内时，得到SNR_ref；当错误概率小于10^-4时，SNR_ref减小一个步长ΔSNR，执行步骤(1.2)～(1.3)直至满足判断条件，得到SNR_ref；当错误概率大于10^-3时，SNR_ref增加一个步长ΔSNR，执行步骤(1.2)～(1.3)直至满足判断条件，得到SNR_ref。

所述的步骤(3)具体为：

(3.1)控制变量其中n为LDPC码的码长，z_l是噪声z的第l个分量，当z_l与调制后的信号符号相反时，H(q_lz_l)＝1，当z_l与调制后的信号符号相同时，H(q_lz_l)＝0；

(3.2)将V的取值范围划分为若干个子区间，SNR_ref下的划分步骤具体为，

(3.2.1)将SNR_ref下的控制变量范围[V_min,ref,V_max,ref]和子区间数目L_ref暂设定为[0,1]和1000，子区间的长度ΔV＝(V_max,ref-V_min,ref)/L_ref；

(3.2.2)运行不包含译码的无约束仿真，获得仿真样本；

(3.2.3)选取第一个噪声数不为零的子区间k对应的V值为V_min,ref，选取最后一个噪声数不为零的子区间k对应的V值为V_max,ref；

(3.2.4)将[V_min,ref,V_max,ref]划分为L_ref＝L₀个子区间，L₀可根据不同的LDPC码的参数设定，并计算ΔV₀＝(V_max,ref-V_min,ref)/L₀；

(3.3)确定可信区间具体为，

(3.3.1)采用包含译码的无约束仿真估计在SNR_ref下噪声落入V的第k个子区间，且导致译码错误的联合概率P_err,k,ref；采用约束仿真估计SNR_ref下在导致译码错误的所有噪声中，落入V的第k个子区间的噪声占总数的条件概率P_k|err,ref；

(3.3.2)选取联合概率P_err,k,ref不为零的一段连续区域；

(3.3.3)在P_err,k,ref不为零的连续区域内，从最左端开始，选取子区间数为δ的区域并暂定为[V^*,V^**]；

(3.3.4)在子区间数为δ的区域内计算错误概率P_err,ref＝P_err,k,ref/P_k|err,ref，与MC的仿真结果进行对比，当错误概率P_err,ref与MC的仿真结果之间的相对误差小于等于10％时，得到区间[V^*,V^**]，当相对误差大于10％时，保持区域内的子区间数为δ，增大区间左右端点的值，重复步骤(3.3.4)以修正区间[V^*,V^**]；

(3.4)辅助函数其中β>1，γ＝2；Ω为[k_min,L]。

所述的步骤(4)具体为：

(4.1)初始化不包含译码的无约束仿真中的控制变量范围[V_min,obj,V_max,obj]，子区间长度ΔV及子区间的数目L_obj，具体为，

(4.1.1)将SNR_obj下的控制变量范围[V_min,obj,V_max,obj]和子区间数目L_obj暂设定为[0,V^**]和1000，子区间的长度ΔV＝(V_max,obj-V_min,obj)/L_obj；

(4.1.2)运行不包含译码的无约束仿真，获得仿真样本；

(4.1.3)选取噪声数不为零且序号k最小的子区间对应的V值为V_min,obj，V_max,obj＝V^**；

(4.1.4)将[V_min,obj,V_max,obj]划分为L_obj个子区间，L_obj＝(V_max,obj-V_min,obj)/ΔV₀；

(4.2)如图2所示，运行包含译码的无约束仿真，估计P_err,k,obj，具体为，

(4.2.1)设定噪声个数i＝0，迭代次数j＝1，迭代过程中需要的噪声数目为N₁和N₁'，第一次迭代中

(4.2.2)从原始高斯分布中抽取一个噪声z₀，若z₀能够导致译码器出现判决错误，令初始噪声等于z₀，i增1；否则，重复步骤(4.2.2)，直至z₀导致译码器出现判决错误，令初始噪声等于z₀，i增1；

(4.2.3)采用Metropolis算法生成一个新的噪声z^*,(i)＝(z₁,…,z_n)，其中n为码长，i增1；

(4.2.4)计算新的噪声z^*,(i)落入的控制变量的子区间号k＝(V(z^*,(i))-V_min,obj)/ΔV；第k个所述子区间内的噪声数目增1；

(4.2.5)将新的噪声加在全零码字上，送入LDPC译码器；若译码器输出端的错误码字数为1，第k个子区间内的错误帧数加1；若译码器输出端的错误码字数为0，第k个子区间内的错误帧数不变；

(4.2.6)若i<N₁，重复步骤(4.2.3)～(4.2.6)直至不满足判断条件；否则，执行步骤(4.2.7)；

(4.2.7)采用如下Berg递归公式计算第j+1次迭代中(1≤k≤L)

其中

(4.2.8)若计算Ω区域内的概率和j增1，执行步骤(4.2.9)；否则，i＝0，j增1，重复所述步骤(4.2.2)～(4.2.8)，直至满足收敛条件，计算执行步骤(4.2.9)；

(4.2.9)采用Metropolis算法从引入了辅助函数之后的偏置分布中抽取一个新的噪声z^*,(i)＝(z₁,…,z_n)，i增1；计算新的噪声落入的控制变量的区间号k＝(V(z^*,(i))-V_min,obj)/ΔV，k内的噪声数增1；

(4.2.10)将噪声z^*,(i)加在全零的LDPC码字上，送入所述译码器进行判决，若新的噪声导致译码器出现判决错误，增1；否则，不变；

(4.2.11)若i≤N₁'，重复所述步骤(4.2.9)～(4.2.11)，直至i＝N₁'，执行步骤(4.2.12)；否则，执行步骤(4.2.12)；

(4.2.12)采用Berg递归公式计算第j+1次迭代中(1≤k≤L)

其中

进一步，

(4.2.13)若输出并计算执行步骤(4.2.14)；否则，i＝0，j增1，重复步骤(4.2.9)～(4.2.13)，直到Ω区域内的P_k,obj估计值收敛，输出并计算P_err|k,obj，执行步骤(4.2.14)。

(4.2.14)计算P_err,k,obj＝P_k,objP_err|k,obj。

如图3所示，步骤(4.2.3)中所述Metropolis算法具体为：

1)提取噪声个数i，迭代次数j，噪声z^*,(i)，表示噪声维数的l＝1；

2)令噪声的第l维分量的下一个可能状态其中，ε₁为最大步长因子，调整ε₁使得每一次迭代过程中被接受的状态数目除以总的状态数目接近0.234，Δz～N(0,σ²)；

3)产生两个服从U(0,1)均匀分布的随机数α₁与α₂，判断

4)若l等于码长n，得到下一个可能的状态执行步骤5)；否则，l增1，重复步骤2)～4)，直至满足条件l＝n，执行步骤5)；

5)计算与

6)判断

7)输出z^*,(i+1)。

所述的步骤(5)具体为：

(5.1)初始化约束仿真中的控制变量范围[V′_min,obj,V′_max,obj]，子区间长度ΔV及子区间的数目L'_obj；具体为，

(5.1.1)将[V′_min,obj,V′_max,obj]暂设定为[0,V_max,obj]即[0,V^**]，子区间长度ΔV＝ΔV₀，运行包含译码的约束仿真，获得仿真样本；

(5.1.2)选取第一个噪声数目不等于零的子区间对应的V值为V′_min,obj，则控制变量范围为[V′_min,obj,V^**]；

(5.1.3)计算L'_obj＝(V′_max,obj-V′_min,obj)/ΔV＝(V^**-V′_min,obj)/ΔV₀；

(5.2)如图4所示，运行约束仿真，估计P_k|err,obj，具体为：

(5.2.1)设定噪声个数i＝0，迭代次数j＝1，迭代过程中需要的总噪声数目为N₂和N'₂，第一次迭代中

(5.2.2)初始噪声等于z₀，i加1；

(5.2.3)采用修改后的Metropolis算法生成一个新的噪声z^*,(i)＝(z₁,…,z_n)，其中n为码长，i加1；

(5.2.4)计算新的噪声z^*,(i)落入的所述控制变量的子区间号k＝(V(z^*,(i))-V′_min,obj)/ΔV，第k个所述子区间内的噪声数目增1；

(5.2.5)若i<N₂，重复所述步骤(5.2.3)～(5.2.5)，直至不满足判断条件，执行步骤(5.2.6)；否则，执行步骤(5.2.6)；

(5.2.6)采用Berg递归公式计算第j+1次迭代中(1≤k≤L'_obj)

其中

(5.2.7)若计算Ω区域内的概率和j增1，执行步骤(5.2.8)；否则，i＝0，j增1，下一次迭代过程中总的仿真样本数N₂＝1.3×N₂，重复所述步骤(5.2.2)～(5.2.7)，直至满足计算P_Ω，执行步骤(5.2.8)；

(5.2.8)采用修改后的Metropolis算法，从引入了辅助函数之后的偏置分布中抽取一个新的噪声z^*,(i)＝(z₁,…,z_n)，i增1；计算噪声z^*,(i)落入的所述控制变量的子区间号k，第k个所述子区间内的噪声数目增1；

(5.2.9)若i<N'₂，重复所述步骤(5.2.8)～(5.2.9)，直至不满足所述判断条件，执行步骤(5.2.10)；否则，执行步骤(5.2.10)；

(5.2.10)采用Berg递归公式计算第j+1次迭代中的(1≤k≤L'_obj)

其中

进一步，

(5.2.11)若输出否则，i＝0，j增1，重复步骤(5.2.8)～(5.2.11)，直到Ω区域内的P_k|err,obj估计值收敛，输出

如图5所示，所述的步骤(5.2.3)中修改后的Metropolis算法具体为：

1)提取噪声个数i，迭代次数j，噪声z^*,(i)，表示噪声维数的l＝1；

2)令噪声的第l维分量的下一个可能状态其中，ε₂为最大步长因子，调整ε₂使得每一次迭代过程中被接受的状态数目除以总的状态数目接近0.234，Δz～N(0,σ²)；

3)产生两个服从U(0,1)均匀分布的随机数α₁与α₂，判断

4)若l等于码长n，得到下一个可能的状态执行步骤5)；否则，l增1，重复步骤2)～4)，直至满足条件l＝n，执行步骤5)；

5)计算与

6)判断E表示能够导致译码器出现判决错误的噪声的集合；

7)输出z^*,(i+1)。

所述的步骤(6)具体为：

在可信区间内计算平均值

具体实施例

本发明以码长为96，码率为1/2的MacKay(96,48)LDPC码为一个特例，介绍本发明的一种快速收敛的估计LDPC码错误概率的方法。

仿真中，发送信号为全零码字，调制方式为二进制相移键控调制，已调信号为{+1,…,+1}。LDPC码的译码器采用置信度传播译码算法，迭代次数为5次。步骤(1)中SNR_ref为4dB，SNR_obj为大于4dB的实数，P_err指的是误帧率(FER)，也可指误比特率(BER)。

步骤(1.1)中SNR_ref＝1dB，步骤(1.3)中ΔSNR＝1dB。步骤(3.3)中的区间[V^*,V^**]＝[0.5,0.55]，δ＝30。步骤(4.1.1)中α＝0.1240。步骤(4.2.1)中，第一次迭代中应产生的总的噪声数目N₁＝5000，N₁'＝10000，图3中的最大步长因子ε₁＝1.5。步骤(5.2.1)中，第一次迭代中应产生的总噪声数目N₂＝10000，N'₂＝100000，图5中的最大步长因子ε₂＝0.5。

如图6所示，图(a)是本发明提供的仿真方法的每个子区间内的噪声数目在Ω区域内趋于平坦的速度，图(b)是DAIS法的在Ω区域内趋于平坦的速度。由于仿真所需的样本数与Ω区域内P_k及P_k|err估计值的收敛速度成正比，收敛速度越快，仿真所需的样本数越少；并且，直方图趋于平坦的速度越快，估计值收敛于真实值的速度越快。从图中可以看出，在Ω区域内，本发明设计的方法的直方图趋于平坦的速度较快，因此估计值收敛于真实值的速度较快，所需的仿真样本数较少。

仿真结果显示，在低信噪比区域，本发明提供的方法得到的FER及BER估计值与MC法的一致；在高信噪比区域，该方法能获得LDPC码的低FER及BER估计值，如图7所示。

表1给出了MC法、DAIS法及快速收敛的错误概率估计法这三种仿真方法所需译码的样本数及仿真增益，仿真增益G₁等于MC法所需译码的样本数与其它两种方法所需译码的样本数之比，仿真增益G₂等于DAIS法所需译码的样本数与其它两种方法所需译码的样本数之比。其中，MC仿真所需的译码样本数是在保证估计值精确度为10％的条件下，由N≥100/FER估计得到。仿真结果表明，本发明的快速收敛的估计LDPC码的错误概率方法所需的仿真样本数少。

表1MC法、DAIS法及快速收敛的错误概率估计法所需译码的样本数及仿真增益

Claims (9)

1.一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述方法包括以下步骤：

(1)利用蒙特卡罗仿真方法确定参考信噪比SNR_ref；

(2)输入待仿真的信噪比SNR_obj，判断输入信噪比SNR_obj是否小于SNR_ref，当SNR_obj≤SNR_ref时，采用蒙特卡罗法估计系统的错误概率，流程结束；当SNR_obj>SNR_ref时，执行步骤(3)；

(3)定义控制变量V，将V的范围划分为若干个子区间，确定可信区间；定义辅助函数，确定仿真预设区域Ω；

(4)估计在SNR_obj下噪声落入V的第k个子区间且导致译码错误的概率P_{err，k，obj}；

(5)估计在SNR_obj下导致译码错误的所有噪声中，落入V的第k个子区间的噪声占总数的条件概率P_k|err,obj；

(6)计算在SNR_obj下的错误概率P_err,obj。

2.根据权利要求1所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述定义控制变量V，将V的范围划分为若干个子区间，确定可信区间；定义辅助函数，确定仿真预设区间Ω步骤具体为：

(3.1)控制变量其中n为LDPC码的码长，z_l是噪声z的第l个分量，当z_l与调制后的信号符号相反时，H(q_lz_l)＝1，当z_l与调制后的信号符号相同时，H(q_lz_l)＝0；

(3.2)将V的取值范围划分为若干个子区间，SNR_ref下的划分步骤具体为，

(3.2.1)将SNR_ref下的控制变量范围[V_min,ref,V_max,ref]和子区间数目L_ref暂设定为[0,1]和1000，子区间的长度ΔV＝(V_max,ref-V_min,ref)/L_ref；

(3.2.2)运行不包含译码的无约束仿真，获得仿真样本；

(3.2.3)当仿真样本数超过10⁷时，选取第一个噪声数不为零的子区间k对应的V值为V_min,ref，选取最后一个噪声数不为零的子区间k对应的V值为V_max,ref；

(3.2.4)将[V_min,ref,V_max,ref]划分为L_ref＝L₀个子区间，L₀可根据不同的LDPC码的参数设定，并计算ΔV₀＝(V_max,ref-V_min,ref)/L₀；

(3.3)确定可信区间具体为，

(3.3.1)采用包含译码的无约束仿真估计在SNR_ref下噪声落入V的第k个子区间，且导致译码错误的联合概率P_err,k,ref；采用约束仿真估计SNR_ref下在导致译码错误的所有噪声中，落入V的第k个子区间的噪声数占总数的条件概率P_k|err,ref；

(3.3.2)选取联合概率P_err,k,ref不为零的一段连续区域；

(3.3.3)在P_err,k,ref不为零的连续区域内，从最左端开始，选取子区间数为δ的区域并暂定为[V^*,V^**]；

(3.3.4)在子区间数为δ的区域内计算错误概率P_err,ref＝P_err,k,ref/P_k|err,ref，与MC的仿真结果进行对比，当错误概率P_err,ref与MC的仿真结果之间的相对误差小于等于10％时，得到可信区间[V^*,V^**]，当相对误差大于10％时，保持区域内的子区间数为δ，增大区间左右端点的值，重复步骤(3.3.4)以修正可信区间[V^*,V^**]；

(3.4)辅助函数其中β>1，γ＝2；Ω为[k_min,L]。

3.根据权利要求2所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述估计在SNR_obj下噪声落入V的第k个子区间且导致译码错误的概率P_err,k,obj具体为：

(4.1)初始化不包含译码的无约束仿真中的控制变量范围[V_min,obj,V_max,obj]，子区间长度ΔV及子区间的数目L_obj，具体为，

(4.1.1)将SNR_obj下的控制变量范围[V_min,obj,V_max,obj]和子区间数目L_obj暂设定为[0,V^**]和1000，子区间的长度ΔV＝(V_max,obj-V_min,obj)/L_obj；

(4.1.2)运行不包含译码的无约束仿真，获得仿真样本；

(4.1.3)选取噪声数不为零且序号k最小的子区间对应的V值为V_min,obj，V_max,obj＝V^**；

(4.1.4)将[V_min,obj,V_max,obj]划分为L_obj个子区间，L_obj＝(V_max,obj-V_min,obj)/ΔV₀；

(4.2)运行包含译码的无约束仿真，估计P_err,k,obj。

4.根据权利要求2所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述估计在SNR_obj下导致译码错误的所有噪声中，落入V的第k个子区间的噪声数占总数的条件概率P_k|err,obj具体为：

(5.1)初始化约束仿真中的控制变量范围[V′_min,obj,V′_max,obj]，子区间长度ΔV及子区间的数目L'_obj；具体为，

(5.1.1)将[V′_min,obj,V′_max,obj]暂设定为[0,V_max,obj]，即[0,V^**]，子区间长度ΔV＝ΔV₀，运行包含译码的约束仿真，获得仿真样本；

(5.1.2)选取第一个噪声数目不等于零的子区间对应的V值为V′_min,obj，则控制变量范围为[V′_min,obj,V^**]；

(5.1.3)计算L'_obj＝(V′_max,obj-V′_min,obj)/ΔV＝(V^**-V′_min,obj)/ΔV₀；

(5.2)运行约束仿真，估计P_k|err,obj。

5.根据权利要求1所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述计算在SNR_obj下的错误概率P_err,obj具体为，

在可信区间内计算平均值

6.根据权利要求3所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，步骤(3.2.2)、(4.1.2)中不包含译码的无约束仿真具体为，

1)设定噪声个数i＝0，迭代次数j＝1，第一次迭代过程中需要的总噪声数目为N₁，噪声落入控制变量第k个子区间的概率其中1≤k≤L；

2)选取一个导致译码器出现判决错误的噪声z₀，i增1；

3)采用Metropolis算法生成一个新的噪声，i增1；

4)计算新的噪声落入的控制变量子区间k，k内的噪声数增1；

5)若i<N₁，重复步骤3)～5)；否则，执行步骤6)；

6)采用Berg递归公式计算

7)若输出否则，i＝0，j增1，下一次迭代过程中的总样本数等于本次迭代中总样本数的1.3倍，重复步骤2)至7)，直至输出

7.根据权利要求2所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，步骤(3.2)中包含译码的无约束仿真具体为，

1)设定噪声个数i＝0，迭代次数j＝1，迭代过程中需要的总噪声数目为N₁和N′₁，第一次迭代中Ω＝[k_min,L]；

2)选取一个导致译码器出现判决错误的噪声z₀，i增1；

3)采用Metropolis算法生成一个新的噪声，i增1；

4)计算新的噪声落入控制变量的子区间号k，k内的噪声数增1；

5)将新的噪声加在全零码字上，送入LDPC译码器；若译码器输出端的错误码字数为1，第k个子区间内的错误帧数增1；若译码器输出端的错误码字数为0，第k个子区间内的错误帧数不变；

6)若i<N₁，重复步骤3)～6)直至不满足判断条件；否则，执行步骤7)；

7)采用Berg递归公式计算

8)若计算Ω区域内的概率和j增1，执行步骤9)；否则，i＝0，j增1，重复步骤2)～8)，直至执行步骤9)；

9)采用Metropolis算法从引入了辅助函数之后的偏置分布中抽取一个新的噪声，i增1；计算新的噪声落入的控制变量的区间号k，k内的噪声数增1；

10)将新的噪声加在全零的LDPC码字上，送入译码器进行判决；若新的噪声导致译码器出现判决错误，k内的错误帧数增1；若未导致错误，k内的错误帧数不变；

11)若i<N′₁，重复步骤9)～11)，直至判断条件成立，执行步骤12)；否则，执行步骤12)；

12)采用Berg递归公式计算进一步计算

13)若输出并计算P_err|k；否则，i＝0，j增1，重复步骤9)～13)，直到Ω区域内的P_k估计值收敛，输出并计算P_err|k；

14)计算P_err,k＝P_kP_err|k。

8.根据权利要求4所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述的步骤(3.3.1)、(5.1)、(5.2)中约束仿真具体为，

1)设定噪声个数i＝0，迭代次数j＝1，迭代过程中需要的总噪声数目为N₂和N'₂，第一次迭代中Ω＝[k_min,L']；

2)选取一个导致译码器出现判决错误的噪声z₀，i增1；

3)采用修改后的Metropolis算法生成一个新的噪声，i增1；

4)计算新的噪声落入控制变量的子区间号k，k内的噪声数增1；

5)若i<N₂，重复步骤3)～5)直至该判断条件不成立；否则，执行步骤6)；若判断条件成立；

6)采用Berg递归公式计算

7)若计算Ω区域内的概率和i＝0，j增1，执行步骤8)；否则，i＝0，j增1，重复步骤2)～7)，直至满足执行步骤8)；

8)采用修改后的Metropolis算法从引入了辅助函数之后的偏置分布中抽取一个新的噪声，i增1；计算新的噪声落入的控制变量的区间号k，k内的噪声数增1；

9)若i<N'₂，重复步骤8)～9)，直至i＝N'₂，执行步骤10)；否则，执行步骤10)；

10)采用Berg递归公式计算进一步计算

11)若输出并计算P_k|err；否则，i＝0，j增1，重复步骤8)～11)，直到Ω区域内的P_k|err估计值收敛，输出并计算P_k|err。

9.根据权利要求8所述的一种快速收敛的估计LDPC码错误概率的方法，其特征在于，所述的步骤3)修改后的Metropolis算法具体为：

1)提取噪声个数i，迭代次数j，噪声z^*,(i)，表示噪声维数的l＝1；

2)令噪声的第l维分量的下一个可能状态其中，ε₂为最大步长因子，调整ε₂使得每一次迭代过程中被接受的状态数目除以总的状态数目接近0.234，Δz～N(0,σ²)；

3)产生两个服从U(0,1)均匀分布的随机数α₁与α₂，判断

$\left\{\begin{array}{cc}{\mathrm{\&alpha;}}_1\&le;\min \&lsqb;1,{\mathrm{\&rho;}}_l\left(z_{b,l}^{*}\right)/{\mathrm{\&rho;}}_l\left(z_{a,l}^{*}\right)\&rsqb;,& z_l^{*,(i+1)}=z_{b,l}^{*}\\ {\mathrm{\&alpha;}}_1>\min \&lsqb;1,{\mathrm{\&rho;}}_l\left(z_{b,l}^{*}\right)/{\mathrm{\&rho;}}_l\left(z_{a,l}^{*}\right)\&rsqb;,& z_l^{*,(i+1)}=z_{a,l}^{*}\end{array}\right.;$

4)若l等于码长n，得到下一个可能的状态执行步骤5)；否则，l增1，重复步骤2)～4)，直至满足条件l＝n，执行步骤5)；

5)计算与

6)判断E表示能够导致译码器出现判决错误的噪声的集合；

7)输出z^*,(i+1)。