基于ENVI和ArcGis的云南省侵蝕模數圖量算方法

劉 剛，劉陽容，楊林山，林 蕓，房林東，田 帥，李紹平

(1.中國電建集團昆明勘測設計研究院有限公司，云南 昆明 650051；2.中國科學院西北生態環境資源研究院，甘肅 蘭州 730000；3.成都市新都區城鄉建設局，四川 成都 610500；4.湘西土家族苗族自治州水利水電勘測設計研究院，湖南 吉首 416000)

采用侵蝕模數圖法來推求設計流域多年平均產沙量的方法已廣泛應用于云南省中小流域水利水電工程泥沙計算中[1]。采用AutoCAD的傳統云南省侵蝕模數圖量算方法需以侵蝕模數圖為底圖，人工手繪侵蝕模數圖圖斑，并進行面積換算，再按照侵蝕類型取對應侵蝕模數，采用面積加權法計算流域平均侵蝕模數，進而推求設計流域多年平均產沙量。上述方法具有工作量大、耗時長、易出錯、量算成果不能重復使用且與ArcGis水文分析成果不能有效銜接等缺點。

目前，在設計工作中使用的云南省侵蝕模數圖為柵格數據(紅、綠、藍3個通道)，模數圖坐標為大地坐標，有75幅之多。設計流域泥沙計算需要水文集水區劃分成果(流域范圍)作為計算邊界，量算時需要單位換算，不便于直接量算。鑒于ArcGis水文集水區劃分在設計中的廣泛運用，且侵蝕模數圖法在云南省使用廣泛，迫切需要提出一種新的方法以實現云南省侵蝕模數圖的快速量算，克服傳統方法的諸多弊端，以提高設計人員的工作效率。

基于上述原因，本文提出了基于ENVI和ArcGis的云南省侵蝕模數圖量算方法。本文利用ArcGis對云南省75幅侵蝕模數圖進行處理，制作帶投影信息的云南省侵蝕模數拼接柵格圖。利用ENVI軟件與ArcGis軟件良好的兼容性，保留侵蝕模數拼接柵格圖的投影信息，使用ENVI圖像分類以及分類后處理功能對帶投影信息的云南省侵蝕模數拼接柵格圖進行矢量化[2]。侵蝕模數圖矢量化后因帶有投影信息，能夠與基于ArcGis的水文集水區劃分成果進行疊圖分析，可推求流域(流域范圍全在云南省)多年平均產沙量。

1 基礎數據制作

基于ENVI和ArcGis的云南省侵蝕模數圖量算方法，需要利用ArcGis與ENVI軟件對設計中常用的侵蝕模數柵格圖進行處理，得到云南省侵蝕模數矢量圖。因設計項目往往分散，為方便工作避免因設計流域不同而重復基礎數據制備環節，本文建議利用ArcGis對云南省75幅侵蝕模數圖進行坐標校正、裁剪、拼接以及投影變化等處理，制作帶投影信息的云南省侵蝕模數柵格拼接圖，再利用ENVI圖像分類、分類后處理等功能制作全云南省侵蝕模數矢量圖。

a)數據制作步驟1。云南省侵蝕模數圖共計75幅，每幅成圖比例尺為1∶20 萬。為保證制作的矢量圖具有很高的空間位置精度，建議利用ArcGis軟件的[Georeferencing]工具對每張侵蝕模數柵格圖進行格網逐點校正，使用[Raster]-[Clip]工具裁去侵蝕模數柵格圖多余邊角，利用[Mosaic To New Raster]工具完成所有柵格數據的拼接。侵蝕模數柵格圖坐標系為大地坐標，使用時需要轉化為投影坐標，可使用[Project Raster]工具對拼接好的柵格數據進行投影轉換，將大地坐標系轉為云南省水文手冊附圖相同的投影坐標系[3]。

