基于多源數據的科爾沁沙地陸地水及地下水儲量變化研究

趙珍珍， 馮建迪

(山東理工大學 建筑工程學院， 山東 淄博 255000)

水是生命的起源，更是人類生存和發展的基本條件之一。作為一種基礎資源，水被廣泛應用于工業、農業、生活、娛樂等人類活動。對于整個生態系統，水資源的地位舉足輕重，影響著全球氣候演變、全球水循環和生物多樣性等基本的生態因素[1-2]。根據不同的空間分布，水資源一般分為兩大類：地表水和地下水。地下水資源具有水質好、蓄水穩定、分布廣泛、易于開發等特點，是重要的淡水資源之一。在人們的生產生活中，地下水也擔任了重要的角色。很多地方，地下水可被直接抽取飲用，作為農村和城市居民生活用水。在農業生產中，地下水是優質的灌溉水源，被大量抽取使用。另外，由于其獲取途徑方便、成本低等因素，地下水也被廣泛應用于工業和礦產開發生產中。

然而，作為一種極其寶貴的淡水資源，地下水正面臨著污染和被過量開采的現狀。目前，水資源匱乏已發展成為一個日益嚴重的全球性生態問題，特別是在干旱和半干旱地區尤為突出[3-5]。根據2012—2016年世界經濟論壇年度全球風險調查，水資源危機是全球十大風險之一，排在前三位。

科爾沁沙地位于中國北方，屬于半干旱和半濕潤的過渡性氣候。年度降水量變化很大，近年來呈現波動下降的態勢。地表水資源日益減少，河流幾近干涸。農業、工業、畜牧業及城市生活用水越來越依賴于地下水。因長期超量開采地下水，地下水位已嚴重下降。

水資源的實時動態監測是實現水資源的保護和管理的重要參考之一。水資源的監測手段有很多。傳統的方法有河流/水庫監測站、地下水位監測井等。傳統的監測站和監測井，易受地理環境影響，測站數量相對稀少、總體分布不均勻，監測能力有限。另外，有些水文數據屬于內部參考資料，出于信息安全等因素，不對外公布。以上因素很大程度上限制了對水資源的監測、管理、保護和研究。

始于2002年的GRACE(gravity recovery and climate experiment，重力恢復和氣候探測任務)衛星重力任務為水資源的監測研究開辟了一條嶄新的道路。GRACE提供的月重力場模型，可用來監測地球表面淺層物質質量遷移。對陸地水儲量變化的監測是GRACE的重要應用之一。另外，結合陸地水文資料或水文模型，GRACE也可用于地下水儲量變化的監測[6]。目前，利用GRACE重力數據監測陸地水儲量變化(地下水儲量變化)的研究很多，研究區域大多是水量豐富的大河流域以及水資源匱乏的干旱和半干旱區。例如，密西西比河流域[7-8]、尼羅河流域[9]、亞馬遜河流域[10-12]和恒河流域[13-14]等。在中國區域，華北平原[15]、新疆天山山區[16]、海河流域[17]、黑河流域[18-19]和長江流域[19]等區域的水資源監測一直是一個熱門的話題。

以上研究多集中在大河流域以及部分干旱和半干旱區，很少涉及科爾沁沙地水資源的研究。為了實現對科爾沁沙地的水資源監測，本文擬聯合GRACE重力數據、水文模型、監測點水井數據及給水度信息，建立科爾沁沙地水資源監測體系，反演該地區陸地水儲量變化和地下水儲量變化，以期為當地水資源的保護和可持續利用提供參考。

1 研究區概況

科爾沁沙地是中國最大的沙地(北緯42°30′—45°30′，東經119°—124°)，平均海拔高度為178.5 m，呈現西高東低，行政區域上包括18個縣(市、旗)，總面積5.06×104km2。科爾沁沙地屬于半干旱農牧交錯區，耕地和草地之間不斷發生互相轉換，土質松散、黏連性差，為土地沙化提供了物質基礎。同時，該地區降水量少，年平均降水量為371.3 mm，而且年際降雨量波動性大，1980年的年均降雨量只有198.7 mm，1991年年均降雨量達到了577.6 mm，經常出現干旱年份，以及偶爾發生洪澇災害；年內降雨不均，降雨主要集中在夏季，占總降水量的69.6%，為豐水期，而春季和秋季總降水比例不足30%，冬季降水量只占百分之一左右。降雨量直接影響該地區的土壤含水量，以及地表水面積和地下水水位變化。

