環渤海低平原水土鹽分與水位埋深的空間變異及協同克立格估值

周在明,張光輝,王金哲,嚴明疆

中國地質科學院水文地質環境地質研究所,河北石家莊 050803

環渤海低平原水土鹽分與水位埋深的空間變異及協同克立格估值

周在明,張光輝,王金哲,嚴明疆

中國地質科學院水文地質環境地質研究所,河北石家莊 050803

針對環渤海低平原水土鹽分和水位埋深,通過選取0～20cm深度內的127個代表性土樣進行土壤全鹽量測定,并對130個水井的水位埋深及128個地下水樣的礦化度進行測定。綜合運用普通克立格(Ordinary Kriging)與協同克立格(CoKriging)方法與GIS技術研究水土鹽分及水位埋深的空間分布狀況。結果表明,土壤鹽分和地下水礦化度的變異系數分別為0.91和 0.87,屬于中等空間變異強度,而水位埋深的變異系數為1.06,屬于強變異強度。單變量空間相關程度均屬于中等,交互變量的空間相關性較強,其空間自相距離(135.6～223.2km)高于單變量(40.6～144.8km)。Kriging與CoKriging分析表明,環渤海低平原土壤鹽分與地下水礦化度分布規律基本一致,而與水位埋深相反。即土壤鹽分和地下水礦化度自內陸平原向東部濱海平原逐漸增加,不同地面高程上存在分布差異,并隨著高程的降低呈增加的趨勢。鹽分含量以唐山—天津—滄州—東營—濱州一線較高,礦化度以天津、滄州部分地區居高,水位在黃河三角洲較淺。通過CoKriging插值的均方根誤差與Kriging相比減少了0.29%～4.78%,而預測值與實測值的相關系數提高了2.03%～20.58%。

環渤海低平原;土壤鹽分;地下水礦化度;水位埋深;協同克立格

水鹽和水位埋深對地下水的利用程度及土壤鹽分有一定的影響(Kahlown et al.,2005),土壤鹽分、鹽漬化對農業可持續生產造成重大威脅(Cetin et al.,2003),土地鹽化程度的增加不但減少了作物產量而且限制了作物種類的選擇(Singh et al.,2009),干旱半干旱地下水淺埋區、濱海低平原區尤為顯著(Il’ichev et al.,2008;Guswa,2002)。對地下水和土壤鹽分含量的準確估計和分布研究在區域和全球水土環境及相關研究中是一項重要的課題(Liu et al.,2005)。對于水土而言其可溶鹽的分布在一定程度上受到地形的影響(Salama et al.,1999),地形特征對地下水及土壤形成與利用具有重要作用(Lin et al.,2007;Wu et al.,2008)。因此在地下水和土壤鹽分含量的插值估算中區域地形因素的影響不容忽視。

協同克立格(CoKriging)插值方法充分利用了待求變量和易測變量的相關性(Odeh et al.,1995),借助輔助變量提高主變量的預測精度,在地下水和土壤鹽分空間插值中得到廣泛的應用,是公認的有效插值方法(Triantafilis et al.,2001;Darwish et al.,2007)。Kim等應用水化學、Kriging與CoKriging方法對首爾以南的五松地區地下水硝酸鹽污染進行了評價,闡明了硝酸鹽的空間分布及與其它物理化學參數的關系,并解釋了 CoKriging在復雜地下水化學數據插值中的優越性(Kim et al.,2009)。Eldeiry等綜合應用地面驗證和遙感數據對美國阿肯色流域土壤鹽分含量的Kriging和CoKriging插值方法進行了比較,認為CoKriging方法對玉米地的土壤鹽分估計精度最高,其次是小麥和苜蓿(Eldeiry et al.,2010)。José等應用CoKriging方法對西班牙東南部塞奎拉河域土地,以土壤水電導率為輔助變量進行了 50cm深度的鹽分預測研究,較好地反應了作物根系層的鹽漬化狀況(José et al.,2010)。在國內相關的研究以西北干旱內陸和黃河三角洲為主,結果表明與普通克立格相比 CoKriging既能夠減少樣點數量又能提高估值精度(趙成義等,2003;姚榮江等,2006)。

