利用GRACE反演昆士蘭州地下水儲(chǔ)量變化

佴晉，魏長(zhǎng)壽，滕雷，趙曰耀

( 1. 山東科技大學(xué) 測(cè)繪與空間信息學(xué)院, 山東 青島 266590；2. 內(nèi)蒙古科技大學(xué) 礦業(yè)與煤炭學(xué)院, 內(nèi)蒙古 包頭 014010 )

0 引 言

水資源對(duì)于整個(gè)生態(tài)系統(tǒng)有著不可替代的作用，其中地下水由于其水質(zhì)好、分布廣、易獲取等特點(diǎn)被廣泛應(yīng)用于日常生活、農(nóng)業(yè)灌溉和工業(yè)生產(chǎn)等領(lǐng)域. 傳統(tǒng)的監(jiān)測(cè)方法只能在局部空間尺度上進(jìn)行有效監(jiān)測(cè)，卻難以在較大空間尺度上實(shí)現(xiàn)地下水位變化監(jiān)測(cè)，不利于對(duì)地下水資源進(jìn)行監(jiān)測(cè)與管理. 有效監(jiān)測(cè)地下水資源的分布和其動(dòng)態(tài)變化并對(duì)其進(jìn)行分析，對(duì)當(dāng)?shù)厣鷳B(tài)環(huán)境與用水管理政策制定具有重要現(xiàn)實(shí)意義，也是應(yīng)對(duì)地下水資源匱乏的必要手段.

重力反演與氣候?qū)嶒?yàn)衛(wèi)星(GRACE)計(jì)劃的實(shí)施[1]為水資源變化監(jiān)測(cè)提供了新的途徑，能監(jiān)測(cè)較大空間尺度上的重力場(chǎng)時(shí)變信息，在短時(shí)間尺度內(nèi)一般認(rèn)為地球重力場(chǎng)變化主要是由于水的質(zhì)量變化引起的. 目前GRACE時(shí)變重力場(chǎng)已被廣泛應(yīng)用于重大空間尺度上陸地水儲(chǔ)量變化研究中. 如LANDERER等[2]在全球范圍內(nèi)，結(jié)合GRACE數(shù)據(jù)和水文學(xué)模型擬合數(shù)據(jù)進(jìn)行比較分析，證實(shí)了GRACE數(shù)據(jù)可以用于陸地水儲(chǔ)量變化監(jiān)測(cè). 此外，GRACE衛(wèi)星提供的月重力場(chǎng)模型結(jié)合陸地水文資料或水文模型也可應(yīng)用于地下水儲(chǔ)量變化的監(jiān)測(cè). 國(guó)內(nèi)外諸多學(xué)者利用GRACE衛(wèi)星月重力場(chǎng)數(shù)據(jù)研究地下水儲(chǔ)量變化取得顯著成果，如CHEN等[3]利用GRACE結(jié)合陸地地面模型分析發(fā)現(xiàn)澳大利亞維多利亞州的地下水一直呈現(xiàn)下降狀態(tài)且降雨異常和年地下水儲(chǔ)量密切相關(guān)；周志才等[4]利用GRACE時(shí)變重力場(chǎng)模型結(jié)合地下水觀測(cè)井?dāng)?shù)據(jù)發(fā)現(xiàn)2003—2009年淮河流域地下水儲(chǔ)量以約0.28 cm/a的趨勢(shì)遞增；謝小偉等[5]利用GRACE衛(wèi)星重力場(chǎng)數(shù)據(jù)和全球陸地?cái)?shù)據(jù)同化系統(tǒng)(GLDAS)水文模型反演陜甘晉高原2005—2014年地下水儲(chǔ)量變化并對(duì)其進(jìn)行時(shí)空分析，得出陜甘晉高原地下水儲(chǔ)量減少主要是由于人為活動(dòng)造成的.

