李丹，梁新強，吳嘉平2*

（1.浙江大學環(huán)境與資源學院，杭州310058；2.浙江大學海洋學院，浙江舟山316000）

當前，中國在水資源領(lǐng)域面臨著巨大的挑戰(zhàn)，逐漸成為限制國民經(jīng)濟、社會發(fā)展的重要因素。其中，飲用水的安全供給關(guān)系到人民群眾的生命健康與安全。據(jù)QIU[1-2]估算，在中國每年因水污染問題造成約6萬人死亡，大約有3億人面臨飲用水短缺問題。而水污染加劇了飲用水供給的壓力，水質(zhì)性缺水在經(jīng)濟迅速發(fā)展的中國越來越嚴重，太湖、滇池就曾爆發(fā)過類似的群體性事件[3-4]，給人民正常生活和身體健康帶來了巨大的威脅。由于江河湖泊及地下水受到污染，水庫已經(jīng)成為中國城市飲用水供給的重要保障[5-6]。近年來，點源污染已得到有效地控制，而非點源污染成了水體污染的主要來源[7-9]。目前，減輕非點源污染的威脅已成為確保水庫水質(zhì)安全的迫切需求[10]。

流域范圍內(nèi)開展營養(yǎng)負荷定量評估是控制并減少水體富營養(yǎng)化污染的有效方法。在世界范圍內(nèi)，流域模型以及情景分析方法都是研究和分析水環(huán)境的常用工具[11-12]。但由于我國農(nóng)業(yè)耕地破碎化程度高，人工干擾強，田間差異性大，因而基于北美、歐洲農(nóng)業(yè)環(huán)境開發(fā)的模型及其參數(shù)不能直接應(yīng)用于中國農(nóng)業(yè)環(huán)境系統(tǒng)[13-15]。因此，參數(shù)的本地化顯得尤其重要，正確有效的參數(shù)化可以大大減少未來模型的調(diào)試與使用成本，為未來模型技術(shù)集成以及業(yè)務(wù)化打下堅實的基礎(chǔ)。在大多數(shù)情況下，決策者制定策略方案時缺乏對方案實施過程中可能產(chǎn)生的結(jié)果的有效預(yù)期，這使得在實際管理中如何獲得更加深入的情景研究結(jié)果變得非常重要，特別是高效準確的情景研究可以增強農(nóng)業(yè)和流域管理的現(xiàn)實意義和促進制定更好的政策和激勵措施。目前，情景分析中多關(guān)注在不同氣候變化、土地利用/覆蓋和土地管理等[16-18]條件下氮、磷營養(yǎng)元素的污染負荷受何影響。但是，關(guān)于政策引導(dǎo)對農(nóng)業(yè)生產(chǎn)、土地利用方式轉(zhuǎn)變和水質(zhì)安全的影響，特別是在流域尺度上的研究還比較少。

改革開放以來，在社會主義市場經(jīng)濟體制的激勵和科技創(chuàng)新的推進下，我國農(nóng)業(yè)產(chǎn)業(yè)結(jié)構(gòu)不斷優(yōu)化。在保障糧食生產(chǎn)的基礎(chǔ)上，畜牧業(yè)、漁業(yè)、經(jīng)濟作物生產(chǎn)量不斷攀升。然而，在飲用水水源地的農(nóng)業(yè)生產(chǎn)問題上仍面臨著許多困難與挑戰(zhàn)。水庫往往地處偏遠，受自然環(huán)境、交通條件的影響，當?shù)鼐用窠?jīng)濟條件較為一般，迫切需要從事各種社會經(jīng)濟活動以改善自身經(jīng)濟條件。因此，如何在保障庫區(qū)人民正常發(fā)展訴求的基礎(chǔ)上，維護水源地的水環(huán)境安全已成為各級地方政府需解決的頭等大事。2013年初，黃浦江死豬事件引發(fā)社會高度關(guān)注。同年5月開始，浙江省大力開展畜牧業(yè)轉(zhuǎn)型升級工作，2014年上半年，對老虎潭庫區(qū)內(nèi)的養(yǎng)雞棚、養(yǎng)豬場進行關(guān)停處理，而隨著茶葉的市場價格不斷攀升，原本經(jīng)營養(yǎng)殖場的農(nóng)民轉(zhuǎn)而進行茶葉種植，極可能帶來更大的面源污染。因此，本文應(yīng)用土壤水文評估模型（soil and water assessment tool,SWAT）分析了老虎潭水庫2010年1月—2015年4月間水體富營養(yǎng)化污染物的主要來源及其作用模式，并通過情景分析方法模擬6種管理措施下庫區(qū)水文、水質(zhì)安全情況及其經(jīng)濟效益，以期為后續(xù)水源地生態(tài)補償機制的建立和提高飲用水安全保障水平提供理論依據(jù)。

1 研究地區(qū)與研究方法

1.1 研究區(qū)概況

如圖1所示，老虎潭水庫位于浙江省西北部（119°48′31″E～119°56′55″E，30°36′37″N～30°43′35″N），湖州市吳興區(qū)埭溪鎮(zhèn)境內(nèi)，水庫攔河壩位于東苕溪支流埭溪上游埭溪鎮(zhèn)張村附近，距離湖州市區(qū)約35 km。水庫集水區(qū)包括湖州吳興區(qū)埭溪鎮(zhèn)、德清縣莫干山鎮(zhèn)以及安吉縣昆銅鄉(xiāng)，總面積110 km2。多年平均氣溫15.8℃，年平均相對濕度81%。流域多年平均降水量1 600.4 mm。庫區(qū)地勢西高東低，地貌類型以上游構(gòu)造剝蝕丘陵及下游河谷平原為主，山脈多呈北西至南東向展布，海拔多在100～500 m之間，最高海拔740 m，山坡坡度一般為30°～50°。上游有莊上、大陳、埭溪3條主要河流，河曲發(fā)育。下游河谷逐漸開闊，兩岸地勢低平，分布少量殘丘，至埭溪后入杭嘉湖平原，向東流經(jīng)東苕溪入太湖。老虎潭水庫總庫容9 966萬m3，向湖州市供水約22萬m3/d，保障市區(qū)近60%的生活飲用水。根據(jù)調(diào)查資料，2013年老虎潭庫區(qū)在冊人口14 321人，常住人口13 110人，主要集中在大沖、鎮(zhèn)水、橫嶺、雙橋、溪北和芳山6個自然村。為保障水庫水質(zhì)安全，在工程建設(shè)過程中當?shù)卣P(guān)停調(diào)整了庫區(qū)內(nèi)23家鄉(xiāng)鎮(zhèn)企業(yè)，減排污染物氨氮（ammonium,NH3-N）5.15 t/a，化學需氧量（chemical oxygen demand,COD）29.34 t/a，廢水55.94萬t/a。老虎潭水庫集水區(qū)內(nèi)有3條主要支流：莊上、大陳、埭溪，所在的流域面積分別為23.22、13.48、67.97 km2。

