基于自适应扫描的图像压缩方法

石翠萍, 张钧萍, 张　晔

(1. 哈尔滨工业大学电子与信息工程学院, 黑龙江 哈尔滨 150001;

2. 齐齐哈尔大学通信与电子工程学院, 黑龙江 齐齐哈尔 161000)



基于自适应扫描的图像压缩方法

石翠萍1,2, 张钧萍1, 张晔1

(1. 哈尔滨工业大学电子与信息工程学院, 黑龙江 哈尔滨 150001;

2. 齐齐哈尔大学通信与电子工程学院, 黑龙江 齐齐哈尔 161000)

摘要：现有的压缩方法通常采用固定的扫描方式。然而,每幅图像都有其自身的特点,采用固定扫描难以使其得到好的压缩性能。针对该问题,提出了一种基于自适应扫描的图像压缩方法。首先,对图像进行小波变换后,根据图像内容决定子带间扫描顺序。其次,对于子带内扫描,根据子带特性设计不同的扫描方法,以尽可能保留图像的纹理特征。最后,采用二叉树编码器对扫描后的系数序列进行编码。实验证明,在相同条件下,该方法明显优于主流的JPEG2000,且优于最新的基于固定扫描的自适应二叉树编码方法。

关键词：图像压缩; 自适应扫描; 图像内容; 二叉树

0引言

绝大多数现有的编码方法都是采用固定扫描方式扫描图像。如最常用的联合图像专家组(joint photographic experts group,JPEG)[1],采用zigzag方式扫描图像。对于嵌入式零树小波(embedded zerotree wavelet,EZW)[2]和基于等级树的集划分(set partitioning in hierarchical trees,SPIHT)[3]方法,均采用morton扫描方式对位平面进行扫描。JPEG2000[4]是一种非常高效的编码方式,然而,其采用的也是固定扫描,对每个位平面,每四行作为一段,对该段按照从上至下、从左至右的方式扫描。这些固定扫描方式,实现相对简单,但都没有考虑到图像的内容。对于图像来说,不同的图像包含不同的信息。如果没有针对图像内容进行分析,则忽视了图像数据自身的信息特性,这往往会造成图像的关键数据丢失,使重建图像的整体效果较差。近年来,一些基于图像内容的压缩方法被提出。文献[5]提出了一种预滤波与高阶空频上下文模型相结合的自适应块变换方法,有效地提高了编码效率。文献[6-7]均从人眼对不同内容的敏感度出发,来设计压缩算法。文献[8]对图像中的不同区域,采用不同的描述方法,以得到最小的图像失真。文献[9-10]分别从基于上下文的量化和基于视觉重要性的角度,对图像进行压缩。文献[11]提出了最新的基于内容的编码方法,能够对建立好的二叉树进行自适应扫描,并得到了较好的压缩效果。这些方法一般是从变换、量化、编码以及人眼视觉的角度,对图像内容进行分析。而且,对于基于上下文的模型,阶数越高,算法复杂度越高。

一些文献也从对图像扫描的角度对压缩算法进行了分析和设计。文献[12]指出了基于内容的预测技术,通常都是针对无损压缩的,且一般采用raster扫描。文献[13]提出了一种新的基于内容的视频编码方法,但和其他视频编码一样,采用的是DCT扫描,并将其与多描述编码相结合。从扫描方式上看,无论是raster扫描还是DCT扫描,从本质上都是一种“线”扫描方式。然而,小波变换后,对各高频子带,不重要的“块”一般比不重要的“线”更为常见。因此,对基于小波变换的压缩,要得到较高的编码性能,必须合理利用这种“块”的关系。

morton扫描可以较好地利用这种“块”的关系,一些经典的扫描方法,如EZW,SPIHT,都是基于morton扫描的。但该扫描方法仍是固定扫描,没有考虑到子带的特性。此外,Hilbert扫描也可以看作是一种“块”扫描,对给定的阶数,Hilbert扫描可以沿着特定的扫描顺序进行扫描。该方法具有较好的局部保持特性,能提供较高的编码性能,但同样地,这种方法也是采用固定的扫描方法。

本文从图像内容出发,结合小波变换后各子带特性不同这一特点,提出了一种新的基于内容自适应扫描的图像压缩方法。

1基于自适应扫描的图像压缩方法