b)數據制作步驟2。ENVI軟件與ArcGis軟件具有良好的兼容性，可保留侵蝕模數拼接柵格圖的投影信息。在ArcGis軟件中，將云南省侵蝕模數拼接柵格圖導出為ENVI可識別的柵格數據形式，并利用ENVI遙感軟件對已經拼接好的云南省侵蝕模數圖進行監督分類[4]。依據侵蝕模數圖類別(無明顯侵蝕、輕度侵蝕、中度侵蝕、強度侵蝕、極強度侵蝕、劇烈侵蝕、凍融侵蝕、壩子、水體、城市)分別建立訓練樣本。使用最小距離分類方法進行分類[5]。其中，訓練樣本的分離性可采用ENVI中ROI Separability Calculation工具進行計算。分離系數范圍為0～2，大于1.9表明樣本之間分離性好；小于1.8需要重新選擇樣本；低于1需合并分類樣本。本次云南省侵蝕模數圖樣本分離系數中，凍融侵蝕與城市類型分離系數為0，兩者區分度較低。可結合凍融侵蝕分布特點，后期在ArcGis軟件中人工區分。

云南省侵蝕模數柵格圖經過監督分類后，分類結果會不可避免的出現細小圖斑，可利用ENVI軟件中Majority/Minority分析、聚類處理(Clump)及過濾處理(Sieve)等功能進行處理[6-7]。經上述處理后，利用ENVI軟件將分類結果文件轉為ArcGis軟件支持的Shapefile(.shp)矢量數據格式，即可得到云南省侵蝕模數矢量(圖1)。

圖1 云南省侵蝕模數矢量圖制作流程示意

2 數據精度評估

為評價上述方法制作的云南省侵蝕模數矢量圖的分類精度，本文借用Kappa系數進行評價。Kappa系數是指一種對遙感圖像分類精度和誤差矩陣進行評價的多元離散方法[8]。Kappa系數計算基于混淆矩陣，為創建混淆矩陣，綜合考慮評估工作量大小，本文利用ArcGis [Create Random Points]工具隨機生成200 個點，并利用云南省侵蝕模數柵格圖人工賦值真實侵蝕類型，使用[Spatial Join]工具對隨機點添加本次制作的云南省侵蝕模數矢量圖侵蝕類型。云南省侵蝕模數矢量圖精度評估點位置分布示意見圖2。依據隨機點屬性表，通過簡單的統計即可構建混淆矩陣(表1)[9]。

圖2 云南省侵蝕模數矢量圖精度評估點位置分布示意

一般Kappa系數計算值在0.00～1.00之間。其中，0.00～0.20表示分類質量差；0.21～0.40表示分類質量一般；0.41～0.60表示分類質量好；0.61～0.80表示分類質量很好；0.81～1.00表示分類質量極好[10]。計算得Kappa系數為1，表明本次云南省侵蝕模數矢量圖分類成果質量極好，與原侵蝕模數柵格圖一致，能夠用于設計工作中。

表1 云南省侵蝕模數矢量圖精度評估混淆矩陣

3 量算方法

基于ArcGis的水文集水區劃分能夠提供流域泥沙計算邊界，計算邊界因帶有空間參考信息，可與步驟2制作的云南省侵蝕模數矢量圖進行疊圖分析，進而推求流域多年平均產沙量。結合水文集水區劃分成果，可使用ArcGis中的 [Clip]工具對步驟2中云南省侵蝕模數矢量圖進行裁剪，可利用[Summary Statistics]工具對裁剪結果文件按侵蝕類型分類匯總，統計各侵蝕類型面積。裁剪前設計流域與云南省侵蝕模數矢量圖應具有相同坐標系，且裁剪后需要更新裁剪成果屬性表的面積字段。

依據SL 190—2007《土壤侵蝕分類標準》，土壤水力侵蝕強度共分為6級，分別為微度、輕度、中度、強烈、極強烈及劇烈。在設計計算中，云南省侵蝕模數圖中城市、水體、壩子無明顯侵蝕類型歸為微度級別，輕度侵蝕類型歸為輕度級別，中度侵蝕類型歸為中度級別，強度侵蝕類型歸為強度級別，極強烈侵蝕類型歸為極強級別，劇烈侵蝕類型歸為劇烈級別。在缺少實測及調查侵蝕模數資料時，不同級別侵蝕模數可參考SL 190—2007《土壤侵蝕分類標準》進行取值。設計中使用的云南省侵蝕模數圖不同侵蝕類型對應侵蝕模數取值見表2。

表2 云南省侵蝕模數圖不同侵蝕類型對應侵蝕模數取值

