朱悅璐，熊俊飛，鄧 煒

(南昌工程學(xué)院 水利與生態(tài)工程學(xué)院，江西 南昌 330099)

土壤含水率是水文、生態(tài)、農(nóng)業(yè)以及巖土工程領(lǐng)域的重要參數(shù)[1-3]，不同學(xué)科對該指標(biāo)研究各有側(cè)重。例如農(nóng)業(yè)工程中，研究人員重點(diǎn)關(guān)注土壤水分變化、判斷作物長勢、評價農(nóng)業(yè)干旱等問題[4]；巖土工程領(lǐng)域則更多偏向于土中水分變化導(dǎo)致的工程災(zāi)害，如地基沉降、邊坡穩(wěn)定等問題[5-6]。可見土壤含水率的測定作為某一研究的“中間環(huán)節(jié)”，對后續(xù)研究結(jié)果的正確性至關(guān)重要。目前測定土壤含水率的方法可分為兩類：其一為單點(diǎn)小尺度測定，其二為大尺度衛(wèi)星遙感測量技術(shù)[7]。

單點(diǎn)小尺度測定法[8-9]對于某點(diǎn)土壤含水率可以較為準(zhǔn)確的測定，但無法對某一區(qū)域內(nèi)土壤含水率的變化進(jìn)行連續(xù)監(jiān)測，只能依靠一點(diǎn)或幾個具有代表性的測點(diǎn)表征整個區(qū)域含水率的變化，若要獲得更精確的數(shù)據(jù)，則需要連續(xù)密集野外作業(yè)，這會導(dǎo)致采樣成本過高且數(shù)據(jù)同步性較差。

而基于衛(wèi)星數(shù)據(jù)反演的現(xiàn)代遙感技術(shù)在一定程度上可以解決上述問題，該方案核心之一是對遙感數(shù)據(jù)進(jìn)行光譜分析。這類研究多集中在通過光譜建立土壤含水率回歸模型[10]、某類土中特殊光譜與含水率對應(yīng)關(guān)系[11]、基于光譜建立土壤含水率反演模型[12]等方面。目前，國內(nèi)已有研究團(tuán)隊(duì)?wèi)?yīng)用無人機(jī)高空攝影對試驗(yàn)區(qū)含水量進(jìn)行監(jiān)測，為防治土壤鹽堿化、農(nóng)業(yè)灌溉方案制定等提供信息指導(dǎo)[13-14]，該技術(shù)本質(zhì)上也是一種小尺度的遙感反演。但應(yīng)該注意到，上述遙感反演分析結(jié)果并不具有普適性，事實(shí)上，目前尚無一種土壤含水率測定手段可以同時滿足不同研究對象。進(jìn)一步，在更大尺度中，例如流域尺度評估農(nóng)業(yè)干旱風(fēng)險問題[15]，受傳感器限制，遙感衛(wèi)星只能在某些特定平坦空曠場合獲得可靠數(shù)據(jù)[16]，因此難于精確指導(dǎo)農(nóng)業(yè)生產(chǎn)，而低空雷達(dá)掃描范圍又不足以支持流域尺度，若要擴(kuò)大監(jiān)測范圍則會提高相應(yīng)成本[17]。此外由于衛(wèi)星用途、運(yùn)行軌道等原因，在現(xiàn)階段，僅為某一課題，通過衛(wèi)星遙感獲取特定區(qū)域長系列逐日土壤含水率數(shù)據(jù)，在我國大多高?；蜓芯繖C(jī)構(gòu)仍不現(xiàn)實(shí)。

基于上述討論，本文擬采用大尺度陸面水文模型VIC為工具，通過氣象數(shù)據(jù)及下墊面數(shù)據(jù)建立流域尺度水文模型，模擬生成研究區(qū)各網(wǎng)格內(nèi)逐日土壤含水量，在此基礎(chǔ)上，從模型內(nèi)和模型外分別對模擬結(jié)果進(jìn)行檢驗(yàn)，以確保滿足工程精度。該方案區(qū)別于傳統(tǒng)方法，第一可以在基于實(shí)測資料基礎(chǔ)上，模擬整體區(qū)域土壤含水率，數(shù)據(jù)結(jié)構(gòu)在空間上具有流域尺度整體性；第二由于構(gòu)建水文模型的氣象資料更容易獲取，且下墊面數(shù)據(jù)在一定時期保持時空一致性，因此可以通過較小的成本生成連續(xù)土壤含水率數(shù)據(jù)。本研究所用方法，可以在流域尺度構(gòu)建具有時空連續(xù)分布特性的土壤含水率數(shù)據(jù)庫，為相關(guān)工程提供數(shù)據(jù)基礎(chǔ)。