圖1 老虎潭庫區(qū)流域位置及氣象和監(jiān)測點Fig.1 Location of Laohutan reservoir watershed,meteorological and monitoring site

1.2 研究數(shù)據(jù)

1.2.1 主要污染源

1.2.1.1 農(nóng)業(yè)面源污染

近年來，隨著城鎮(zhèn)化的推進，人口逐漸向大中城市聚攏。由于勞動力外流，庫區(qū)內(nèi)農(nóng)業(yè)耕作活動由原有的精耕細作轉(zhuǎn)向較低勞動強度的方式，粗放的耕種模式勢必會造成新的環(huán)境壓力，其中化肥農(nóng)藥的過量施用表現(xiàn)最為突出。如表1所示，庫區(qū)內(nèi)水稻田施用氮肥量（以純氮計）高達每年500～600 kg/hm2，約為國家環(huán)保部發(fā)布的《化肥使用環(huán)境安全技術(shù)導(dǎo)則》中推薦施氮量150～180 kg/hm2的3.3倍。毛竹林復(fù)合肥施用量約每年750～1 125 kg/hm2，略高于國家推薦用量750 kg/hm2。茶園施氮量約 800～1 000 kg/hm2，是推薦施氮量225～450 kg/hm2的2.8倍左右。同時，庫區(qū)內(nèi)廣泛使用草甘膦除草劑，致使土壤表層缺少植被覆蓋，加重了水土及營養(yǎng)元素的流失。

1.2.1.2 農(nóng)村生活污染

農(nóng)村生活污染主要來源于居民生活污水、糞便污水、生活垃圾及散養(yǎng)畜禽等。農(nóng)村生活污染源成為影響水環(huán)境的重要因素，特別在水源地、庫區(qū)尤為明顯。根據(jù)《第一次全國污染源普查城鎮(zhèn)生活源產(chǎn)排污系數(shù)手冊》，在二區(qū)二類城鎮(zhèn)中，農(nóng)村人均生活污水量為175 L/d，人均COD、總氮（total nitrogen,TN）、總磷（total phosphorus,TP）、氨氮排放量分別為73、12.9、1.05、9.1 g/d。各村鎮(zhèn)的農(nóng)村生活污染排放量見附表1（http://www.zjujournals.com/agr/CN/article/showSupportInfo.do?id=10547）。根據(jù)調(diào)查核算，老虎潭水庫流域內(nèi)農(nóng)村生活污染源總氮排放量達到53.94 t/a，氨氮為38.05 t/a，總磷排放量為4.39 t/a，COD排放量為305.27 t/a。

表1 老虎潭庫區(qū)農(nóng)作物施肥情況Table 1 Fertilization management of Laohutan reservoir watershed in SWAT operations kg/(hm2·a)

1.2.1.3 畜禽養(yǎng)殖污染

2013年9月，我們對老虎潭水庫流域內(nèi)的畜禽養(yǎng)殖情況進行了詳細調(diào)查統(tǒng)計，結(jié)果顯示，庫區(qū)內(nèi)共有127家養(yǎng)殖戶，各村鎮(zhèn)主要養(yǎng)殖種類有豬、雞和鴨等。根據(jù)《第一次全國污染源普查畜禽養(yǎng)殖業(yè)源產(chǎn)排污系數(shù)手冊》中華東地區(qū)畜禽養(yǎng)殖產(chǎn)排污系數(shù)，結(jié)合調(diào)查的畜禽養(yǎng)殖數(shù)量，核算老虎潭水庫流域內(nèi)各畜禽養(yǎng)殖污染物排放量和產(chǎn)排污系數(shù)，結(jié)果見表2，各村鎮(zhèn)畜禽養(yǎng)殖污染負荷的核算結(jié)果見附表2（http://www.zjujournals.com/agr/CN/article/showSupportInfo.do?id=10547）。根據(jù)調(diào)查結(jié)果核算，老虎潭水庫流域畜禽養(yǎng)殖污染源中的總氮排放量達14.76 t/a，氨氮為3.56 t/a，總磷排放量為2.51 t/a，COD排放量為260.54 t/a。

表2 畜禽養(yǎng)殖污染物排污系數(shù)一覽表Table 2 Pollutant discharge coefficient of livestock and poultry breedingg/（頭·d）

1.2.2 土地利用

土地利用分類基于多光譜衛(wèi)星數(shù)據(jù)，包括Landsat TM/ETM+和SPOT-5獲得的該地區(qū)在2008—2014年的土地利用情況。老虎潭水庫集水區(qū)內(nèi)總土地面積約11 026 hm2。2008年，土地利用以林地為主，面積約8 881 hm2，占總土地面積的80.55%，其中大部分為竹林（7 655 hm2），還有少量松竹混合林（1 226 hm2）；其次為耕地和水體，面積分別為595、589 hm2，分別占總土地面積的5.40%和5.34%。茶園面積447 hm2（4.06%），園地面積218 hm2（1.98%），村莊251 hm2（2.27%），道路44 hm2（0.40%）。從2008年到2014年，流域內(nèi)的土地利用變化主要由林地和部分園地轉(zhuǎn)化成茶園。茶園面積擴張明顯，2014年茶園面積為606 hm2，較2008年增加了159 hm2，總增幅35.48%，年均增幅4.43%。水庫集水區(qū)內(nèi)土地利用類型見表3，土地利用現(xiàn)狀及其變化見圖2。

1.2.3 土壤類型

本文所使用的土壤數(shù)據(jù)來源于《全國第二次土壤普查》，將研究區(qū)內(nèi)1∶5萬土壤圖經(jīng)數(shù)字化所得。水庫壩址上游流域內(nèi)土壤類型主要有黃泥沙土、黃泥土2大類。其中，黃泥沙土一般分布于埭溪流域，黃泥土多分布在莊上及大陳流域。另有黃泥沙田、泥砂田、灰黃泥土等在內(nèi)的水稻土分布于溪谷兩側(cè)的耕地、林地中。土壤圖和詳細土壤屬性分別見附圖1和附表 3（http://www.zjujournals.com/agr/CN/article/showSupportInfo.do?id=10547）。

