本發(fā)明涉及土壤學(xué)與水土保持科學(xué)，特別是涉及一種基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法。

背景技術(shù)：

1、輸沙過(guò)程是土壤學(xué)和水土保持科學(xué)的重要研究?jī)?nèi)容。流域侵蝕輸沙很大程度上影響著土壤、營(yíng)養(yǎng)物質(zhì)以及污染物等的遷移，造成水土流失區(qū)域土地肥力降低、農(nóng)業(yè)生產(chǎn)受到污染、區(qū)域生態(tài)環(huán)境被破壞等嚴(yán)重后果。因此，輸沙過(guò)程的精確化模擬具有重要的科學(xué)價(jià)值和應(yīng)用價(jià)值。

2、輸沙率和含沙量是輸沙過(guò)程中的兩個(gè)重要指標(biāo)，因此，輸沙過(guò)程的精細(xì)化模擬需要進(jìn)行輸沙率和含沙量的精細(xì)化模擬。流域作為天然的地理實(shí)體，成為水沙治理等工程措施的重點(diǎn)實(shí)驗(yàn)單元。前人以流域?yàn)閱卧獙?duì)流域侵蝕產(chǎn)沙機(jī)制展開(kāi)了建模與研究，研究中將流域內(nèi)部當(dāng)作“黑箱”處理，僅以某幾個(gè)點(diǎn)位的結(jié)果或均值代替整體的產(chǎn)流產(chǎn)沙結(jié)果，忽略了流域內(nèi)部的泥沙輸移過(guò)程，缺乏對(duì)泥沙輸移過(guò)程模擬的研究。

3、近年來(lái)，隨著無(wú)人機(jī)攝影測(cè)量、雷達(dá)測(cè)量等技術(shù)的發(fā)展，獲取高精度地形數(shù)據(jù)越來(lái)越方便快捷，有不少學(xué)者開(kāi)始使用高精度地形數(shù)據(jù)計(jì)算地形變化，以像元為單位，使用流向算法對(duì)流域進(jìn)行徑流輸沙路徑模擬，模擬輸沙率空間分布。但尚未有學(xué)者使用流向算法對(duì)流域內(nèi)部的含沙量分布進(jìn)行精細(xì)化模擬，因此亟需結(jié)合輸沙率、集水口徑流量、流域內(nèi)部產(chǎn)流量等水文數(shù)據(jù)，計(jì)算含沙量空間分布。

技術(shù)實(shí)現(xiàn)思路

1、為了解決以上技術(shù)問(wèn)題，本發(fā)明提供一種基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，包括以下步驟：

2、s1、獲取研究區(qū)域時(shí)序觀測(cè)的多個(gè)測(cè)量時(shí)期的原始dem和水文數(shù)據(jù)，水文數(shù)據(jù)包括平均降雨強(qiáng)度、降雨歷時(shí)、平均土壤容重、徑流時(shí)長(zhǎng)以及集水口實(shí)測(cè)徑流量；

3、s2、對(duì)原始dem數(shù)據(jù)進(jìn)行填洼處理，得到填洼后的dem數(shù)據(jù)；

4、s3、針對(duì)于每個(gè)測(cè)量時(shí)期，將相鄰期的原始dem數(shù)據(jù)相減，得到dem數(shù)據(jù)的地形變化量，進(jìn)一步乘以平均土壤容重轉(zhuǎn)化得到地形變化質(zhì)量；

5、s4、根據(jù)填洼后的dem數(shù)據(jù)，采用質(zhì)量守恒原理約束下的多流向算法，計(jì)算泥沙流向和對(duì)應(yīng)方向的輸沙量，并根據(jù)徑流時(shí)長(zhǎng)數(shù)據(jù)和輸沙量數(shù)據(jù)，得到輸沙率數(shù)據(jù)；

6、s5、根據(jù)集水口實(shí)測(cè)徑流量和預(yù)估總流量，計(jì)算得到流域內(nèi)的水流下滲率；