1 水文模型構(gòu)建

1.1 研究流域

本文選擇渭河流域作為研究區(qū)，地理坐標(biāo)為東經(jīng)106°20′～110°37′，北緯33°40′～37°18′，其地勢東西高于南北。按照0.5°×0.5°經(jīng)緯度將流域劃分為75個網(wǎng)格，全流域包括21個氣象站點(diǎn)和5個水文站點(diǎn)，如圖1所示。

1.2 氣象數(shù)據(jù)

氣象數(shù)據(jù)輸入時段為1980—2015年，包括日最高氣溫、最低氣溫、降水和風(fēng)速。分布式水文模型運(yùn)算中，各個網(wǎng)格的數(shù)據(jù)相互獨(dú)立，考慮到距離作為降水空間分布的最重要因素，本研究采用距離權(quán)重方案，將21個氣象站點(diǎn)數(shù)據(jù)插值到各網(wǎng)格中心，其公式為

(1)

式中n為格網(wǎng)內(nèi)站點(diǎn)數(shù)；di為第i個站點(diǎn)至每個網(wǎng)格中心點(diǎn)的距離；Pi為第i個站點(diǎn)的實(shí)測氣象數(shù)據(jù)，即為每個網(wǎng)格插值后所生成的氣象數(shù)據(jù)。

1.3 植被參數(shù)

我國的地類代碼為三級分類，共有6大類67個小類。該分類方法與VIC模型采用的Maryland全球1 km陸面覆蓋類型分類不同，因此在構(gòu)建模型下邊界條件時，本研究具體的處理方法是將土地利用類型圖中多個子項(xiàng)合并成VIC模型中一種植被分類，使其與Maryland發(fā)展的陸面覆蓋類型參數(shù)庫相同。

應(yīng)用Arcgis提取圖2中每個網(wǎng)格的信息，并對覆蓋類型進(jìn)行統(tǒng)計(jì)，其中i為陸面類型分類號，Vi表示第i類植被在網(wǎng)格內(nèi)所占的面積，滿足公式

(2)

應(yīng)用式(2)對研究區(qū)75個網(wǎng)格內(nèi)植被信息逐一計(jì)算，進(jìn)行數(shù)字化處理，每種信息對應(yīng)不同的顏色和相應(yīng)的數(shù)據(jù)以供模型讀??；其制作過程以80s植被分布為例，如圖3所示。

上述文件中，不同植被類型對應(yīng)不同根區(qū)深度，各植被12個月的葉面積指數(shù)則通過VIC模型從植被參數(shù)庫文件調(diào)取，該套數(shù)據(jù)可以反映流域真實(shí)植被覆蓋情況。

1.4 土壤參數(shù)

本研究選取的VIC模型網(wǎng)格土壤參數(shù)來源于NOAA(The National Oceanic and Atmospheric Administration)提供的全球5′土壤質(zhì)地分類描述，共有12類，分別為砂土、壤質(zhì)沙土、沙壤土、粉質(zhì)壤土、粉土、壤土、沙質(zhì)粘土壤、粉質(zhì)粘土壤、粘土壤、沙質(zhì)粘土、粉質(zhì)粘土和粘土。

圖1 水文模型網(wǎng)格劃分

圖2 研究區(qū)植被參數(shù)

圖3 構(gòu)建數(shù)字化植被類型圖

在模型運(yùn)行過程中，各土壤類型分別包括53種參數(shù)，整體可分為兩類，一類是與土壤特性有關(guān)的參數(shù)，共有47個，另一類是需要對比實(shí)測輸出結(jié)果進(jìn)行調(diào)整的參數(shù)，共有6個。第1類參數(shù)一旦研究區(qū)域選定，當(dāng)土壤結(jié)構(gòu)已知時，即標(biāo)定為常數(shù)。第2類參數(shù)在研究區(qū)域選定后，要根據(jù)流域土壤實(shí)際情況進(jìn)行調(diào)整，調(diào)整結(jié)果的好壞，直接影響最終模擬精度，這一工作即為模型參數(shù)率定。

1.5 參數(shù)率定