表3 2008年和2014年老虎潭庫區(qū)土地利用類型代碼及面積Table 3 Soil texture and land use distribution in Laohutan reservoir watershed in the years of 2008 and 2014

圖2 老虎潭庫區(qū)土地利用變化(2008—2014年)Fig.2 Land use and land cover change in Laohutan reservoir watershed during 2008—2014

1.2.4 水文水質(zhì)

本文所使用的水文、水質(zhì)數(shù)據(jù)來源于老虎潭水庫管理局。水質(zhì)采樣點有4個（圖1），分別位于3條支流匯入水庫的末端以及水庫壩前取水口。采樣時間從2010年1月—2015年4月，采樣周期為每月1次。水質(zhì)指標包括：pH值、溶解氧、高錳酸鉀指數(shù)、生化需氧量、總氮、總磷、氨氮含量、氟化物含量、部分重金屬含量。

1.2.5 其他數(shù)據(jù)

GIS數(shù)據(jù)庫中的數(shù)字高程模型（digital elevation model,DEM）、河網(wǎng)數(shù)據(jù)（1∶2 500）來源于中國國家測繪局。氣象數(shù)據(jù)來源于庫區(qū)內(nèi)的5個氣象站點，主要包括每日最高和最低氣溫、降雨量、蒸發(fā)量、相對濕度4類氣象要素。結(jié)合杭州氣象監(jiān)測點獲得的每日風速以及光照強度構(gòu)成模型驅(qū)動數(shù)據(jù)庫。

1.3 研究方法

1.3.1 SWAT模型簡介

SWAT模型是由美國農(nóng)業(yè)部農(nóng)業(yè)研究中心開發(fā)的流域尺度模型。該模型可以在不同土壤類型、土地利用和管理條件下模擬流域內(nèi)產(chǎn)水、產(chǎn)沙以及農(nóng)藥化肥等污染物負荷變化[19-20]。SWAT模型以其完善的模型結(jié)構(gòu)和高效的計算功能，在世界范圍內(nèi)得到了廣泛的應(yīng)用[11,21]。在老虎潭庫區(qū)所在的苕溪以及太湖流域都得到了成功的應(yīng)用[22-27]。ZHANG等[22]模擬徑流的納什效率系數(shù)（Nash-Sutcliffe coefficient of efficiency,NSCE）指數(shù)在0.70～0.83之間。LAI等[26]模擬徑流、總氮、總磷的NSCE指數(shù)在0.58～0.91、0.53～0.78、0.55～0.72之間。ZHAO等[27]模擬徑流的NSCE指數(shù)在0.77～0.81之間。這些結(jié)果都說明SWAT模型率定后能較好地模擬當?shù)貜搅骷八|(zhì)。本文所使用的模型版本為SWAT 2012（rev.627）。

1.3.2 模型構(gòu)建

1.3.2.1 敏感性分析

SWAT模型的參數(shù)眾多，在校準和驗證模型時如果對每個參數(shù)都進行調(diào)節(jié)將會費時費力，且對模型的應(yīng)用和發(fā)展極為不利。但是在模型的研究及應(yīng)用過程中，人們發(fā)現(xiàn)眾多參數(shù)中有些對模型模擬結(jié)果的影響很小，而一些敏感性參數(shù)對模擬結(jié)果具有舉足輕重的作用，所以測算出模型各參數(shù)的敏感性后再進行校準和驗證將會事半功倍。本研究根據(jù)以往工作經(jīng)驗，參考SWAT模型手冊及國內(nèi)外相關(guān)文獻研究結(jié)果[28-31]，選取25個對徑流、氮、磷模擬有較大影響的參數(shù)，詳見表5和表6。

本文采用全局敏感性分析（global sensitivity analysis）[32]以及單次敏感性分析（one-at-a-time sensitivity analysis,OAT）[33]2種方法計算參數(shù)的敏感性。在全局敏感性分析中通過拉丁超立方抽樣（Latin hypercube sampling）產(chǎn)生的參數(shù)值和目標函數(shù)值做回歸分析，得到式（1）多元回歸方程，然后利用t檢驗確定每個參數(shù)的顯著性水平。

式中：g為目標函數(shù)；α、β為回歸系數(shù)；bi為參數(shù)值。全局敏感性的大小使用2個評價指標進行衡量：t檢驗和P值。t代表參數(shù)的敏感性，其絕對值越大，表示該參數(shù)的敏感性越大；P代表敏感性的顯著性水平，其值越接近0，表示該參數(shù)越顯著，當P＜0.05時表示該參數(shù)對結(jié)果影響極為顯著[25]。

1.3.2.2 模型的率定和不確定性分析

模型率定基于SUFI-2算法。采用拉丁超立方抽樣的參數(shù)值帶入模型中進行模擬，并與實際觀測值進行比對，獲得新的參數(shù)值調(diào)整范圍。重復(fù)迭代多次直至獲得理想的參數(shù)值，將其作為模型率定結(jié)果。每次模擬次數(shù)為500次[34]，從中選出滿足預(yù)先設(shè)置的閾值（徑流：NSCE≥0.5；營養(yǎng)元素：NSCE≥0.35）的模擬結(jié)果組合，并求出95%置信度下模型模擬的不確定性范圍（95%prediction uncertainty,95PPU）[30]。

1.3.2.3 模型驗證評估

采用納什效率系數(shù)（NSCE）、徑流模擬的偏差百分比（BIAS）和決定系數(shù)（R2）對SWAT模型運行效果進行評估。

式中：Qo為實地觀察數(shù)據(jù)；Qs為模型模擬數(shù)據(jù)為觀測數(shù)據(jù)的平均值；-Qs為模擬數(shù)據(jù)的平均值；n為觀測數(shù)據(jù)的總數(shù)。

根據(jù)MORIASI等[35-36]的研究，當模擬月徑流時NSCE＞0.5、R2＞0.70、|BIAS|＜25%，認為模型模擬得到了滿意的結(jié)果。當模擬營養(yǎng)元素時，NSCE＞0.35、R2＞0.30、|BIAS|＜30%即可滿足模型應(yīng)用要求。SWAT模型對總氮、總磷的模擬難度大于對徑流和泥沙的模擬。