同時，科爾沁沙地人口數量穩步上升，從1980年的2.37×106人，增長到2016年的3.13×106人，35 a的時間里人口增加了32.1%。伴隨著人口的快速增長，為滿足生活需求，耕地面積和牲畜存欄量也在不斷增加。科爾沁沙地經濟水平不斷提高，GDP從1999年的1.35×1010元，經過16年的發展，到2016年末達到了1.95×1011元，增長了14倍；農牧民人均純收入從1999年的2 028元，增長到2014年的9 932元，增加到原來的近5倍；產業結構由最初的以第一產業為主，逐漸轉變為第二產業占據主要地位，第三產業次之，第一產業的比重最小。科爾沁沙地特殊的地理位置、氣候特征和人為因素，共同決定了該地區地表水和地下水的不斷減少、沙地范圍大并呈現出波動性變化等生態環境問題。

2 反演方法

利用GRACE時變重力場反演地表面質量變化的計算模型可以表示如下[20]：

〔ΔClmcos(mλ)+ΔSlmsin(mλ)〕

(1)

由于GRACE任務的飛行特點、共軌星間跟蹤觀測模式、軌道誤差、大氣和海洋模型誤差等因素，由原始觀測數據解算的地球重力場模型球諧系數存在條帶誤差和高頻誤差。為了削弱這些誤差的影響，本文采用去相關濾波與高斯濾波(平滑半徑為300 km)的組合濾波方法，經過去相關濾波和高斯濾波以后，公式(1)變為：

(2)

式中：Wl——l階高斯平滑核函數。該函數可以通過一組遞推公式表示：

W0=1

(3)

(4)

(5)

其中，

(6)

式中：r——“半寬”常數，即高斯平滑半徑；R——地球的平均半徑。

經過去相關濾波和高斯平滑濾波數據處理，GRACE數據中的干擾因素已被剔除，剩余的信息為地球表面質量遷移帶來的重力異常。在較短的時間尺度內，地球表面質量遷移一般是指陸地水儲量變化。陸地水儲量變化可進一步轉換為等效水高的形式，公式為

(7)

式中：ρ——水的密度(1 000 kg/m3)；H——等效水高。

本文采用Landerer和Swenson(2012)提出的尺度因子法，其基本原理是利用最小二乘技術，以濾波前GLDAS-NOAH得到的水儲量信號與濾波后的水儲量信號之差的最小為目標，求解得到尺度因子[21]為0.7。一個區域的陸地水儲量變化包括：土壤含水量、地表徑流、積雪和冰川以及地下水的變化，可以用公式(8)表示[22]。將土壤含水量、地表徑流、積雪和冰川的變化從GRACE計算的陸地水儲量變化中扣除，即可得到地下水儲量變化，從而實現利用GRACE重力數據監測地下水儲量的目的[23]。

ΔTWS=ΔSMS+ΔSMES+ΔRESS+ΔGWS

(8)

其中：ΔTWS——陸地水儲量變化； ΔSMS——土壤含水量變化； ΔSWES——積雪和冰川變化； ΔRESS——地表徑流變化； ΔGWS——地下水儲量變化。

鑒于科爾沁沙地的實際情況和已有的資料，在計算科爾沁沙地的地下水儲量變化時，本文不考慮積雪和冰川以及地表徑流的變化，僅計算了占地表水主導地位的土壤含水量變化。那么，公式(8)簡化為：

ΔTWS=ΔSMS+ΔGWS

(9)

該簡化的水資源平衡公式已被應用于多個地區的地下水儲量變化的研究中[18,24-25]。

3 結果與分析

3.1 基于GRACE數據反演陸地水儲量變化