參數(shù)率定是保證模擬結(jié)果正確的重要環(huán)節(jié)。VIC模型需要率定的參數(shù)有6個，分別為B，Ds，Dm，Ws，d2，d3。以總量精度誤差BLAS%、納什效率系數(shù)NES作為判斷函數(shù)，其中BLAS代表模擬和實(shí)測值在數(shù)量上的擬合程度，NES代表模擬值和實(shí)測值在趨勢上的擬合程度。前者絕對值越小，模型模擬精度越高，后者越接近于1，模型模擬精度越高。其表達(dá)式如下：

(3)

(4)

其中Qsim.t和Qobs.t分別為徑流序列的模擬值和觀測值。

1.6 匯流模型

采用Dag Lohmann研發(fā)的匯流程序構(gòu)建匯流模型，共有7類文件，其中流向文件是匯流模型的核心，它包括匯流方向、網(wǎng)格大小、經(jīng)緯坐標(biāo)、流域形狀等信息，其余6類為可選項(xiàng)。

該程序分別利用單位線法和圣維南方程計(jì)算坡面匯流和河道匯流，整個匯流過程假定為線性，網(wǎng)格內(nèi)的水分流入相鄰8個網(wǎng)格中的一個，這一過程由上游至下游，直至匯入流域出口斷面。

對于坡面匯流，不論任何土壤，都假定由降雨引起的快速產(chǎn)流和慢速產(chǎn)流兩部分組成，其表達(dá)式為

(5)

其中Qs(t)為快速產(chǎn)流，QF(t)為慢速產(chǎn)流，

對于河道匯流，假定由地表產(chǎn)流和地下產(chǎn)流兩部分組成，并應(yīng)用圣維南方程計(jì)算，其表達(dá)式為

(6)

式中C為波速，D為擴(kuò)散系數(shù)。參數(shù)C，D決定了水流在河道間流動的性質(zhì)，該值可通過實(shí)測獲取。

1.7 模型運(yùn)行

上述氣象數(shù)據(jù)、植被數(shù)據(jù)、下墊面數(shù)據(jù)以及匯流程序整體上構(gòu)成了VIC模型的運(yùn)行環(huán)境，利用本文準(zhǔn)備的前處理數(shù)據(jù)，即可驅(qū)動VIC模型模擬長系列日土壤含水量。

2 結(jié)果與討論

2.1 模型參數(shù)率定結(jié)果

本文選取渭河上游林家村站、中游魏家堡站、下游華縣站、以及兩大支流張家山、狀頭控制站作為率定對象，通過調(diào)整參數(shù)，以式(3)～(4)為標(biāo)準(zhǔn)對比模擬值和實(shí)測值。率定結(jié)果表明，林家村站NES=0.88，總量精度誤差BLAS=4.4%，魏家堡站NES=0.85，BLAS=0.62%，該結(jié)果可以達(dá)到一般工程要求，如圖4所示。其余各站點(diǎn)參數(shù)率定結(jié)果及指標(biāo)類似，結(jié)果見表1。

圖4 各站率定結(jié)果

表1 研究流域最終參數(shù)方案

應(yīng)用上述全流域優(yōu)化后的方案驅(qū)動VIC模型，對全流域經(jīng)行模擬，整體模擬結(jié)果NES=0.81，BLAS=6.9%，這一結(jié)果滿足工程要求，在全流域尺度下，模型通過率定。

2.2 模型參數(shù)驗(yàn)證結(jié)果

為驗(yàn)證上述參數(shù)的合理性，選擇緊接率定期后一半時期作為驗(yàn)證期，將表2各子流域參數(shù)代入VIC模型。渭河干流上中下游、兩大支流驗(yàn)證期NES、BLAS檢驗(yàn)結(jié)果分別為：(0.84，5.3%)，(0.85，6.2%)，(0.78，7.9%)，(0.86，5.0%)，(0.87，4.3%)，各流域模擬效果如圖5(a)～(e)所示。若以各小流域驗(yàn)證結(jié)果平均值計(jì)算全流域模擬結(jié)果，此時NES=0.84，BLAS=5.7%，該值比率定期效果更好，與實(shí)際經(jīng)驗(yàn)不符。

圖5 各流域驗(yàn)證期實(shí)測模擬對比

進(jìn)一步的，以華縣+狀頭站徑流之和作為渭河流域出口斷面對全流域進(jìn)行驗(yàn)證，模擬結(jié)果為：NES=0.81，BLAS=7.3%。如圖5(f)所示，這一驗(yàn)證期結(jié)果與率定期結(jié)果自洽，且滿足一般工程精度，因此認(rèn)為模型在全流域通過驗(yàn)證。

