汪邦穩(wěn),方少文,楊勤科
(1.北京師范大學(xué)地理學(xué)與遙感科學(xué)學(xué)院,北京100875;2.江西省水土保持科學(xué)研究所,江西南昌330029;3.西北大學(xué)城市與環(huán)境學(xué)院,陜西西安710127)
水土流失是全球頭號(hào)環(huán)境問題,嚴(yán)重威脅著區(qū)域生態(tài)安全和糧食安全[1]。水土流失評(píng)價(jià)是水土保持科學(xué)研究的重要內(nèi)容,準(zhǔn)確評(píng)價(jià)區(qū)域水土流失可以為水土保持規(guī)劃和綜合治理提供科學(xué)依據(jù)。長(zhǎng)期以來,我國(guó)的水土流失研究是在比較微觀的尺度如坡面、小流域進(jìn)行的,而對(duì)全國(guó)和區(qū)域水土流失、區(qū)域地表過程和全球變化進(jìn)行研究,則需宏觀尺度水土流失定量評(píng)價(jià)方法的支持,以認(rèn)識(shí)土壤侵蝕對(duì)環(huán)境要素的影響、全球變化與區(qū)域土壤侵蝕的關(guān)系[2]。本研究以模型和“3S”技術(shù)為基礎(chǔ),探討贛南地區(qū)水土流失及其影響因子的評(píng)價(jià)方法,為準(zhǔn)確評(píng)價(jià)該地區(qū)的水土流失和科學(xué)進(jìn)行水土保持規(guī)劃提供技術(shù)、數(shù)據(jù)支持。
贛州市位于贛江上游、江西省南部,面積3.94萬 km2,轄1區(qū)2市15縣,2008年末總?cè)丝?88.95萬人,生產(chǎn)總值達(dá)834.77億元,現(xiàn)有耕地面積30.4萬hm2;地貌以山地丘陵為主,占區(qū)域總面積的83%,平均海拔300—500 m,坡度16°~45°;地帶性土壤主要有紅壤、黃壤和黃棕壤,其中紅壤分布最廣、面積最大;屬典型的亞熱帶濕潤(rùn)季風(fēng)氣候,年均降水量1 586.9 mm,年均蒸發(fā)量1 416.3 mm,降水年內(nèi)、年際分配不均,多年平均氣溫18.9℃;全市森林覆蓋率為74.2%,林相單一,屬人工次生林,以針葉林為主,林分結(jié)構(gòu)不合理,主要樹種為馬尾松和濕地松。贛州市是江西省最重要的糧油林果生產(chǎn)區(qū),但由于長(zhǎng)期不合理的經(jīng)濟(jì)活動(dòng),該地區(qū)的原生植被已遭到嚴(yán)重破壞,土壤侵蝕嚴(yán)重,生態(tài)環(huán)境亟待整治。
我國(guó)的區(qū)域水土流失評(píng)價(jià)工作起步較晚,由于缺乏基礎(chǔ)數(shù)據(jù)等原因,目前尚未有成熟完善的區(qū)域水土流失評(píng)價(jià)模型。近年來有學(xué)者利用美國(guó)的通用土壤流失方程(USLE)對(duì)我國(guó)的一些區(qū)域進(jìn)行了水土流失評(píng)價(jià),取得了較好的效果。但是,由于USLE模型是在美國(guó)平緩坡度的條件下建立的,其環(huán)境條件與我國(guó)不同,因此該模型在我國(guó)應(yīng)用具有一定的局限性。鑒于此,劉寶元等[3]基于USLE的建模經(jīng)驗(yàn),充分考慮到我國(guó)水土流失的環(huán)境特征,建立了中國(guó)土壤流失預(yù)報(bào)模型(China Soil Loss Erosion,簡(jiǎn)稱 CSLE),即
式中:A為多年平均土壤流失量,t/(hm2·a);R為降雨侵蝕力,MJ·mm/(hm2·a);K為土壤可蝕性因子,t·h/(MJ·mm·a);S為坡度因子,無量綱;L為坡長(zhǎng)因子,無量綱;B為水土保持生物措施因子,無量綱;E為水土保持工程措施因子,無量綱;T為水土保持耕作措施因子,無量綱。
該模型的特點(diǎn),一是根據(jù)我國(guó)水土保持的實(shí)際情況,將USLE中的作物和水土保持措施兩大因子變?yōu)樗帘3秩蟠胧┮蜃?,即變?yōu)樯铩⒐こ毯透鞔胧┮蜃?二是模型的結(jié)構(gòu)相對(duì)簡(jiǎn)單,便于推廣應(yīng)用。鑒于此,本研究選擇該模型作為贛南地區(qū)水土流失的評(píng)價(jià)模型。
CSLE中降雨侵蝕力(R)是方程中一個(gè)最基本的因子,它是降雨量、降雨強(qiáng)度、雨型和雨滴動(dòng)能的函數(shù)。國(guó)內(nèi)外的研究表明,坡面土壤流失量與降雨能量E及某一時(shí)段最大降雨強(qiáng)度It的乘積EIt有著密切的關(guān)系。但是,由于次降雨過程資料難以獲得,因此一般利用氣象站常規(guī)觀測(cè)統(tǒng)計(jì)資料如日降雨量、月降雨量、年降雨量建立降雨侵蝕力簡(jiǎn)易算法。我國(guó)學(xué)者在這方面做了大量的研究工作,建立了全國(guó)及各地區(qū)R值計(jì)算模型并分析了其時(shí)空分布特征。
3.1.1R值的計(jì)算方法
由于日降雨量資料在贛南地區(qū)可以獲得,因此研究中采用章文波[4]基于全國(guó)日降雨量數(shù)據(jù)構(gòu)建的R值簡(jiǎn)易算法計(jì)算區(qū)域R值,計(jì)算式為
式中:R為降雨侵蝕力,MJ·mm/(hm2·a);α、β為模型參數(shù);k為某半月內(nèi)侵蝕性降雨日數(shù),d,k=1,2,…,m;Pk為半月內(nèi)第k天的日降雨量,mm,利用紅壤區(qū)的研究成果,贛南地區(qū)侵蝕性降雨標(biāo)準(zhǔn)為10 mm;為多年平均降雨侵蝕力,MJ·mm/(hm2·a);n為年數(shù);Pd10為一年內(nèi)侵蝕性降雨日雨量的平均值(即一年中≥10 mm日降雨量的總和與相應(yīng)日數(shù)的比值),mm;Py10是侵蝕性降雨年總量的多年平均值(即≥10 mm日雨量年累加值的多年平均值),mm。
3.1.2R值柵格數(shù)據(jù)獲取
采用式(2)計(jì)算獲取的各站點(diǎn)的降雨侵蝕力R值,借助ArcMap軟件,采用Kriging內(nèi)插方法對(duì)降雨侵蝕力點(diǎn)數(shù)據(jù)進(jìn)行空間表面插值,得到降雨侵蝕力柵格數(shù)據(jù)。
K為土壤可蝕性因子,是表征土壤性質(zhì)對(duì)侵蝕敏感程度的指標(biāo),是進(jìn)行水土流失定量評(píng)價(jià)的重要依據(jù)。USLE中對(duì)土壤可蝕性因子的定義為:標(biāo)準(zhǔn)小區(qū)春天順坡耕翻一次,裸土狀態(tài)下單位降雨侵蝕力的年土壤流失量。對(duì)土壤可蝕性的研究主要是用來評(píng)價(jià)土壤性質(zhì)和水土流失的關(guān)系,國(guó)外早在20世紀(jì)20年代就開始進(jìn)行,我國(guó)土壤可蝕性研究大體上可分為土壤可蝕性評(píng)價(jià)指標(biāo)確定、USLE中K值的測(cè)定和估算、土壤可蝕性指標(biāo)的時(shí)空變化及不確定性分析3個(gè)階段。張科利在分析布設(shè)在我國(guó)不同水土流失區(qū)徑流小區(qū)實(shí)測(cè)資料的基礎(chǔ)上,對(duì)我國(guó)土壤可蝕性K值進(jìn)行了系統(tǒng)研究,觀測(cè)計(jì)算出了一組土壤可蝕性K值,并對(duì)無實(shí)測(cè)資料地區(qū)土壤可蝕性估算問題進(jìn)行了探討。
3.2.1K值的計(jì)算方法
在區(qū)域尺度上土壤可蝕性難以用小區(qū)觀測(cè)獲得,因此國(guó)內(nèi)外學(xué)者構(gòu)建了土壤可蝕性估算模型,盡管關(guān)于土壤可蝕性值估算的研究很多,但具有代表性的成果為EPIC模型。本研究根據(jù)現(xiàn)有的數(shù)據(jù)資料,選用式(3)EPIC模型計(jì)算贛南地區(qū)土壤可蝕性K值,然后利用張科利[5]的修訂公式[式(4)]對(duì)獲得的數(shù)據(jù)進(jìn)行修正,最后得到贛南地區(qū)的土壤可蝕性因子K值。上兩式中:SAN、SIL、CLA分別為土壤中砂粒、粉粒、黏粒的含量,%;C為土壤有機(jī)碳含量,%;SN1=1-SAN/100;Kepic為EPIC模型計(jì)算的土壤可蝕性,t·h/(MJ·mm·a)。
3.2.2 土壤可蝕性因子K柵格數(shù)據(jù)獲取
利用江西省1∶50萬土壤圖,借助ArcMap軟件,查閱《江西紅壤》、《中國(guó)土種志》獲取贛南各土壤類型的粉粒、極細(xì)砂粒、砂粒、有機(jī)質(zhì)、有機(jī)碳含量,依據(jù)式(3)獲取EPIC_K值矢量數(shù)據(jù),依據(jù)式(4)獲取CSLE_EPIC_K值矢量數(shù)據(jù),然后進(jìn)行數(shù)據(jù)轉(zhuǎn)換獲取CSLE_EPIC_K柵格數(shù)據(jù),詳見圖1。
圖1 土壤可蝕性因子K值計(jì)算流程
3.3.1LS值計(jì)算方法
在USLE和RUSLE(修正的通用土壤流失方程)中,地形對(duì)水土流失的影響用坡度坡長(zhǎng)因子表示。坡度因子(slope gradient factor)定義為某一坡度土壤流失量與坡度為5.13°、其他條件都一致的坡面產(chǎn)生的土壤流失量之比率。坡長(zhǎng)因子(slope length factor)在USLE中定義為某一坡面土壤流失量與坡長(zhǎng)為22.1 m、其他條件都一致的坡面產(chǎn)生的土壤流失量之比率。坡度與坡長(zhǎng)因子的乘積稱為地形因子,反映了坡度坡長(zhǎng)對(duì)土壤流失的綜合影響。
坡度坡長(zhǎng)因子難以通過小區(qū)觀測(cè)獲取,目前主要以DEM為基礎(chǔ)按一定的算法提取,在RUSLE中,坡度因子、坡長(zhǎng)因子的算法如下
上二式中:L為坡長(zhǎng)因子,無量綱;λ為由DEM提取的坡長(zhǎng),m;22.1為22.1 m標(biāo)準(zhǔn)小區(qū)斜坡長(zhǎng);α為指數(shù);S為坡度因子,無量綱;θ為由DEM提取的坡度值,(°)。
由于RUSLE是在緩坡條件下的天然徑流小區(qū)觀測(cè)資料基礎(chǔ)上建立的,其坡度坡長(zhǎng)因子的計(jì)算也是針對(duì)這種緩坡地形,直接采用上面的算法計(jì)算陡坡地區(qū)的坡度坡長(zhǎng)因子值存在較大的誤差,因此劉寶元對(duì)RUSLE中坡度坡長(zhǎng)因子用于陡坡時(shí)的計(jì)算方法進(jìn)行了修正,修正后的坡度因子計(jì)算式為[6]
坡長(zhǎng)因子算法為
式中:m為坡長(zhǎng)指數(shù),根據(jù)坡度的不同取不同的值;其他符號(hào)意義同上。
3.3.2 基礎(chǔ)數(shù)據(jù)處理方法
由于國(guó)家基礎(chǔ)地理信息中心生產(chǎn)的DEM上面有一些殘留的平三角,不能完全滿足水土流失評(píng)價(jià)工作的需求,因此本研究基于1∶5萬地形圖,在ANUDEM專業(yè)軟件支持下插值生成水文地貌關(guān)系正確的DEM。ANUDEM是Hutchinson教授在長(zhǎng)期研究成果的基礎(chǔ)上開發(fā)的,其開發(fā)應(yīng)用歷史已超過20年,作為專業(yè)化DEM生產(chǎn)軟件,它已受到許多國(guó)家的廣泛關(guān)注。把等高線、高程點(diǎn)、河流等信息輸入到軟件中,設(shè)置各參數(shù)值,運(yùn)行輸出高質(zhì)量的DEM。
3.3.3LS值獲取方法
用AML語言編程,在ArcGIS/Info下進(jìn)行LS因子值提取,主要計(jì)算過程包括:填洼、流向和柵格坡長(zhǎng)的計(jì)算,局地山頂點(diǎn)和坡度變化點(diǎn)的提取,坡度和坡長(zhǎng)計(jì)算,坡度和坡長(zhǎng)因子值計(jì)算。計(jì)算流程如圖2所示。
圖2 坡度坡長(zhǎng)因子計(jì)算流程
B為水土保持生物措施因子,是水土流失評(píng)價(jià)的基本參數(shù)。B因子值是指在一定條件下,耕作農(nóng)地上的土壤流失量與同等條件下適時(shí)翻耕的連續(xù)休閑對(duì)照地上的土壤流失量之比。它是一個(gè)無量綱參數(shù),其值變化于0~1之間。
從國(guó)內(nèi)外的研究來看,計(jì)算B值的方法之一是根據(jù)作物小區(qū)與裸露小區(qū)多年平均土壤流失量比較計(jì)算得到。目前國(guó)內(nèi)B因子的研究大多局限于一定作物,數(shù)據(jù)一般從實(shí)測(cè)徑流小區(qū)資料計(jì)算得到,或是根據(jù)作物農(nóng)作期劃分以及所在區(qū)域的降雨侵蝕力季節(jié)分布等計(jì)算,所得數(shù)據(jù)僅能滿足徑流小區(qū)或小流域研究的需要,而不能滿足區(qū)域B值研究的需要。
區(qū)域尺度水土保持生物措施因子一般是通過遙感圖像提取的。為了滿足區(qū)域水土流失在遙感監(jiān)測(cè)中的應(yīng)用,有學(xué)者將區(qū)域B因子研究轉(zhuǎn)為土地利用以及不同的土地利用類型的植被蓋度研究。根據(jù)土地利用圖和植被蓋度圖,結(jié)合前人研究成果,賦予不同土地利用類型、不同植被蓋度不同的B值,得到區(qū)域B值柵格數(shù)據(jù)。該方法在陜西省水土流失普查試點(diǎn)研究中應(yīng)用,效果比較理想。因此,本研究利用遙感解譯獲取土地利用和提取植被蓋度,在此基礎(chǔ)上獲取B值。
3.4.1 基礎(chǔ)數(shù)據(jù)處理方法
(1)遙感解譯。目前,遙感解譯比較成熟的方法有三種,即目視解譯、計(jì)算機(jī)屏幕解譯和計(jì)算機(jī)輔助自動(dòng)分類??紤]到工作量和多期解譯的系統(tǒng)性,研究中統(tǒng)一采用監(jiān)督分類和非監(jiān)督分類相結(jié)合的方法,同時(shí)結(jié)合DEM、Google Earth、河流行政區(qū)劃圖以及專家意見等進(jìn)行解譯,具體的數(shù)據(jù)處理流程如圖3所示。
圖3 B值基礎(chǔ)數(shù)據(jù)處理流程
(2)NDVI提取。植被蓋度信息的提取,是區(qū)域水土流失監(jiān)測(cè)和評(píng)價(jià)的基礎(chǔ)。目前學(xué)者提出的植被指數(shù)多達(dá)幾十種,不同的植被指數(shù)有不同的優(yōu)缺點(diǎn),適用范圍也不盡相同,但NDVI是目前應(yīng)用最為廣泛的植被指數(shù)。遙感圖像處理軟件中都具有植被指數(shù)提取功能,如在ERDAS中可以直接提取TM、MSS等有關(guān)的植被指數(shù),在這些軟件下可以快捷地獲取所需要的植被指數(shù)。
植被蓋度可以通過NDVI推算獲得,其計(jì)算式為
式中:f為植被蓋度;NDVI為所求像元的NDVI值,無量綱;NDVImin、NDVImax分別為研究區(qū)內(nèi)NDVI的最小值與最大值,無量綱。
3.4.2 植被因子B值獲取方法
鑒于贛南地區(qū)歷史和自然的原因,林下流現(xiàn)象普遍存在且十分嚴(yán)重,“遠(yuǎn)看青山在,近看水土流”是該地區(qū)水土流失的真實(shí)寫照。江西省水保所對(duì)贛南的考察發(fā)現(xiàn),海拔500 m以上的林地植被結(jié)構(gòu)很好,基本上不存在林下流現(xiàn)象,而位于海拔200 m以下的馬尾松和濕地松林水土流失普遍存在。林下水土流失的主要原因是林下無灌木和草被層,林地對(duì)水土流失的控制作用十分有限,個(gè)別情況下甚至沒有作用。因此,B值獲取的方法主要是根據(jù)已有的研究成果,結(jié)合植被蓋度給贛南地區(qū)不同的土地利用進(jìn)行賦值。對(duì)林地B值的計(jì)算,既考慮了植被蓋度又考慮了海拔的影響。贛南地區(qū)不同土地利用類型和不同植被蓋度下的B值見表1。
表1 贛南地區(qū)不同土地利用類型和不同植被蓋度下的B值
國(guó)內(nèi)外學(xué)者對(duì)水土保持措施因子(E)進(jìn)行了大量的研究,但由于研究方法和途徑的不同,所以目前還沒有一個(gè)有效的E值計(jì)算模型。本研究參考影響水土流失的其他因子的定義方法,將有工程措施和無工程措施水土流失的比值定義為水土保持工程措施因子。對(duì)于面積達(dá)3.94萬km2的贛南地區(qū)來說,要收集該區(qū)域內(nèi)所有的含有空間信息的工程措施需要很長(zhǎng)的時(shí)間,而且過去的工程措施情況也無法獲得,為此我們通過統(tǒng)計(jì)贛南各縣(區(qū)、市)歷年的數(shù)據(jù),結(jié)合野外調(diào)查和專家咨詢,以縣(區(qū)、市)為單位收集了贛南地區(qū)的梯田、水平竹節(jié)溝、谷坊、淤地壩等工程措施資料,參考現(xiàn)有的研究成果[7],嘗試計(jì)算了區(qū)域水土保持措施中的工程措施因子值。計(jì)算公式為
式中:St為梯田面積,km2;Sglt為谷坊、攔沙壩和塘壩的控制面積,km2;Sz為水平竹節(jié)溝控制面積,km2;S為土地面積,km2;α為梯田的減沙系數(shù);β為谷坊、攔沙壩、塘壩的減沙系數(shù);γ為水平竹節(jié)溝減沙系數(shù)。
實(shí)際工作中,因水土保持工程措施的類型不同,其控制面積也會(huì)不同,因此根據(jù)現(xiàn)場(chǎng)調(diào)查,結(jié)合《水利部水土保持工程建設(shè)標(biāo)準(zhǔn)》來確定單位水土保持工程措施的控制面積。確定的結(jié)果為:谷坊控制面積1.5 hm2,塘壩控制面積30 hm2,水平竹節(jié)溝控制面積0.4 hm2/km(水平竹節(jié)溝每間隔4 m挖一條帶),攔沙壩控制面積200 hm2。有關(guān)各種水土保持工程措施的減沙效益已經(jīng)有大量的研究成果,本研究根據(jù)對(duì)水保工作者的咨詢,同時(shí)參考前人的研究成果,確定梯田的減沙系數(shù)α為76.3%,谷坊、攔沙壩、塘壩的減沙系數(shù)β為90%,水平竹節(jié)溝減沙系數(shù)γ為65%。
利用式(10)實(shí)現(xiàn)對(duì)贛南各縣(區(qū)、市)水土保持工程措施因子E值的獲取,借助ArcGIS軟件,把獲取的E值作為屬性數(shù)據(jù)轉(zhuǎn)化成贛南水土保持工程措施因子E值的柵格數(shù)據(jù)。
水土流失是造成生態(tài)環(huán)境惡化的重要影響因素,科學(xué)、快速、準(zhǔn)確地評(píng)價(jià)水土流失可以為區(qū)域水土保持規(guī)劃和治理提供科學(xué)依據(jù)。但是,由于區(qū)域差異性的存在,每個(gè)地方的水土流失影響主導(dǎo)因素各有區(qū)別,因此本研究以模型和“3S”技術(shù)為基礎(chǔ),在弄清各個(gè)因子含義和算法的基礎(chǔ)上,探討了適合贛南地區(qū)侵蝕環(huán)境的水土流失評(píng)價(jià)模型及其因子值的獲取途徑。同時(shí),對(duì)贛南地區(qū)林下水土流失區(qū)植被因子B值以及贛南水土保持工程措施因子E值的獲取方法進(jìn)行了探討,為科學(xué)、快速、準(zhǔn)確地評(píng)價(jià)該地區(qū)的水土流失提供了一定的技術(shù)支持。
[1]李銳.中國(guó)水土流失基礎(chǔ)研究的機(jī)遇與挑戰(zhàn)[J].自然雜志,2008,30(1):6-11.
[2]楊勤科,李銳.我國(guó)區(qū)域土壤侵蝕與環(huán)境研究述評(píng)[J].中國(guó)人口·資源與環(huán)境,2006,16(6):90-94.
[3]劉寶元,畢小剛,付素華.北京土壤流失方程[M].北京:科學(xué)出版社,2010.
[4]章文波,謝云.利用日雨量計(jì)算降雨侵蝕力的方法研究[J].地理科學(xué),2002,22(6):705-711.
[5]張科利,彭文英,楊紅麗.中國(guó)土壤可蝕性值及其估算[J].土壤學(xué)報(bào),2007,44(1):7-13.
[6]Liu B Y,Risse L M.Slope gradient effects on soil loss for steep slopes[J].Transactions of the ASAE,1994,37(6):1 835-1 840.
[7]謝紅霞,李銳,楊勤科,等.退耕還林(草)和降雨變化對(duì)延河流域土壤侵蝕影響[J].中國(guó)農(nóng)業(yè)科學(xué),2009,42(2):569-576.