摘要：水文分析實戰(zhàn)教程河流平均比降計算本章導(dǎo)讀河流比降指的是河流任意兩端點間的高程差與兩點間的水平距離之比值，簡單的概括為單位河長內(nèi)的落差，其單位一般都是。讀者可以參考水文分析實戰(zhàn)教程河流提取與河網(wǎng)分級的做法去提取。

本章導(dǎo)讀：河流比降指的是河流任意兩端點間的高程差與兩點間的水平距離之比值，簡單的概括為單位河長內(nèi)的落差，其單位一般都是‰ 。 河流（縱）比降對地質(zhì)災(zāi)害、水運、水電、礦砂和新構(gòu)造運動研究有實用價值。比降是水文學(xué)上重要的一項指標(biāo)。在絕大部分的GIS通用平臺上，都沒有實現(xiàn)平均比降的計算，原因在于計算比降的條件比較苛刻，不符合通用GIS的要求。本章同樣是一個實操章節(jié)，筆者將會利用ArcGIS軟件，結(jié)合其ArcPy編寫一個平均比降計算的工具，在通用GIS軟件上實現(xiàn)河流平均比降計算。 BY 李遠(yuǎn)祥

河流比降是很多水文研究者研究的內(nèi)容，其水文意義非常的大。這里研究的河流比降指的是河流的縱比降，是河流任意兩端點間的高程差與兩點間的水平距離之比值，即單位河長的落差。

正是因為直接表達(dá)的是單位河長的落差，所以在水文學(xué)上其最直接的影響就是流速，帶來更多的水動力引發(fā)的能量轉(zhuǎn)移，如勢能與動能的轉(zhuǎn)化（涉及到水電站的建設(shè)）；河水對河床和河岸的沖刷，影響泥沙量；大面積持續(xù)降雨時，洪峰到來的時間和強度等。所以，比降在水文學(xué)上從來就不缺少話題。

河流的比降分為床面比降和水面比降。床面比降，用以表示河床縱斷面地形的變化；水面比降即河流中任意兩端點間的瞬時水面高程差與其相應(yīng)水平距離之比，用以表明河流全程或分段的水面坡度，故又稱水力坡度，通常說的河流比降就是河流水面比降。

在恒定均勻流情況下，水面坡度恰好等于河底坡度。所以，在沒有大量斷面數(shù)據(jù)支撐的情況下，根據(jù)地形和河流線型去研究的河流比降都是床面比降。

筆者在這里研究的都是床面比降，因為只考慮地形與河流本身的幾何形狀，沒有任何的站點數(shù)據(jù)支撐。但在非洪水情況下，這個床面比降與水面比降基本上是吻合的。

在河流相距較短的兩個點上計算河流比降，意義不是十分大。一般需要計算的是整個河流或者某一河段的平均比降，這對于空間跨度較大的河流來說才有研究價值。

本章研究的內(nèi)容也就是河流（段）的平均比降。研究河流平均比降的方法有幾種，比較流行的是【約翰斯通-克羅斯】法。其原理是根據(jù)比降大致均勻的原則，將河段劃分為若干個子河段，然后去計算整體的平均比降。其公式如下：

公式中個參數(shù)如下
S---- 河段平均河底比降
Li----第i個子河段長度
Si----第i個子河段河底比降
N---- 子河段數(shù)量

另一種方法是【勞倫斯】法，其需要經(jīng)驗系數(shù)和流速，相對來說比較復(fù)雜，需要結(jié)合測站的一些數(shù)據(jù)。由于筆者沒有這些數(shù)據(jù)，就不作論述。

本章考慮的是在ArcGIS中實現(xiàn)【約翰斯通-克羅斯】法的計算。

利用GIS手段來研究水文，其核心方式是使用DEM數(shù)據(jù)去模擬水文過程。GIS的空間疊加能力非常的管用，非常簡單的就能獲取到河流線對應(yīng)的DEM像元值，也就是河流對應(yīng)的高程值，這樣用于計算比降是非常有用的。最起碼解決了河流高程點獲取的問題。

河流沿線高程的獲取，不能直接使用原始的DEM數(shù)據(jù)。原因有兩點：

