基于分段抽取软判决加权Walsh Hadamard变换的卷积码识别算法

张岱，张玉，杨晓静，樊斌斌（.电子工程学院脉冲功率激光技术国家重点实验室，安徽合肥007；.电子工程学院504教研室，安徽合肥007；.电子工程学院40教研室，安徽合肥007）

基于分段抽取软判决加权Walsh Hadamard变换的卷积码识别算法

张岱1，张玉2，杨晓静3，樊斌斌2
（1.电子工程学院脉冲功率激光技术国家重点实验室，安徽合肥230037；2.电子工程学院504教研室，安徽合肥230037；3.电子工程学院402教研室，安徽合肥230037）

针对低信噪比环境下卷积码识别研究存在的不足，提出一种基于分段抽取软判决加权Walsh Hadamard变换（WHT）的卷积码识别算法。该算法利用接收比特的解调软判决信息求取软判决频次序列，并构造加权Walsh Hadamard矩阵，从而识别得到基本校验序列。利用基本编码矩阵构造规则，最终实现记忆长度及基本生成矩阵的识别。算法结合分段抽取思想，降低了所需运行存储量。仿真实验表明，该算法可在低信噪比环境下对不同码率卷积码进行有效识别，具有较好的容错性，且对大约束度的卷积码性能提高更为显著。

信息处理技术；卷积码识别；低信噪比；软判决信息；加权Walsh Hadamard变换；分段抽取

0　引言

信道编码识别技术在非合作信号处理、智能通信等领域具有广泛应用[1]。其中，在非合作信号处理领域的应用通常是在低信噪比环境下（对于处理的基带信号，信噪比通常在5 dB以下，解调后得到的编码序列误码率在10-1～10-2量级）进行，因此研究一种容错性高的卷积码识别算法将具有重要意义。

目前，对卷积码的识别研究主要有基于欧几里德的识别算法[2]、基于 Boyer-Moore的快速合冲法[3]、基于Walsh Hadamard变换（WHT）的识别算法[4]、基于校验序列的识别算法[5-6]、基于矩阵分析列相关性的识别算法[7-8]和基于校验方程组成立概率的识别算法[9]等。其中基于欧几里德的识别算法适用于1/n码率卷积码，运算量小但容错性能差；基于Boyer-Moore的快速合冲法运算量小且具备一定容错性，但该方法只适用于1/2码率卷积码的识别；基于WHT的识别算法同样适用于1/n码率卷积码，具有较好的容错性，但算法运行所需存储量巨大；基于校验序列的识别算法可对不同码率卷积码进行识别，但由于其基于对硬判决序列编码矩阵分析的原理，容错性一般；基于矩阵分析列相关性的识别算法，通过对编码矩阵的迭代分解，设置判决门限寻找满足列相关性的编码参数及校验矩阵，该方法同样只利用了解调硬判决序列，即便经历多次迭代后识别容错性仍有待提高；基于校验方程组成立概率的识别算法将解调软判决信息与方程成立概率结合，容错性高，但只适用于对（2，1，m）卷积码的识别。

综上所述，现有的卷积码识别方法主要利用解调硬判决信息（即解调过程中经判决得到的编码序列）进行识别，识别容错性有待提高。对此，本文在经典WHT算法的基础上，提出一种基于分段抽取软判决加权WHT（S-SWWHT）的卷积码识别算法。与经典WHT识别算法相比，该算法充分利用了解调软判决信息（即解调过程中软判决所利用的符号可靠度信息），并根据信噪比环境设定Walsh Hadamard矩阵加权系数，提高了识别容错性。同时通过分段抽取，降低了对运行存储量的要求，在识别得到基本校验序列之后，进一步识别得到基本生成矩阵。仿真结果表明，该算法可在低信噪比环境下实现对多种码率卷积码的有效识别，具有很好的应用价值。

1　卷积码识别问题描述

对于k/n码率卷积码，其基本校验多项式矩阵h（x）可表示如下：

式中：K为码字约束度。提取基本校验多项式矩阵行向量的系数，得基本校验序列为

将卷积码码字序列按一定规格构建编码矩阵C，其和基本校验序列存在如下校验关系：

而对于原卷积码的基本生成矩阵G，其与基本校验序列之间存在着如下校验关系：