不同的图像具有不同的内容。图像的内容可以从不同角度进行描述。例如,对于图像的特征,一般包括如颜色、纹理、形状,拓扑关系等。对于具体应用,图像内容又可定义为与该应用紧密相关的特征。如对于生物特征识别,掌纹、人脸、虹膜、指纹均可作为特征。本文针对图像压缩,从图像能量的角度来定义图像内容。原因为:不同的图像具有不同的方向信息,而方向信息的丰富程度可以从子带的能量分布体现出来。如果能够先扫描能量较大的子带,则在相同的码率下,那些信息更丰富的系数能够先扫描,这对提高图像编码性能是极为有利的。总之,图像内容的定义与应用紧密相关。对于图像压缩,将图像能量定义为图像内容,结合子带特性确定扫描方法,理论上应能提高编码性能,得到更好的重建图像质量。

根据上面的分析,本文提出了一种新的基于自适应扫描的图像压缩方法,算法整体框图如图1所示。其编码过程包括正交小波变换、基于内容的自适应扫描,以及对扫描后的系数序列进行二叉树编码等部分。下面对编码过程中的主要模块进行详细的阐述。

图1　提出方法的总体框图

1.1基于内容的自适应扫描

对于一幅小波变换图像,扫描和编码通常是紧密联系在一起的。先扫描到的系数可以先编码,也就是说,对于那些先被扫描的系数,对应的码流应在整个码流的前面。相应的,这部分码流可以先被解码并显示。因此,在相同比特率下,不同的扫描方式会生成不同质量的重建图像。这说明了一个重要的问题:对图像编码来说,扫描方式非常重要。如果能先扫描那些对图像重建有更大贡献的系数,那么重建图像的质量必然会得到改善,尤其是在低码率的情况下。

根据前面的分析,对小波变换图像来说,“块”扫描比“线”扫描更为适合。然而,对每个小波子带,应采用何种“块”扫描方式？此外,针对不同子带重构时,对重建图像的贡献也不同的问题,应如何设计子带间的扫描方式？这两个问题是本节需要解决的问题。

本文提出的基于内容自适应的扫描算法步骤如下:

设原始图像为X,

步骤 1用A表示小波变换,经L级小波变换后,小波图像可以表示为

这里,N表示小波子带总数。

步骤 2对于每个小波子带(XAT)(i)(i=1,…,N),用l表示该子带对应的小波层数(l=1,2,…,L),d表示子带方向(d=1,2,3,4)。其中,“1”表示最低频子带,“2”表示水平方向子带,“3”表示对角方向子带,“4”表示垂直方向子带。则该子带能量为

c(i,j)表示当前子带中,位于(i,j)点处的系数值。M和N分别表示当前子带的行数和列数。

步骤 3对所有子带,按能量El,d降序的顺序,确定子带间的扫描顺序。从数学角度,可表示为

式中，permu表示子带顺序的重排。

步骤 4确定各子带内的扫描顺序,方法如下

对于每个小波子带(XAT)(i)(i=1,…,N):

(1) 若该子带方向是“1”或“2”,则采用horizontal_z扫描方式;

(2) 若该子带方向是“4”,则采用vertical_z扫描方式;

(3) 若该子带方向是“3”,扫描方式需由该层方向为“2”和“4”的子带共同决定。

若El,2≥El,4,则该子带采用horizontal_z扫描方式。

若El,2

最后,根据上面确定的子带间扫描顺序和子带内扫描顺序,对变换图像XAT进行扫描,生成一维的系数序列。

可见,对于一幅图像,其扫描顺序和扫描方法是由图像内容决定的。首先,根据子带能量确定子带间扫描顺序,目的是使码流整体上能够按照子带重要性进行组织。其次,对于不同的子带,根据子带特性设计扫描方法,是为了尽可能保留图像的纹理信息。通过这种方式,使那些对重建图像贡献较大的系数尽量排在前面,以提高整体压缩性能。提出的horizontal_z扫描和vertical_z扫描如图2所示。

图2　子带内采用的扫描方法