本文采用德克薩斯大學CSR(center of space research)空間研究中心提供的200307—201012 GRACE Release-05 Level-2數據，對其進行去相關濾波與高斯濾波(平滑半徑為300 km)的組合濾波方法處理后，反演了科爾沁沙地陸地水儲量變化時間序列(圖1)。

圖1 基于GRACE數據和CPC水文模型的科爾沁沙地陸地水儲量變化(200307-201012)

由GRACE時變重力場數據反演得到科爾沁沙地200307—201012的陸地水儲量下降速率為-13.2±2.6 mm/a(下降速率即為等效水高的變化率，下同)。同時，本文采用全球水文模型CPC (climate prediction center)進行對比分析，該模型是根據全球觀測的降水分布而建立，采用的數據包括全球降水數據、溫度、潮濕度、水平風速、地表大氣壓和太陽輻射數據等。CPC提供的產品包括地表積雪分布及其厚度、地表以下4層的土壤含水量、土壤溫度等。作為驗證數據，CPC水文模型得到的陸地水儲量變化與GRACE的結果進行比較。從圖1可以看出：①GRACE數據和CPC水文模型反演的科爾沁沙地陸地水儲量結果整體上均呈現持續減小趨勢，具有較好的一致性；②二者具有相似的季節變化特性，相關系數為0.72。需要注意的是，CPC水文模型計算的水儲量變化振幅略小于GRACE的結果。這是因為GRACE的水儲量變化是所有陸地水變化的總和，相對于CPC的結果來說，還包含了其他的地表質量變化(如地表生物總含水量變化和固體地球質量變化)。

3.2 基于GRACE和GLDAS水文模型反演地下水儲量變化

本文構建了兩種科爾沁沙地地下水儲量變化反演方法：一是聯合GRACE數據和水文模型反演地下水儲量變化；二是利用監測點水井數據和給水度數據確定地下水儲量變化。最后，將這兩種方法計算的科爾沁沙地地下水儲量變化結果進行比較，相互驗證。

本文利用GLDAS水文模型計算了科爾沁沙地的土壤含水量變化(如圖2所示)。從圖2中可以看出，科爾沁沙地的土壤含水量存在明顯的季節變化特性，在夏季達到最大值，在春冬季達到最小值。科爾沁沙地位于中國北方，屬于溫帶大陸性氣候特性。夏季炎熱，受東南季風的影響，降水較多，占全年降水量的70%～80%；春冬季多偏北風和西北風，空氣干燥，降水稀少。

圖2 基于GLDAS水文模型的科爾沁沙地土壤含水量變化(200307-201012)

圖3記錄了2003年到2010年科爾沁沙地的降雨量月均值序列。通過對比圖2和圖3可知，科爾沁沙地的土壤含水量與降雨量是密切相關的，二者存在明顯的正相關關系。根據GRACE數據計算的陸地水儲量變化減去GLDAS水文模型計算的土壤含水量變化，得到了科爾沁沙地的地下水儲量變化(如圖4所示)。從圖4中可以看出，2003年7月到2010年12月期間，科爾沁沙地的地下水儲量總體上呈現下降趨勢，總體下降速率為-13.5±1.9 mm/a。其中，2003年7月到2005年7月，地下水儲量是增加的，速率為30.2±4.1 mm/a；從2005年1月開始，地下水儲量急劇下降，到2009年12月，地下水儲量下降速率為-26.7±2.7 mm/a；2009年7月到2010年12月，地下水儲量依然處于下降態勢，但下降速率明顯減緩，為-6.0±1.9 mm/a。

圖3 2003-2010年科爾沁沙地月平均降水量變化

圖4 基于GRACE數據和GLDAS水文模型的200307-201012科爾沁沙地地下水儲量變化

3.3 基于監測點水井數據和給水度反演地下水儲量變化

3.3.1 基于監測點水井數據反演地下水位變化 本文使用了《中國地質環境監測地下水位年鑒》記錄的地下水數據，研究了科爾沁沙地的地下水位變化。根據年鑒記錄，內蒙古自治區通遼市科爾沁沙地區域共分布著6個國家級地下水監測點。將每個監測點的2005—2010年地下水月平均水位提取出來，計算得到6個監測點的地下水位變化序列。再將6個監測點的地下水位變化序列求平均值，得到了科爾沁沙地2005—2012年的月平均地下水變化序列(如圖5所示)。從圖5可以看出，科爾沁沙地地下水位一直處于下降的趨勢。

