水-气界面CO2通量监测的静态箱法与薄边界层模型估算法比较*

姚 骁，李 哲，郭劲松，林初学，陈永柏，李 翀

(1：中国科学院水库水环境重点实验室，重庆 400714) (2：中国科学院重庆绿色智能技术研究院，重庆 400714) (3：重庆大学城市建设与环境工程学院，重庆 400045) (4：中国长江三峡集团公司，北京 100038)

水-气界面CO2通量监测的静态箱法与薄边界层模型估算法比较*

姚 骁3，李 哲1，2，4**，郭劲松1，2，林初学4，陈永柏4，李 翀4

(1：中国科学院水库水环境重点实验室，重庆 400714) (2：中国科学院重庆绿色智能技术研究院，重庆 400714) (3：重庆大学城市建设与环境工程学院，重庆 400045) (4：中国长江三峡集团公司，北京 100038)

模型估算法与静态箱法是水-气界面气体通量监测的主要方法，因原理不同监测结果通常存在一定差异.目前对引起上述差异的主要环境因素仍不清晰.本研究使用自行设计的静态箱对三峡支流澎溪河水-气界面CO2通量进行监测，并与同步开展的CO2通量薄边界层模型估算法结果相比较，探讨该水域引起这两种监测方法结果产生差异的主要环境因素.结果表明，瞬时风速、水汽温差及水深均会对静态箱法及模型估算法的监测结果产生影响.风速越强、水汽温差越大、水深越大，这两种方法监测结果的差异就越小；而水域面积对两种方法的差异没有影响.比较发现，两种方法所获通量数据呈显著正相关，但静态箱法所获通量数据离散性显著高于薄边界层模型估算法.从方法的稳定性角度，在峡谷河道型水库水体温室气体监测中薄边界层模型估算法可能更为适宜.

静态箱法；薄边界层模型估算法；CO2通量；环境因素

水-气界面气体交换是水生生态系统中碳、氮等生源要素同大气进行物质交换的重要途径.通常，对水-气界面间气体交换通量的监测方法有微气象法、模型估算法、箱法、遥感反演法等[1].其中，模型估算法和箱法因其简单、灵活、易操作等特点在野外现场监测中最为常用.在原理上，模型估算法根据空气和水体内气体浓度梯度差运用Fick定律来估算通量，气体在两介质中的浓度及其传质系数的确定是该方法的核心[2-5].箱法原是测量陆地生态系统温室气体排放通量应用最广泛的方法，1990s中期被引用到水-气界面间气体通量监测中[6].考虑到水-气界面箱法监测下垫面并不稳定，在监测中通常使箱体浮于水面，定长时间段内监测箱体内气体浓度变化梯度以获取水-气界面气体交换通量速率值.

采用模型估算法和箱法对天然水域水-气界面间气体交换进行监测的相关研究不胜枚举，但由于模型估算法与箱法自身差异，它们间的监测结果也存在显著区别.在实践中，Duchemin等将模型估算法与静态箱法通量监测结果进行比较，发现后者通量值普遍高于前者，甚至高出一个数量级，认为模型估算法可能会低估水-气界面实际通量，箱法所造成的人工微环境又可能使得监测结果偏高[7].Vachon等通过对箱内外表层水体扰动强度的监测，认为上述差异随风速的增加而降低，当风速较大时(>4m/s)，两种方法的监测结果趋近一致[8].但Soumis等认为静态箱内热传输(箱内气温与表层水体温度差异造成)将使白天静态箱法监测结果较模型估算法的偏高，而晚上则偏低[9].尽管如此，不同方法所适用的环境特征目前并不明晰；对于特定水域，如具有显著流动性的河道型水库，两种方法监测结果是否存在明确数学关系使其能够相互适配也没有明确回答.

三峡水库是典型的深水峡谷型水库，其水库水体温室气体效应近年来备受关注.同国际上(南美或加拿大)已经开展相关研究的水库不同，三峡水库的峡谷河道型特征使其保持了明显的流动性，水-气界面气体交换过程不仅受到风速、温度等的影响，且同所在水域的水深、水域面积等密切相关.在该水域实际运用模型估算法和箱法监测水-气界面温室气体(CO2、CH4等)交换通量，需充分考虑前述因素的影响，而目前仍鲜有报道.本文以三峡水库典型支流澎溪河为研究对象，采用前述两种监测方法，对其回水区不同典型区段的水-气界面CO2通量进行监测，结合同期相关环境参量(现场风速、温度)和水域特征(水深、水域面积)，对两种方法的监测结果及其干扰因素进行初步比较分析，以期为进一步完善该类型水库水-气界面温室气体监测方法提供前期研究积累.