1.3.3 情景分析

情景分析已經(jīng)被廣泛地應(yīng)用于指導(dǎo)自然資源保護和環(huán)境管理規(guī)劃[37]。通過情景分析方法，可以考慮未來可能發(fā)生的事件過程、各種結(jié)果及其影響，進而幫助決策者做出更加明智的選擇。

預(yù)測未來本身存在著很大的難度，為了實現(xiàn)對客觀世界的再現(xiàn)，必須對其進行高度概化，選取重要的不定因素進行討論。在水體富營養(yǎng)化污染中，污染源復(fù)雜多變，河川徑流的不同水文模式?jīng)Q定流域內(nèi)營養(yǎng)元素的遷移、轉(zhuǎn)化、衰減過程，影響最終進入水環(huán)境的污染物總量[11]。因此，為了控制并減少營養(yǎng)負荷，需開展流域范圍內(nèi)定量評估工作，建立可能的營養(yǎng)負荷削減措施方案，并通過模型量化各方案的影響程度。不同的行政管理措施會對水環(huán)境產(chǎn)生不同的影響。模型率定時的不確定性作為假設(shè)情景的不確定性范圍。本文考慮3個不確定因素，即施肥量、畜禽養(yǎng)殖是否取締、茶園是否開墾，設(shè)置了6種情景模式（表4）。其中：情景S1（現(xiàn)有施肥量、保持畜禽養(yǎng)殖、未開墾茶園）是庫區(qū)原有狀態(tài)；情景S2與情景S1不同之處在于取締了畜禽養(yǎng)殖，這種狀態(tài)是庫區(qū)實施產(chǎn)業(yè)升級、保護水源地政策所期望實現(xiàn)的情景；情景S3與情景S2所不同之處在于開墾了茶園，這種狀態(tài)是政策實施所造成的真實結(jié)果；情景S4是在情景S3的基礎(chǔ)上改變了施肥量，是在現(xiàn)有政策基礎(chǔ)上提出的改進方案；情景S5是假設(shè)在情景S1的原始狀態(tài)下不取締畜禽養(yǎng)殖，僅僅是推行更優(yōu)的施肥方案；情景S6是在情景S5的基礎(chǔ)上只取締畜禽養(yǎng)殖而不進行茶園開墾的環(huán)境最優(yōu)方案。

2 結(jié)果

2.1 徑流模擬

模型模擬時間從2008年1月1日至2014年12月31日，前2年（2008—2009）僅作為模型的預(yù)熱期，不輸出結(jié)果，率定期3年（2010—2012），驗證期2年（2013—2014）。使用埭溪監(jiān)測點每日及每月徑流量實測數(shù)據(jù)進行模型率定及驗證。

表4 老虎潭水庫流域營養(yǎng)元素污染情景分析Table 4 Scenarios for analyzing effects of fertilization,livestock,land use and land cover change(LUCC)on nutrient pollution

2.1.1 參數(shù)敏感性分析結(jié)果

選取SWAT中與徑流相關(guān)的10個參數(shù)進行敏感性分析，如表5所示。徑流模擬中敏感性最高的參數(shù)是徑流曲線系數(shù)（curve number,CN2）、土壤蒸發(fā)補償系數(shù)（soil evaporation compensation factor,ESCO）、淺層地下水再蒸發(fā)系數(shù)（groundwater revaporization coefficient,Gw_revap）和土壤飽和水力傳導(dǎo)度（saturated hydraulic conductivity,Sol_K）。這4個參數(shù)的P值均小于0.05，顯著影響著徑流模擬結(jié)果，并且t檢驗絕對值都大于2，對徑流影響幅度大。前人的研究也得到了類似的結(jié)果[38-40]。

2.1.2 模型率定結(jié)果

根據(jù)參數(shù)敏感性分析結(jié)果，選取敏感性前10的參數(shù)，以埭溪監(jiān)測點2010—2012年每日及每月徑流量監(jiān)測數(shù)據(jù)作為參照數(shù)據(jù)進行率定，得到參數(shù)最優(yōu)值（表5）。模型率定期每日徑流量驗證參數(shù)NSCE、R2和BIAS分別是0.71、0.84和－0.63%；每月徑流量NSCE、R2和BIAS分別是0.85、0.88和－2.32%，表明率定結(jié)果總體令人滿意。

表5 徑流模擬參數(shù)敏感性分析及率定結(jié)果Table 5 Sensitivity analysis results for hydrologic parameters

2.1.3 模型驗證結(jié)果

采用表5中的參數(shù)率定結(jié)果，運行SWAT模型，以埭溪監(jiān)測點2013—2014年每日及每月徑流數(shù)據(jù)作為驗證參照數(shù)據(jù)。每日及每月徑流模擬的偏差百分比（BIAS）均控制在20%以內(nèi)（分別為－17.57%、－12.42%），模型的納什效率系數(shù)（NSCE）分別達0.75、0.77，決定系數(shù)（R2）分別達0.79、0.80，精度滿足模擬要求。說明模擬結(jié)果能較好地反映老虎潭庫區(qū)徑流量變化規(guī)律。

圖3和圖4展示了老虎潭流域埭溪監(jiān)測點每日及每月徑流量的觀測值和模擬值，其中，在每月徑流量模擬結(jié)果中包含了納什系數(shù)等于0.5作為似然函數(shù)閾值下徑流模擬95%的置信范圍（95PPU）[30]。徑流模擬結(jié)果同觀測結(jié)果具有相同的變化趨勢。在某些洪水季節(jié)，峰值流量的模擬值略低于真實觀測值，例如2010年7月、2012年8月和2013年6月。與此同時，在模擬枯水季節(jié)徑流量（小于1 m3/s）時極易造成過量模擬。這主要是因為枯水季節(jié)徑流主要依靠地下水來維持，而地下水流向和流量的復(fù)雜程度增加了模型模擬的難度。

圖3 埭溪監(jiān)測點每日徑流量模擬情況Fig.3 Observed and simulated daily discharge in Daixi

圖4 埭溪監(jiān)測點每月徑流量模擬情況Fig.4 Observed and simulated monthly discharge in Daixi

2.2 營養(yǎng)元素模擬