如果河流數(shù)據(jù)與DEM數(shù)據(jù)不是同一份數(shù)據(jù)，那么像元點與河流線的重合程度有待商榷。

如果河流數(shù)據(jù)與DEM數(shù)據(jù)吻合（或者河流數(shù)據(jù)直接從DEM數(shù)據(jù)提取出來的），其高程的獲取所依賴的DEM必須是填洼后的無凹陷點DEM。

河流線屬于一維數(shù)據(jù)，DEM在平面上屬于二維數(shù)據(jù)，其相交之后會產(chǎn)生一個像元點的連續(xù)像元。如果了解過D8單流向算法的話，就可以很容易理解，如果直接從原始DEM中獲取像元值，其值有可能是一個匯點（凹陷點）。然而使用填洼后的DEM，可以確保該點與周邊一點的高程持平。

下圖是一個測試數(shù)據(jù)，黃色部分為河流線，點為河流的節(jié)點。包含兩份DEM，一份是原始的DEM，一份是填洼過的DEM(Fill_DEM) 。河流是使用了前面章節(jié)介紹到的河流提取方法提取出來的線狀河流。讀者可以參考《ArcGIS水文分析實戰(zhàn)教程（6）河流提取與河網(wǎng)分級》的做法去提取。
測試的主要內(nèi)容是分兩種情況打印河流各節(jié)點的高程值。一是獲取原始DEM的高程值，而是獲取填洼后的DEM的高程值。通過高程值進行對比選用不同數(shù)據(jù)引起的差異。

第一張圖是原始DEM的沿河高程值，第二張圖是填洼后的DEM的沿河高程值?？梢钥吹绞褂迷糄EM的時候，有可能出現(xiàn)下游高程比上游高程還要高的情況，也就是說這個上游地勢低的點剛好就是一個匯（凹陷點）。而填洼后的DEM則會將這個點給填平，確保不會斷流，從而使得下游高程最少不會高于上游高程。對比其他高程值，是沒有變化的，變化的部分就只是出現(xiàn)【匯】的點。

如果原始DEM的精度足夠高的話，這種情況一般很少碰到。但為了穩(wěn)妥起見，筆者還是建議對原始DEM做填洼，生成無凹陷點的DEM。讀者可以參考之前章節(jié) 《ArcGIS水文分析實戰(zhàn)教程（4）地形預(yù)處理》

如果所使用的河流是從DEM中提取出來的，那么在矢量河流提取的時候使用【柵格河流矢量化】時是否需要使用【簡化折線】參數(shù)，這個參數(shù)是默認(rèn)選上的，如下圖所示：

不選擇簡化折線，河流線會基本上就是按照流量線走向，河流折點密度會比較高；選擇簡化的話，河流線的折點會被抽希一部分。這個怎么取舍？如果在地形變化比較大的情況下，不建議簡化；如果地勢趨平，那就直接簡化線。這種方式會直接影響到最終的比降數(shù)值。

下圖黃線部分是簡化節(jié)點后的河流線，紅線是沒有簡化節(jié)點的河流線，其形狀基本上是重合的，而且在關(guān)鍵的河流節(jié)點上，也是完全重合。

也就是說，簡化與不簡化，其參與計算的河段數(shù)量是不一樣的。上圖數(shù)據(jù)中簡化后的子河段數(shù)為32，不簡化子河段的數(shù)量是92 。如果人為計算的話，那么可定是后者的計算量比較大了。

ArcGIS 除了前面章節(jié)提及到的水文分析工具及其相關(guān)工具之外，還有非常龐大的工具集對空間計算進行支撐。由于沒有斷面數(shù)據(jù)的支撐，這里只考慮通過地形和河流線型進行計算，也就是說實現(xiàn)【約翰斯通-克羅斯】法。

根據(jù)【約翰斯通-克羅斯】法的計算公式，可以推導(dǎo)出其分析流程

通過河流線提取河流折點，利用這些折點與DEM數(shù)據(jù)疊加，獲取高程值

利用河流折點，獲取其坐標(biāo)值，根據(jù)兩點距離公式，計算每一段子河段的長度

根據(jù)每個點的高程值，計算出相鄰兩點之間的高差

根據(jù)高差與子河段長度，計算出各子河段的河底比降