針對整個環渤海低平原水土鹽分問題,進行大尺度區域布點,結合地形因素對水土鹽分的空間分布格局進行研究的尚不多見。為此,本文以環渤海低平原為研究對象,以研究區高程為輔助變量應用協同克立格插值方法對水土鹽分及水位埋深的空間分布進行研究,為提高水土鹽分含量的插值估算精度,為地下水和土地資源管理利用和農業生產的合理布局提供科學依據。

1 材料與方法

環渤海低平原位于華北東部,包括河北、山東內陸低平原和濱海低平原兩大類型,含天津、滄州等80多個縣市。由黃河、海河、灤河沖積而成,地勢低平,大部分海拔50 m以下,濱海區10 m左右?？偯娣e9.6×104km2;總耕地面積6000多萬畝,是我國重要的糧棉、果蔬產區。

環渤海低平原屬于歐亞大陸東岸暖溫帶半干旱季風氣候區,冬春寒冷干燥,夏季炎熱多雨。多年平均氣溫12.2℃,平均年降水量500～600mm,全區年內降水分配不均,主要集中在 6—9月,占全年降水量的60%～80%,年均蒸發量900～1400mm,干燥度達1.5左右。

1.1 數據來源

根據該區的土地利用類型和地下水狀況,全區共測量水井130個,布設水樣點128個、土樣點127個(圖1),取樣點間距約 20km。取樣范圍為(N36°03′～N39°35′;E114°36′～E119°28′),取 樣 面 積8.97×104km2。取樣時間為2010年4—5月,每個點采用GPS定位,土樣深度為0～20cm。

水位埋深用皮尺與測繩測定,水樣礦化度采用重量法測定。土樣全鹽量測定過程如下,稱取過2mm篩的風干土試樣50～100 g,按土水比1:5配制浸出液。采用蒸干法,添加15%雙氧水溶液與2%碳酸鈉溶液進行測定,詳細操作方法參見《土壤農業化學分析方法》(魯如坤等,1999)。

以1:5萬地形圖(等高距20 m)與SRTM高程數據為基礎,用實測高程點進行高程驗證(圖2),最終確定高程點1120個,平均間距約6.7km,覆蓋研究區主要地貌單元。利用均方根誤差描述高程精度,Zi為高程數據真值,zi為計算值,n為誤差個數,計算得RMSE=6.52,說明所得到的高程空間插值結果滿足數據分析要求。

圖1 環渤海低平原取樣點分布圖Fig.1 Distribution of soil and groundwater samples in the low plain around the Bohai Sea

圖2 環渤海低平原高程分布圖Fig.2 Digital Elevation Model of the low plain around the Bohai Sea

1.2 普通克立格與協同克立格方法

在地統計學中,半方差函數的一些重要參數如塊金值、基臺值和變程等可以用來表示區域化變量在一定尺度上的空間變異和相關程度,它是研究土壤特性空間變異的關鍵,也是克立格插值的精度要素(Yang et al.,2008;姜勇等,2005)。在本征平穩假設下,半方差計算公式為(式1)

式中,r(h)為步長h的半方差函數,N(h)是間距為h的計算對數,Z(xi)和Z(xi+h)分別是區域化變量Z(xi)和Z(xi+h)在空間位置xi和xi+h處的實測值。

兩個隨機變量的協同區域化可以用交互半方差函數來表示(式2)

式中,rij(h)是兩個變量的交互半方差值,N(h)是具有相同間距h的變量Zi(x)和Zj(x)的離散點的數目。如果兩個變量是正相關的,那么變量Zi從xa到xa+h的增加或減少會引起Zj的增加或減少,交互半方差就是正值。

如果變異函數和相關分析的結果表明某一屬性的空間相關性存在,則可以利用普通克立格進行插值(式 3),

式中,z*(x0)是待估點處的估計值,z(xi)是實測值,λi是分配給每個實測值的權重且∑λi=1。n是參與點估值的實測值的數目。

協同克立格是普通克立格的擴展形式,它要用到兩個或兩個以上的變量,其中一個是主變量,其它的作為輔助變量,將主變量的自相關性和主輔變量的交互相關性結合起來用于無偏最優估值中。其公式為(式4),

式中,z*(x0)是待估點x0處的估計值,z1(xi)和z2(xj)分別是主變量z1和輔助變量z2的實測值,λi和λj分別是分配給主變量z1和輔助變量z2的實測值的權重,且。n和p是參與x0點估值的主變量z1和輔助變量z2的實測值數目。

