匯水面積

地形因子是土壤侵蝕模型中的主要構成部分，也是目前研究的重點和難點。土壤侵蝕模型中的地形因子包括坡長因子(L) 和坡度因子(S)，坡長因子反映了土壤侵蝕量與坡長之間的量化關係，坡度因子則反映了坡度對侵蝕的影響。通用土壤流失方程(USLE）及其修訂版（RUSLE）中的L 因子表示標準化到22.13 m 坡長上的土壤侵蝕量。

坡長為從地表徑流源點到坡度減小至有沉積出現地方之間的距離，或到一個明顯的渠道之間的水平距離。坡長可以從野外直接測量，也可以通過等高線量算確定。RUSLE 中一般用匯水渠道表示坡長結束的地方，由於地形圖上一般沒有標明匯水渠道，所以通過等高線確定的坡長一般過長。隨著數字製圖技術的出現，利用DEM（Digital Ele-vation Model）計算坡長得到廣泛應用。其中基於徑流路徑柵格累積提取坡長的演演算法思路為：在計算水流流向的基礎上，定義局部高點作為坡長累計計算起點，根據水流來向和流向關係，由高到低，通過不斷尋求徑流結束點的方式，利用多重循環和迭代方法，完成對累計坡長的計算，該演演算法中以Hichy和Remortel 的研究為代表。國內學者以Remortel的AML程序為基礎進行了改進和應用。

事實上，在真實的二維空間中，地表徑流及其產生的土壤流失，不是決定於到分水嶺或農田邊界的距離，而是決定於單位等高線長度上的徑流面積（簡稱為單位匯水面積）。由於受水流匯聚或分散的影響，單位匯水面積與坡長存在較大的差別。所以，土壤侵蝕模型中的地形因子計算僅考慮坡面縱向形態的影響還不全面，還要考慮坡面平面形態的影響，也就是要考慮徑流匯聚程度的影響。現有研究中對於基於徑流路徑的LS 因子值與基於單位匯水面積的LS 因子值的差異研究尚不多見。為此，本文嘗試對兩者的差異及其區域特徵進行探討。

本研究所涉及區域位於黑龍江九三農墾鶴山農場鶴北流域中的2 號小流域，該流域地理位置介於125°16′12″E~125°18′7″E，49°0′0.8″N~49°1′5.1″ N，流域面積為3 km。為研究方便，本文所討論2號小流域為通過徑流演演算法得到的純自然流域，與受諸如機耕路、壟作方向等人類活動而使自然徑流流向發生改變，並最終形成的流域形態稍有不同。該區是典型東北漫川漫崗地帶，呈坡長坡緩的地形特徵，海拔高度一般在400 m 以下，相對高差≤100 m，坡度一般為1°~3°，大坡度在3°~6°之間。該區域基本是在1949年後開始大規模開墾，儘管時間相對較短，但到目前為止侵蝕相當嚴重，據相關研究，該區的侵蝕速率達到了14.5 t/(hm·a)，遠超容許土壤流失量。研究區氣候屬寒溫帶大陸性半濕潤氣候，氣溫冷熱相差懸殊，夏季最熱月份在7 月，平均氣溫20.8℃，最高氣溫可達37℃；冬季最冷在1 月份，平均氣溫為-22.5℃，最低氣溫可達-43.7℃。初春溫差較大，年均氣溫0.4℃左右。年降水量在500~550 mm 間，降水年際變化大，分佈不均，集中在7~9 月份，佔到全年的64%。本區漫崗基本都已開墾成耕地，耕種作物主要為大豆和小麥（Triticum aestivum），大豆（Glycine max）實行壟作耕種。由於難以做到完全等高起壟，所以往往造成順（坡）壟溝狀面蝕和斷壟溝狀侵蝕。

研究的總體思路為，在構建柵格DEM 的基礎上，分別計算基於單位匯水面積的LS 因子值和基於徑流路徑柵格累積的LS 因子值，並對兩類演演算法計算值進行比較分析。本研究所用柵格DEM 是在1:1 萬地形圖基礎上矢量化所得，首先利用矢量化數據構建不規則三角網TIN，進而將其轉化為柵格DEM。通過量算等高線間的間距，最終確定柵格DEM 的解析度為2 m。依據構建轉化得到的柵格DEM，分別計算基於單位匯水面積和基於徑流路徑柵格累積兩類不同演演算法LS 因子值。

為對比分析兩類不同演演算法差異，考慮到數據量，利用Hawths Tools 工具在計算得到的LS 因子柵格圖上隨機生成300 個點，並提取其LS 因子值。由於RUSLE極少使用超過305 m的坡長，加之本研究所用基於徑流路徑的坡長演演算法是根據RUSLE 所得，因此，本研究中對於坡長超過305 m的隨機點不予考慮，最終實際參與對比分析的隨機點為291個。流域隨機分佈291 點的坡度和坡長的最小值、最大值、平均值和標準偏差如表1 所示。

本研究基於單位匯水面積的坡長因子計算，採用的是Desmet 提出的演演算法，該演演算法以Foster提出的不規則坡面坡長因子計算公式推演而來。地形因子中的坡度因子則分別採用Wischmeier、McCool、Govers、Nearing等人的研究成果。上述坡度坡長因子的具體演演算法請參閱相應文獻，在此不再贅述。對於基於單位匯水面積的地形因子而言，由於坡長因子和坡度因子分別來自於不同的研究成果，為表述方便，將LS 因子演演算法的名稱，以LS 因子中所採用的坡度因子演演算法進行命名。研究中面積計算涉及到的流向演演算法，統一採用Quinn 等提出的多流向演演算法（Multiple Flow）。

基於徑流路徑柵格累計演演算法是在Remortel的AML（arc macro language）程序基礎上對相關參數修改得到。Remortel 演演算法提取的坡長反應了侵蝕-沉積發生的空間位置，由於其簡單易算而得到廣泛的應用。目前Remortel 的AML 程序已經發展到第四版，不同版本的AML 程序演演算法在坡長指數因子值的選取和坡度因子演演算法上有所不同。以其最新版第四版為例，代碼反映出該版本m 取值較之前版本更為精細，m 值根據McCool等研究並經局部內插得到。該版演演算法中的坡度因子，則是根據McCool 等於1987 年提出的S 因子計算公式得來。第四版AML 程序中坡長指數（m）取細溝和細溝間侵蝕之比較低情況下的數值，而這僅適用於諸如草地和其他有覆蓋的緊實土壤，這與本研究區細溝和淺溝較為發育不符，因此將m修正為RUSLE 中採用的Foster 等的計算方法。

4.1 Remortel 修正演演算法與原版演演算法比較

通過與McCool 特定坡度坡長條件下的LS 因子值的相關分析，發現Remortel 修改演演算法與McCool 的LS 因子值的相關係數達到了0.999，而兩者線性回歸方程的斜率也近乎等於1，可以看出兩者近乎呈1：1 的關係。與之相比，Remortel 第四版AML 原始程序演演算法與McCool 之LS 因子值的相關性儘管依然相對較高，達到了0.985，但從回歸方程的斜率看則明顯要小於相等坡度坡長條件下的Mc-Cool 之值。究其原因，與Remortel 第4 版AML 程序演演算法中坡長指數取值有關，其計算中採用的坡長指數（m）為細溝和細溝間侵蝕之比較低情形下的數值，而本研究採用的是在細溝和細溝間侵蝕比率中等情形下的坡長指數作為參照。通過與實測徑流泥沙數據的對比，坡長指數採用RUSLE 中Foster 等的方法計算結果與實測值結果吻合較好。總之，無論從相關係數來看，還是從線性回歸方程的斜率看，參數修正後的Remortel 演演算法計算值明顯要好於第四版原始演演算法計算值。

4.2 不同演演算法流域值比較

對流域整體而言，以單位匯水面積為基礎得到的4 種演演算法LS 因子值存在明顯差異，特別是Govers 演演算法與其他3 種演演算法在最大值、平均值及標準偏差間相差較大。與之相比，McCool、Wisch-meier 和Nearing 演演算法計算值在流域值整體方面相差不大。從表2可以看出，無論地形因子的最大值、平均值還是標準偏差，根據Govers演演算法計算得到的流域LS因子值在4種演演算法中都是最大的，其流域最大值達到85.88，而McCool、Wischmeier、Nearing 演演算法的流域最大值僅有不到Govers 演演算法計算值的一半，如根據McCool 演演算法計算得到的流域最大值是34.78，Nearing 的最大值在三者中相對較大，也只有35.24。地形因子流域最大值在4 種演演算法中最小的為Wischmeier 演演算法，只有29.20。對於流域平均值而言，也呈現與流域最大值相類似的特徵，所不同的是McCool 演演算法計算值稍大於Nearing 演演算法計算值。標準偏差方面，也是Govers 演演算法最大，其他3 種演演算法相對較小。總之，對流域整體LS 因子計算相關參數而言，4 種演演算法中，Govers演演算法除去最小值外都遠遠高於其他3 種演演算法計算值，對於McCool、Wischmeier 和Nearing 演演算法而言，流域整體計算值都相差不大。與前述基於單位匯水面積LS 因子演演算法計算值相比，基於徑流路徑柵格累計的Remortel 修改演演算法LS 因子值，在流域最大值、平均值以及標準偏差上都要小於前者，僅僅只有最小值要高於Gov-ers 演演算法的最小值。

4.3 兩類演演算法LS 因子值比較

4 種演演算法中的點基本都位於1：1 線以上，僅有Wischmeier 和Nearing 演演算法中，有少數點位於1：1 線以下，但數量相當有限。這說明291 點的4 種演演算法計算值與參照值相比，基本上都要大於參照值。從回歸趨勢線的斜率來看，Nearing演演算法與參照值最為接近，其回歸趨勢線的斜率僅有1.485，其次依次為McCool演演算法及Wischmei-er演演算法，偏離最大的為Govers演演算法計算值，回歸線斜率達到了3.786。對於不同演演算法的相關性而言，Mc-Cool、Nearing及Wischmeier與參照值的相關性相差不大，最好的為Nearing 演演算法計算值，相關係數達到了0.753，其次為McCool和Wischmeier演演算法，相關性最差的為Govers演演算法，相關係數只有0.630。從以上分析可以看出，在其他因子相同的條件下，由基於單位匯水面積地形因子計算得到的土壤侵蝕量要大於由柵格累計演演算法得到的土壤侵蝕量。以往的徑流小區實驗基本都是圍繞坡度和坡長展開，坡度、坡長與侵蝕量間的量化關係得到了實驗數據的支持，而截至目前單位匯水面積與侵蝕量關係的實驗研究尚未見相關報道。坡度和坡長參數，尤其是坡長的測量或計算有一定的不確定性和複雜性，儘管理論推導上看，單位匯水面積代替坡長具有一定的合理性，但從上面的分析可以看出，要確切的認識單位匯水面積得到的地形參數與侵蝕量的關係，還需進一步從試驗得到的實測數據來分析。

變化在不同演演算法間相差不大。與之形成鮮明對照的是，L 因子隨面積的變化則差異明顯。此處需要說明的是，由於在計算基於單位匯水面積的LS 因子時，除涉及面積參數外，還有修正因子，即與坡向相關的參數，故而其隨著面積的增加並非呈現規則的遞增趨勢。在5 種演演算法中，尤其以Gover 的坡長因子演演算法變化最大，在相同的匯水面積條件下，Gover 的L 因子演演算法計算值要比其他四種演演算法計算值大幾倍甚至十幾倍之多。這說明在研究區坡長坡緩的地形條件下，不同LS 因子演演算法對於坡長因子的影響更為明顯。

4.4 坡長指數m的影響

在4種演演算法LS因子計算中，除去Govers的L因子演演算法中坡長因子m 為定值外，其他演演算法都為變值。由於Desmet不規則坡面L因子演演算法涉及坡度、坡向及面積等諸多因子，為簡化相關參數具體分析坡長指數的影響，本文以基於徑流路徑的坡長為自變數，結合研究流域的坡度分佈情況，將1.57%和4.13%及12.67%作為固定坡度（此處因3.97%和4.13%坡度條件下計算得到的坡長因子相差不大，故只採用4.13%的坡度進行分析），分析不同演演算法條件下，坡長因子L隨坡長增加的變化情形。

在坡長大於約22 m 的坡面上，無論是RUSLE 還是USLE 中的坡長因子計算值，在等坡長條件下，坡長因子值都隨著坡度的增加而增加；在坡長小於約22 m 的坡面上，則正好相反，即在等坡長條件下，坡長因子值都隨著坡度的增加反而減小。這說明在大於22 m的坡面上，等坡長條件下土壤侵蝕隨著坡度的增加而增強，而在小於22 m 的坡面上，則正好相反。RUSLE 和USLE模型計算可以看出，坡長與侵蝕間關係的複雜性。RUSLE 中的m取值有所不同，但在有限的坡長條件下，兩者的差異相對有限。只不過，在大於22 m 的坡面上，1.57%坡度下USLE 演演算法坡長因子值大於等坡長條件下RUSLE 中的坡長因子值，但在4.13%和12.67%坡度條件下，正好相反，RUSLE 中的坡長因子值大於等坡長條件下USLE 演演算法得到的坡長因子值。對於坡長指數為定值的Govers 演演算法而言，由於m 為定值，所以不同坡度條件下由其計算得到的坡長因子不存在隨坡度變化的情形，而且由於m值較大，由其計算得到的坡長因子值在大於22 m 坡長的坡面上，都要遠遠大於USLE 和RUSLE 中L因子的計算值。

對於坡長因子的差異，一般可以通過直接比較坡長因子指數的大小來表示。McCool 等通過對不同小區數據的總結得到，坡長指數一般變化於0~0.9 之間，多集中於0.27~0.68 間。國內學者在不同地區針對特定坡度和坡長小區，研究了地形因子與侵蝕量的量化關係。縱觀國內相關研究，全國不同地區的坡長指數變化於0.14~0.46 之間，在相同坡度級別條件下，中國學者得到的m 值普遍低於美國所採用的坡長指數值，而且國內外的研究都得出坡長指數隨坡度增加而增大的趨勢。

而Govers 採用的坡長指數為1.45，要大於絕大部分國外有關的研究數值。由此可見，對於坡面侵蝕而言，應用Govers 演演算法的計算值要遠遠大於實際值，正如Govers 所言，該演演算法更為適用於細溝侵蝕佔比更高的坡面侵蝕預測。

現有不同區域坡長與侵蝕量的試驗研究，基本是在特定的坡度條件下進行，這也就意味著研究得出的坡長因子指數為定值，而如果坡度都採用定值的話，坡長因子則意味著土壤侵蝕只與坡長有關。從侵蝕機理來講，坡度和坡長以及坡型坡向等地形因子共同決定了降水在彙集流動過程中能量轉化能力。所以，如果就某一具體區域單純從地形角度來講，侵蝕機理分析可以看出，坡長因子指數採用與坡度有關的變值更為合理。當然，如果就不同區域而言，坡長指數m 還要受到諸如土壤特性、地表植被類型和田間管理措施的影響，說明坡長指數存在著區域異質性和差異性。

信息技術的發展，使得運用DEM 數據進行流域土壤侵蝕量預測成為現實。在具體應用中，計算是基於單個像元數據進行，而在整個流域中由於坡度存在非均一性，特別是對於面積有限而地形高差懸殊的流域來講更是如此，這就使得應用定值坡長指數m 的坡長因子計算受到限制。所以在對流域土壤流失量進行計算時，使用可變坡長指數m 值計算坡長因子更為合理。

1）在研究區坡長坡緩的地形條件下，基於單位匯水面積的地形因子計算值普遍要大於基於徑流路徑柵格累積的地形因子計算值。這也就意味著在其他因子相同的條件下，由前者計算得出的土壤侵蝕量要高於後者計算得到的侵蝕量。在現有坡長測量或計算具有不確定性和複雜性條件下，展開單位匯水面積與侵蝕量關係的試驗研究尤為必要。

2）坡度因子和坡長因子分別來看，研究區不同演演算法計算得到的坡度因子差別不大，與之相比，坡長因子則差異明顯。這說明在研究區坡長坡緩的地形條件下，坡長因子對LS因子演演算法響應更為敏感。

3）從不同演演算法坡長因子計算值來看，坡長因子指數為變值的坡長因子計算值間的差別不大，而對於坡長因子指數相對較大並為定值的Govers演演算法而言，其計算值則要遠遠大於其他演演算法的計算值。侵蝕機理分析認為，在應用DEM 進行流域尺度地形因子計算時，坡長指數採用與坡度有關的變值更為合理。

排水系統是現代城市中重要的市政基礎設施，隨著合流制排水體制帶來的水體污染問題越來越受到人們的重視，新建或改建的排水系統多採用分流排水體制。其中雨水管網系統建設投資大，關係到民生安全，其設計工作不容忽視。匯水區域的地形決定了雨水管網的布置形式及每段管線的服務面積、管徑選擇、雨水口的形式和布置等一系列問題。在《室外排水設計規範》中並沒有對雨水設計流量公式中的匯水面積給出計算方法。通常的劃分原則是當匯水區有適當的地形坡度時，依照雨水匯入低側的原則，按照地面雨水徑流的水流方向劃分匯水面積，將雨水干管布置在地形低處或溪谷線上。當地形平坦時，則根據就近排除的原則，把匯水面積按周圍管渠的布置用等分角線劃分，干管布置在排水流域的中間，儘可能擴大重力流排除雨。也有研究採用其他的劃分方式，如水的範圍Thiessen 多邊形法、面積管長比法用中以目視估計和手工劃分為主，難以準確地體現地形對匯流的影響或直接忽略了這種影響，劃分精度和效率較低，容易受人為因素影響。為準確體現地形對匯流的影響，筆者提出了在 GIS 環境下基於數字高程模型( DEM) 的匯水面積劃分方法。DEM 數據是地理信息系統的基礎數據，自產生以來在很多領域都得到了廣泛應用 也是目前用於流域地形分析的主要數據。通過 DEM 可提取大量的地表形態信息，如流域柵格單元的坡向、坡度及單元格之間的關係等，基於 DEM 提取流域的水文特徵也發展了多種較成熟的演演算法。但目前基於 DEM 的水文特徵提取主要用於天然流域中河網及相應匯水區的分析研究，用於城市環境下的情況較少，雨水管網作為人工設施，每個管段的匯水區形狀及面積除了受地形影響外，同時也受城市規劃布局的約束，尤其是道路布局的影響。筆者利用 DEM 提取流域水文特徵的原理和方法，結合城市雨水系統的布置形式與工作特點，提出了一種新的劃分方法。

1.1 DEM 數據中窪地的處理

窪地是指低於周圍柵格的區域，分為偽窪地和自然窪地。偽窪地在 DEM 數據中非常普遍，主要來自輸入數據的錯誤、不合適的插值方法和柵格大小等方面。自然窪地則是實際中存在的窪地 較小的如地面坑窪，較大的如湖泊、蓄洪設施等。無論哪種類型的窪地，在流向分析時都會造成水流在窪地彙集而無法流出的現象，影響匯水面積分析的準確性，需要作填窪處理。但不同類型的窪地要分別對待，將由於數據採集誤差產生的偽窪地和蓄水量不足以影響整體分析的小型自然窪地進行填窪處理，對於大型自然窪地，需要在填窪操作中覆蓋窪地範圍圖層( 例如湖泊的匯水邊界) ，對圖層下的柵格不進行填窪處理，避免產生與實際情況不符的結果。

1.2 排水路線

道路通常是街區內地面徑流的集中地，也是雨水管線的定線基礎，城市規劃中一般是將雨水管渠鋪設在路面以下，在沒有鋪設雨水管渠的路段，路面本身也承擔著彙集和輸送雨水徑流的作用。因此，在分析匯流面積時將道路和管渠共同組成的雨水排水網路作為研究對象，稱之為排水路線。排水路線在空間上與道路中心線一致。排水路線具有坡向，坡向為排水路線上雨水的實際流動方向，當道路鋪設有排水管渠時，以管渠坡向作為排水路線坡向，無管渠時則以道路縱坡向為排水路線坡向。最終的分析結果即為每段排水路線的匯水面積。

1.3 道路匯水面積

因為道路邊溝通常低於相鄰街區的地面標高，雨水在匯入道路之後、發生溢流之前不會再次流出路面，作用類似於管渠，所以要對路面柵格高程加以修改以模擬這一作用。道路數據通常是採集道路中心線而產生的線性矢量數據，本身沒有路面寬度信息，路面寬度作為一項屬性存儲在關聯的屬性表內。為了在DEM 中描繪路面寬度，以排水路線為中心線，單側路面寬度為緩存半徑作緩衝分析，緩衝區的柵格即為路面柵格。首先將排水路線下的路面中心柵格降低一定高程值，兩側緩存區內的柵格高程相應降低，在中心柵格和邊緣柵格間形成一個連續的斜面，這樣就準確地描繪了道路匯水面積，並使路面雨水彙集到一處。

1.4 基於DEM的流向計算

在對窪地和道路處理完之後，進一步作流向分析。計算水流流向的演演算法有多種 其中D8 流向分析演演算法是較早提出並得到廣泛應用的一種實用演演算法，該演演算法是通過比較中心單元格與相鄰8 個單元格間的高程大小與落差，將高程下降最大的方向視為該單元格的流向，產生流向柵格分析結果。從技D8 演演算法來分析匯水區的地表徑流方向。在流向處理時同樣需要剔除影響匯流的大型自然窪地，方法與填窪處理時相同。

1.5 流向柵格的修正

由於在道路匯水面積處理過程中只是將沿排水路線的中心柵格降低了同一高度，沿排水路線的地形依然呈高低起伏，相應位於排水路線下方的流向單元柵格還不能準確地反映排水路線上水流的流動方向，沿線流向不一致，所以需要修正這些流向柵格值使水流沿排水路線坡向流動。

排水路線坡向為道路數據或管渠數據採集時的數字化方向，例如，數字化過程中以A 端為起點、B端為終點，則該段排水路線由A 坡向B。所以數字化時要以排水路線的高端為起點、低端為終點，使數字化的方向與路面上或管段中的水流方向一致，並按坡向分段進行數字化。

流向修正僅僅對排水路線下的單元柵格的流向作了修正，而兩側緩衝區內流向並沒有修正，緩衝區內徑流依舊匯入到排水路線下中心柵格部分。

1.6 匯水區的生成

在對排水路線下流向柵格修正之後，在流向柵格圖層上就形成了連續的類似天然河網的排水路線，根據流向柵格和排水路線即可分析獲得每段排水路線兩側的街區匯水面積。分析結果符合地表徑流沿最陡方向流動的自然現象，同時各路段路面上的降水也都在本段彙集。

通過以上分析可獲得各段排水路線兩側街區的匯水面積，而在實際情況中，無論是道路還是管渠排水，都不只是承擔其兩側街區的匯水量，而是同時擔負著輸送上游來水的作用，所以需要分析各段排水路線的上游匯水面積，才能準確計算各段雨水徑流量。分析上游匯水面積需要明確各段排水路線間的連接情況，才能在此基礎上根據排水路線的坡向分析得到總匯流面積。

在管網系統龐大、管段眾多、布置複雜的情況下，依靠人工去判斷某段排水路線上游的匯水面積是一項費時費力的工作。為了使這項工作實現自動化，定義排水路線的連接情況，即上游路線的末端與下游路線的起端為同一點，這樣在排水路線的交匯處，上游排水路線上的徑流就可自動流入與其末端相連的下游排水路線，由此分析得到與各段路線的起點連通的所有排水路線，將它們各自的匯流面積組合即為該排水路線上游的匯水面積，這些信息將自動存在屬性數據表格中，如此就可以輕鬆求解某管段的總匯水面積大小及組成。

利用目前GIS 領域中基於DEM 的水文特徵提取演演算法的成果，採用了一種新的劃分匯水面積的方法: 首先確定排水路線，然後分析流向並獲得各段排水路線相應的匯水面積，並考慮排水路線的連接情況，實現總匯水面積的自動化計算與劃分。分析過程自動化程度高、結果可靠，比傳統方法準確、快速。整個過程在GIS 環境中實現，有效利用了GIS 強大的資料庫管理功能和圖形顯示功能，對獲得的劃分結果可以方便地提取匯水區特徵參數，如面積、最長匯流路徑、特徵寬度等，和其他基於計算機的優化及

分析方法之間的數據交換功能也更容易實現，成果表現更為直觀明了，有助於城市雨水管網的規劃與設計工作更加精確、直觀、快速，具有很強的實用性。