昆士蘭州是澳大利亞的第二大州，水資源處于長(zhǎng)期短缺狀態(tài)，昆士蘭州已成為澳大利亞第一個(gè)將循環(huán)水作為飲用水的州，監(jiān)測(cè)地下水儲(chǔ)量變化對(duì)該地區(qū)的水資源管理以及水資源可持續(xù)利用具有重要意義. 本文利用2003-01—2015-12共142個(gè)月的GRACE時(shí)變重力場(chǎng)模型反演澳大利亞昆士蘭州地下水儲(chǔ)量變化，并結(jié)合正向建模方法恢復(fù)研究區(qū)域信號(hào)，總結(jié)其在空間與時(shí)間上的變化趨勢(shì). 結(jié)合地下水井實(shí)測(cè)數(shù)據(jù)進(jìn)行對(duì)比分析，研究降雨對(duì)于昆士蘭州地下水儲(chǔ)量變化的影響.

1 數(shù)據(jù)與方法

1.1 GRACE數(shù)據(jù)

本次反演采用的數(shù)據(jù)是美國(guó)德克薩斯大學(xué)空間研究中心(CSR)發(fā)布的最新版本GRACE Level-2的RL06月重力場(chǎng)模型數(shù)據(jù)，由于高階項(xiàng)信號(hào)基本被噪音覆蓋，所以采用最高階數(shù)為60階[6]. 低階項(xiàng)數(shù)據(jù)誤差較大，采用衛(wèi)星激光測(cè)距(SLR)解算的值替換月重力場(chǎng)模型中的C20項(xiàng)數(shù)據(jù)[7]，并采用Swenson估算的一階重力位系數(shù)來替代原有的一階重力位系數(shù)[8]，時(shí)變重力場(chǎng)球諧系數(shù)方法反演陸地水儲(chǔ)量變化模型可描述為[9]

式中： Δh為等效水柱高； α 為 地球平均半徑； ρa(bǔ)ve為地球平均密度； ρwat為 地球水密度；Pˉlm為規(guī)格化的締合勒讓德函數(shù)；kl為一階負(fù)荷Love數(shù)； ΔClm和 ΔSlm表示地球重力場(chǎng)球諧系數(shù)相對(duì)于基準(zhǔn)重力場(chǎng)的變化量.

由于GRACE受衛(wèi)星軌道誤差、傳感器誤差以及模型誤差等因素影響，其重力場(chǎng)模型高階項(xiàng)含有較多噪聲和南北條帶噪聲[10]. 為了削弱這兩種誤差的影響，本文采用P4M15去相關(guān)濾波和半徑為300 km的高斯濾波組合濾波的方法對(duì)研究區(qū)域的GRACE數(shù)據(jù)進(jìn)行濾波處理[11]. 經(jīng)過去相關(guān)濾波和高斯濾波后，在GRACE數(shù)據(jù)處理過程中，對(duì)其球諧系數(shù)進(jìn)行高階項(xiàng)截?cái)嗖?duì)其使用高斯平滑濾波，導(dǎo)致初始信號(hào)損失形成泄露誤差. 為了減小泄露誤差的影響，本文利用正向建模方法對(duì)信號(hào)進(jìn)行恢復(fù)，基本原理是通過多次迭代，使得數(shù)據(jù)逐漸向真值逼近[12].

1.2 水文數(shù)據(jù)

1) GLDAS水文模型

GLDAS由NASA與NCEP共同建立，包括NOAH、VIC、CLM、MOSAIC四種陸地表面模式. 本文采用2003-01—2015-12的GLDAS-Noah水文模型數(shù)據(jù)，其時(shí)間分辨率為1個(gè)月. 利用該數(shù)據(jù)的積雪數(shù)據(jù)以及土壤水(0～200 mm)數(shù)據(jù)，將其球諧展開至與GRACE相同階次，采用相同的濾波方法進(jìn)行數(shù)據(jù)處理并對(duì)其進(jìn)行去平均處理計(jì)算區(qū)域地表水儲(chǔ)量變化.

2)降雨數(shù)據(jù)