1.3 數據處理

采用SPSS 16.0軟件進行水土全鹽量及水位埋深的統計分析,變異函數及協同變異函數模型的擬合采用地學統計軟件 GS+7.0,克立格與協同克立格空間插值應用ArcGIS 9.3軟件。

2 結果分析

2.1 描述性統計特征

根據測定結果對地下水礦化度、水位埋深和土壤鹽分含量進行統計(表1)。環渤海低平原區耕層土壤鹽分含量均值為1.08 g/kg屬于輕度鹽漬化土(王遵親等,1993),地下水礦化度平均值為2.06 g/L屬于微咸水,水位埋深平均為16.61 m。由表1可見,環渤海低平原區的地下水埋深屬于強變異強度(變異系數為1.06),地下水礦化度與土壤鹽分的空間變異強度屬于中等,但也較大(變異系數分別為0.87和0.91)。這主要是由于研究區獨特的水文、地形結構以及復雜的人類活動所致。由于地下水過量開采出現區域性地下水位變化不均,內陸低平原淺層地下水位下降幅度較大,東部濱海平原區淺層地下水下降幅度較小,地下水超采是地下水位埋深空間變異加大的主要原因。從內陸到濱海地下水礦化度逐漸增大,水質變差,再加上人為排灌的影響,其共同作用使得平原區地下水礦化度空間變異相對較大。而研究區范圍較大且存在地勢差異、土地耕種方式差異、田塊灌溉制度差異等因素的共同作用導致耕層土壤鹽分變異性較大。

表1 環渤海低平原土壤鹽分、地下水礦化度和水位埋深統計特征值Table 1 Descriptive statistics of soil salinity,TDS and GWL in the low plain around the Bohai Sea

2.2 高程與觀測項目的相關性

根據相應取樣點的高程分別與土壤鹽分、地下水礦化度和水位埋深進行 person相關分析發現,其相關性分別滿足 0.05顯著性水平下的顯著負相關,相關系數?0.253;0.01顯著性水平下的極顯著負相關,相關系數?0.394;和 0.01顯著性水平下極顯著的正相關,相關系數為0.322。這說明土壤鹽分含量、地下水礦化度、水位埋深在不同高程上存在分布差異,并隨著高程的增大,土壤鹽分含量和地下水礦化度呈減少的趨勢,水位埋深則呈現增加的趨勢。在研究區內高程與土壤鹽分、地下水礦化度、水位埋深受區域化現象或空間過程的影響,屬于協同區域化變量,且應用偏度、峰度聯合計算法進行檢驗發現,通過對數轉換后 P＞0.05,因此區域高程屬于對數正態分布類型。

2.3 普通克立格與協同克立格插值

根據半方差函數理論,計算得到半方差函數擬合模型,其決定系數在0.58～0.96之間均達到顯著性水平(表2)?？梢?土壤鹽分、礦化度和水位埋深半方差和交互半方差函數可分別應用球狀模型和指數模型進行擬合。由塊金值與基臺值的比值C0/Sill可知,土壤鹽分、礦化度和水位埋深單變量空間相關程度均屬于中等,其值分別為41.78%、25.68%和40.00%,對交互變量而言,除土壤鹽分的空間相關性屬于中等外,礦化度和水位埋深的空間相關性均較強。此外在空間自相關距離上交互變量(135.6km、200.2km 和 223.2km)也高于單變量(40.6km、144.8km和90.0km)?？傮w上看,土壤鹽分、地下水礦化度和水位埋深的交互變量半方差函數理論模型較之單變量都有一定程度的改變。

分別根據半方差和交互半方差函數模型對研究區土壤鹽分、礦化度和水位埋深數據進行插值,獲得其空間分布圖(圖3)。本研究中 Kriging 和CoKriging是基于環渤海低平原區80多個縣市的實測數據做出的,在這種范疇廣、數據點較少的條件下,我們制作分布圖時,盡量通過參數的選擇使均方根標準誤(RMSSE)接近于 1,以保證預期誤差的變異盡可能地小,從而使得到的分布圖是我們現有數據資料的最佳和最優效的預期結果。

表2 土壤鹽分、地下水礦化度和水位埋深的半方差、交互半方差函數模型Table 2 Models of semivariogram and cross-semivariogram for soil salinity,TDS and GWL