匯總這些數(shù)據(jù)計算出整個河段的平均河底比降

在這5各步驟中，只有第1步與ArcGIS工具有關(guān)。河流的折點提取，可以使用【數(shù)據(jù)管理】--【要素】--【要素折點轉(zhuǎn)點】工具來實現(xiàn)，如下圖所示：

至于利用折點與DEM數(shù)據(jù)疊加，獲取其高程值，可以使用【數(shù)據(jù)管理】--【柵格】--【柵格屬性】--【獲取像元值】工具來獲取。這個工具每次只能實現(xiàn)一個坐標(biāo)點與DEM數(shù)據(jù)疊加查詢。如下圖

所以，為了方便計算，筆者采用ArcPy去寫一個腳本工具，實現(xiàn)參數(shù)的靈活輸入。其中heightList就是用來保存高程值的數(shù)組。

    # 獲取河流的節(jié)點，以用作計算高差和河長
    for row in arcpy.da.SearchCursor(outLine1, ["SHAPE@XY"]):
        # 從點圖層數(shù)據(jù)獲取所有的點的XY
        x, y = row[0]
        #從地形數(shù)據(jù)中獲取其高程的象元值
        pointStr = str(x) +" "+str(y)
        result = arcpy.GetCellValue_management(dem,pointStr)
        riverPoints.append(row[0])
        heightList.append(str(result))

第2-5步基本上就是一些數(shù)學(xué)的運算，利用python就很好實現(xiàn)了。為了讓廣大水文研究者用起來，筆者將這段代碼公開，并且制作成為工具。在代碼中筆者已經(jīng)預(yù)留了ArcGIS的工具參數(shù)，方便制作成為腳本工具。

河流平均比降【約翰斯通-克羅斯】法的實現(xiàn)代碼如下

#############################################################
#                                                           #
#   該工具利用【約翰斯通-克羅斯】方法計算河底平均比降           #
#   在使用的時候必須先選擇一條河流，工具會判斷輸入的河流數(shù)量     #
#   大于一條河流時不予以計算                                  #
#   用于計算的DEM為無凹陷點DEM                                #
#   BY 李遠(yuǎn)祥                                                #
#                                                            #
##############################################################


# -*- coding: gbk -*-
import arcpy
import math
inputLine =arcpy.GetParameterAsText(0) #河流圖層

dem = arcpy.GetParameterAsText(1) #無凹陷點DEM
riverLen = 7873.232038

arcpy.AddMessage(u"..................")
arcpy.AddMessage(u"......開始執(zhí)行...........")

# 記錄高程值的數(shù)據(jù)
heightList = []

iputCout=0
mmCursor = arcpy.da.SearchCursor(inputLine, ["SHAPE@LENGTH"])
for row1 in mmCursor:
    riverLen = row1[0]      #順便獲取河流長度
    iputCout=iputCout+1

arcpy.AddMessage(u"......獲取了河流長度......")
#如果選中的是一個線要素，則執(zhí)行計算，否者不計算多條河流
if iputCout==1:
    ##先將線轉(zhuǎn)為折點
    outLine1 = arcpy.FeatureVerticesToPoints_management(inputLine, "in_memory/outline")
    arcpy.AddMessage(u"......河流線已經(jīng)轉(zhuǎn)換為點集.......")

    #河流節(jié)點
    riverPoints = []

    # 獲取河流的節(jié)點，以用作計算高差和河長
    for row in arcpy.da.SearchCursor(outLine1, ["SHAPE@XY"]):
        # 從點圖層數(shù)據(jù)獲取所有的點的XY
        x, y = row[0]
        #從地形數(shù)據(jù)中獲取其高程的象元值
        pointStr = str(x) +" "+str(y)
        result = arcpy.GetCellValue_management(dem,pointStr)
        riverPoints.append(row[0])
        heightList.append(str(result))

    arcpy.AddMessage(u"......已經(jīng)獲取了折點的高程值......")
    #子河段長度

    subRiverList = []
    hDiffList = []
    # 計算每一個子河段的長度
    for i in range(0,len(riverPoints)-1):
        x,y = riverPoints[i]
        x1,y1 = riverPoints[i+1]
        #計算河段長度
        dis =math.sqrt(pow((x1-x),2)+pow((y1-y),2))
        subRiverList.append(dis)

    arcpy.AddMessage(u"......河流子河段長度計算完成完成.......")

    # 計算從高程從低到高部分的差距
    listLen = len(heightList)
    for j in range(0,listLen-1):
        h1 = float(heightList[j])
        h2 = float(heightList[j+1])
        xheight = 0
        xheight = h1-h2
        hDiffList.append(xheight)

    arcpy.AddMessage(u"......各折點的高程差計算完成......")

    # 計算每個河段的比降
    subbijiang = []
    for i in range(len(hDiffList)):
        Si= hDiffList[i]/subRiverList[i]
        subbijiang.append(Si)

    arcpy.AddMessage(u"......子河段比降計算完成......")
    # 計算最終的比降
    # 計算分子
    sumUp = 0
    for i in range(len(hDiffList)-1):
        sumUp = sumUp + subRiverList[i]*pow(subbijiang[i],0.5)

    S = pow((sumUp/riverLen),2)

    arcpy.AddMessage(u"......該河流的總長為： " + str(riverLen)+u"   平均比降為: " +str(S))
    arcpy.AddMessage(u"......計算完畢.......")