至此，卷积码的识别问题可描述为：将接收序列按一定规格构建编码矩阵，求解（6）式的解向量，此时的解向量即为基本校验序列，基本生成矩阵进而可由（7）式原理求解得到。

2　基于S-SWWHT的卷积码识别

为了方便理解，首先介绍经典WHT识别算法，进而引入本文研究内容。

2.1经典WHT识别算法

由第1节内容，卷积码识别问题最终转化为有限域内求解线性方程组的问题。由于在非合作信号处理领域，信号接收常在低信噪比环境下进行，经解调后的编码序列将存在较多误码，（6）式进一步演化为对码字校验方程组的求解：式中：R为接收编码序列构建的N×l的矩阵（N为校验方程组中方程的个数，l=n（K+1）为基本校验序列长度）；h为待求的基本校验序列；E为接收序列错误图样。由此引入经典WHT识别算法，WHT通过对以R中行向量为地址的单元进行累加，得到频次序列再对 f进行WHT运算：

2.2求解软判决频次序列

由2.1节，经典WHT识别算法仅仅利用解调硬判决信息进行识别，在低信噪比环境下，解调后的数据误码率较高，该方法很难再准确识别。为此引入软判决WHT（SWHT）识别算法，对非合作方而言，增加了对接收信息的利用率。

在信号解调过程中，通常利用到每个比特的对数似然比（LLR）信息，借以完成输出符号的判决，即软判决解调。对大多数通信方或非合作方而言，解调和信道解码二者统一进行，因此可将解调信息存储并加以利用。由（9）式可知，在基于WHT的识别方法中，接收编码数据首先需转化为接收频次序列再进行处理。对此，本文利用接收软判决信息，将编码矩阵行向量的可靠度信息作为软判决频次计算接收数据的频次序列。

以加性高斯白噪声（AWGN）信道下的二进制相移键控（BPSK）调制为例，设高斯白噪声 w～N（0，σ2），当信道发送ci，输出为yi时，该比特的LLR可表示如下：

|LLR（yi）|在解调时可表示输出比特ri的可靠度信息，记为L（ri）.将单个比特的可靠度信息引入整个序列中，对于接收编码矩阵R，其行向量Ri的可靠度可表示为

然而，当直接选用LLR作为比特可靠度信息进行处理时，序列中正确码字的可靠度将远高于错误码字的可靠度，进而减弱了误码对序列可靠度信息的影响，不利于WHT运算的进行。对此，本文通过引入模糊隶属度的概念对可靠度信息作重新设定，将可靠度值隶属于一较小范围内。

对于BPSK调制，有：

则在信号层对信噪比进行估计之后[11]，根据不同信噪比环境下的序列中比特可靠度分布，选取一两边型拋物型F分布如下：

式中：st=E（|LLR（yi）|）=2/σ2.由此改进后的比特可靠度为L（ri）=A（LLR（yi）），与之对应，接收编码矩阵行向量Ri的可靠度为

对以R中行向量为地址的单元，将其可靠度进行累加，得到软频次序列

2.3构建加权Walsh Hadamard矩阵

式中：Φ（zα）=α（正态分布函数，可查表得到）；m= weight（ibin）.针对解向量抽头数不同导致代入校验方程组时方程成立概率不同的现象，本文采用加权WHT（WWHT）的方法，对抽头数较大的向量地址所对应的Walsh Hadamard矩阵列向量进行补偿。

对于解向量ibin，若其为正确解向量，则代入方程组（8）式中，方程成立与不成立个数差值的估计值为

若其为错误解向量，则代入方程组（8）式中，方程成立与不成立个数差值的估计值为

式中：λ∈[-1，1].根据不同重量的解向量代入含错线性方程组中的方程成立情况，对 Walsh Hadamard矩阵的列向量设置加权系数如下：

式中：k为Walsh Hadamard矩阵的列数。结合（17）式可对Walsh Hadamard矩阵列加权系数求解，最终完成加权Walsh Hadamard矩阵的构建，即

2.4卷积码的识别

利用2.2节的软频次序列f′以及2.3节得到的加权Walsh Hadamard矩阵，运行SWWHT算法得到Walsh谱如下：

再根据待识别卷积码码率的不同，分别设计1/n及（n-1）/n码率卷积码的基本生成矩阵识别算法。