对提出的基于内容的自适应扫描方法,下面给出一个具体例子。原始图像如图3(a)所示。假设小波分解层数为3层,小波变换后,先计算各子带的能量,结果见表1。根据表1,可以确定子带间的扫描顺序,分别为LL3,LH3,HL3,HH3,LH2,HL2,HH2,LH1,HL1,HH1。对于每个子带,扫描方式由子带特性决定。对LL3,HL3,HL2和HL1,采用horizontal_z扫描方式,这种扫描在保证“块”扫描的同时,扫描是向着水平方向移动的。同样,对LH3,LH2和LH1,采用vertical_z扫描方式,保证在给定码率下,尽可能多的保留垂直方向信息。对HH3,HH2和HH1,扫描方式由对应小波层的水平子带和垂直子带共同决定。由表1可知,对于每级子带,垂直信息均多于水平信息,故对这3个对角子带,均采用vertical_z扫描方式。整个扫描过程如图3(b)所示。扫描后,生成一维系数序列。为了比较不同的“块”扫描方法的有效性,图3(c)~图3(e)分别给出了morton扫描、Hilbert扫描和自适应扫描后,一维系数序列的能量分布图示。可以看出,相比于其他两种方法,提出的自适应扫描方法能够尽可能先扫描重要系数,这为后续的编码打下了良好的基础。

表1　小波变换图像中各子带对应的能量(×109)

图3　自适应扫描过程和几种扫描方法的结果

1.2二叉树编码

对图像压缩来说,目前采用嵌入式编码方法的居多,因为这种编码方法支持图像的渐进式传输,并能在码流的任意位置解码。

文献[11]提出了一种最新的基于二叉树的嵌入式编码方法,即自适应二叉树编码方法(binary tree coding adaptively,BTCA)。该方法不需要复杂的步骤,如基于上下文的模型、率失真最优化等,大大降低了算法复杂度。其基本过程如下:首先,根据扫描后的一维序列建立二叉树。然后,对该二叉树进行编码,编码是按照从下至上,从左至右的顺序进行。在编码过程中利用了重要系数的兄弟通常也是重要的这一思想,在一定程度上提高了编码性能。函数BTCA() 详细给出了二叉树编码的实现过程。

由于二叉树编码方法相对简单,且效率较高,故本文对变换图像进行基于内容的自适应扫描后,采用二叉树编码方法对生成的一维系数序列进行编码。

Functioncode=BTCA(Tk)

初始化设整个二叉树的树高为D,令d=D

While (d>1 )

{

Letct={ }.If Γ(i)≥Tk-1

若Γ(i)的兄弟为Γ(i+1),且Γ(i+1)不重要,则

ct=Bin_Tree_Enc(Γ,i+1,Tk);

否则

若Γ(i)的兄弟为Γ(i-1),且Γ(i-1)不重要,则

ct=Bin_Tree_Enc(Γ,i-1,Tk);

code={code,ct}Z;

d=d-1Z。

输出在阈值为Tk时,对应位平面输出的码流。

}

在上述BTCA算法中,对给定节点编码,采用了Bin_Tree_Enc算法。该算法过程如下：

Functioncode=Bin_Tree_Enc(Γ,i,Tk)

若Γ(i)在更大的阈值被编码,即Γ(i)≥Tk-1,则

若Γ(i)不在树的最底层,编码Γ(i)的两个孩子;否则编码Γ(i)的符号位。

若Γ(i)有重要父母,且Γ(i)的兄弟是不重要的,则

若Γ(i)不在树的最底层,编码Γ(i)的两个孩子;否则编码Γ(i)的符号位。

若Γ(i)≥Tk

若Γ(i)不在树的最底层,编码Γ(i)的两个孩子,并且在码流前加‘1’;否则编码Γ(i)的符号位,并在码流前加“1”

否则,输出‘0’。

输出以Γ(i)为树根的子树的码流。

1.3额外比特开销

对于提出的基于自适应扫描的压缩方法,为了使接收端能够正确解码,子带间扫描顺序以及对角子带的扫描顺序,应作为边信息,同压缩码流一同传输至接收端。

对J层小波分解,仅需3J+1个数来表示子带间的扫描顺序。在本文中,子带间扫描顺序是按LLJ, HLi, HHi, LHi(i=J,J-1,…,1)顺序来记录的。此外,对于子带内的扫描方式,仅需J个数来表示所有对角子带的扫描方式,记录方式是按照HHi(i=J,J-1,…,1)的顺序进行。这里,用“0”表示horizontal_z扫描,用“1”表示vertical_z扫描。因此,对一幅给定图像,采用提出的自适应扫描方法,总共需要的边信息为(3J+1)+J=4J+1个整数。