由圖3可見,兩種插值方法所獲得的分布圖在整體趨勢和具體斑塊形狀上基本相似,表明協同克立格和普通克立格都能夠較好地反應研究區土壤鹽分、礦化度和水位埋深的空間分布情況。環渤海低平原土壤鹽分含量,自內陸平原向東部濱海平原逐漸增加,鹽分含量較高的地區出現在唐山—天津—滄州—東營—濱州一線。鹽分含量小于1 g/kg的非鹽化土以內陸平原為主,分布在保定—衡水—邢臺—邯鄲一線。這與地下水礦化度的分布規律基本吻合,地下水礦化度大于5 g/L咸水主要分布在東部濱海天津、滄州地區,1～5 g/L的微咸水分布于廣大內陸平原區,其鹽分含量也是逐漸向濱海地區增大,而水位埋深與之相反,在黃河三角洲較淺。

從地形因子看,土壤鹽分和地下水礦化度高的地區地勢平坦,海拔較低,水位較淺;與之相反的是海拔與地勢起伏相對較大的內陸區,土壤鹽分和地下水礦化度不大。從成土母質上看,濱海平原的土體一部分是由河流入海沖積成的三角洲如黃河三角洲,一部分是海積平原。質地多為粉砂、細砂組合,粗粉砂含量多達 60%以上。這種砂壤土體構型使得水鹽遷移量高于內陸平原的粘壤土和亞砂土的組合。從農業活動看,河北平原區的灌溉水源以抽取深層地下淡水為主,山東平原區以引黃灌溉為主,并建有完善的水利工程,有合理的耕作制度,因此農田長期處于脫鹽狀態,鹽分累積小。濱海平原區農田主要依靠自然降雨,長期的地面蒸發使得鹽化加重。

圖3 環渤海低平原土壤鹽分(a、a’)、礦化度(b、b’)和水位埋深(c、c’)空間分布(a、b、c 普通克立格插值;a’、b’、c’協同克立格插值)Fig.3 Spatial distribution of soil salinity (a,a’),TDS (b,b’)and GWL (c,c’)in the low plain around the Baohai Sea(a,b ,c-Ordinary Kriging;a’,b’,c’-Cokriging)

2.4 兩種插值方法預測精度比較

用評價方法(協同克立格CK)的均方根誤差相對于參考方法(普通克立格 OK)的均方根誤差減少的百分數(RRMSE)表示預測精度的提高程度(式6)。

用RR表示評價方法相對于參考方法相關系數的提高程度(式7)。

式中,RMSEOK和ROK分別是參考方法的預測均方根誤差及預測值與實測值間的相關系數,RMSEOK和RCK分別表示評價方法的預測均方根誤差及預測值與實測值間的相關系數。

由表3可見,土壤鹽分、礦化度和水位埋深通過協同克立格插值的均方根誤差與普通克立格插值相比分別減少了0.29%、2.18%和4.78%;而預測值與實測值的相關系數分別提高了 20.58%、2.03%和11.31%。表明在相同的取樣條件下,協同克立格插值由于融合了更豐富的空間信息,預測的精度要高于普通克立格插值,可優化樣點的插值精度。

表3 土壤鹽分、地下水礦化度和水位埋深的普通克立格與協同克立格精度比較Table 3 Comparison of prediction accuracy of soil salinity,TDS and GWL between Ordinary Kriging (OK)and CoKriging(CK)

3 結論

環渤海低平原區土壤鹽分和地下水礦化度屬于中等空間變異強度(Cv值分別為0.91和0.87)而水位埋深屬于強變異強度(Cv值為1.06)。半方差和交互半方差函數可分別應用球狀和指數模型進行擬合。單變量空間相關程度均屬于中等,交互變量的空間相關性較強,空間自相距離(135.6～223.2km)明顯高于單變量(40.6～144.8km)。

克立格與協同克立格空間插值表明,環渤海低平原土壤鹽分與地下水礦化度分布規律基本一致,而與水位埋深相反。即土壤鹽分和地下水礦化度自內陸平原向東部濱海平原逐漸增加,不同地面高程上存在分布差異,并隨著高程的降低呈增加的趨勢。鹽分含量以唐山—天津—滄州—東營—濱州一線較高,礦化度以天津、滄州部分地區居高,水位在黃河三角洲較淺。由于高程與上述觀測項目存在很好的相關性,通過協同克立格插值的均方根誤差與普通克立格插值相比減少了 0.29%～4.78%,而預測值與實測值的相關系數提高了2.03%～20.58%。