对于常用的（n-1）/n码率卷积码，由（4）式知其基本校验序列唯一，且一般有码字约束度K≤7.由于待求解的线性方程组（8）式基于多项式方程组的原理构建，可假设一较大参数K′作为解向量的最高阶次，利用预估参数n，根据本文算法对基本校验序列h识别如下：

将其转化为基本校验多项式矩阵如下：严格的变换关系，因此可由h（x）直接变换得到[13]。对于（n，n-1，m）非系统卷积码基本生成多项式矩阵的识别原理如下：

解多项式方程组（25）式时将得到多解，这些解均是基本生成多项式矩阵行向量的线性组合，可通过对解向量进行初等变换使其化到最简，利用基本编码生成多项式矩阵的构造原理[14]，最终完成原卷积码基本生成多项式矩阵g（x）的识别。

对于1/n码率卷积码，可将接收编码数据分成n路，按顺序选取其中两路，按照1/2码率卷积码的方法进行识别，在得到n-1组参数之后，对于相邻两组参数中表示同一位置多项式的可通过对X的加权（X为x的函数），使得gt（x）=最终确保gcd（g1（x），g2（x），…，gn（x））=1，完成识别[4]。

2.5基于分段抽取降低运行存储

需要注意的是，在算法运行过程中，需构建2n（K+1）×2n（K+1）规格的Walsh Hadamard矩阵并调度运算，即至少需要s=4n（K+1）的运行内存，当接收编码参数为n=6，K=7时，存储量需求s达到1028bit量级，远远超过一般计算机的内存量级（109bit）[15]。对此，本文基于分段抽取的思想，继续对本文识别算法进行改进，降低运行存储。

设分段长度为d，对于2l长的软频次序列f′，以k2d+1，k=0，1，…，2l-d为起点取连续2d个元素，记为f′d，利用SWWHT算法求的Walsh谱，其最大值点的二进制序列即为（8）式解向量的低d位hd；再以k，k=0，1，…，2d为起点，每隔2d个元素取一个元素，组成2l-d长序列，记为，利用SWWHT算法求的 Walsh谱，其最大值点的二进制序列即为（8）式解的高l-d位hl-d[16].将hd与hl-d连接，最终得到接收编码序列的基本校验序列为h= [hdhl-d].对于分段后d、l-d仍较大的情况，可继续对hd与hl-d执行分段抽取运算，直至算法所需运行存储s小于计算机存储性能为止。再利用2.4节原理完成对卷积码基本生成多项式矩阵的识别。

2.6算法识别流程

综合以上识别步骤，本文算法的识别流程可表示如图1所示。

2.7算法运算量及运行存储分析

由于分段抽取是基于对解向量低、高位分别处理的原理，其与整体WHT运算所需运算量一致，因此本文所提出的基于S-SWWHT的卷积码识别算法

对于（n，n-1，m）系统卷积码其基本生成多项式矩阵g（x）与基本校验多项式矩阵h（x）之间存在相对于经典WHT算法，运算量增加在软频次序列的计算和Walsh Hadamard矩阵列向量加权系数的求取上。

图1　基于S-SWWHT的卷积码识别流程Fig.1　Flow chart of convolutional code recognition based on S-SWWHT

在软频次序列的求取时，需将比特LLR信息代入模糊隶属度函数A（x）并进行叠加运算，在低信噪比情况下，约需进行N·l次乘（除）法和2N·l次加（减）法运算；在Walsh Hadamard矩阵列向量加权系数的求解上，不同抽头数时的方程成立概率Pw，i可作为先验信息存储。由此，在计算矩阵加权系数时，需进行l次开方运算、2l+2（l-1）次乘（除）法和2l次加（减）法运算。

上述分析表明，基于S-SWWHT的卷积码识别算法相对于经典WHT算法在运算量上增加了l次开方运算、2l+（2N+2）l-2次乘（除）法和（2N+ 2）l次加（减）法运算。然而本文算法在运行存储需求上，较经典WHT算法有明显降低，利用分段抽取的原理可有效地将WHT算法的运行存储降低至计算机性能范围内，对于不同长度的基本校验序列，进行分段抽取后的运行存储量变化情况可如图2所示（例如，对编码参数为n=6，K=7的情况运行SWWHT算法，当做两次分段抽取之后，运行存储量s 由1029bit降低至108bit）.

图2　分段抽取后算法运行存储量变化Fig.2　Change in an storage space of the algorithm after subsection extraction