7、s6、根據(jù)水流下滲率，采用多流向算法，計(jì)算得到每個(gè)像元的徑流數(shù)據(jù)；

8、s7、根據(jù)步驟s4得到的輸沙率數(shù)據(jù)和步驟s6得到的徑流數(shù)據(jù)，計(jì)算得到含沙量空間分布數(shù)據(jù)。

9、本發(fā)明進(jìn)一步限定的技術(shù)方案是：

10、進(jìn)一步的，步驟s1中，原始dem數(shù)據(jù)是m×n的二維數(shù)組，一個(gè)測(cè)量時(shí)期的dem數(shù)據(jù)包括多個(gè)中心像元，一個(gè)中心像元對(duì)應(yīng)8個(gè)與該中心像元相鄰的鄰域像元，領(lǐng)域像元如下式所示：

11、

12、其中，第ri行第rj列的元素dem(ri,rj)表示中心像元處的高程值；使用bk表示中心像元(ri,rj)的鄰域像元(rp,rq)，也即k連通區(qū)域，k=8。

13、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s2中，識(shí)別原始dem數(shù)據(jù)中的洼地，若某一像元的高程低于周圍8個(gè)像元，則被識(shí)別為洼地；對(duì)于每個(gè)洼地，計(jì)算每個(gè)像元的洼地深度，并通過(guò)將洼地點(diǎn)的高程值增加到周圍鄰域的最高高程值來(lái)填充洼地；

14、每個(gè)像素的洼地深度通過(guò)將周圍鄰域的最高值減去該像素的高程值來(lái)計(jì)算：

15、?，

16、且不是空洞，

17、其中，d(ri,rj)表示中心像元(ri,rj)處的洼地深度，dem(ri,rj)表示原始的中心像元(ri,rj)的高程值，dem(rp,rq)表示原始的鄰域像元(rp,rq)的高程值，bk表示鄰域像元bk的集合，(m,n)分別表示原始dem數(shù)據(jù)的二維數(shù)組的數(shù)值。

18、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s2中，使用洼地深度的值來(lái)提高洼地的高程，如下式所示：

19、

20、其中，hfilled(ri,rj)表示中心像元(ri,rj)處填充后的高程值，h(ri,rj)表示中心像元(ri,rj)處填充前的原始高程值，d(ri,rj)表示中心像元(ri,rj)處的洼地深度。

21、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s3中，使用柵格計(jì)算，令前一測(cè)量時(shí)期的原始dem數(shù)據(jù)減去后一測(cè)量時(shí)期的原始dem數(shù)據(jù)，計(jì)算相鄰期原始dem數(shù)據(jù)之間的地形變化量；根據(jù)地形變化量，通過(guò)預(yù)先獲取的dem分辨率和平均土壤容重，得到地形變化質(zhì)量，如下式所示：

22、

23、其中，表示第n個(gè)測(cè)量時(shí)期的第r個(gè)dem像元處的地形變化質(zhì)量；表示第n個(gè)測(cè)量時(shí)期的第r個(gè)dem像元處的地形變化量，通過(guò)兩個(gè)測(cè)量時(shí)期的原始dem數(shù)據(jù)相減獲取；x表示預(yù)先獲取的dem分辨率；表示平均土壤容重，p表示孔隙率。

24、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s4中，采用質(zhì)量守恒原理約束下的多流向算法進(jìn)行路徑分配，得到中心像元向鄰域像元的輸沙量，并根據(jù)徑流時(shí)長(zhǎng)數(shù)據(jù)和輸沙量數(shù)據(jù)，得到輸沙率數(shù)據(jù)，如下式所示：

25、