1 材料与方法

1.1 模型估算法计算通量的原理

水-气界面CO2交换通量主要受以下几个因素的影响：1) 表层水体CO2的分压；2) 气体传质系数，而气体传质系数又受流速、风速、温度等因素影响.根据Fick定律，对于淡水水体，水-气界面CO2交换通量(正为释放、负为吸收)可由公式(1)计算得出[10]：

FCO2=kx(Cwater-Cair)

(1)

式中，FCO2为水-气界面CO2扩散通量(mmol/(m2·h))，kx为CO2传质系数(cm/h)，Cwater为CO2在水中的浓度(mmol/L)，Cair为现场温度及压力下CO2气体在水中的饱和浓度(mmol/L).

对于传质系数kx的估算，目前通用的2种模型分别是薄边界层模型(TBL)和表面更新模型(SRM).TBL模型是假定气体转移是由水表面的薄边界层控制的水-气界面的浓度与大气中气体浓度形成溶解平衡；而SRM模型是假定水面漩涡可取代水表面薄层，且取代速度决定于水的被搅动程度.目前世界范围内对kx系数的确定绝大多数采用的是1989年Jähne等建立的数学经验公式[11]：

(2)

式中，k600为六氟化硫(SF6)气体的交换系数(cm/h).对于湖泊、库生态系统，选用Cole等在1998年建立的经验公式[12]：

(3)

Sc为t℃下CO2的Schmidt常数，对淡水而言，可按下式进行计算[13]：

Sc(CO2)=1911.1-118.11t+3.4527t2-0.04132t3

(4)

U10为水面上方10m风速(m/s).通常现场监测所得的水体上方风速Uz可用下式进行换算[14]：

(5)

式中，z为测量风速时的高度(m)，Uz为z高度风速大小(m/s)，Cd10为10m时的阻力系数(取0.0013)，κ为Von Karman常数(取0.41).

对水中CO2浓度的确定有水化学平衡法和顶空平衡法.水化学平衡法是通过测定水样中DIC浓度或碱度、pH值及水温求得水中CO2浓度[15]；而顶空平衡法是在装有水样和初始惰性气体的密封玻璃瓶中，通过剧烈摇晃使水体内气体浓度与上方空气中气体浓度达到平衡，测量上方空气中气体的浓度值而得到平衡前水体待测气体的分压[16].水化学平衡法可在野外现场测得，避免了水样的保存带来的误差，且可用于历史数据序列的反演分析，故本研究采用水化学平衡法来计算水中CO2的浓度.

(6)

(7)

(8)

通常K1、K2可由公式(9、10)计算得出[17]，其中T为采样时的水温(℃)：

pK1=-6320.813/TK-19.569224lnTK+126.34048

(9)

pK1=-5143.692/TK-14.613358lnTK+90.18333

(10)

图1 模型估算法(水化学平衡)计算示意图Fig.1 The process of model calculated method(water chemical equilibrium)

上述模型实施技术路线见图1.对于kx的确定，现场采用YSI®ProODO溶解氧仪测定水温(精度0.1℃)，采用SMAR®AR-836风速仪测量水面上方2m的瞬时风速(U2)(精确度为0.1m/s)，并通过公式(5)换算成水面上方10m的风速(U10，m/s)；对于Cwater的确定，现场采用YSI®63测定水体pH值(精度为0.01)、采用HACH®的微量滴定器(精度为1.25μl)及标准硫酸溶液现场滴定碱度；水体中p(CO2)根据pH、碱度及亨利常数计算而得；而Cair则获取自静态箱法测试中的起始浓度值[18].现场气温、大气压来自手持式数字大气压计.

1.2 测定CO2气体通量的静态箱法