過去幾十年，農(nóng)業(yè)生態(tài)系統(tǒng)中的氮元素流失對水生環(huán)境造成了嚴重影響[41]。在流域尺度下模擬氮元素的遷移轉(zhuǎn)化對水質(zhì)保護、生態(tài)環(huán)境可持續(xù)發(fā)展具有重要的意義。模型模擬了老虎潭庫區(qū)2010年1月—2015年4月的總氮流失情況，在參數(shù)敏感性分析的基礎(chǔ)上，以壩前取水點（Laohutan dam,LHT）2010年1月—2012年12月3年的水質(zhì)數(shù)據(jù)對模型進行率定校正，得到相關(guān)參數(shù)最優(yōu)值（表5）。將模型模擬結(jié)果與流域4個實測點在2013年1月—2015年4月數(shù)據(jù)進行對比驗證。

2.2.1 參數(shù)敏感性分析結(jié)果

選取SWAT中與泥沙、營養(yǎng)元素流失相關(guān)的15個參數(shù)進行敏感性分析。結(jié)果（表6）顯示，在氮磷模擬中，氮土壤流失比例系數(shù)（nitrate percolation coefficient,NPERCO）、泥沙流失線性因子（linear parameter of maximum concentration of sediment,SPCON）、磷土壤分區(qū)系數(shù)（phosphorus soil partitioning coefficient,PHOSKD）、有機質(zhì)礦化速率（residue decomposition coefficient,RSDCO）和反硝化速率（denitrification exponential rate coefficient,CDN）是敏感性最高的5個參數(shù)。其中：NPERCO和SPCON的P值小于或等于0.05，是極為敏感的2個參數(shù)，剩下的3個參數(shù)P值在0.05～0.20之間，屬于較為敏感的參數(shù)。

2.2.2 模型率定結(jié)果

以壩前取水點（LHT）2010年1月—2012年12月間每月總氮負荷監(jiān)測數(shù)據(jù)作為參照數(shù)據(jù)進行率定，得到表6所示的參數(shù)最優(yōu)值。率定期驗證參數(shù)NSCE、R2和BIAS分別是0.80、0.89和－4.65%，表明率定結(jié)果總體令人滿意。

表6 氮磷模擬參數(shù)敏感性分析及率定結(jié)果Table 6 Sensitivity analysis results for nutritive parameters

2.2.3 模型驗證結(jié)果

圖5展示了4個監(jiān)測點總氮負荷量的觀測值和模擬值?？偟昱欧拍M值與觀測值基本保持一致，除上方外4個監(jiān)測點的NSCE在0.53～0.58之間，決定系數(shù)（R2）在 0.69～0.74之間，|BIAS|＜15%。然而，模型在2014年之后對TN負荷量的估算在4個監(jiān)測點都有低估，其中上方流域最為明顯，埭溪流域其次，這2個流域是原有畜禽養(yǎng)殖主要分布區(qū)，同時也是新增茶園主要所在流域。主要原因可能是受到滯后效應(yīng)[42]的影響，即土地利用方式的變化對土壤氮周轉(zhuǎn)與流失過程以及地下水氮水平的作用不會在短期內(nèi)就呈現(xiàn)[43]，可能存在幾年到幾十年的滯后。另一方面，新生茶園（開墾期1～3年）地表覆蓋度低，水土流失、氮磷流失嚴重，模型未能準確模擬其氮磷流失負荷，造成了一定程度的低估，但是隨著茶樹逐漸成熟，模型對其模擬能力將會恢復(fù)。

2013年庫區(qū)水質(zhì)較2011年和2012年有了明顯的改觀，總氮質(zhì)量濃度從之前的2.24 mg/L下降到了1.69 mg/L，庫區(qū)取水口水質(zhì)已經(jīng)大體恢復(fù)到了2009年的水平。然而2014年水質(zhì)惡化，總氮質(zhì)量濃度重新超過劣Ⅴ類水質(zhì)標準，達到2.07 mg/L。庫區(qū)水體處于富營養(yǎng)化狀態(tài)，仍然有爆發(fā)硅藻、藍藻等環(huán)境危機事件的可能性。

圖5 老虎潭水庫流域每日總氮負荷模擬值和實測值Fig.5 Simulated and observed daily total nitrogen(TN)loads at Laohutan reservoir watershed

2.3 最佳管理模式

圖6 情景分析總氮來源組成及經(jīng)濟效益Fig.6 Scenario analysis of total nitrogen sources and economic effectiveness

情景分析結(jié)果（圖6）顯示，6種情景中S1產(chǎn)生的入庫量最高，達到（181.6±20.5）t/a，經(jīng)濟效益較高，年產(chǎn)值為（1.72±0.09）億元，總氮污染物中有74.4%源自農(nóng)業(yè)非點源污染[（135.1±20.3）t/a]，其次是人類生活污染[（17.3±2.0）t/a，9.5%]、畜禽養(yǎng)殖[（14.7±1.5）t/a，8.1%]和其他來源[（14.5±2.9）t/a，8.0%]。從中可以看出，庫區(qū)內(nèi)農(nóng)業(yè)非點源污染是最大的污染源，同時庫區(qū)內(nèi)畜禽養(yǎng)殖為當?shù)鼐用駧砹讼喈斂捎^的經(jīng)濟收入。情景S2取締畜禽養(yǎng)殖對其他污染源并無實質(zhì)控制效果，入庫總氮較情景S1減少8.1%，達到（166.9±21.0）t/a，水質(zhì)保護效果一般；而當?shù)鼐用窠?jīng)濟收入受到較大 影響 ，降 到（0.93±0.06）億 元/年 ，減少了45.9%。情景S3在情景S2的基礎(chǔ)上增加茶園開墾，結(jié)果導(dǎo)致農(nóng)業(yè)非點源污染略有增加[（138.7±20.8）t/a]，占總污染的81.4%，致使總氮入庫總量略有增加（170.5±21.5）t/a；經(jīng)濟收入略有增加[（1.11±0.06）億元/年]。情景 S4、S5、S6均采用國家推薦施肥管理模式，可以有效減少農(nóng)業(yè)生產(chǎn)帶來的非點源污染，入庫總氮較情景S1削減幅度巨大（51.0%～62.0%）；同時，由于在生態(tài)種植下稻田的經(jīng)濟效益明顯提升，在一定程度上提升了農(nóng)民的經(jīng)濟收入。其中：情景S4是在現(xiàn)有管理模式上提出的改進措施，推廣生態(tài)施肥既可以減少污染物總量，又可以切實增加居民收入[（1.42±0.08）億元/年]。情景S5是假設(shè)在實施畜禽養(yǎng)殖產(chǎn)業(yè)升級前就實施科學施肥管理的效果，由于沒有新開墾茶園，農(nóng)業(yè)面源污染量較情景S4略少，但是每年會產(chǎn)生（14.5±2.9）t/a的畜禽養(yǎng)殖總氮污染，所造成的總污染比情景S4更多；此外，情景S5獲得了最高的經(jīng)濟效益[（2.04±0.10）億元/年]。情景S6是假設(shè)嚴格控制所有因素（取締畜禽養(yǎng)殖、不開墾新茶園、科學施肥）得到的環(huán)境最優(yōu)型管理模式，總氮入庫量為（87.4±7.9）t/a，較情景 S5[（101.9±10.1）t/a]減少14.2%，但是經(jīng)濟效益較情景S5減少38.5%[（1.26±0.07）億元/年]。綜合經(jīng)濟效益和環(huán)境保護，情景S5應(yīng)該是最佳管理模式，情景S4和S6為較優(yōu)管理模式。綜合情景分析可知，推廣科學施肥，減少化學肥料使用，多施用有機肥，充分利用自然條件，減少人為對自然的干擾，可以獲得優(yōu)質(zhì)農(nóng)產(chǎn)品，并保障當?shù)厮h(huán)境安全。