依據各侵蝕級別面積，結合表1中不同侵蝕級別對應的平均侵蝕模數，可依據式(1)計算設計流域平均侵蝕模數。依據計算得到的流域平均侵蝕模數采用式(2)計算得到設計流域多年平均產沙量。

Mave=∑Mi·Ai

(1)

SA=Mave·F

(2)

式中Mave——流域平均侵蝕模數，t/(km2·a)；Mi——不同侵蝕類型侵蝕模數，t/(km2·a)；Ai——不同侵蝕類型面積，km2；F——設計流域面積，km2；SA——設計流域多年平均產沙量，t/a。

4 實例應用

4.1項目區概況

金龍河位于云南省中西部，屬瀾滄江流域，是瀾滄江一級支流黑惠江的上游段，在匯入劍湖之前稱為金龍河，出劍湖口之后始稱黑惠江。金龍河發源于麗江市古城區九河鄉白漢場羅鳳山，先后流經九河、三河、梅園、金龍、金星等地，于金華鎮柳營村附近匯入劍湖。金龍河流域面積674.1 km2，河長40.6 km，河道比降3.7‰。下文將以此流域多年平均產沙量計算為例說明本方法的實現過程。

4.2量算過程

本文以金龍河流域范圍ASTER GDEM 2 DEM為水文分析地形數據[11]，利用ArcGis水文分析工具進行水文集水區劃分，初步確定金龍河流域邊界。在地勢平坦地區，利用ASTER GDEM 2數據進行水文分析，其集水區劃分成果可能存在較大偏差。本文采用該流域附近1∶5萬掃描地形圖為底圖，對水文集水區劃分成果進行校正。

使用校正后的金龍河流域邊界，依本文前述量算方法對云南省侵蝕模數矢量圖進行裁剪，分類統計各侵蝕類型面積。金龍河流域侵蝕模數圖類型有壩子、城市、水體、無明顯侵蝕、輕度侵蝕、中度侵蝕、強度侵蝕。對比金龍河流域侵蝕模數矢量圖斑(圖3)與原侵蝕模數柵格圖可知，金龍河流域侵蝕模數矢量圖圖斑邊界與柵格圖圖斑邊界幾乎完全一致，表明基于ENVI和ArcGis的云南省侵蝕模數圖量算方法是可行的。本次金龍河流域共涉及140個圖斑，采用本法能夠大大縮短侵蝕模數圖量算時間，提高設計人員的工作效率。

a)金龍河流域侵蝕模數矢量

b)金龍河流域侵蝕模數柵格矢量邊界對比圖3 金龍河流域侵蝕模數矢量、柵格對比情況

按照表1侵蝕級別合并侵蝕類型，經統計，微度級別面積占比72.6%，輕度級別面積占比19.3%，中度級面積別占比8.0%，強度級別面積占比為0.1%。依據侵蝕級別面積匯總成果，按式(1)計算流域平均侵蝕模數。在實際工作中，通常擬定3個計算方案(平均侵蝕模數取高值、中值及低值)，并結合相關資料及工程經驗選取合適的計算方案。本次考慮到金龍河流域主要以微度侵蝕為主，且隨著流域內小流域綜合治理工作的深入開展，流域水土流失將減少。因此，本次金龍河不同侵蝕級別侵蝕模數取低值方案，即微度級別取200 t/(km2·a)，輕度級別取500 t/(km2·a)，中度級別取2 500 t/(km2·a)，強度級別取5 000 t/(km2·a)。依據計算，金龍河流域多年平均侵蝕模數為448 t/(km2·a)，流域多年平均產沙量為30.2 萬t。

5 結語

本方法對于單幅侵蝕模數圖量算同樣適用，只是數據制作步驟1更加簡單。從方便工作角度出發，建議制作全云南省侵蝕模數矢量圖。本文采用的侵蝕模數圖量算方法具有量算速度快、精度高、計算簡單的優點。可有效減少設計人員工作量，對設計人員實際工作有較大幫助。同時，在使用設計流域裁剪侵蝕模數矢量圖后，可將裁剪結果與侵蝕模數柵格圖進行對比，根據項目需要對局部不一致的地方進行修正。