以第1.1节中的图3(a)为例,根据表1,子带间的扫描顺序为1,3,4,2,6,7,5,9,10,8,对角子带内的扫描顺序为1,1,1。因此,边信息的总花销为10+3=13个字节。由于原始图像“fishingboat”的大小为512×512,故边信息所占的比例为13/(512×512)≈0.004 9%。根据上面的分析,边信息的开销仅取决于小波分解层数J,与其他因素无关。也就是说,当图像更大,或图像的位深超过8bit时,边信息所占的比例更小。而且,如果对这些边信息采用熵编码,则边信息开销更是进一步减少。根据该分析,提出的基于内容的自适应扫描方法,用于传输边信息所需的开销是极少的,几乎可以忽略不计。

从数据传输上讲,边信息应放在最前面,同压缩码流一起传输至解码端。原因在于,本文提出的扫描方法是基于图像内容的,也就是说,对于不同的图像,编码时扫描的顺序也不同。为了保证解码的正确性,扫描顺序应作为边信息,先传输至解码端。由于边信息开销极小,通常最多为十几个字节,只要网络条件不是极差,差到连边信息都无法传输完毕,解码端就能根据该信息,“还原”出编码端的扫描过程,从而对压缩码流正确解码。此外,本文方法对网络要求极低,这是因为采用的二叉树编码方法为嵌入式编码,解码可以在任意位置停止。在边信息能传输完毕的前提下,即使遇到网络中断的情况,也能对已接收的压缩码流进行解码,生成重建图像。

2图像质量评估测度

绝大多数压缩方法均采用峰值信噪比(peak signal to noise ratio,PSNR)来评估重建图像的质量,然而,PSNR已被证实并不与人眼感知相一致[14]。因为实际中,有时具有较高PSNR的重建图像,其视觉效果并不好。为了更好地评估所提算法,本文除采用PSNR作为客观评估方法外,还采用基于人眼视觉系统(human vision system,HVS)特性的结构相似性度量方法(structural similarity index measure,SSIM)[15]作为客观评估方法,以综合评定重建图像的质量。

设X和Y分别表示原始图像和重建图像,则重建图像的PSNR可表示为

(1)

式中，b表示原始图像的位深。

SSIM可以用来衡量两幅图像的结构相似度，其可由下式计算：

(2)

式中,μX，μY分别表示X和Y的均值；σX，σY分别表示X和Y的标准差；σXY表示X和Y的协方差；c1,c2,c3,k1,k2均为常数。

3实验及结果分析

文献[11]提出了一种最新的基于二叉树的压缩方法,并将其与SPIHT、基于集合划分的嵌入式块编码(set partitioning embedded block,SPECK)相比,得出了在相同的条件下,提出的“morton扫描+BTCA方法”要优于两者的结论。故在本文实验中,只需将所提算法与文献[11]中效果最好的“morton扫描+BTCA方法”对比即可。为了更充分地验证本文算法,结合前面对Hilbert扫描的分析,本文还与“Hilbert扫描+BTCA方法”对比。此外,目前应用最广泛的、具有较高编码性能的JPEG2000方法也被用作对比方法。

实验中采用双正交9/7小波变换,小波分解层数均为5层。测试图像大小均为512×512。

为了证实本文方法的视觉效果,先采用Lena图像进行实验,实验中采用的码率为0.0313bpp,实验结果如图4所示。图4(b)~图4(e)分别是利用文献[11]的方法,基于Hilbert扫描的方法、本文方法以及JPEG2000方法压缩Lena图像,得到重建图像的局部放大图,对应图4(a)中上面白框的区域。同样的,图4(f)~图4(i)也是采用上述4种方法,得到的重建图像的局部放大图,对应的是图4(a)中下面白框的区域。可以看出,在视觉质量上,本文方法得到的结果更清晰,保留的纹理更多,这是由于提出的扫描方法是基于图像内容进行设计的,且在扫描过程中考虑到了子带特性。

为了充分验证本文算法,实验中采用一组常用的自然图像作为测试图像。在多个码率下,分别采用上述4种方法进行压缩,并用PSNR和SSIM作为评价重建图像的质量的标准。PSNR结果见表2和表3,SSIM结果见表4和表5。在表2~表5中,用Fishboat表示测试图像Fishingboat,用Cam-man表示测试图像Cameraman,用J2K表示JPEG2000。