3 討論

老虎潭庫區(qū)內(nèi)河流均屬山溪性河流，具有明顯的流域邊界，匯水流向相對簡單，是我國華東、華南水庫的典型形式。河床縱坡降較大，源短流急，水位、水量受降水強度影響變化幅度較大，非常適合SWAT模型模擬。由敏感性分析結(jié)果可知，水文模擬時重點調(diào)整的參數(shù)是SCS徑流曲線系數(shù)（CN2）、土壤蒸發(fā)補償系數(shù)（ESCO）、淺層地下水再蒸發(fā)系數(shù)（Gw_revap）和土壤飽和水力傳導(dǎo)度（Sol_K）。氮、磷模擬時重點調(diào)整的參數(shù)是氮土壤流失比例系數(shù)（NPERCO）、泥沙流失線性因子（SPCON）、磷土壤分區(qū)系數(shù)（PHOSKD）、有機質(zhì)礦化速率（RSDCO）和反硝化速率（CDN）。CN2參數(shù)嚴重影響著地面匯流過程，而ESCO和Gw_revap的高敏感性來源于老虎潭流域的地理特性，該地緯度較低（30°N），累計全年太陽輻射強度較大（4 000～5 000 MJ/m2），水分蒸發(fā)、植物蒸騰作用對土壤含水率影響較大，對流域水循環(huán)起著至關(guān)重要的影響。Sol_K參數(shù)大小表征土壤中水分向下滲透速度的快慢，與土壤的機械組成有直接關(guān)系，對水分在土壤中的垂直分布有著至關(guān)重要的影響。除此之外，土壤表層可利用水量（Sol_AWC）和淺層地下水再蒸發(fā)閾值（REVAPMN）對徑流模擬結(jié)果也有一定影響。多項研究的徑流參數(shù)優(yōu)化結(jié)果與以上結(jié)果類似[28,38,44]。NPERCO代表氮元素從土壤中滲透到地表水中的能力，NPERCO值越大，代表更多的溶解態(tài)氮滲透到地表水或者側(cè)流當中，說明大量的溶解態(tài)氮通過地表徑流和側(cè)流的形式流失。事實上，相關(guān)研究[45-46]在模擬氮流失時僅僅對這個參數(shù)進行了率定。SPCON是控制流域內(nèi)泥沙流失強度最敏感的參數(shù)[47]，總氮的流失很大一部分是泥沙吸附態(tài)氮，氮磷流失量同泥沙的攜帶能力強烈相關(guān)[48]。PHOSKD指表層土壤（10 mm）中可溶性磷的濃度和地表徑流中可溶性磷的濃度之比，這個參數(shù)直接控制徑流中磷的流失量。RSDCO和CDN控制著氮循環(huán)中固化以及反硝化的速率，是氮循環(huán)中至關(guān)重要的參數(shù)。

農(nóng)業(yè)耕種模式和土地利用變化對氮的流失具有顯著影響。自然植被通?？梢杂行p少流失，而農(nóng)業(yè)活動則會加劇流失。在天然竹林地中氮的流失主要以泥沙結(jié)合態(tài)為主。而在耕作土地中氮的流失主要以可溶態(tài)氮為主，因此其流失量主要受徑流強度影響。雨季（6—9月）的營養(yǎng)鹽輸出大于旱季（12月—次年3月）。降雨影響地表水質(zhì)有2種途徑，第1種是造成面源污染的沖刷作用，第2種是對原有水體的稀釋作用，它們對氮磷遷移起正反2方面作用[49]。由于施肥時間不同，這2種作用對老虎潭庫區(qū)的主要影響不同。每年9—10月是竹林施肥主要時間段，每年3—4月是茶園主要施肥期，這段時間內(nèi)河流中的污染物濃度同降雨量呈正相關(guān)，表明這期間的主要污染模式為農(nóng)田面源污染。而在5—7月土壤中所施用的氮磷含量下降，地表植被覆蓋度增加，雨水沖刷所帶入的污染物數(shù)量減少，而降雨量可以更有效地對水體進行稀釋，從而降低污染物的濃度。

從情景分析可以看出，關(guān)停畜禽養(yǎng)殖可以快捷、便利地減少總氮的入庫量，但是由于該項措施所引起的后續(xù)效應(yīng)（新建茶園）對污染物控制不利，同時在一定程度上損害了當?shù)剞r(nóng)戶的經(jīng)濟利益，地方政府對農(nóng)業(yè)產(chǎn)業(yè)升級的行政干預(yù)并不能達到完美控制減排并保障經(jīng)濟發(fā)展的目的。在現(xiàn)有條件下依照情景S4進行，從農(nóng)業(yè)面源污染著手控制污染來源，雖然政策實施難度比較大，但是可以得到最優(yōu)的污染控制效果并且取得一定的經(jīng)濟效益。如引導(dǎo)庫區(qū)農(nóng)民采用環(huán)境友好型耕種技術(shù)，即減少化肥施用，多施用有機肥，可減少面源污染，逐步降低老虎潭水庫水體中氮、磷含量，確保水庫水源的安全。情景S4是現(xiàn)階段最容易實現(xiàn)的最佳管理模式。