3　仿真实验与结果分析

由于（n，1，m）卷积码的识别可在（2，1，m）卷积码识别的基础上进行，本节针对常用的（n，n-1，m）卷积码进行仿真实验，其中在（2，1，m）卷积码的识别上与经典WHT识别算法[4]及基于校验方程成立概率的识别算法[9]作性能对比，在（n，n-1，m）卷积码（n≥3）的识别上与基于矩阵分析列相关性的识别算法[8]及基于校验序列的识别算法[6]作性能对比。

对生成的卷积码编码数据，仿真实验使用BPSK调制，信道模型为AWGN模型，实验中假定已准确估计码率k/n.下面通过设置仿真的不同信噪比情况，对本文算法进行有效性验证及性能分析。

3.1有效性验证

以常用的（2，1，5）和（3，2，4）非系统卷积码为例验证本文算法的有效性。仿真实验采用BPSK调制，信道模型为AWGN模型。

例1 （2，1，5）卷积码的基本生成多项式矩阵用8进制数表示为g（53 25），即

设定信噪比为1 dB，数据量为1 000的情况下，估计基本校验序列长度为14，采用S-SWWHT算法与传统WHT算法分别对（2，1，5）卷积码进行识别，结果如图3所示。

图3　S-SWWHT与传统WHT识别结果对比Fig.3　The recognition results of S-SWWHT and traditionalWHT algorithms

首先验证S-SWWHT算法的有效性，由图3可得，S-SWWHT算法识别得到最优解向量为2 459，将其表示为二进制：100110011011.另一归一化Walsh谱值较接近的解向量表示成二进制形式为10011001101100，将其高阶的0抽头去除后，可得二者为同一解向量。这是由于估计校验序列长度大于实际校验序列长度，从而导致存在两个最优解向量。系数抽取得到基本校验多项式矩阵如下：

（2，1，m）卷积码存在如下性质：h1，1（x）=g1，2（x），h1，2（x）=g1，1（x）.结合（27）式求得基本生成多项式矩阵，与原卷积码相符。

同时，通过观察可得WHT算法识别得到多个较优解向量，且其中含有错误解向量。因此，S-SWWHT算法的识别效果优于传统WHT算法，且在估计解向量长度较大时，S-SWWHT算法识别更为准确。

例2 （3，2，4）非系统卷积码的基本生成多项式矩阵用8进制数表示为：g（23 35 0；0 5 13），即

在信噪比为2.5 dB，数据量为5 000的情况下，估计基本校验序列长度为24.按照传统WHT算法，需构建规格为224×224的Walsh-Hadamard矩阵，这将超出一般计算机的运行存储性能。为此，将采用S-SWWHT算法对（3，2，4）卷积码进行识别。

图4展示了对（3，2，4）卷积码的识别结果，可得解向量的低12位表示成二进制形式为111101111111，高12位为001010101110.进而得到基本校验序列：111101111111001010101110，按码长3抽取后得到基本校验矩阵为 h2=[11110011 10110101 11111010]T，表示成多项式矩阵形式如下：

对（3，2，4）卷积码含错编码序列，根据（2）式，设基本生成多项式矩阵的最高阶数（编码记忆长度）为m，建立F2（x）上的多项式方程如下：

对多项式方程（30）式分析可得，其最高阶次为m+7，未知项系数有3（m+1）个，由此，将（30）式按阶次提取系数可构建F2上方程个数为m+8、未知数个数为3（m+1）的线性方程组。又基本生成多项式矩阵的行数即为式（30）最大线性无关解的个数，即解空间维数。对线性方程组维数进行分析可得如下关系：

解得卷积码编码记忆长度为m=4.求解由多项式方程（30）式系数提取构造的线性方程组，得到基础解向量为 p1=110011011100000，p2= 110011110101101.表示成多项式矩阵形式如下：

图4　S-SWWHT算法对（3，2，4）卷积码的识别结果Fig.4　Recognition results of（3，2，4）convolutional code based on S-SWWHT algorithm

对p（x）作行初等化简，使其行最高阶次最小（即使基本生成多项式矩阵结构最简），得到：

（33）式与原卷积码基本生成多项式矩阵相同，验证了本文算法在低信噪比环境下对高码率非系统卷积码盲识别的有效性。

3.2识别性能对比