上述相同參數(shù)所得出不同模擬結(jié)果的現(xiàn)象并不矛盾，這是由于采用不同計(jì)算方案所造成的，若直接采用各個小流域模擬結(jié)果平均值作為全流域模擬結(jié)果，由于沒有考慮各小流域徑流量、面積等要素在全流域的權(quán)重，因此計(jì)算結(jié)果不如對出口斷面進(jìn)行檢驗(yàn)效果準(zhǔn)確。

2.3 土壤含水量模擬結(jié)果

模型通過率定和驗(yàn)證后，即可用本文所生成的全套參數(shù)驅(qū)動VIC模型模擬研究區(qū)土壤含水率，該含水率數(shù)據(jù)為1985—2015年各網(wǎng)格上、中、下3層土逐日土壤含水量，共12 784組，保存在輸出文件中，將這些數(shù)據(jù)單獨(dú)提取，編制為Excel表格，即形成滑坡體所在網(wǎng)格含水量數(shù)據(jù)庫。

為進(jìn)一步驗(yàn)證該數(shù)據(jù)庫的正確性，選用1981—1991年寶雞土壤監(jiān)測站的實(shí)測含水量數(shù)據(jù)進(jìn)行驗(yàn)證。該站在流域內(nèi)所處位置如圖6所示。

該監(jiān)測站在上述時段內(nèi)，每月18日實(shí)測一次，10年共有120個土壤含水量數(shù)據(jù)，以含水量數(shù)據(jù)庫中同時段月平均數(shù)據(jù)與實(shí)測數(shù)據(jù)進(jìn)行對比驗(yàn)證，其結(jié)果如圖7所示。

從圖7可以看出，網(wǎng)格模擬含水量和實(shí)測含水量在月尺度下，變化趨勢基本相同，但二者數(shù)值相差較大，網(wǎng)格的平均含水量深度在300～700 mm之間，而站點(diǎn)實(shí)測含水量數(shù)值在150～400 mm之間，模擬值在整體上高于實(shí)測值。歸因分析認(rèn)為，這種系統(tǒng)性偏差主要在于模擬土層和實(shí)測土層的最大深度不同，受現(xiàn)場檢測手段、監(jiān)測目的限制，站點(diǎn)最大監(jiān)測深度為500 mm左右，而本文構(gòu)建的水文模型3層土平均厚度約為3.5 m，二者相差7倍左右。

2.4 模擬結(jié)果討論

考慮到月尺度下網(wǎng)格平均含水量可能會掩蓋某些關(guān)鍵信息，例如極端降雨或者干旱所生成的偶然極值，因此選擇相同時間點(diǎn)的模擬值與實(shí)測值逐一對比以判斷模擬效果，如圖8所示。由圖中可以看出，相同時刻含水量模擬值和實(shí)測值分布與月平均值對比具有類似趨勢，即整體上模擬值大于實(shí)測值，但二者散點(diǎn)圖分布更為緊密，圖中已無法直接判斷有無系統(tǒng)性偏差；這表明同時刻逐一對比效果好于平均值對比。

圖7 土壤含水量實(shí)測模擬對比 圖8 土壤含水量實(shí)測模擬逐一散點(diǎn)對比

將120組“實(shí)測-模擬”含水量數(shù)據(jù)應(yīng)用SPSS軟件進(jìn)行統(tǒng)計(jì)分析，取最大值、最小值、均值、中值以及標(biāo)準(zhǔn)差等5個量作為統(tǒng)計(jì)檢驗(yàn)指標(biāo)，其結(jié)果如表2所示；從極值和均值的角度來看，模擬結(jié)果和實(shí)測結(jié)果相差不大，模擬效果較好。

表2 含水量統(tǒng)計(jì)分析表

對上述結(jié)果做進(jìn)一步相關(guān)性分析，其結(jié)果表明，模擬值序列和實(shí)測值序列相關(guān)性尚好，相關(guān)系數(shù)R2=0.749 2。各實(shí)測值散點(diǎn)與對應(yīng)模擬值散點(diǎn)差值絕對值基本在±100 mm之間，平均誤差約為15.2%。差值最大絕對值為108.1 mm，最大誤差為32%。考慮到模型的不確定性，該極端偶然誤差可以接受。如圖9所示。

圖9 土壤含水量相關(guān)性分析及差值表