本文采用自行设计的改良型静态箱对水-气界面CO2通量进行监测.采样箱主体由塑料材质构成(430mm×330mm×140mm)，箱内四周贴有塑料泡沫，增加箱体的稳定性及增大浮力(图2).根据《淡水水库温室气体监测导则》建议，采样前先使采样箱内充满空气，然后倒置在水面上，采样前10s开启箱内风扇，使箱内空气混合均匀，确保取样稳定.每隔5min用针筒抽取采样箱内气体100ml，并注入铝箔采气袋后保存.每个采样点共取5次，整个采气过程为20min.集气袋取回实验室后48h内完成CO2浓度的分析.气体采用安捷伦7820型气相色谱仪.气体采用十通阀进样，CO2先经TDX-01柱分离后再通过甲烷转化炉转化，最后用FID检测器检验.通过测定气样浓度变化率计算水-气界面的CO2交换通量，计算公式为[19]：

(11)

式中，Flux为CO2气体扩散通量(mg/(m2·h))；Ls为时间-浓度关系图中的斜率(μatm/min或ppm/s)；P为监测时的大气环境压力(kPa)；F1为ppm到μg/m3的转换系数；F2为分钟到小时的转换系数，为60；V为浮箱内套入的空气体积(m3)；P0为标准大气压，为101.325 kPa；A为箱体所覆盖的水面面积(m2)；T为监测时的箱内温度(℃).计算后使用摩尔质量换算成CO2的摩尔浓度通量，监测中气体浓度变化的线性相关性系数R2均大于0.9才能满足要求.2011年3月高阳平湖水-气界面CO2通量静态箱现场取样测试及拟合图见图3.

图2 水-气界面悬浮静态箱Fig.2 Floating chamber applied in the study

图3 高阳平湖采样点静态箱取样测试与拟合结果Fig.3 Measurement and fitting results from static-chamber method in Lake Gaoyang

1.3 研究区域

澎溪河(亦称：小江)，地处四川盆地东部边缘，流域面积5173km2，干流全长182.4km，是三峡库区北岸中段流域面积最大的支流(图4a).三峡水库畜水后，澎溪河流域包含了峡谷、消落带区以及湖库等各种特征的水域，其在三峡水库具有代表性.2010年5月至2011年5月对该流域上游开县至下游河口设置了8个水-气界面温室气体监测断面(图4b)，每月10日对上述区域进行样品采集及实验分析工作，各断面采样点位于河道深弘线处.采样中确保在同一个采样点上模型估算法和静态箱法同步进行，使它们之间通量数据配对可比.采样期间，同步采集相关环境参量指标.

图4 澎溪河流域水系(a)及采样点(b)示意图 Fig.4 Sketch map of backwater area(a) and sampling spot(b) in the Pengxi River

1.4 采样点所代表的水域面积的确定

为明确采样点，根据前期研究成果[20]，本研究采用环境因素控制法来对澎溪河各采样点所能代表的水域面积进行划分.该方法认为水域生态系统中温室气体产汇过程受水质理化特征影响显著，不同采样点间水质理化特征的空间差异是导致采样点间温室气体通量存在差别的原因.故可采用超标倍数赋权的方法，将不同采样点间环境因素权值归一化，确定不同采样点水域温室气体通量水平所占权重，进而计算每个采样点所代表的水域面积.初步研究发现，澎溪河流域CO2通量的主要影响因素是pH、Chl.a[21-22]，故以它们的均值为标准，计算采样点的水面面积控制权重.

1.5 数据处理方案

化学实验分析及计算所得数据全部录入SPSS®或Origin®进行统计分析，采用非参检验K个独立样本法来辨别不同数据组的差异性，采用变异系数CV值(标准差同均值的比值)表征数据变化幅度，用Spearman相关性分析说明数据间变化的线性相关性.

为了表征使用基于水化学平衡的TBL估算法监测所得的CO2通量与对应的使用静态箱法监测所得的CO2通量之间的差异性，引入比值α，具体计算式为：

α=F静态箱/FTBL

(12)

式中，F静态箱为静态箱法监测的CO2通量(以下简称“静态箱通量”)(mmol/(m2·h))，FTBL为TBL估算法监测的CO2通量(以下简称“TBL通量”)(mmol/(m2·h)).

2 结果与分析

2.1 TBL估算法与静态箱法间的总体差别

静态箱通量同TBL通量具有显著的正相关关系，Spearman相关性系数为0.783(P≤0.01)，所获数据的线性拟合结果见图5，两种方法所获通量结果具有较好的一致性.研究期间，TBL通量均值为0.59±0.85mmol/(m2·h)，静态箱通量均值为1.70±2.54mmol/(m2·h).总体上，静态箱通量高于TBL通量，静态箱通量与TBL通量比值α的范围为0.24～23.96，均值为4.02±0.33，中位值为2.97，频次分布见图6.可以看到α分布较集中，主要分布在1～4之间，当α大于8之后只有个别特例出现.