圖5 2005-2010年科爾沁沙地地下水位月平均變化量

3.3.2 基于監測點水井數據和給水度反演地下水儲量變化 本文根據研究區域內估算的潛水含水層的給水度結果和相關的給水度室內試驗數據，利用Kringing模型算法，估計了科爾沁沙地的給水度數值，估值的均方根誤差優于0.001。考慮到區域均值的不確定性，科爾沁沙地給水度的取值區間為：(0.09～0.1)。將監測點水井數據得到的地下水位變化乘以給水度，即可得到地下水儲量變化。圖6給出了200501—201012科爾沁沙地的地下水儲量變化，包含了給水度取0.09和0.10的計算結果。根據地下水消耗的速率，將科爾沁沙地的地下水儲量變化分為兩個階段：200501—200912，地下水下降速度較快，下降速率為-29.1 mm/a(給水度取0.1)和-26.04 mm/a(給水度取0.09)；200907—201012，地下水下降速度有所減緩，下降速率為-6.1 mm/a給水度取0.1)和-5.5 mm/a(給水度取0.09)。

3.4 科爾沁沙地地下水儲量變化結果比較分析

基于GRACE數據和GLDAS水文模型得到的地下水儲量變化，與監測點水井數據得到的地下水儲量變化進行比較(如表1所示)。由于《中國地質環境監測地下水位年鑒》記錄的地下水信息始于2005年，所以對比時間段選取200501—201012。從表1可以看出，這兩種方法得到的地下水儲量變化在此期間是一致的。兩種方法相互驗證，增加了科爾沁沙地地下水儲量變化結果的可靠性。

表1 不同階段GRACE和水井數據監測的地下水變化結果比較

綜合分析GRACE數據、GLDAS水文模型、監測點水井數據和給水度信息得到的科爾沁沙地地下水儲量變化可知，2003—2010年科爾沁沙地的地下水儲量呈現消耗的趨勢，消耗速率為-13.5±1.9 mm/a，其中200501—200912月，地下水儲量消耗最快，2009年7月以后，地下水儲量有所恢復，下降速率減緩。

3.5 科爾沁沙地水儲量變化成因分析

結合氣象數據和當地農業模式分析科爾沁沙地水儲量變化的成因。首先，分析了科爾沁沙地降水量與陸地水儲量變化的關系。圖7給出了2003—2010年科爾沁沙地的年降水量數據，其中，年降水量的最大值出現在2010年，最小值出現在2009年。根據相關研究報道，2009年中國北方發生了嚴重的干旱[26]。同時結合圖2可以看出，由GRACE數據得到的2009年科爾沁沙地陸地水儲量下降最為明顯。隨著降水量的增加，從2010年開始陸地水儲量有所恢復。此外，2006—2009年的降水量均處于較低的水平，平均約為310 mm，是導致科爾沁沙地陸地水儲量多年下降態勢的主要原因之一。

然后，從科爾沁沙地的農業生產模式方面分析了當地地下水儲量消耗的原因。科爾沁沙地的農作物種類繁多，包括玉米、向日葵、果樹、各種蔬菜等。其中，玉米為當地的主要糧食作物，被大面積廣泛種植。研究表明，當地玉米的種植面積占農作物總種植面積的70%[27-28]。玉米是一年生單子葉草本植物，秸稈很高，枝葉茂盛，生長過程中需水量很大。

圖7 2003-2010年科爾沁沙地年降水量變化趨勢

