神東礦區2005—2018年水土流失動態變化研究

王麗云，王小燕，王 義，岳本江，屈 創

(1.黃河上中游管理局，陜西 西安 710021； 2.神華神東煤炭集團有限公司，陜西 神木 719300)

神東礦區位于晉陜蒙接壤地區，是我國最重要的煤炭生產和能源基地之一[1]。面對大規模開采與極脆弱生態環境保護之間的矛盾，礦區不斷創新生產方式，探索生態治理與環境保護的新理念，提出了以控制外圍風沙侵蝕為前提，內外圍結合治理，促進礦區整體恢復和改善礦區生態環境的思路，走出了一條資源保護型開采與生態環境治理相協調的綠色礦業之路[2]。為全面體現神東礦區2005—2018年水土流失綜合治理成果，客觀評價礦區水土保持工作成效，本研究采用中國土壤流失方程(CSLE)和風力侵蝕模型，定量計算礦區2005、2010、2015、2018年4期土壤侵蝕模數，并對不同時期不同侵蝕強度土地面積進行對比分析，希望能為礦區水土流失綜合治理提供數據支撐。

1 研究區概況

神東礦區分布于陜西省榆林市神木市北部、府谷縣西部，內蒙古自治區鄂爾多斯市東勝區及伊金霍洛旗的南部、準格爾旗的西南部，地理位置為109°51′～110°46′E、38°52′～39°41′N，包括大柳塔礦、活雞兔礦、哈拉溝礦等14個礦井，總面積約1 064 km2，黃河一級支流窟野河(烏蘭木倫河)縱貫礦區中部。地處黃土高原北部與毛烏素沙地交錯區，具有典型的干旱荒漠高原氣候特征，多年平均降水量僅362 mm，年水面蒸發量為2 297.4～2 838.7 mm[1]。區內大部分為典型的風成沙丘及沙漠灘地，尤其是北部、東南部為黃土丘陵溝壑區，梁峁起伏，溝壑縱橫，地表支離破碎，風水蝕交互作用導致水土流失嚴重。

2 監測方法

2.1 數據來源與預處理

本研究主要涉及遙感影像、地形數據、氣象資料、土壤資料等基礎數據。

(1)遙感影像。采用2005、2010年的SPOT5影像和2015、2018年的ZY-3數據進行2005、2010、2015、2018年礦區土地利用和水土保持措施遙感解譯。采用監測年前3年空間分辨率為250 m的MOD13Q1數據計算植被覆蓋度。采用監測年逐日AMSR-E Level 2A亮溫數據產品(空間分辨率為25 km)計算表土濕度數據。遙感影像使用前進行輻射校正、幾何校正、配準、影像融合、鑲嵌與裁切等預處理，融合后的影像既提高了多光譜影像的空間分辨率，又保留了其多光譜特性。遙感數據統一采用CGCS2000坐標系、Albers投影。

(2)地形數據。基于礦區1∶10 000比例尺的數字高程模型DEM，用符素華等[3]開發的土壤侵蝕模型地形因子計算工具計算坡長、坡度因子等。

(3)其他數據。收集礦區8個站點2000年以來的逐日降水量資料和逐日整點風速資料，用于計算降雨侵蝕力因子和風力因子。植被覆蓋度修訂系數采用第一次全國水利普查水土保持情況普查成果數據。

2.2 土壤侵蝕模型

2.2.1 水力侵蝕模型

在水力侵蝕區，采用中國土壤流失方程(CSLE)計算土壤侵蝕模數[4-5]。計算公式為

A=R·K·L·S·B·E·T

(1)

式中：A為土壤侵蝕模數，t/(hm2·a)；R為降雨侵蝕力因子，MJ·mm/(hm2· h·a)；K為土壤可蝕性因子，t·hm2·h/(hm2·MJ·mm)；L為坡長因子；S為坡度因子；B為生物措施因子；E為工程措施因子；T為耕作措施因子。

(1)降雨侵蝕力因子。將8個站點多年平均1～24個半月降雨侵蝕力矢量化，利用普通克里金空間插值方法，生成空間分辨率為2 m的24個半月降雨侵蝕力柵格數據。

(2)土壤可蝕性因子。采用第一次全國水利普查水土保持情況普查土壤可蝕性因子成果作為礦區土壤可蝕性因子值。

(3)坡長、坡度因子。基于區域DEM數據，采用北京師范大學符素華等[3]開發的地形因子計算工具計算坡長、坡度因子值，并對林草地采用緩坡公式修正。

(4)生物措施因子。基于MOD13Q1 NDVI數據采用像元二分法計算植被覆蓋度，利用植被覆蓋度修訂系數，獲得每年24個半月30 m空間分辨率的植被覆蓋度，再結合24個半月降雨侵蝕力因子比例與土地利用類型計算園地、林地、草地的生物措施因子；參照《區域水土流失動態監測技術規定(試行)》對非園地、林地和草地生物措施因子查表賦值。

(5)工程措施因子。根據野外調查和遙感影像解譯，礦區水土保持工程措施主要為土坎水平梯田，工程措施因子賦值為0.084。

(6)耕作措施因子。礦區在全國輪作區中涉及兩個一級區、三個二級區。一級區北部中高原半干旱喜涼作物一熟區下的二級區后山壩上晉北高原山地半干旱喜涼作物一熟區和隴中青東寧中南黃土丘陵半干旱喜涼作物一熟區，耕作措施因子取值0.488；一級區北部低高原易旱喜溫一熟區下的二級區黃土高原東部易旱喜溫作物一熟區，耕作措施因子取值0.417。

2.2.2 風力侵蝕模型