降雨數(shù)據(jù)為全球降水氣候?qū)W項(xiàng)目(GPCP)數(shù)據(jù).采用2003-01—2015-12的月降雨數(shù)據(jù)，進(jìn)行平均化處理得到昆士蘭州月平均降雨數(shù)據(jù).

3)地下水井監(jiān)測(cè)資料

地下水實(shí)測(cè)數(shù)據(jù)根據(jù)聯(lián)合國(guó)教科文組織國(guó)際水文計(jì)劃地下水門戶網(wǎng)站(https://groundwaterportal.net/project/ggmn)下載. 實(shí)測(cè)數(shù)據(jù)為昆士蘭州監(jiān)測(cè)井觀測(cè)的地下水位埋深數(shù)據(jù)，觀測(cè)時(shí)隔為1天，水位單位為m. 經(jīng)過篩選，最終從20 444口監(jiān)測(cè)井中挑選出128個(gè)數(shù)據(jù)分布較好的地下水監(jiān)測(cè)井，監(jiān)測(cè)水井分布如圖1所示，對(duì)監(jiān)測(cè)水井?dāng)?shù)據(jù)進(jìn)行去平均處理得出每月的變化值.

圖1 昆士蘭州地下水監(jiān)測(cè)井點(diǎn)位分布圖

1.3 地下水儲(chǔ)量計(jì)算

在短時(shí)間尺度內(nèi)一般認(rèn)為地球重力場(chǎng)變化主要是由于水的質(zhì)量、地表生物量和冰雪質(zhì)量的變化引起的，但研究區(qū)域內(nèi)地表生物量相對(duì)穩(wěn)定，所以在研究過程中忽略生物變化量的影響[13]. 根據(jù)GRACE月重力場(chǎng)模型反演得到的陸地水儲(chǔ)量變化是包含地表徑流、土壤水、雪水當(dāng)量、冠層水和地下水變化等變化的總變化，而地表徑流和冠層水基本處于平衡狀態(tài)，變化可以忽略不計(jì)，所以昆士蘭地區(qū)的陸地水的總變化減去土壤水和雪水當(dāng)量的變化可以得到地下水儲(chǔ)量的變化，其中土壤水和雪水當(dāng)量可以由GLDAS獲取.

式中： ΔGGW為 地下水儲(chǔ) 量 變 化量； ΔTTWS為GRACE反演得到的陸地水水儲(chǔ)量總變化量； ΔSSM和 ΔSSWE為由GLDAS水文模型得到的土壤水和雪水當(dāng)量的變化.

2 結(jié)果與分析

2.1 水儲(chǔ)量時(shí)間變化分析

圖2為昆士蘭州2003—2015年間陸地水、地表水以及地下水儲(chǔ)量變化的時(shí)間序列圖. 對(duì)比圖2中的GRACE和GLDAS時(shí)間序列可知，11 a間GRACE和GLDAS反演得出的時(shí)間序列具有明顯的季節(jié)性和周期性變化，且兩者間變化趨勢(shì)基本一致. 兩者時(shí)間序列間的相關(guān)系數(shù)高達(dá)0.832，其中每年的夏秋多表現(xiàn)為上升的趨勢(shì)，而冬季和春季多表現(xiàn)為下降的趨勢(shì). 這是因?yàn)槔ナ刻m州的夏秋兩季比較溫潤(rùn)，雨季為1～3月，水資源得到有效的補(bǔ)充；春季和冬季昆士蘭州的降雨量較少，且農(nóng)作物需要大量灌溉，導(dǎo)致水儲(chǔ)量減少.

圖2 水儲(chǔ)量變化時(shí)間序列

由圖2可知，2003—2006年地下水儲(chǔ)量緩慢下降，2006—2007年地下水儲(chǔ)量緩慢上升，2007—2013年地下水儲(chǔ)量明顯上升，2013—2016年地下水儲(chǔ)量又略微減少. 利用線性擬合方法擬合昆士蘭州地下水儲(chǔ)量變化，2003—2015年昆士蘭州的地下水儲(chǔ)量整體變化速率為1.3±0.09 cm/a，其中2007-07—2012-11地下水儲(chǔ)量變化呈最快上升趨勢(shì)，有資料顯示昆士蘭州于2007年啟動(dòng)污水凈化計(jì)劃并宣布將循環(huán)水作為飲用水. 地下水最大虧損出現(xiàn)在2004年1月，為-14.1 cm；地下水最大盈余出現(xiàn)在2012年2月，達(dá)到12.8 cm.