下面通过蒙特卡洛实验，分两种情况，对不同卷积码识别算法的识别性能进行对比。

图5表示接收数据量为1 000时，对（2，1，5）卷积码的识别性能对比。由图5可以看出，在识别码字约束度较大的卷积码时，本文算法误码适应能力明显优于其他算法。设置准确识别的识别率门限为90%，则本文算法在识别性能上，较基于方程成立概率的识别算法提高约0.5 dB，较经典WHT算法提高约1.2 dB.

图5　（2，1，5）卷积码识别性能对比Fig.5　Performance comparison of the recognition of（2，1，5）convolutional code

图6表示接收数据量为5 000时，对（3，2，4）非系统卷积码的识别性能对比图。由于该码的编码约束度为24，从中可以看出，识别该码较识别（2，1，5）卷积码容错能力有所下降，然而，与其余高码率卷积码的识别算法相比，本文算法识别容错能力明显占优。设置准确识别的识别率门限为90%，则本文算法在识别性能上，较基于列相关性的识别算法提高约3 dB，较基于校验序列的识别算法提高约3.5 dB.

图6　（3，2，4）卷积码识别性能对比图Fig.6　Performance comparison of the recognition of（3，2，4）convolutional code

4　结论

本文提出了一种基于S-SWWHT的卷积码识别算法，该算法在预估码率的情况下首先构建Walsh Hadamard矩阵，通过引入软判决频次序列以及矩阵加权系数，实现Walsh谱的求解，进而完成对基本校验序列的识别。同时，结合分段抽取降低算法运行存储。最后利用基本编码矩阵构造规则，实现对卷积码记忆长度及基本生成多项式矩阵的识别。该算法充分利用每个比特的解调软判决信息，可在低信噪比环境下实现对不同码率卷积码的识别，具有较好的容错性，在非合作信号处理中的信息截获及合作通信中的智能通信等领域具有重要应用意义。算法不足之处在于，该算法对编码参数的盲识别部分仍有待进一步研究完善。

（References）

[1]解辉，黄知涛，王丰华.信道编码盲识别技术研究进展[J].电子学报，2013，41（6）：1166-1176. XIE Hui，HUANG Zhi-tao，WANG Feng-hua.Research progress of blind recognition of channel coding[J].Acta Electronica Sinica，2013，41（6）：1166-1176.（in Chinese）

[2]Wang FH，Huang Z T.A method for blind recognition of convolutional code based on euclidean algorithm[C]∥IEEE International Conference on Wireless Communications Networking and Mobile Computing.Shanghai，China：IEEE，2007：1414-1417.

[3]Lu P Z，Zou Y.Fast computationsofGröbner basesand blind recognitions of convolutional codes[J].Lecture Notes in Computer Science，2007，4547：303-317.

[4]刘健，王晓君，周希元.基于Walsh-Hadamard变换的卷积码盲识别[J].电子与信息学报，2010，32（4）：884-888. LIU Jian，WANG Xiao-jun，ZHOU Xi-yuan.Blind recognition of convolutional coding based on Walsh-Hadamard transform[J]. Journal of Electronic&Information Technology，2010，32（4）：884-888.（in Chinese）

[5]刘建成，杨晓静.基于求解校验序列的（n，1，m）卷积码盲识别[J].电子与信息学报，2012，34（10）：2363-2368. LIU Jian-cheng，YANG Xiao-jing.Blind recognition of（n，1，m）convolutional code based on solving check-sequence[J].Journal of Electronic&Information Technology，2010，34（10）：2363-2368.（in Chinese）

[6]杨晓静，刘建成，张玉.基于求解校验序列的（n，k，m）卷积码盲识别[J].宇航学报，2013，34（4）：568-573. YANG Xiao-jing，LIU Jian-cheng，ZHANG Yu.Blind reco-gnition of（n，k，m）convolutional codes based on solving check-sequence [J].Journal of Astronautics，2013，34（4）：568-573.（in Chinese）

[7]Marazin M，Gautier R，Burel G.Dual codemethod for blindidentification of convolutional encoder for cognitive radioreceiver design [C]∥GLOBECOM Workshops 2009.Honolulu，Hawaii，US：IEEE，2009：1-6.

[8]Marazin M，Gautier R，Burel G.Blind recovery of k/n rate convolutional encoders in a noisy environment[J].EURASIP Journal on Wireless Communications and Networking，2011，168（1）.1-9.