在風力侵蝕地區，根據土地利用類型，選用耕地、草(灌)地、沙地(漠)等風力侵蝕模型，計算風力侵蝕模數[6-7]。耕地、草(灌)地、沙地(漠)風力侵蝕模型計算公式分別為

(2)

(3)

(4)

式中：Qfa為每半個月內耕地風力侵蝕模數，t/(hm2·a)；W為每半個月內表土濕度因子，介于0～1之間；Tj為每半個月內各風速等級的累計時間，min；Z0為地表粗糙度，cm；j為風速等級序號；Uj為第j等級的平均風速，m/s；Qfg為每半個月內草(灌)地風力侵蝕模數，t/(hm2·a)；V為植被覆蓋度，%；Qfs為每半個月內沙地(漠)風力侵蝕模數，t/(hm2·a)。

(1)風力因子。風力因子是指不同等級風速和風向對風力侵蝕的潛在能力，與臨界侵蝕風速和各等級風速的累積時間有關。風力因子計算采用收集的礦區逐日24 h整點風速和風向資料，按1 m/s間隔統計全年大于等于臨界侵蝕風速的各等級風速的累積時間，利用普通克里金方法，插值生成24個半月各風速等級對應的累計時間。耕地的臨界侵蝕風速取值為5 m/s，草(灌)地、沙地(漠)的臨界侵蝕風速按照李智廣等[7]在我國風力侵蝕抽樣調查方法研究中不同植被覆蓋度的臨界侵蝕風速取值。

(2)表土濕度因子。表土濕度因子表征了土壤表層0～2.5 cm深度范圍內含水率抗拒土壤風力侵蝕的潛在能力。基于AMSR-E Level 2A亮溫數據計算每天的表土濕度因子。

(3)植被覆蓋度。計算方法同水力侵蝕植被覆蓋度。

(4)地表粗糙度因子。地表粗糙度因子指因植被、微地形和農田耕作導致的零風速位置的高度。在風力侵蝕區，除耕地因呈斑塊狀分布，地表粗糙度在邊界發生突變外，其他諸如草(灌)地和沙地(漠)等的邊界地表粗糙度都呈逐漸過渡特征。根據不同土地利用類型，采用李智廣等[7]在我國風力侵蝕抽樣調查方法研究中不同下墊面粗糙度的賦值。

3 結果與分析

3.1 水土流失狀況

神東礦區各時期風力侵蝕和水力侵蝕面積見表1。2005年礦區平均土壤侵蝕模數為3 322.53 t/(km2·a)，水土流失面積943.36 km2，占礦區總面積的88.66%，侵蝕強度以輕、中度為主。2010年礦區平均土壤侵蝕模數為2 497.58 t/(km2·a)，水土流失面積855.47 km2，占礦區總面積的80.40%，侵蝕強度以輕、中度為主。2015年礦區平均土壤侵蝕模數為1 843.98 t/(km2·a)，水土流失面積730.66 km2，占礦區總面積的68.67%，侵蝕強度以輕度為主。2018年礦區平均土壤侵蝕模數為852.17 t/(km2·a)，水土流失面積587.61 km2，占礦區總面積的55.23%，侵蝕強度以輕度為主。

表1 神東礦區各時期不同侵蝕強度土地面積統計

3.2 水土流失動態變化

對比4期土壤侵蝕數據(表2)，2005—2010年神東礦區水土流失面積減少了87.89 km2，中度及以上侵蝕面積均有所減少，其中：中度侵蝕面積減少最明顯，減少了89.75 km2；其次為強烈侵蝕，面積減少了64.76 km2。2010—2015年水土流失面積減少了124.81 km2，中度、強烈侵蝕面積減少明顯，分別減少了117.17、72.41 km2。2015—2018年水土流失面積減少了143.05 km2，中度侵蝕面積減少最多，減少了139.97 km2。

表2 神東礦區各時期不同侵蝕強度土地面積變化統計

總的來看，2005—2018年神東礦區水土流失面積明顯減少，共減少了355.75 km2，中度侵蝕面積減少最多，其次為強烈侵蝕，分別減少了346.89、170.96 km2。對比分析礦區不同時期水土流失面積動態變化可知，2005—2018年礦區水土流失狀況為高強度侵蝕向低強度侵蝕轉移，水土流失面積明顯減少，水土流失強度明顯降低，生態環境顯著改善。

4 結 語

2005—2018年，神東礦區土壤侵蝕模數由3 322.53 t/(km2·a)降低至852.17 t/(km2·a)，水土流失面積減少了355.75 km2，相當于礦區總面積的1/3還多。取得如此成效主要得益于礦區始終堅持“生態礦區、環保礦井、清潔煤炭”的建設理念。礦區不僅實施了傳統水土保持工程措施及生態水保林、經濟林、灌草等林草植被措施，還創造性地利用柳桿障蔽生態錨固坡等技術防治水土流失。截至2018年底，礦區已累計治理水土流失面積339 km2，水土流失和土地沙化得到有效控制，植被覆蓋度大幅度提升，由開發初的3%～11%提高到60%以上[8]，生態環境得到顯著改善。

本研究采用資料收集、遙感監測、模型計算(風力侵蝕、水力侵蝕模型)等方法，定量推算了神東礦區2005、2010、2015、2018年4期土壤侵蝕模數，并對不同侵蝕強度土地面積進行了對比分析，有效評價了礦區土壤侵蝕狀況。但是，神東礦區侵蝕營力存在水力、風力及重力侵蝕疊加情況，重力侵蝕主要是井工采煤引起的滑坡和沉陷，如何率定重力侵蝕模數，更準確地掌握礦區水土流失狀況，還需進一步研究。