2.2 水儲(chǔ)量空間變化分析

圖3～5為陸地水、地表水、地下水儲(chǔ)量變化趨勢(shì)的空間分布圖，三者空間分布基本一致. 從圖中可以發(fā)現(xiàn)陸地水儲(chǔ)量的整體變化呈東增西減的趨勢(shì)，西北和西南虧損嚴(yán)重，但整體呈上升趨勢(shì)；地表水的整體變化多呈現(xiàn)明顯上升趨勢(shì)；而地下水儲(chǔ)量變化趨勢(shì)與陸地水儲(chǔ)量變化趨勢(shì)大致相同，整體都呈現(xiàn)上升趨勢(shì)，且其在空間上也表現(xiàn)出東增西減的變化趨勢(shì).

圖3 陸地水儲(chǔ)量變化趨勢(shì)

圖4 地表水儲(chǔ)量變化趨勢(shì)

由圖5可以看出，整個(gè)東部沿海地下水大多呈增長(zhǎng)趨勢(shì)，每年約以0.1～3.78 cm/a的速度增加；整個(gè)昆士蘭州中西部地區(qū)的地下水大多呈現(xiàn)虧損狀態(tài)，每年約以0.84～0.2 cm/a的速度減??；昆士蘭中部和西北部地區(qū)的地下水消耗最為嚴(yán)重，約以2～4 cm/a的速度減小. 地下水儲(chǔ)量變化的差異可能與不同地區(qū)的降水情況、人類活動(dòng)以及農(nóng)作物種植情況等因素相關(guān)[14].

圖5 地下水儲(chǔ)量變化

3 驗(yàn)證與討論

3.1 降雨數(shù)據(jù)對(duì)比

文中采用全球降水氣候?qū)W項(xiàng)目提供的GPCP數(shù)據(jù)對(duì)昆士蘭州GRACE反演的地下水進(jìn)行對(duì)比分析.對(duì)2003—2015年間的GPCP數(shù)據(jù)進(jìn)行處理得到2003—2015年年均降雨異常. 圖6(a)為2003—2015年昆士蘭州地下水儲(chǔ)量變化時(shí)間序列及同時(shí)間的月降雨數(shù)據(jù)，圖中折線為反演得到的地下水儲(chǔ)量變化時(shí)間序列，黑色柱狀為對(duì)GPCP數(shù)據(jù)處理得到的昆士蘭地區(qū)的月降雨數(shù)據(jù). 由圖6(a)可以發(fā)現(xiàn)降雨量的最大值多出現(xiàn)在1～3月，而地下水的峰值相對(duì)降水具有3～5個(gè)月的延遲，這是由于降水下滲對(duì)地下水進(jìn)行補(bǔ)給需要一定的時(shí)間. 在每年的12月至次年3月降雨較多，對(duì)應(yīng)的地下水呈上升狀態(tài)；每年的6～10月降雨較少，地下水儲(chǔ)量整體呈下降趨勢(shì)；2003-06—09地下水出現(xiàn)驟降現(xiàn)象，在此期間降雨也急劇減少；2009年的1月和2月由于降雨突然增多導(dǎo)致短時(shí)間內(nèi)地下水發(fā)生快速上升的現(xiàn)象，已有資料顯示2009年2月昆士蘭州由于連續(xù)的強(qiáng)降雨，昆士蘭州遭到洪水侵襲；同樣在2011年的雨季，由于強(qiáng)降雨影響導(dǎo)致昆士蘭州發(fā)生洪災(zāi)，同期的地下水儲(chǔ)量呈快速上漲的趨勢(shì). 圖6(b)為2003—2015年年均降雨異常，由圖6(b)可以發(fā)現(xiàn)2003—2006年、2012—2015年年降水量低于年均值，與昆士蘭州干旱較為吻合.2006—2012年年降水偏多明顯高于年均值，特別是2010年，昆士蘭州年降雨異常值達(dá)到最大的393 mm，地下水立馬呈現(xiàn)快速上升的趨勢(shì)；從2006年起降雨開始增多，地下水儲(chǔ)量下降的趨勢(shì)得到有效緩解；2007—2011年地下水呈快速增長(zhǎng)趨勢(shì)，同期降雨也較常年偏多. 由此可見，昆士蘭州地下水儲(chǔ)量變化趨勢(shì)與該地區(qū)降雨量有較好的一致性，降雨是引起地下水變化的主要原因.