在科爾沁沙地，年平均降水量僅為350 mm，發生在玉米生長周期內的降水占全年降水的70%左右。玉米生長需水量很大，僅依靠自然降雨是不能滿足玉米對水分需求的。本文引用張銅會等[27]的研究結果(如表2所示)，展示了科爾沁沙地玉米耗水量及用水來源[29]。玉米的生長周期分為4個階段：播種—拔節、拔節—抽穗、抽穗—灌漿、灌漿—成熟。從表2可以看出，在不同的生長階段，玉米對水的需求是不同的。此外，從該表中還可以看出，玉米總的耗水量已超過當地降雨量，因此，需要抽取地下水進行灌溉。灌溉水量占到了玉米總耗水量的37%左右。特別是在播種—拔節和拔節—抽穗階段，灌溉水量所占比例高達60%左右。

表2 科爾沁沙地玉米不同生長階段的需水特點[27]

結合圖7可知，2006—2009年間，科爾沁沙地的降水量嚴重不足，低于平均值。為了抵御干旱，確保農作物的產量，更多的地下水被抽取到地面，用于灌溉。與此同時，地表蒸發量也隨之更加。水資源大量被消耗，地下水得不到足夠的補給，導致地下水位不斷下降。

圖8是利用2005—2010年《中國地質環境監測地下水位年鑒》地下水位數據計算的科爾沁沙地地下水位季節變化。從圖中可以看出，每年的地下水位從3月份開始下降，一直到8月份下降至最低點。8月份以后，地下水得以補給，水位上升。根據前文描述的當地農業模式和降水特點，地下水消耗和恢復的過程可以理解為：春季降雨稀少，當地居民大量抽取地下水以對抗春旱，導致地下水從3月份開始下降；夏季農作物生長的過程中，水的需求量是巨大的，降雨量不能滿足農作物的生長需求，因此大量地下水被不間斷的抽取到地面，導致地下水位持續性下降；秋季作物開始成熟，對水的需求也在逐漸減少，地下水得以恢復；秋季以后直到冬季結束，該地區無農作物生長，地下水在冬末春初時，恢復到最大值。另外，從圖8中還可以發現，雖然地下水存在消耗和恢復的周期性，但是每年的地下水位都不能恢復到最初的水位，地下水位在不斷的下降。

圖8 基于監測點水井數據的2005-2010年科爾沁沙地地下水位季節變化

4 討論與結論

本文采用200307—201012的GRACE時變重力場模型Release-05 Level-2，采用去相關濾波與高斯平滑濾波(平滑半徑300 km)結合的濾波方法，反演了科爾沁沙地陸地水儲量變化，將其與CPC水文模型反演結果進行了對比分析；結合GLDAS水文模型計算了土壤含水量變化，并分析了科爾沁沙地地下水處理變化的時空特征；利用監測點水井數據和給水度信息確定地下水儲量變化，并與GRACE反演結果進行了初步驗證。研究結果表明：①科爾沁沙地陸地水儲量整體呈現減少趨勢，減少速率為-13.2±2.6 mm/a，與CPC水文模型反演結果具有較好的一致性；②科爾沁沙地地下水儲量整體上呈現下降趨勢，下降速率為-13.5±1.9 mm/a，與監測點水井數據反演結果具有較好的一致性。

結合氣象數據和當地農業模式分析了科爾沁沙地水儲量變化的成因：首先，降水量變化與陸地水儲量變化有直接關系，2003-2010年期間，年降水量最小值出現在2009年，此時對應的科爾沁沙地的陸地水儲量下降最為明顯，年降水量最大值出現在2010年，此時期的陸地水儲量也開始有所恢復；第二，根據科爾沁沙地玉米不同生長階段需水量以及對應階段降水量的特點，發現干旱和農業灌溉是導致地表水減少和地下水位降低的重要原因之一。

基于GRACE衛星時變重力場數據反演長時間序列、大尺度的陸地水儲量變化和地下水儲量變化，該方法有效地彌補了傳統水文監測中因觀測站點分布不均勻、觀測范圍有限等因素導致資料獲取不充分、數據不均勻的缺陷。隨著后續衛星重力探測計劃的相繼實施、時變重力場反演方法的改進和水文模型的日趨完善，將進一步提高衛星重力觀測結果中地球物理信號的時空分辨率，以及信號分離和提取的精度，為研究固體地球物質遷移以及精化和驗證全球和局部水文模型提供可靠的地球物理信息。