4 結(jié)論

4.1 本文應(yīng)用SWAT模型計算得出2010—2014年老虎潭水庫流域年均入庫水量約1萬m3、總氮180～200 t。其中：總氮74.4%來自農(nóng)業(yè)生產(chǎn)，9.5%來自生活污染，8.1%來自畜禽養(yǎng)殖，8.0%來自其他污染源。

4.2 通過模型模擬及情景分析，發(fā)現(xiàn)老虎潭水庫庫區(qū)內(nèi)耕地過量施肥、土壤侵蝕和水土流失是非點源污染產(chǎn)生的主要原因。因此，推廣合理施肥、加強水土保持措施、控制非點源污染是對庫區(qū)水環(huán)境最優(yōu)的管理方式。

4.3 地方政府應(yīng)謹慎運用行政手段干預(yù)農(nóng)業(yè)產(chǎn)業(yè)升級，應(yīng)當在充分調(diào)研和科學預(yù)測的基礎(chǔ)上加強引導(dǎo)和調(diào)控，從而促進污染減排，同時保障生態(tài)效益和經(jīng)濟效益，實現(xiàn)社會福利最大化。

[1] QIU J.China to spend billions cleaning up groundwater.Science,2011,334(6057):745.

[2] QIU J.China pledges to get wealthier with less water.Nature,2009.

[3] 沈建軍,李柏山,徐海萍.太湖水污染原因分析及治理措施.環(huán)境科學導(dǎo)刊,2009,28(2):27-29.SHEN J J,LI B S,XU H P.Water pollution analysis and treatment measures of Taihu Lake.Environmental Science Survey,2009,28(2):27-29.(in Chinese with English abstract)

[4] 丁莞韻.中國水污染事件記實.環(huán)境保護,2007(7B):83-85.DING W Y.China’s water pollution event documentary.Environmental Protection,2007(7B):83-85.(in Chinese)

[5] XIN X,LI K,FINLAYSON B,et al.Evaluation,prediction,and protection of water quality in Danjiangkou Reservoir,China.Water Science and Engineering,2015,8(1):30-39.

[6] ZHANG H,HUANG T,LIU T.Sediment enzyme activities and microbial community diversity in an oligotrophic drinking water reservoir,eastern China.PLoS One,2013,8(10):e78571.

[7] ONGLEY E D,ZHANG X L,YU T.Current status of agricultural and rural non-point source pollution assessment in China.Environmental Pollution,2010,158(5):1159-1168.

[8] JIANG J,LI S,HU J,et al.Amodeling approach to evaluating the impacts of policy-induced land management practices on non-point source pollution:A case study of the Liuxi River watershed,China.Agricultural Water Management,2014,131:1-16.

[9] HAN H,ALLAN J D.Uneven rise in N inputs to the Lake Michigan Basin over the 20th century corresponds to agricultural and societal transitions.Biogeochemistry,2012,109(1/2/3):175-187.

[10]NIELSEN A,TROLLE D,ME W,et al.Assessing ways to combat eutrophication in a Chinese drinking water reservoir using SWAT.Marine and Freshwater Research,2013,64(5):475.

[11]HASHEMI F,OLESEN J E,DALGAARD T,et al.Review of scenario analyses to reduce agricultural nitrogen and phosphorus loading to the aquatic environment.Science of the Total Environment,2016,573:608-626.

[12]ISIK S,KALIN L,SCHOONOVER J E,et al.Modeling effects of changing land use/cover on daily streamflow:An artificial neural network and curve number based hybrid approach.Journal of Hydrology,2013,485:103-112.

[13]XIE X H,CUI Y L.Development and test of SWAT for modeling hydrological processes in irrigation districts with paddy rice.Journal of Hydrology,2011,396(1/2):61-71.

[14]ZHOU F,SHANG Z Y,ZENG Z Z,et al.New model for capturing the variations of fertilizer-induced emission factors of N2O.Global Biogeochemical Cycles,2015,29(6):885-897.

[15]ZHENG J,LI G Y,HAN Z Z,et al.Hydrological cycle simulation of an irrigation district based on a SWAT model.Mathematical and Computer Modelling,2010,51(11/12):1312-1318.

[16]WILSON C O,WENG Q H.Simulating the impacts of future land use and climate changes on surface water quality in the Des Plaines River watershed,Chicago Metropolitan Statistical Area,Illinois.Science of the Total Environment,2011,409(20):4387-4405.

[17]LIN K,Lü F,CHEN L,et al.Xin’anjiang model combined with curve number to simulate the effect of land use change on environmental flow.Journal of Hydrology,2014,519:3142-3152.

[18]SHAMIR E,MEGDAL S B,CARRILLO C,et al.Climate change and water resources management in the Upper Santa Cruz River,Arizona.Journal of Hydrology,2015,521:18-33.

[19]ARNOLD J G,SRINIVASAN R,MUTTIAH R S,et al.Large area hydrologic modeling and assessment partⅠ:Modeldevelopment.JournaloftheAmerican Water Resources Association,1998,34(1):73-89.

[20]NEITSCH S L,ARNOLD J G,KINIRY J R,et al.Soil&Water assessment tool theoretical documentation,version 2009//Texas Water Resources Institute Technical Report No.406.USA:TexasA&M University System,2011.

[21]LI Z L,XU Z X,LI Z J.Performance of WASMOD and SWAT on hydrological simulation in Yingluoxia watershed in northwest of China.Hydrological Processes,2011,25(13):2001-2008.

[22]ZHANG Q L,CHEN Y X,JILANI G,et al.Model AVSWAT apropos of simulating non-point source pollution in Taihu Lake basin.Journal of Hazardous Materials,2010,174(1/2/3):824-830.

[23]金婧靚,王飛兒,戴露瑩,等.苕溪流域非點源污染特征及其影響因子.應(yīng)用生態(tài)學報,2011,22(8):2119-2125.JING J L,WANG F E,DAI LY,et al.Characteristics of nonpoint source pollution in Tiaoxi watershed and related affecting factors.Chinese Journal of Applied Ecology,2011,22(8):2119-2125.(in Chinese with English abstract)

[24]李麗嬌,薛麗娟,張奇.基于SWAT的西苕溪流域降雨-徑流關(guān)系及水量平衡分析.水土保持通報,2008(5):81-85.LI L J,XUE L J,ZHANG Q.SWAT-based rainfall-runoff relationship and water balance analysis in Xitiaoxi catchment.Bulletin of Soil and Water Conservation,2008,28(5):81-85.(in Chinese with English abstract)