为分析两种方法所获通量数据的离散性，以TBL通量作为参考值，按其通量分为<0mmol/(m2·h)、0～1mmol/(m2·h)和>1mmol/(m2·h) 3个区间，分别计算上述不同区间内静态箱通量、TBL通量数据序列的变异系数(CV值)和它们CV值之间的比值(图7).分析发现，随着TBL通量增加，TBL通量CV值呈现先增加后减小的特征.TBL通量的CV值在0～1mmol/(m2·h)最大；在>1mmol/(m2·h)的区间范围内最小.静态箱通量CV值变化亦呈现出先增加后减少的特征，静态箱通量CV值在0～1mmol/(m2·h)最大；在<0mmol/(m2·h)的区间范围内最小.静态箱通量数据序列离散性在所有区间内均显著大于TBL通量数据序列(二者CV值的比值大于1，图7)，二者差异随着TBL通量增加而增大.

图5 TBL估算法与静态箱法的通量比较Fig.5 Comparison of flux between TBL and static chamber

图6 α值频次分布图Fig.6 Histogram of α

2.2 主要环境要素对两种方法监测结果的影响

图7 TBL估算法与静态箱法的CV值以及CV值之比Fig.7 The radio of TBL and static chamber’s CV

经过换算之后的10 m风速(U10)范围在0～5m/s，将α值按U10风速的大小分为<1m/s、1～3m/s和>3m/s 3组，通过非参数检验K个独立样本法进行显著性检验(图8)，发现风速对α值影响显著(P<0.05)，随着风速增大，α值呈显著下降的趋势.

图8 风速及温差对两种通量监测方法的影响Fig.8 Influence of wind speed and water temperature on observation of flux

当U10<1m/s时，α值为4.82±2.95，变化范围为0.24～20.61，α值的分布较为分散；U10在1～3m/s时，α值为4.17±2.12，最大不超过13.2；而当U10>3m/s时，α平均值为2.32±0.90，变化范围为0.68～3.73，分布也较为集中.

研究期间，现场表层水汽温差变化范围为-13.73～5.27℃，将其按<-2℃、-2～0℃和>0℃分为3组，检验结果表明水汽温差会对α造成较显著的影响(P<0.05)，α值随着水汽温差的升高而变小，其值从4.91倍逐渐降低到3.25倍，并且α值的分布趋于集中.

将澎溪河流域各采样点现场水深按<10m、10～30m、>30m分为3组，分析结果表明，水深会对α值造成显著的影响(P<0.01)，α值随水深的增加而呈减小的趋势.当水深<10m时，α均值为4.93±5.39，变化范围为0.46～20.61，同前述结果相似，该区间内α值分布亦极不均匀.当水深在10～30m之间时，α值有所减小，但均值仍为4.36±1.92，变化范围为1.44～9.50.当水深>30m时，α值明显变小且更集中，数据范围仅为0.54～4.96，且大部分α值在1～3倍之间，均值为2.71±1.17(图9a).以各采样点所代表的水域面积为基础，将总体样本划分为<500m2、500～1000m2及>1000m23组.分析结果表明，不同水域面积内α值的差异性并不显著(P>0.05，图9b).

图9 水深(a)及水域面积(b)对α值的影响Fig.9 Influence of water depth(a) and water area(b) on observation of flux

2.3 讨论

采用静态箱法和TBL估算法对水-气界面温室气体通量开展监测是目前水体温室气体监测研究普遍采用的监测方法.但目前尚无统一、明确的标准方法对水-气界面静态箱法、TBL估算法分别进行标定与误差分析，对它们监测结果间进行数据弥合还有很多困难.Soumis等认为白天水汽温差较小时静态箱法会高估CO2的排放量，而夜间水汽温差较大时静态箱法会低估CO2的排放量[9].Vachon等认为风速造成的水表扰动减弱了由于静态箱所产生的人为干扰的影响，风速增加，静态箱法所获通量数据可能更接近于真实值[8].Duchemin等指出在深水区，TBL估算法会高估风速的影响，而在浅水区当风速较低时，TBL估算法又会低估通量的大小[7].本研究所选择的峡谷河道型水库水域具有独特性.在相对狭窄的峡谷河段风场峡谷效应显著，表现为风力变化大，瞬时风速可能极大，易出现乱流涡旋风和升降气流.同时，峡谷河道型水库水域局部时期亦可能具有显著流动性，水域表层水体扰动对两种方法的影响明显.