為考察上述含水量誤差的來源，按照渭河流域?qū)嶋H氣候條件，將該數(shù)據(jù)按月份劃分為旱季和雨季，其中1至3月為旱季，6至8月為雨季。分段考察兩組數(shù)據(jù)的相關(guān)性，如圖10所示。分析結(jié)果表明，含水率實(shí)測值與模擬值相關(guān)系數(shù)旱季R2=0.612 5，雨季R2=0.927 9；雨季模擬結(jié)果遠(yuǎn)好于旱季，因此可以認(rèn)為，造成含水量模擬值與實(shí)測值誤差的主要原因是旱季土壤水分模擬不準(zhǔn)確。事實(shí)上大量前人研究結(jié)果表明，VIC模型在濕潤地區(qū)模擬效果要好于干旱地區(qū)，本文所述土壤含水率模擬結(jié)果亦是該結(jié)論的一個特例。

圖10 土壤含水率旱季相關(guān)性分析，雨季相關(guān)性分析

至此，本研究以上述方案生成的土壤含水量數(shù)據(jù)，在流域尺度、區(qū)域尺度和網(wǎng)格尺度上皆通過檢驗(yàn)，在精度上可以滿足工程要求。網(wǎng)格內(nèi)土壤含水量模擬值與實(shí)測值的檢驗(yàn)結(jié)果表明：(1)VIC模型模擬的土壤含水量在長時間范圍內(nèi)，變化趨勢具有代表性；(2)在具體對應(yīng)的時間節(jié)點(diǎn)上，模擬效果要好于月平均趨勢對比；(3)在季節(jié)性對比中，雨季模擬效果要優(yōu)于旱季模擬效果。上述結(jié)論與VIC模型區(qū)域模擬的傳統(tǒng)經(jīng)驗(yàn)相符，可用于指導(dǎo)實(shí)際工程。

3 結(jié)束語

本文以實(shí)測氣象數(shù)據(jù)和下墊面數(shù)據(jù)為基礎(chǔ)，應(yīng)用VIC模型模擬生成研究區(qū)長系列日土壤含水率數(shù)據(jù)。這一工作的核心在于將一種較為容易獲取的資料轉(zhuǎn)為現(xiàn)階段某些區(qū)域較難獲取的數(shù)據(jù)，可為監(jiān)測土壤含水率困難的偏遠(yuǎn)地區(qū)提供一種有效的解決途徑。為保證水文模型模擬在大尺度條件下的準(zhǔn)確性，通過模型內(nèi)參數(shù)率定、參數(shù)驗(yàn)證，模型外實(shí)測站點(diǎn)對比分析等手段對結(jié)果進(jìn)行檢驗(yàn)，認(rèn)為本研究方案生成的土壤含水率可滿足一般工程精度。本文主要結(jié)論如下：

模型參數(shù)率定結(jié)果表明，渭河干流上中下游，和兩大支流涇河北洛河NES均大于80%，BLAS小于10%，全流域整體率定結(jié)果NES=0.81，BLAS=6.9%，模型在子流域尺度和全流域尺度均通過率定。

模型參數(shù)驗(yàn)證結(jié)果表明，各子流域驗(yàn)證期NES均大于75%，BLAS小于10%，全流域整體驗(yàn)證結(jié)果NES=0.81，BLAS=7.3%，該結(jié)果比率定期稍差，但符合水文模型模擬的一般規(guī)律，模型在子流域尺度和全流域尺度通過驗(yàn)證。

120組站點(diǎn)實(shí)測土壤含水率與其所在網(wǎng)格土壤含水率月平均值對比結(jié)果表明，二者具有相同趨勢但數(shù)值相差較大，整體上模擬值大于實(shí)測值，二者相差7倍左右，歸因分析認(rèn)為，這種系統(tǒng)性的偏差來自于土層深度不同以及網(wǎng)格月平均值所帶來的誤差。

實(shí)測值與模擬值在相同時刻的逐一對比結(jié)果表明，序列相關(guān)系數(shù)R2=0.7 492，平均誤差為15.2%，最大偶然誤差為32%，這一結(jié)果要遠(yuǎn)好于平均值對比，可滿足一般工程精度。分季節(jié)逐一對比結(jié)果表明，含水率相關(guān)系數(shù)旱季R2=0.612 5，雨季R2=0.927 9，可以認(rèn)為，造成模擬整體偏差的另一原因在于VIC模型對干旱區(qū)域模擬的不準(zhǔn)確，這一結(jié)論與傳統(tǒng)的經(jīng)驗(yàn)相符合。