姜勇,莊秋麗,梁文舉,施春健,歐偉.2005.空間變異在土壤性質長期定位觀測及取樣中的應用[J].土壤通報,36(4):531-535.

魯如坤.1999.土壤農業化學分析方法[M].北京:中國農業科學技術出版社.

王遵親,祝壽泉,俞仁培.1993.中國鹽漬土[M].北京:科學出版社.

姚榮江,楊勁松,劉廣明.2006.土壤鹽分和含水量的空間變異性及其 CoKriging估值——以黃河三角洲地區典型地塊為例[J].水土保持學報,20(5):133-138.

趙成義,王玉潮,李子良,李國振.2003.田塊尺度下土壤水分和鹽分的空間變異性[J].干旱區研究,20(4):252-256.

CETIN M,KIRDA C.2003.Spatial and temporal changes of soil salinity in a cotton field irrigated with low-quality water[J].Journal of Hydrology,272:238-249.

DARWISH K H M,KOTB M M,ALI R.2007.Mapping soil salinity using collocated cokriging in Bayariya,Oasis,Egypt[C]//Proceedings of the 5th International Symposium on Spatial Data Quality,ITC Enschede,The Netherlands.

ELDEIRY A A,GARCIA L A.2010.Comparison of ordinary Kriging,regression Kriging,and Cokriging techniques to estimate soil salinity using Landsat images[J].Journal of Irrigation and Drainage Engineering,136(6):355-364.

GUSWA A J.2002.Models of soil moisture dynamics inecohydrology:a comparative study[J].Water Resources Research,38(9):1-15.

IL’ICHEV A T,TSYPKIN G G,PRITCHARD D,RICHARDSON C N.2008.Instability of the salinity profile during the evaporation of saline groundwater[J].The Journal of Fluid Mechanics,614:87-104.

JIANG Yong,ZHUANG Qiu-li,LIANG Wen-ju,SHI Chun-jian,OU Wei.2005.Application of sptial variability in long-term site-specific observatory study of soil properties and sampling strategy[J].Chinese Journal of Soil Science,36(4):531-535(in Chinese with English abstract).

JOSé M P,FERNANDO V,JOSé L R.2010.Spatial evaluation of soil salininty using the WET sensor in the irrigated area of Seuura river lowland[J].Journal of Plant Nutrition and Soil Science,doi:10.1002/jpln.200900221.

KAHLOWN M A,ASHRAF M,ZIAUL H.2005.Effect of shallow groundwater table on crop water requirements and crop yields[J].Agricultural Water Management,76:24-35.

KIM K H,YUN S T,CHOI B Y,CHAE G T,JOO Y S,KIM K,KIM H S.2009.Hydrochemical and multivariate statistical interpretations of spatial controls of nitrate concentrations in a shallow alluvial aquifer around oxbow lakes (Osong area,central Korea)[J].Journal of Contaminant Hydrology,107:114-127.

LIN Y S,LIN Y W,WANG Y,CHEN Y G,HSU M L,CHIANG S H,CHEN Z S.2007.Relationship between topography and spatial variations in groundwater and soil morphology within the Taoyuan-Hukou Tableland,Northwestern Taiwan[J].Geomorphology,90(1-2):36-54.

LIU X Y,PETERSON J,ZHANG Z Y,CHANDRA S.2005.Improving soil salinity prediction with high resolution DEM derived from LIDAR data[J].International Archives of the Photogrammetry,Remote Sensing and Spatial Information Sciences,36(7):41-43.

LU Ru-kun.1999.Analytical methods of soil agricultural chemistry[M].Beijing:China Agricultural Science Press(in Chinese).

ODEH I O A,MCBRATNEY A B,CHITTLEBOROUGH D J.1995.Further results on prediction of soil properties from terrain attributes:Heterotopic cokriging and regression-kriging[J].Geoderma,67(3-4):215-226.

SALAMA R B,OTTO C J,FITZPATRICK R W.1999.Contributions of groundwater conditions to soil and water salinization[J].Hydrogeology Journal,7:46-64.