[9]于沛东，李静，彭华.一种利用软判决的信道编码识别新算法[J].电子学报，2013，41（2）：301-306. YU Pei-dong，LI Jing，PENG Hua.A novel algorithm forchannel coding recognition using soft-decision[J].Acta Electronica Sinica，2013，41（2）：301-306.（in Chinese）

[10]伍文君，黄芝平，唐贵林，等.含错扰码序列的快速恢复[J].兵工学报，2009，30（8）：1134-1138. WUWen-jun，HUANG Zhi-ping，TANG Gui-lin，etal.Fastrecovery of interfered scrambling codesequence[J].Acta Armamentarii，2009，30（8）：1134-1138.（in Chinese）

[11]许华，王爱粉，杨晓宇.常规数字通信信号信噪比估计综述[J].信号处理，2013，29（6）：723-733. XU Hua，WANG Ai-fen，YANG Xiao-yu.Survey of the SNR estimation of conventional digital communication signals[J].Journal of Signal Processing，2013，29（6）：723-733.（in Chinese）

[12]黄开枝，陈松.基于加权WHT的软判决序列快速估计算法[J].电子与信息学报，2013，35（1）：74-79. HUANG Kai-zhi，CHEN Song.A soft fastestimationmethod of PN sequence based on weightingWalsh-Hadamard transform（WHT）[J].Journal of Electronic&Information Technology，2013，35 （1）：74-79.（in Chinese）

[13]杨勇.快速盲识别（n，k，m）系统卷积码的解方程方法[J].计算机应用，2013，33（S2）：77-79. YANG Yong.Fast blind recognition algorithm of（n，k，m）system convolutional coding based solving equation[J].Journal of Computer Applications，2013，33（S2）：77-79.（in Chinese）

[14]Katsiotis A，Kalouptsidis N.On（n，n-1）punctured convolutional codes and their trellismodules[J].IEEE Transactions on Communications，2011，59（5）：1213-1217.

[15]Lukusa T，Ouahada K，Ndjionggue A R，etal.Experimental verication of frequencymappingswith Hadamard Tran-sform for power line communications channel[C]∥2014 18th IEEE International Symposium on Power Line Com-munications and Its Applications.Glasgow，UK：IEEE，2014：190-195.

[16]Agaian S，Sarukhanyan H，Astola J.Gold Hadamard transform matrices[C]∥Proceedings of the 2006 International TICSP Workshop on Spectral Methods and Multirate Signal Processing. Florence，Italy：TICSP，2006：107-112.

An Algorithm for Convolutional Codes Recognition Based on Sectionally Extracting Soft-decision Weighted Walsh Hadamard Transform

ZHANG Dai1，ZHANG Yu2，YANG Xiao-jing3，FAN Bin-bin2
（1.State Key Laboratory of Pulse Power Laser Technology，Electronic Engineering Institute，Hefei230037，China；2.504 Staff Room，Electronic Engineering Institute，Hefei230037，Anhui，China；3.402 Staff Room，Electronic Engineering Institute，Hefei230037，Anhui，China）

An algorithm for convolutional code recognition based on sectionally extracting soft-decision weighted Walsh Hadamard transform（WHT）is proposed for the issue of convolutional code recognition in the low signal-to-noise ratio（SNR）environment.The proposed algorithm is used to obtain a soft-decision frequency sequence by using received bits'soft-decision information，and construct a weighted Walsh Hadamard matrix to recognize the basic check sequence.The code memory length and basic generated matrix are recognized based on the construction rule of basic encodingmatrixt.The algorithm has the reduced run storage space.Simulation results show that the algorithm could recognize the convolutional codeswith different coding rates in the low SNR environment，and has good fault-tolerant performance，which is improved greatly especially in recognition of long restraint-length convolutional codes.

information processing technology；convolutional codes recognition；low signal-to-noise ratio；soft-decision information；weighted Walsh Hadamard transform；section extraction

TN911.7

A

1000-1093（2015）12-2298-08

10.3969/j.issn.1000-1093.2015.12.012

2014-10-08

国家自然科学基金项目（61201379）；安徽省自然科学基金项目（1208085QF103）

张岱（1993—），男，博士研究生。E-mail：daishen1208@126.com；张玉（1962—），男，教授，硕士生导师。E-mail：zy_57325@126.com