从前述分析上可以看出，尽管两个方法所获通量数据呈显著正相关，但静态箱法所获通量数据离散性(变异系数)显著高于TBL估算法，而TBL估算法所获通量结果总体上相对稳定.且α分布主要集中在1～4，只有极个别值在8以上.

在所遴选的4个环境参量中，瞬时风速、水汽温差与采样点水深均对α值影响显著，总体趋势是随着上述环境参量增加，两种方法差异逐渐减少，离散性亦逐渐降低，说明在风速加大、水汽温差变大、水深增加的环境条件下，静态箱法和TBL估算法所获结果趋于一致.但在风速相对较小、水汽温差相对较低的环境条件下，α值离散性较大，不确定度较高，两种方法可比性相对较差.依据前述观点[7-9]，其原因可能同静态箱自身密闭效应对通量监测产生干扰有关，该干扰容易减弱气体的空间变化，造成很多不确定因素，包括扰动水体改变CO2浓度梯度、空气压力梯度、物质流动等；而TBL估算法本身并不受这些干扰的影响.从方法的稳定性角度来说，在峡谷河道型水库水体温室气体监测中，TBL估算法可能更为适宜，图7中TBL估算法与静态箱法各自的CV值趋势以及它们CV值之比的趋势可以支撑这个观点.但由于TBL估算法自身也受到监测参量的影响，一些参量监测(如pH值、风速等)对TBL估算法结果影响显著，特别是当pH小于7.5时，就不适合用TBL估算法来计算CO2通量了，上述影响不可忽略.

[1] 赵 炎，曾 源，吴炳方等.水库水气界面温室气体通量监测方法综述.水科学进展，2011，22(1)：135-146.

[2] Urabe J， Iwata T， Yagami Yetal. Within-lake and watershed determinants of carbon dioxide in surfacewater： A comparative analysis of a variety of lakes in the Japanese Islands.LimnologyandOceanography， 2011， 56(1)：49-60.

[3] Jonsson A， Åberg J， Jonsson M. Variations in pCO2during summer in the surface water of an unproductive lake in northern Sweden.TellusB， 2007， 59(5)：797-803.

[4] Gelbrecht J， Faitm， Dittrich Metal. Use of GC and equilibrium calculations of CO2saturation index to indicate whether freshwater bodier in north-eastern Germany are net sources or sinks for atmospheric CO2.Fresenius’JournalofAnalyticalChemistry， 1998， 361(1)：47-53.

[5] Huotari J， Ojala A， Peltomaa Eetal. Temporal variations in surface water CO2concentration in a boreal humic lake based on high-frequency measurements.BoraelEnvironmentResearch， 2009， 14：48-60.

[6] Anthony WH， Hutchinson GL， Livingston GP. Chamber measurement of soil-atmosphere gas exchange： linear vs. diffusion-based flux models.SolidScienceSocietyofAmericaJournal， 1995， 59(5)：1308-1310.

[7] Duchemin E， Lucotte M， Canuel R. Comparison of static chamber and thin boundary layer equation methods for measuring greenhouse gas emissions from large water bodies.EnvironmentalScienceandTechnology， 1999， 33(2)：350-357.

[8] Vachon D， Prairie YT， Cole JJ. The relationship between near-surface turbulence and gas transfer velocity in freshwater systems and its implications for floating chamber measuiements of gas exchange.LimnologyandOceanography， 2010， 55(4)：1723-1732.

[9] Soumis N， Canuel R， Lucotte M. Evaluation of two current approaches for the measurement of carbon dioxide diffusive fluxes from lentic ecosystems.EnvironmentalScienceandTechnology， 2008， 42(8)：2964-2969.

[10] Trfmbly A， Varflvy I， Rochm Cetal. Greenhouse gas emissions fluxes and processes hydroelectric reservoirs and natural environments. New York： Springer， 2005：725-732.

[11] Jähne B， Libner P， Fischer Retal. Investigating the transfer processes across the free aqueous viscous boundary layer by the controlled flux method.TellusB， 1989， 41(2)：177-195.