图4　码率为0.0313bpp时,不同方法得到的重建图像局部放大图

dB

表3　在不同码率条件下,本文方法与其他方法的PSNR性能比较2　dB

对表2和表3进行分析,可以看出,对于所有的测试图像,在给定的比特率下,采用提出的基于自适应扫描的压缩方法得到的重建图像的PSNR,均高于其他3种基于扫描的方法。原因在于,提出的方法能够先扫描包含信息较多的系数,故在均方误差意义下,能得到更好的编码性能。

对表4和表5进行分析,可以看出,对于这些测试图像,几乎在所有给定码率下,采用提出的压缩方法得到的SSIM,均高于其他三种方法。极个别情况下,对特定的图像,如peppers,在码率为0.125bpp时,采用JPEG2000得到的重建图像的SSIM更高。值得注意的是,JPEG2000的高性能是以高复杂度为代价的,如率失真最优化,基于上下文的模型等,都极大的增加了JPEG2000的复杂度。尽管如此,在绝大多数情况下,本文方法得到的重建图像的SSIM仍然最高。从另一个角度,由于每幅图像都具有不同的特性,包括纹理、复杂度等,相比于采用一副图像来评估算法,更好的方法是采用多幅图像进行测试,然后取结果的平均值。根据表4和表5,在多个码率下,对多幅测试图像的SSIM结果取平均,可以看到,在所有给定码率下,本文方法仍然优于其他3种方法。

表4　在不同码率条件下,本文方法与其他方法的SSIM性能比较1

表5　在不同码率条件下,本文方法与其他方法的SSIM性能比较2

根据上述对实验结果的分析,可以得出结论,相比于其他3种基于扫描的方法,提出的基于内容自适应扫描的压缩方法,能够提供更好的编码性能,以及更好的视觉效果,证实了提出方法的有效性。

从应用价值上讲,本文提出的压缩方法,子带间扫描顺序以及对角子带内的扫描顺序均通过子带能量大小比较即可确定。对扫描后的系数序列进行二叉树编码时,主要的计算也是两个数的大小比较。因此,本文提出的压缩方法,主要通过简单的大小比较即可实现,不需要诸如最佳率失真截断、基于上下文模型等复杂的操作,故复杂度非常低,这使该算法极适合于硬件实现,并用于实际的压缩系统中。

从性能提高的角度,若在该算法的基础上,采用熵编码,则编码性能还会进一步提高。熵编码同样可以采用硬件实现。

从复杂度的角度,本文方法的复杂度远远低于JPEG2000的复杂度。与最新的基于固定扫描的BTCA方法相比,增加的子带间和子带内自适应扫描步骤,通过少数几个比较大小的步骤即可实现。在硬件实现中,几个数的比较时间是极短的,几乎可以忽略不计。

4结论

针对现有压缩算法通常采用固定扫描,没有考虑到图像自身特性的问题,在对图像内容进行分析的基础上,结合小波子带特性,提出一种新的基于内容自适应扫描的图像压缩方法。该方法将子带间扫描和子带内扫描分别考虑,并采用二叉树编码器对扫描后的系数序列进行编码。实验结果证明了本文方法明显优于最主流的JPEG2000,且优于目前最新的基于固定扫描的自适应二叉树编码方法。提出的扫描方法是无冗余的,与其他固定扫描方法相比,只是改变了系数扫描的顺序。当然,本文方法也有一定的不足,如子带内的扫描是沿着水平或垂直方向进行,对图像来讲,若能根据其内容的纹理特征,计算最佳扫描方向,并沿着该方向扫描子带,从理论上应能得到更好的性能。但该方法无疑会增加算法的复杂度。如何在提高算法编码性能的基础上,尽可能降低算法复杂度,是下一步要研究的重点。

参考文献：

[1] Wallace G K. The JPEG still picture compression standard[J].IEEETrans.onConsumerElectronics,1992,34(1),30-44.

[2] Shapiro J M. Embedded image coding using zerotrees of wavelet coefficients[J].IEEETrans.onSignalProcessing,1993,41(12),3445-3462.

[3] Said A, Pearlman W A. A new, fast, and efficient image codec based on set partitioning in hierarchical trees[J].IEEETrans.onCircuitsandSystemsforVideoTechnology,1996,6(3),243-250.