26、其中，表示第n測(cè)量時(shí)期中心像元向第k個(gè)鄰域像元r輸送泥沙量的輸沙率；k表示與第r個(gè)中心像元相鄰的鄰域像元的序號(hào)；tansknr表示坡度，上標(biāo)為控制流向分散或集中的指數(shù)，為1時(shí)表示多流向，為無(wú)窮大時(shí)表示單流向；lknr表示第n測(cè)量時(shí)期第k號(hào)鄰域像元的等高線長(zhǎng)度加權(quán)因子；表示第n測(cè)量時(shí)期第r個(gè)中心像元處的地形變化質(zhì)量；qnr表示鄰域像元向第r個(gè)中心像元輸送的總泥沙量，；tknr表示第n測(cè)量時(shí)期第k個(gè)鄰域像元向第r個(gè)中心像元輸送的泥沙量。

27、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s4中，當(dāng)?shù)趎測(cè)量時(shí)期中與第r個(gè)中心像元相鄰的第k號(hào)鄰域像元高于第r個(gè)中心像元時(shí)，lknr=0；當(dāng)?shù)趎測(cè)量時(shí)期中與第r個(gè)中心像元相鄰的第k號(hào)鄰域像元低于第r個(gè)中心像元且位于水平或垂直方向時(shí)，；當(dāng)?shù)趎測(cè)量時(shí)期中與第r個(gè)中心像元相鄰的第k號(hào)鄰域像元低于第r個(gè)中心像元且位于對(duì)角線方向時(shí)，。

28、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s5中，假設(shè)各個(gè)位置下滲率一致，流域內(nèi)的水流下滲率表達(dá)式如下：

29、

30、其中，rn表示第n個(gè)測(cè)量時(shí)期的流域內(nèi)水流下滲率；表示集水口實(shí)測(cè)徑流量；表示流域內(nèi)預(yù)估總流量，計(jì)算方法如下式所示：

31、

32、其中，r表示降雨平均強(qiáng)度，tr表示降雨歷時(shí)，sc表示流域面積。

33、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s6中，每個(gè)像元的徑流數(shù)據(jù)計(jì)算公式如下：

34、

35、其中，qsknr表示第n個(gè)測(cè)量時(shí)期第r個(gè)中心像元向第k個(gè)與該中心像元相鄰的鄰域像元輸送的徑流量；rn表示第n個(gè)測(cè)量時(shí)期的流域內(nèi)水流下滲率；qsnr表示鄰域像元向第r個(gè)中心像元輸送的總徑流量，；tsknr表示第n個(gè)測(cè)量時(shí)期第k個(gè)鄰域像元向第r個(gè)中心像元輸送的徑流量；qsr表示第n個(gè)測(cè)量時(shí)期第r個(gè)中心像元上降雨產(chǎn)生的流量。

36、如前所述的基于時(shí)序dem的流域徑流含沙量空間分布計(jì)算方法，步驟s7中，含沙量空間分布數(shù)據(jù)計(jì)算公式如下：

37、

38、其中，cnr表示第n個(gè)測(cè)量時(shí)期第r個(gè)中心像元上的含沙量；表示第n個(gè)測(cè)量時(shí)期中心像元向第k個(gè)鄰域像元r輸送泥沙量的輸沙率；qsknr表示第n個(gè)測(cè)量時(shí)期中心像元向第k個(gè)鄰域像元r輸送的徑流量。

39、本發(fā)明的有益效果是：

40、本發(fā)明中，基于質(zhì)量守恒原理，采用多流向算法，對(duì)泥沙輸移路徑進(jìn)行模擬，得到輸沙率空間分布；再通過(guò)該多流向算法對(duì)水流輸移路徑進(jìn)行模擬，得到徑流量空間分布，通過(guò)計(jì)算得到含沙量空間分布，能夠提高對(duì)含沙量的監(jiān)測(cè)精度，用于流域水土流失的監(jiān)測(cè)，以評(píng)價(jià)區(qū)域的生態(tài)環(huán)境，為生態(tài)環(huán)境治理提供重要參考。