else:
    arcpy.AddMessage(u"......您選中了 "+str(iputCout) +  u"  條河流......")
    arcpy.AddMessage(u"......工具已經(jīng)停止計算,沒有算出您想要的結(jié)果......")
    arcpy.AddMessage(u"......因為您目前選中了多條河流。只能選中一條河流進行計算!......")

從代碼中可以看到，這個工具針對的是一條河流進行計算的。如果想指定河流的某一段進行計算，那么就需要將河流轉(zhuǎn)為折點，選擇對應(yīng)的折點，生成新的河段才能參與計算。這個過程需要用到【點集轉(zhuǎn)線】工具進行重構(gòu)要分析的河段。

可以將上面的python代碼制作成為腳本工具，方便在ArcGIS中進行可視化的調(diào)用。工具的制作可以參考筆者另一篇博文《ArcMap自定義腳本工具制作》，里面有詳細(xì)的論述。

該工具用到的ArcGIS工具涉及到ArcInfo級別 和 Spatial 擴展模塊，請先確保ArcGIS軟件具備相應(yīng)的授權(quán)。

制作好工具之后可以雙擊調(diào)用該工具了，【約翰斯通-克羅斯】法需要用到無凹陷點DEM和河流線型數(shù)據(jù)。

從上面代碼可以看到，工具會判斷選中的河流數(shù)量，如果數(shù)量大于1，就不作計算，確保計算的是某一河段。ArcGIS的交互性非常的優(yōu)秀，只要選中某一條河流，就會將選中的部分作為參數(shù)輸入。

對應(yīng)工具的流程，這里第一步需要選中要分析的一條河流。如下圖

接著雙擊調(diào)用工具，兩個參數(shù)分別是河流圖層和無凹陷點DEM。一定要注意，這個DEM是填洼過的DEM，不是原始DEM數(shù)據(jù)。如下圖

點擊確定就可以執(zhí)行計算結(jié)果，結(jié)果會在最終的提示框中找到，如下圖所示：

至此，復(fù)雜的河流平均比降就計算完畢了！

水文分析是一個復(fù)雜的過程，需要有深厚的水文學(xué)基礎(chǔ)作為背景。而GIS是一個通用工具，提供非常龐大的分析功能。二者要結(jié)合起來，不但需要學(xué)習(xí)兩個專業(yè)的知識，還需要具備知識整合的能力。

代碼和制作ArcGIS工具的方法已經(jīng)在文中貼出來了。讀者可以自行參考制作。當(dāng)然想省事的話，筆者也已經(jīng)上傳制作好的工具，可以點擊跳轉(zhuǎn)到 下載頁面 。為了加強交流，防止伸手，多動手操作，筆者已經(jīng)將下載分值設(shè)置為最高分10分。 積分不夠的讀者請多動手制作 ^_^ 
如果讀者對水文分析感興趣的話，可以持續(xù)關(guān)注CSDN的GIS制圖樂園，以及微信公眾號【GIS制圖樂園】。BY 李遠(yuǎn)祥