[25]LIANG X Q,WANG Z B,ZHANG Y X,et al.No-tillage effects on N and P exports across a rice-planted watershed.Environmental Science and Pollution Research,2016,23(9):8598-8609.

[26]LAI G Y,YU G,GUI F.Preliminary study on assessment of nutrient transport in the Taihu Basin based on SWAT modeling.Science in China,Series D:Earth Sciences,2006,49(Suppl.1):135-145.

[27]ZHAO G J,H?RMANN G,FOHRER N,et al.Development and application of a nitrogen simulation model in a data scarce catchment in South China.Agricultural Water Management,2011,98(4):619-631.

[28]NIRAULA R,MEIXNER T,NORMAN L M.Determining the importance of model calibration for forecasting absolute/relative changes in streamflow from LULC and climate changes.Journal of Hydrology,2015,522:439-451.

[29]NIRAULAR,NORMAN LM,MEIXNER T,et al.Multi-gauge calibration for modeling the semi-arid Santa Cruz Watershed in Arizona-Mexico Border area using SWAT.Air,Soil and Water Research,2012,5:41-57.

[30]YANG J,REICHERT P,ABBASPOUR K C,et al.Comparing uncertainty analysis techniques for a SWAT application to the Chaohe Basin in China.Journal of Hydrology,2008,358(1/2):1-23.

[31]張平,劉云慧,肖禾,等.基于SWAT模型的北京密云水庫沿湖區(qū)氮磷流失風險分區(qū).中國農(nóng)業(yè)大學學報,2011(3):53-59.ZHANG P,LIU Y H,XIAO H,et al.Identification of the risk area of nitrogen and phosphorus loss at lakeshore of Miyun reservoir by using SWAT model.Journal of China Agricultural University,2011,16(3):53-59.(in Chinese with English abstract)

[32]VAN GRIENSVEN A,MEIXNER T,GRUNWALD S,et al.A global sensitivity analysis tool for the parameters of multivariable catchment models.Journal of Hydrology,2006,324(1/2/3/4):10-23.

[33]VAN GRIENSVEN A,FRANCOS A,BAUWENS W.Sensitivity analysis and auto-calibration of an integral dynamic model for river water quality.Water Science and Technology,2002,45(9):325-332.

[34]ABBASPOUR K C.SWAT-CUP:SWAT calibration and uncertainty programs:Auser manual.Duebendorf,Switzerland:Swiss Federal Institute of Aquatic Science and Technology,Eawag,2011.

[35]MORIASI D N,GITAU M W,PAI N,et al.Hydrologic and water quality models:Performance measures and evaluation criteria.Transactions of the ASABE,2015,58(6):1763-1785.

[36]MORIASI D N,ARNOLD J G,VAN LIEW M W,et al.Model evaluation guidelines for systematic quantification of accuracy in watershed simulations.Transactions of the ASABE,2007,50(3):885-900.

[37]曹茜,于德永,孫云,等.土地利用/覆蓋變化與氣候變化定量關(guān)系研究進展.自然資源學報,2015,30(5):880-890.CAO Q,YU D Y,SUN Y,et al.Research progress on the quantitative relationship between land-use/land-cover change and climate change.Journal of Natural Resources,2015,30(5):880-890.(in Chinese with English abstract)

[38]HOLVOET K,VAN GRIENSVEN A,SEUNTJENS P,et al.Sensitivity analysis for hydrology and pesticide supply towards the river in SWAT.Physics and Chemistry of the Earth,Parts A/B/C,2005,30(8/9/10):518-526.

[39]RAJIB M A,MERWADE V,YU Z Q.Multi-objective calibration of a hydrologic model using spatially distributed remotely sensed/in-situ soil moisture.Journal of Hydrology,2016,536:192-207.

[40]WELLEN C,ARHONDITSIS G B,LONG T,et al.Quantifying the uncertainty of nonpoint source attribution in distributed water quality models:A Bayesian assessment of SWAT’s sediment export predictions.Journal of Hydrology,2014,519(Part D):3353-3368.

[41]GAO M,QIU J,LI C,et al.Modeling nitrogen loading from a watershed consisting of cropland and livestock farms in China using Manure-DNDC.Agriculture,Ecosystems&Environment,2014,185:88-98.

[42]LIU J,DIETZ T,CARPENTER S R,et al.Complexity of coupled human and natural systems.Science,2007,317(5844):1513-1516.

[43]BüRGI M,?STLUND L,MLADENOFF D J.Legacy effects of human land use:Ecosystems as time-lagged systems.Ecosystems,2017,20(1):94-103.

[44]LI Z,LIU W Z,ZHANG X C,et al.Impacts of land use change and climate variability on hydrology in an agricultural catchment on the Loess Plateau of China.Journal of Hydrology,2009,377(1/2):35-42.

[45]BEHERA S,PANDA R K.Evaluation of management alternatives for an agricultural watershed in a sub-humid subtropical region using a physical process based model.Agriculture,Ecosystems&Environment,2006,113(1/2/3/4):62-72.

[46]SCHILLING K E,WOLTER C F.Modeling nitrate-nitrogen load reduction strategies for the Des Moines River,Iowa using SWAT.Environmental Management,2009,44(4):671-682.

[47]SHAN N,RUAN X H,XU J,et al.Estimating the optimal width of buffer strip for nonpoint source pollution control in the Three Gorges Reservoir Area,China.Ecological Modelling,2014,276:51-63.

[48]王曉燕,王一峋,王曉峰,等.密云水庫小流域土地利用方式與氮磷流失規(guī)律.環(huán)境科學研究,2003,16(1):30-33.WANG X Y,WANG Y X,WANG X F,et al.The character of nutrient loss and land use in a small watershed of Miyun Reservoir.Research of Environmental Sciences,2003,16(1):30-33.(in Chinese with English abstract)

[49]李瑞玲,張永春,劉莊,等.太湖緩坡丘陵地區(qū)雨強對農(nóng)業(yè)非點源污染物隨地表徑流遷移的影響.環(huán)境科學,2010,31(5):1220-1226.LI R L,ZHANGYC,LIU Z,et al.Rainfall intensity effects on nutrients transport in surface runoff from farmland in gentle slope hilly area of Taihu Lake basin.Environmental Science,2010,31(5):1220-1226.(in Chinese with English abstract)