[4] JPEG2000imagecodingsystem[M].ISO/IEC Std,15444-1, 2000.

[5] Tu C J, Tran T D. Context-based entropy coding of block transform coefficients for image compression[J].IEEETrans.onImageProcessing, 2002, 11(11): 1271-1283.

[6] Fahmy G, Jr J B, Panchanathan S. Texture characterization for joint compression and classification based on human perception in the wavelet domain[J].IEEETrans.onImageProcessing, 2006, 15(6): 1389-1396.

[7] Liu D, Sun X Y, Wu F, et al. Image compression with edge-based inpainting[J].IEEETrans.onCircuitsandSystemsforVideoTechnology, 2007, 17(10): 1273-1287.

[8] Shirani S. Content-based multiple description image coding[J].IEEETrans.onMultimedia, 2006, 8(2): 411-419.

[9] Cagnazzo M, Antonini M, Barlaud M. Mutual information-based context quantization[J].SignalProcessing:ImageCommunication, 2010, 25(1): 64-74.

[10] Deng C W, Lin W S, Cai J F. Content-based image compression for arbitrary-resolution display devices[J].IEEETrans.onMultimedia, 2012, 14(4): 1127-1139.

[11] Huang K K, Dai D Q. A new on-board image codec based on binary tree with adaptive scanning order in scan-based mode[J].IEEETrans.onGeoscienceandRemoteSensing,2012,50(10)：3737-3750.

[12] Memon N, Neuhoff D L, Shende S. An analysis of some common scanning techniques for lossless image coding[J].IEEETrans.onImageProcessing, 2000, 9(11): 1837-1848.

[13] Sawant S, Adjeroh D A. Balanced multiple description coding for 3D DCT video[J].IEEETrans.onBroadcasting, 2011, 57(4): 765-776.

[14] Mathews J, Nair M S, Jo L. Modified BTC algorithm for gray scale images using max-min quantizer[C]∥Proc.oftheInternationalConferenceonAutomation,Computing,Communication,ControlandCompressedSensing, 2013: 377-382.

[15] Wang Z, Bovik A C, Sheikh H R, et al. Image quality assessment: from error visibility to structural similarity[J].IEEETrans.onImageProcessing, 2004, 13(4): 600-612.

石翠萍(1980-),女,讲师,博士研究生,主要研究方向为图像压缩、图像处理。

E-mail:scp1980@126.com

E-mail:zhangjp@hit.edu.cn

张晔(1960-),男,教授,博士研究生导师,主要研究方向为图像处理。

E-mail:zhye@hit.edu.cn

网络优先出版地址：http://www.cnki.net/kcms/detail/11.2422.TN.20150901.0928.002.html

Image compression method based on adaptive scanning

SHI Cui-ping1,2, ZHANG Jun-ping1, ZHANG Ye1

(1.SchoolofElectronicsandInformationEngineering,HarbinInstituteofTechnology,Harbin150001,China;

2.SchoolofCommunicationandElectronicEngineering,QiqiharUniversity,Qiqihar161000,China)

Abstract:Existing compression schemes usually scan an image by a fixed way. However, each image has its own characteristics, which makes it difficult to obtain a good coding performance with a fixed scanning manner. To solve this problem, an adaptive scanning-based compression method is proposed. Firstly, the wavelet transform is performed. Then the scanning order among the subbands is determined by the content of the image. Secondly, for the scanning of a subband, in order to preserving the texture information as much as possible, the scanning method is designed with the considering of the characteristic of the subband. Finally, the binary tree code is exploited to encode the one-dimensional coefficient sequence. Experiments show that, under the same conditions, the proposed compression method is superior to JPEG2000 and the state-of-the-art compression methods based on binary tree coding adaptively (BTCA).

Keywords:image compression; adaptive scanning; image content; binary tree

通讯作者张钧萍(1970-),，女,教授,博士研究生导师,主要研究方向为图像处理。

作者简介：

中图分类号：TP 751.1

文献标志码：A

DOI:10.3969/j.issn.1001-506X.2016.01.30

基金项目：国家自然科学基金(61271348)；黑龙江省高校科技成果产业化前期研发培育项目(1254CGZH04)；齐齐哈尔大学青年教师科研 启动项目(2010W-M31)资助课题

收稿日期：2014-10-24；修回日期：2015-06-04；网络优先出版日期：2015-09-01。