[12] Cole JJ， Caraco NF. Atmospheric exchange of carbon dioxide in a low-wind oligotrophic lake measured by addition of SF6.LimnologyandOceanography， 1998， 43(4)：647-656.

[13] Roehm CL， Prairie YT， Giorgio PA. The pCO2dynamics in lakes in the boreal region of northern Québec， Canada.GlobalBiogeochemCycles， 2009， 23：GB3013.

[14] Sibgh S， Bhatti TS， Kothari DP. Wind power estimation using artificial neural network.JournalofEnergyEngineering， 2007， 133(1)：46-52.

[15] Worrall F， Lancaster A. The release of CO2from riverwaters-the contribution of excess CO2from groundwater.Biogeochemistry， 2005， 76：299-317.

[16] Sobek S， Tranvik LJ， Cole JJ. Temperature independence of carbon dioxide supersaturation in global lakes.GlobalBiogeochemCycles， 2005， 19：GB2003.

[17] Solomon S， Qin D， Manning Metal. Historical overview of climate change science. Cambridge， UK and USA： Intergovernmental Panel on Climate Change(IPCC)， 2007：95-127.

[18] Matthe CJD， Stlouisv L， Hesslein RH. Comparison of three techniques used to measure diffusive gas exchange from sheltered aquatic surface.EnvironmentalScienceandTechnology， 2003， 37(4)：772-780.

[19] Zhu RB， Liu YS， Xu Hetal. Carbon dioxide and methane fluxes in the littoral zones of two lakes， east Antarctica.AtmosphericEnvironment， 2010， 44(3)：204-311.

[20] 李 哲，白 镭，蒋 滔等.三峡澎溪河水域CO2与CH4年总通量估算.水科学进展，2013，24(4)：551-559.

[21] 郭劲松，蒋 滔，李 哲等.三峡水库澎溪河春季水华期p(CO2)及影响因素分析.水科学进展，2011，22(6)：829-838.

[22] 李 哲，白 镭，郭劲松等.三峡水库两条支流水-气界面CO2、CH4通量比较初探.环境科学，2013，34(3)：1008-1016.

Comparison between closed static chamber method and thin boundary layer method on monitoring air-water CO2diffusion flux

YAO Xiao3， LI Zhe1，2，4， GUO Jinsong1，2， LIN Chuxue4， CHEN Yongbai4& LI Chong4

(1：KeyLaboratoryofReservoirAquaticEnvironment，ChineseAcademyofSciences，Chongqing400714，P.R.China)(2：ChongqingInstituteofGreenandIntelligenceTechnology，ChineseAcademyofSciences，Chongqing400714，P.R.China)(3：CollegeofUrbanConstructionandEnvironmentalEngineering，ChongqingUniversity，Chongqing400045，P.R.China)(4：ChinaThreeGorgesCorporation，Beijing100038，P.R.China)

Closed static chamber method and thin boundary layer method are frequently applied in monitoring air-water gas diffusive fluxes. However， there are differences in results between the methods due to theoretical principles which are still unclear. A type of closed static chamber was applied in Pengxi River， Three Gorges Reservoir to monitor air-water CO2flux and compared with the results from thin boundary layer method that simultaneously monitored to discuss the causes for the differences. Results indicated that wind speed， air-water temperature difference， and depth of the water were the impact factors related to the difference between the methods above. It was found that the smaller difference between the methods above was created by the higher wind speed， or the greater air-water temperature difference， or the greater water depth. However， water surface area showed no significant impact on the difference between the methods. Through a comparison of CO2fluxes data sets from the methods， it could be noted that there was a significant positive correlation on the monitoring results of the CO2fluxes.But monitoring results from the closed static chamber method had higher variance compared to thin boundary layer method.From the viewpoint of the stability of the methods， the thin boundary layer method might be more feasible in the river-channel based reservoir area.

Closed static chamber method； thin boundary layer method； CO2flux； environmental factor

*国家自然科学基金项目(51179125)、中国科学院西部行动计划项目(KZCX2-XB3-14)和重庆市自然科学基金重点项目(CSTC2012JJB20004)联合资助.2014-03-26收稿；2014-07-26收修改稿.姚骁(1990～)，男，硕士研究生；E-mail：yaoxiao1990@163.com.

**通信作者；E-mail：zheli81@sina.com.