SINGH R B,CHAUHAN CPS,MINHAS P S.2009.Water production functions of wheat (Triticum aestivum L.)irrigated with saline and alkali waters using double-line source sprinkler system[J].Agricultural Water Management,96:736-744.

TRIANTAFILIS J,ODEH I O A,MCBRATNEY A B.2001.Five geostatistical models to predict soil salinity from electromagnetic induction data across irrigated cotton[J].Soil Science Society of America Journal,65:869-878.

WANG Zun-qin,ZHU Shou-quan,YU Ren-pei.1993.Saline Soil in China[M].Beijing:Science Press(in Chinese).

WU W,FAN Y,WANG Z Y,LIU H B.2008.Assessing effects of digital elevation model resolutions on soil-landscape correlations in a hilly area[J].Agriculture,Ecosystems and Environment,126:209-216.

YANG J,HUANG Z C,CHEN T B,LEI M,ZHENG Y M,ZHENG G D,SONG B,LIU Y Q,ZHANG C S.2008.Predicting the probability distribution of Pb-increased lands in sewage-irrigated region:A case study in Beijing,China[J].Geoderma,147:192-196.

YAO Rong-jiang,YANG Jin-song,LIU Guang-ming.2006.Spatial variability of soil salinity and moisture and their estimations by CoKriging method – A case study in characteristic field of Yellow River Delta[J].Journal of Soil and Water Conservation,20(5):133-138(in Chinese with English abstract).

ZHAO Cheng-yi,WANG Yu-chao,LI Zi-liang,LI Guo-zhen.2003.Study on the spatial variability of soil moisture content and salt content in the field scale[J].Arid Zone Research,20(4):252-256(in Chinese with English abstract).

Spatial Variability of Soil Salinity,Total Dissolved Solid and Groundwater Depth Based on Cokriging in the Low Plain around the Bohai Sea

ZHOU Zai-ming,ZHANG Guang-hui,WANG Jin-zhe,YAN Ming-jiang
Institute of Hydrogeology and Environmental Geology,Chinese Academy of Geological Sciences,Shijiazhuang,Hebei050803

Taking into account the problem of soil salinity,total dissolved solid (TDS)of groundwater and groundwater level (GWL)existient in the low plain around the Bohai Sea,the authors collected soil samples at the depths of 0～20cm from 127 sites in the plain,determined the soil salinity,and measured shallow groundwater depth in 130 sites andTDS ofwater samples in 128 sites.Classical statistical and geostatistical methods combined with GIS technique were used to analyze spatial variability of soil salinity,TDS and GWL.The results show that spatial variability of soil salinity and TDS belongs to the moderate degree while GWL belongs to the strong degree with Cv being 0.91,Cv being 0.87,and Cv being 1.06.The semivariogram and cross-semivariogram of soil salinity,TDS and GWL were fitted by the spherical and the exponential model respectively.Spatial correlation of single variables belongs to the moderate degree,while cross variable belongs to the strong degree,with spatial correlation distance 135.6km to 223.2km longer than 40.6km to 144.8km.Ordinary Kriging and CoKriging maps show that soil salinity and TDS spatial distribution are almost the same,increasing from the inland plain to the east coastal plain,spatial distribution difference exists in different altitudes,and the soil salinity and TDS rise with the decrease of the altitude;nevertheless,things are just the opposite for GWL.Soil salinity is high in Tangshan-Tianjin-Cangzhou-Dongying-Binzhou,TDS is high in some areas of Tianjin and Cangzhou,and GWL is low in the Yellow River Delta.In comparison with the ordinary Kriging with the same sampling numbers,the root-mean-square error produced by CoKriging decreases by 0.29%～4.78%,while the correlation coefficient between the predicted value and the measured value increases by 2.03%～20.58%.

low plain around the Bohai Sea;soil salinity;total dissolved solid;groundwater level;CoKriging

P641.51;P322.3

A

10.3975/cagsb.2011.04.14

本文由國家科技支撐計劃項目(編號:2007BAD69B02;2009BADA3B05)和河北省科技廳重點基礎研究項目(編號:08966711D)聯合資助。

2011-03-16;改回日期:2011-04-15。責任編輯:魏樂軍。

周在明,男,1980年生。博士研究生。主要從事土壤水鹽管理的研究工作。E-mail:tougaozhou@163.com。