圖6 2003—2015年間昆士蘭州降雨數(shù)據(jù)

3.2 水井實(shí)測(cè)數(shù)據(jù)對(duì)比

昆士蘭地區(qū)地下水位監(jiān)測(cè)井位置分布如圖1所示. 對(duì)挑選出的128口地下水位監(jiān)測(cè)井進(jìn)行去平均處理得到2003—2015年昆士蘭州地下水位變化的時(shí)間序列，對(duì)GRACE反演得到的地下水儲(chǔ)量變化進(jìn)行驗(yàn)證. 若要將二者進(jìn)行等效對(duì)比，需要將由監(jiān)測(cè)井得到的地下水位乘以研究區(qū)域的給水度[15]，由于資料的缺失并未獲得研究區(qū)域的相關(guān)資料，因此僅對(duì)二者的變化趨勢(shì)進(jìn)行分析. 如圖7所示，兩者變化趨勢(shì)基本一致，反演結(jié)果與地下水位數(shù)據(jù)間的相關(guān)系數(shù)高達(dá)0.841. 其中地下水位在2003—2007年呈下降趨勢(shì)，2007—2011年表現(xiàn)出明顯的上升趨勢(shì)，與GRACE反演結(jié)果相符合，驗(yàn)證了GRACE反演結(jié)果的正確性.

圖7 地下水儲(chǔ)量變化與地下水位變化

4 結(jié) 論

本文采用2003—2015年共13年的GRACE數(shù)據(jù)反演了昆士蘭州的地下水儲(chǔ)量變化，利用正向建模法恢復(fù)研究區(qū)域信號(hào)，并與GPCP降雨數(shù)據(jù)以及地下水井實(shí)測(cè)數(shù)據(jù)進(jìn)行綜合比較分析，研究結(jié)論如下：

1)時(shí)間上，昆士蘭州地下水在2003—2015年大約以1.3±0.09 cm/a的速率增加，且地下水儲(chǔ)量呈現(xiàn)出明顯的季節(jié)性變化，夏秋兩季降雨增多地下水呈上升趨勢(shì)，冬季和春季降雨減少地下水則呈現(xiàn)下降趨勢(shì).

2)空間上，地下水儲(chǔ)量呈現(xiàn)東增西減的顯著空間差異，其中昆士蘭州東北部由于雨水充沛導(dǎo)致地下水上升速度最快，以3.78 cm/a的速度增長(zhǎng)；中部和西北部地區(qū)的地下水消耗最為嚴(yán)重，約以2～4 cm/a的速度減小.

3)水井水位變化與反演得到的地下水儲(chǔ)量變化趨勢(shì)基本一致，初步證明了反演結(jié)果的正確性. 與GPCP降雨數(shù)據(jù)對(duì)比發(fā)現(xiàn)，地下水儲(chǔ)量變化與降雨密切相關(guān)，且年地下水儲(chǔ)量變化與降雨異常密切相關(guān)，降雨是引起地下水變化的主要因素；2007—2011年降雨較常年偏多，地下水呈快速增長(zhǎng)趨勢(shì)；2012—2015年降雨量較往年偏少，因此，2012—2015年地下水儲(chǔ)量呈下降趨勢(shì).