CN110717232A - 一种考虑梯田淤地坝影响的子流域汇流模拟方法 - Google Patents

Abstract

本发明公开了一种考虑梯田淤地坝影响的子流域汇流模拟方法，涉及水文模型模拟技术领域；包括如下步骤：采用子流域套等高带的方式作为模型基本计算单元，基于DEM将子流域根据高程划分成数目不等的等高带并统计相关参数，提取沟道栅格，设置沟道长度阈值，剔除伪沟道；对等高带内的梯田和淤地坝进行概化，并计算梯田对坡面径流的拦截比例，并统计淤地坝相关参数；基于运动波方程模拟各等高带坡面、单条沟道和河道汇流过程。该方法能够模拟计算单元内淤地坝建设和梯田建设对汇流过程的影响。

Description

一种考虑梯田淤地坝影响的子流域汇流模拟方法

技术领域

本发明涉及分布式水文模型子流域汇流模拟技术领域，尤其涉及一种考虑梯田淤地坝影响的子流域汇流模拟方法。

背景技术

分布式水文模型是研究水循环过程和机理的有效手段，也是解决许多水文实际问题的有效工具。随着研究问题的复杂程度增加，模型功能也越来越复杂。在黄土高原山区，梯田分布广泛，沟道纵横，为了有效利用雨水增加粮食产量，该地区在沟道上修建了大量的淤地坝，用于拦蓄山区洪水沉积泥沙，淤地造田。梯田淤地坝的建设对流域汇流过程影响巨大。梯田和淤地坝区域平坦，面积相对较大，能够有效拦截上游径流，增加局部区域含水率，改变了天然水循环过程，减小汇流流速，使得泥沙沉积。为了进一步提高子流域汇流过程模拟精度，评估梯田和淤地坝建设对水循环过程的影响，有必要在模型模拟中增加对梯田和淤地坝过程的模拟。

发明内容

本发明的目的在于提供一种考虑梯田淤地坝影响的子流域汇流模拟方法，从而解决现有技术中存在的前述问题。

为了实现上述目的，本发明采用的技术方案如下：

一种考虑梯田淤地坝影响的子流域汇流模拟方法，包括如下步骤：

S1，采用子流域套等高带的方式作为模型基本计算单元，确定子流域范围及等高带范围；包括坡度、长度、宽度等等；

基于DEM对子流域和等高带内的栅格进行分析统计，并获取子流域、等高带相关参数；其中，子流域是地表水文过程相对独立的区域单元，等高带是子流域内按高程大小分割的数目不等的长方形概化单元。

S2，基于DEM提取等高带内的模拟沟道，并统计模拟沟道的数量及平均长度；根据子流域内各栅格最终汇入的河道沟道栅格类型，统计沟道和河道汇流栅格个数，并计算沟道河道径流分配系数；

具体实现步骤为：减小阈值提取模拟河网，其中非河道栅格作为备选沟道栅格，并计算各直接汇入河道栅格的备选沟道的最长汇流长度；

设定沟道长度阈值，剔除沟道长度低于阈值的伪河道，统计概化沟道的数量以及评价长度。根据子流域内各栅格最终汇入的河道沟道栅格类型，统计子流域内栅格汇流属性，并统计各等高带内坡面沟道河道径流分配系数；

S3，基于统计资料或遥感资料，分析确定各等高带计算单元内的梯田面积；根据各等高带内的梯田面积确定当前等高带内的梯田对上一等高带的坡面径流拦截比例，即认为该部分坡面径流完全进入当前等高带的梯田内，参与当前等高带梯田的产汇流过程；

S4，设定淤地坝只存在于沟道上，且将每个沟道上的所有淤地坝概化为一个等效淤地坝，并假定每条沟道上淤地坝大小一样且位于沟道出口处；统计淤地坝相关参数；

淤地坝作为一个小型水库拦蓄上游沟道汇流径流量，在淤地坝区域内蒸发下渗等调蓄后排入河道；淤地坝所需的参数包括每个等高带内的淤地坝面积、最大洼地储留深、沟道径流拦截比例、淤地坝透水系数。

S5，基于运动波方程模拟各等高带坡面汇流过程，其中各等高带上游流入量为上一等高带坡面径流乘以对应等高带坡面径流分配系数减去被当前等高带梯田拦截的部分，侧向汇入量等于当前等高带非梯田部分产流量加上梯田调蓄后的产流量；梯田调蓄过程主要采用水量平衡法进行模拟；

S6，基于运动波方程模拟单条沟道汇流过程，其中侧向汇入量为各等高带坡面径流乘以对应等高带沟道径流分配系数除以沟道数量。对模拟计算的单条沟道汇流量进行淤地坝调蓄模拟，即淤地坝先拦蓄一部分沟道径流，加入淤地坝储量中，再根据淤地坝水量平衡模拟出流量，从而实现淤地坝对沟道径流量的调节作用。

S7，采用运动波方程对子流域河道进行汇流模拟，此时需要考虑上游子流域河道汇入量、当前子流域侧向进入量(包括各等高带坡面径流直接进入量，以及经淤地坝调蓄后的沟道径流进入量)。重复步骤S4-S6，直到所有时间、所有子流域汇流模拟结束。

优选地，步骤S1中：

模型基本单元是子流域套等高带，其中，子流域是地表水文过程相对独立的区域单元，主要用于沟道、河道汇流运算；等高带是子流域内部的细分，用于反映高程变化的差异影响，主要用于产流、坡面汇流运算；

确定等高带范围包括：首先确定等高带个数，再依据等高带的个数确定等高带范围。

优选地，子流域范围确定后，可通过两种方法确定等高带的个数；

方法一：给定最大等高带个数和最小等高带面积阈值参数计算获得；首先采用子流域面积除以最小等高带面积阈值确定等高带个数，再同最大等高带个数求最小值得到最终的等高带个数；

方法二：给定等高带平均面积参数计算获得；采用子流域面积除以等高带面积，向上取整得出最终的等高带个数；

根据等高带个数确定等高带范围的方法，包括对子流域内所有栅格按高程降序排列，获得最高和最低高程值，采用等高带个数将这个高程区间进行等分，按高程从高到底依次从1开始自然数编号；遍历子流域内栅格，根据栅格高程位于的高程区间分别赋值对应的区间编号，从而确定各等高带范围，即相同编号的栅格位于同一个等高带。

优选地，根据划分的等高带范围，以DEM为基础数据提取等高带的坡度、长和宽参数，其中，等高带坡度是等高带范围内所有栅格坡度的平均值，等高带宽是沿着坡面流向的长度，等高带长是垂直坡面流向的长度，长和宽计算公式如下：

wid_c＝(H_max-H_min)/sin(slope)

len_c＝A_c/wid

其中，wid_c为等高带宽，len_c为等高带长，H_max是等高带范围内的最大高程，H_min是等高带范围内的最小高程，slope是等高带平均坡度，A_c是等高带面积。

优选地，步骤S2中具体包括：

S21，由于提取的沟道存在着长度过小的情况，认为这些短的沟道属于提取出的伪沟道，需要加以剔除，因此直接由阈值提取出的沟道为备选沟道；

S22，遍历子流域内的河道栅格，统计流入河道栅格的备选沟道栅格作为备选沟道的出口；对每个备选沟道出口栅格，统计所有流入该栅格的上游备选沟道栅格，确定该备选沟道的最大汇流路径长度；

S23，设定一个沟道最小长度阈值，并认为长度小于该阈值的备选沟道为伪沟道；将伪河道栅格属性赋值为坡面栅格，统计确定该子流域内沟道条数N_g，并对所有沟道长度求平均值得到沟道的平均概化长度len_g；

S24，对流域内每个河道和沟道栅格进行溯源遍历确定各栅格的汇流属性；

如果沟道栅格位于最低等高带，则给所有汇入该沟道栅格的栅格赋值最低等高带沟道汇流属性(F1)，同理，对所有汇入最低等高带河道栅格的栅格赋值最低等高带河道汇流属性(F2)，对所有汇入其他等高带沟道栅格的栅格赋值其他等高带沟道汇流属性(F3)，对所有汇入其他等高带河道栅格的栅格赋值其他等高带河道汇流属性(F4)。

S25，统计沟道和河道汇流栅格个数，并计算等高带坡面、沟道、河道径流分配系数；

优选地，最低等高带径流分配系数计算公式如下：

r_p＝0

r_g＝N_F1/(N_F1+N_F2)

r_r＝N_F2/(N_F1+N_F2)

其他等高带径流分配系数计算公式如下：

r_p＝(N_F1+N_F2)/(N_F1+N_F2+N_F3+N_F4)

r_g＝N_F3/(N_F1+N_F2+N_F3+N_F4)

r_r＝N_F4/(N_F1+N_F2+N_F3+N_F4)

其中，r_p，r_g，r_r表示坡面、沟道、河道的径流分配系数；N_F1为该等高带内最低等高带沟道汇流属性栅格个数；N_F2为该等高带内最低等高带河道汇流属性栅格个数；N_F3为该等高带内其他等高带沟道汇流属性栅格个数；N_F4为该等高带内其他等高带河道汇流属性栅格个数。

优选地，步骤S3中：

等高带内的梯田面积确定可有以下两种方法确定：

(1)收集地市统计年鉴中的梯田面积，并根据土地利用类型，按山区旱地斑块进行空间展布，并同等高带计算单元进行相交处理，得到等高带内梯田的面积；

(2)通过遥感手段，基于DEM、遥感影像资料，根据一定的标准，解译得出流域范围内的梯田面积；

当前等高带梯田对上一等高带坡面径流拦截比例由以下公式计算：

r_tj＝(2.2409-0.4075*A_sj/A_C1)*A_s/A_C0

其中，r_tj为当前等高带梯田对上一等高带坡面径流拦截比例，A_sj为当前等高带内梯田面积，A_Cj为当前等高带面积，A_C(j-1)为上一等高带面积，j表示等高带编号；如果当前等高带是最上面的等高带，则r_t1＝0。

优选地，步骤S4中：

统计淤地坝相关参数包括每个等高带内的淤地坝面积、最大洼地储留深、沟道径流拦截比例、淤地坝透水系数；

每个等高带内的淤地坝面积可根据DEM和遥感影像，基于一定的准则解译出子流域范围内的淤地坝范围，并同各等高带范围进行叠合处理，从而得出各等高带内的淤地坝面积；

淤地坝最大洼地储留深可根据参考文献设置；

淤地坝沟道径流拦截比例用于表征淤地坝对沟道径流的拦截程度，范围在0-1之间，是反映淤地坝位置和溢流信息的一个参数，如果实际情况下淤地坝大多位于沟道下游，则该比例系数较大，如果实际情况下淤地坝大多位于沟道上游，则该比例系数较小，如果淤地坝有溢流设施，则该比例相对较小；

淤地坝透水系数用于表示通过坝体渗透的水量占储水量的比例；

将沟道汇流过程均等分为N_g个，且根据实际情况确定子流域内是否有淤地坝；如果有淤地坝，则将流域内的所有淤地坝概化为N_g个大小一样的等效淤地坝，且假定位于各沟道出口，从而起到淤地坝拦蓄沟道径流的作用。

优选地，步骤S5中：

运动波方程模拟等高带的汇流过程计算公式为：

ql_i＝q_i/len_i；

S_fi＝S_0i；

其中,ql_i为汇流单宽流入量，q_i为参与汇流过程侧向径流量，len_i为汇流长度，A_i为汇流过水断面面积，Q_i为汇流断面出口流量，n_i为曼宁糙率系数，R_i为水力半径，S_0i为坡降，S_fi为摩擦坡降，t为时间，x为沿汇流方向的坐标，i表示坡面(＝pj,j为等高带编号)、沟道(＝g)、河道(＝r)不同汇流情况；

对每个等高带采用以下方法模拟梯田的产流过程：

H_j＝(Q_p(j-1)*r_tj*dt/A_sj+P_j-Inf_j)；

其中，H_j为当前等高带梯田内堆积的水深，H_max为梯田最大可接受的储留深，Q_p(j-1)为上一等高带流入当前等高带的流量(如果j为第一等高带，该值为0)，dt为模拟时间间隔，P_j为当前时段降水量，Inf_j为当前时段梯田部分入渗量，R_tj为当前等高带梯田部分产流深，其余参数意义同前面。

优选地，步骤S6中：

对于沟道汇流而言，只在子流域层面进行汇流演算，且认为子流域内有N_g条长度为len_g的沟道，对每条沟道分别进行汇流演算得到各沟道汇流量Q_g，再求和得到总的沟道汇流量Q_g*N_g。沟道的断面概化为倒等腰梯形，相关参数根据调查或文献获得。沟道汇流上游入流量＝0，参与侧向径流的q_g等于各等高带计算的坡面汇流量Q_pj*r_g之和，坡降S₀₂为外部给定的参数。

值得注意的是，淤地坝调蓄模拟过程如下：首先，将根据各等高带淤地坝面积和洼地储留深，计算得出当前时段淤地坝上洼地储留的水量以及最大储水量，计算得出当前淤地坝可容纳的最大水量＝最大储水量-时段初储水量-时段降水+时段蒸发+时段入渗量；其次，根据淤地坝拦截比例和当前时段所有沟道汇入径流量，计算得出当前时段拦截沟道径流进入淤地坝区域的水量，并采用当前淤地坝可容最大水量进行修正，计算得到淤地坝储水量；第三，根据淤地坝储水量以及透水系数，计算得出当前时段渗透淤地坝排入下游的水量；第四，计算当前时段淤地坝净储水量＝时段拦截水量-时段淤地坝透水量；第五，修正所有沟道汇流量＝修正前流量-当前时段淤地坝净储水量/时段时长；最后，计算时段末淤地坝储水量＝时段初储水量+当前时段淤地坝净储水量，并根据各等高第面积，计算各等高第上淤地坝的时段末洼储水深。在计算时段末洼地储留深的时候，需要根据当前时段淤地坝净储水量的符号以及时段初各等高带洼地储留深的值进行分类处理，如果净储水量为负，则对时段初各等高带洼地储留深按比例缩减即可，如果净储水量为正，则将所有储水量平摊到各等高带上。

优选地，步骤S7中：

对河道采用运动波方程进行河道汇流模拟；对于河道而言，也只在子流域层面进行汇流演算，只有1个河段，长度为根据DEM提取的len_r,断面也概化为倒等腰梯形，相关参数根据河流上水文实测断面进行概化得出；河道上游入流量等于上游子流域出口径流量，参与侧向径流的q_r为各等高带计算的坡面汇流量Q_pj*r_r以及经过淤地坝调蓄后的沟道汇流量Q′_g*N_g总和。

本发明的有益效果是：

本发明中提供一种考虑梯田淤地坝影响的子流域汇流模拟方法，该方法能够模拟计算单元内淤地坝建设对汇流过程的影响；能够模拟计算单元内梯田建设对汇流过程的影响。

附图说明

图1是本发明应用实例计算流程图；

图2是本发明应用实例等高带、备选沟道、河道栅格以及淤地坝位置信息；

图3是本发明应用实例沟道选择示意图；

图4是本发明应用实例汇流属性示意图；

图5是本发明应用实例各等高带径流分配系数；

图6是本发明应用实例坡面沟道河道汇流模拟示意图A；

图7是本发明应用实例坡面沟道河道汇流模拟示意图B。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施方式仅仅用以解释本发明，并不用于限定本发明。

实施例

本发明实施例提供了一种考虑梯田淤地坝影响的子流域汇流模拟方法，如图1所示，包括如下步骤：

S1，采用子流域套等高带作为模型基本计算单元。本实例以1个子流域为例进行说明整个汇流过程。其中，假定该子流域划分为3个等高带，并按栅格高程排序分割得到，等高带1高程最高，等高带3高程最低，如图2所示。在应用实例提取中，河道提取阈值取的是5km²。

S2，设置0.5km²阈值提取备选沟道。查找直接汇入河道栅格的备选沟道栅格作为出口栅格，并对其进行溯源遍历搜索，查找各出口栅格以上最长的备选沟道汇流长度。设定最小长度阈值，剔除小于该长度的伪沟道。本应用实例中共有10条备选沟道流入河道，根据各栅格汇流路径长度，求10条备选沟道的长度，如图3所示。其中点栅格为对应备选沟道最长的汇流路径。设置200m作为备选沟道最小长度阈值，可以发现图3中红圈中的备选沟道长度小于该阈值(180m)，作为伪沟道剔除(见图4中沟道栅格点)。因此本应用实例共9条沟道。为了计算机处理方便，这里直接采用沟道栅格个数(155)乘以栅格边长(90m)作为所有沟道长度，再除以沟道数量，得到每条概化沟道的平均长度1550m。

对流域内不同等高带内的沟道河道栅格进行溯源遍历，确定各栅格的汇流属性。首先将子流域内所有栅格赋值为F2；其次对等高带3中的沟道栅格赋值F1，并逐个沟道栅格进行溯源遍历，将汇入其中的所有坡面栅格赋值为F1；第三将等高带1和2中的河道栅格赋值F4，并逐河道栅格进行溯源遍历，将汇入其中的所有坡面栅格赋值为F4；最后将等高带1和2中的沟道栅格赋值F3，并逐沟道栅格进行溯源遍历，将汇入其中的所有坡面栅格赋值为F3，结果如图4所示。统计各等高带中F1-F4栅格的个数，并根据公式计算坡面、沟道、河道径流分配系数，结果如图5所示。

S3，根据收集的资料确定各等高带内的梯田面积。从地市统计年鉴中收集对应年份的梯田面积，一般情况下，收集的该类面积会大于子流域套等高带的面积，需要进行空间展布。一般认为地市内部的梯田分布比较均匀，这里将收集的梯田统计资料面积按地市内的山区旱地空间分布进行展布，得到单位面积上的梯田面积。然后，累加等高带范围内的梯田面积即采用对应山区旱地面积乘上其单位面积的梯田面积。如果没有相关统计数据，则可以通过遥感手段，基于DEM和遥感影像资料，采用一定的准则，解译提取流域内的梯田范围，并同等高带进行叠合相交处理，从而得出等高带内的梯田面积。

采用以下公式计算某个等高带对上一等高带的坡面径流拦截比例：

r_tj＝(2.2409-0.4075*A_sj/A_C1)*A_s/A_C0

其中，r_tj为当前等高带梯田对上一等高带坡面径流拦截比例，A_sj为当前等高带内梯田面积，A_Cj为当前等高带面积，A_C(j-1)为上一等高带面积，j表示等高带编号。如果当前等高带是最上面的等高带，则r_t1＝0。

S4，采用DEM和遥感影像解译出子流域范围内的淤地坝范围，假设分布情况如图2所示。采用GIS叠合操作，统计等高带范围内淤地坝面积；通过文献查找淤地坝最大洼地储留深；根据淤地坝控制区域和沟道控制区域的比例求得淤地坝拦截比例系数。

S5，对每个等高带采用以下方法模拟梯田的产流过程：

H_j＝(Q_p(j-1)*r_tj*dt/A_sj+P_j-Inf_j)；

其中，H_j为当前等高带梯田内堆积的水深，H_max为梯田最大可接受的储留深，Q_p(j-1)为上一等高带流入当前等高带的流量(如果j为第一等高带，该值为0)，dt为模拟时间间隔，P_j为当前时段降水量，Inf_j为当前时段梯田部分入渗量，R_tj为当前等高带梯田部分产流深，其余参数意义同前面。

采用运动波方程对各等高带坡汇流过程进行模拟。在本应用实例中，先对等高带1坡面汇流进行模拟计算，其中上游入流量为0，侧向入流量q_p1为等高带1非梯田区域产流量和梯田区域产流量R_t1的面积加权值，模拟计算得到等高带1出口坡面径流量Q_p1，然后按0.18,0.80和0.02的比例将Q_p1进行匹分，得到等高带2坡面汇流的上游流入量、直接进沟道量和直接进河道的量；其次对等高带2坡面汇流进行模拟计算，其中上等高带入流量为0.18*Q_p1*(1-r_t2),侧向入流量q_p2为等高带2非梯田区域产流量和梯田区域产流量R_t2的面积加权值，模拟计算得到等高带2出口坡面径流量Q_p2，然后按0.64,0.33和0.03的比例将Q_p2进行匹分，得到等高带3坡面汇流的上游流入量、直接进沟道量和直接进河道的量；最后对等高带3坡面汇流进行模拟计算，其中上等高带入流量为0.64*Q_p2*(1-r_t3),侧向入流量q_p3为等高带3非梯田区域产流量和梯田区域产流量R_t3的面积加权值，模拟计算得到等高带3出口坡面径流量Q_p3，然后按0.53和0.47的比例进行匹分得到直接进沟道量和直接进河道的量。汇流示意图如图6所示。

S6，所有等高带坡面汇流完成后进行沟道汇流模拟，仍采用运动波方式进行,沟道上游入流量为0，侧向入流量q_g＝(0.8Q_p1+0.33Q_p2+0.53Q_p3)/N_g，模拟计算得到单个沟道出口径流量Q_g，则沟道总径流量为Q_g*N_g。

采用淤地坝调蓄模块对沟道总径流量Q_g*N_g进行修正，水量平衡计算原理如图7所示。首先，根据各等高带产流模块计算得到的降水、蒸发、下渗量计算整个淤地坝区域的对应水量，计算淤地坝当前储水量，并根据最大储水量计算得到当前淤地坝最大可容纳水量；其次，根据淤地坝拦截比例计算得到拦截上游沟道径流量I_u，并同最大可容纳水量进行对比，得出超淤地坝库容的溢流量O₁；第三，根据淤地坝储水量以及透水系数，计算得出当前时段渗透淤地坝排入下游的水量O₂；第四，计算当前淤地坝净储水量ΔS＝P+I_u-E-Inf-O₁-O₂；第五，修正沟道汇流量Q′_g＝Q_g-ΔS/dt；最后，根据ΔS符号对各等高带淤地坝上的洼地储留深进行修正，如果净储水量为负，则对时段初各等高带洼地储留深按比例缩减即可，如果净储水量为正，则将所有储水量平摊到各等高带上。汇流示意图如图6所示。

S7，对河道采用运动波方程进行河道汇流模拟。河道上游入流量为上游子流域出流量，侧向径流入流q_r＝0.02Q_p1+0.03Q_p2+0.47Q_p3+Q′_g*N_g，模拟得出当前子流域河道径流量Q_r，汇流示意图如图6和图7所示。重复步骤S5到S7，直到所有时间、所有子流域汇流模拟结束。

Claims (10)

1.一种考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，包括如下步骤：

S1，采用子流域套等高带的方式作为模型基本计算单元，确定子流域范围及等高带范围；

基于DEM对子流域和等高带内的栅格进行分析统计，并获取子流域、等高带相关参数；

S2，基于DEM提取等高带内的模拟沟道，并统计模拟沟道的数量及平均长度；根据子流域内各栅格最终汇入的河道沟道栅格类型，统计沟道和河道汇流栅格个数，并计算沟道河道径流分配系数；

S3，基于统计资料或遥感资料，分析确定各等高带计算单元内的梯田面积；根据各等高带内的梯田面积确定当前等高带内的梯田对上一等高带的坡面径流拦截比例，即认为该部分坡面径流完全进入当前等高带的梯田内，参与当前等高带梯田的产汇流过程；

S4，设定淤地坝只存在于沟道上，且将每个沟道上的所有淤地坝概化为一个等效淤地坝，并假定每条沟道上淤地坝大小一样且位于沟道出口处；统计淤地坝相关参数；

淤地坝作为一个小型水库拦蓄上游沟道汇流径流量，在淤地坝区域内调蓄后排入河道；

S5，基于运动波方程模拟各等高带坡面汇流过程，其中各等高带上游流入量为上一等高带坡面径流乘以对应等高带坡面径流分配系数减去被当前等高带梯田拦截的部分，侧向汇入量等于当前等高带非梯田部分产流量加上梯田调蓄后的产流量；

S6，基于运动波方程模拟单条沟道汇流过程，其中侧向汇入量为各等高带坡面径流乘以对应等高带沟道径流分配系数除以沟道数量；

S7，采用运动波方程对子流域河道进行汇流模拟，重复步骤S4-S6，直到所有时间、所有子流域汇流模拟结束。

2.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S1中：

模型基本单元是子流域套等高带，其中，子流域是地表水文过程相对独立的区域单元，主要用于沟道、河道汇流运算；等高带是子流域内部的细分，用于反映高程变化的差异影响，主要用于产流、坡面汇流运算；

确定等高带范围包括：首先确定等高带个数，再依据等高带的个数确定等高带范围。

3.根据权利要求2所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，子流域范围确定后，可通过两种方法确定等高带的个数；

方法一：给定最大等高带个数和最小等高带面积阈值参数计算获得；首先采用子流域面积除以最小等高带面积阈值确定等高带个数，再同最大等高带个数求最小值得到最终的等高带个数；

方法二：给定等高带平均面积参数计算获得；采用子流域面积除以等高带面积，向上取整得出最终的等高带个数；

根据等高带个数确定等高带范围的方法，包括对子流域内所有栅格按高程降序排列，获得最高和最低高程值，采用等高带个数将这个高程区间进行等分，按高程从高到底依次从1开始自然数编号；遍历子流域内栅格，根据栅格高程位于的高程区间分别赋值对应的区间编号，从而确定各等高带范围，即相同编号的栅格位于同一个等高带。

4.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，根据划分的等高带范围，以DEM为基础数据提取等高带的坡度、长和宽参数，其中，等高带坡度是等高带范围内所有栅格坡度的平均值，等高带宽是沿着坡面流向的长度，等高带长是垂直坡面流向的长度，长和宽计算公式如下：

wid_c＝(H_max-H_min)/sin(slope)

len_c＝A_c/wid

其中，wid_c为等高带宽，len_c为等高带长，H_max是等高带范围内的最大高程，H_min是等高带范围内的最小高程，slope是等高带平均坡度，A_c是等高带面积。

5.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S2中具体包括：

S21，由于提取的沟道存在着长度过小的情况，认为这些短的沟道属于提取出的伪沟道，需要加以剔除，因此直接由阈值提取出的沟道为备选沟道；

S22，遍历子流域内的河道栅格，统计流入河道栅格的备选沟道栅格作为备选沟道的出口；对每个备选沟道出口栅格，统计所有流入该栅格的上游备选沟道栅格，确定该备选沟道的最大汇流路径长度；

S23，设定一个沟道最小长度阈值，并认为长度小于该阈值的备选沟道为伪沟道；将伪河道栅格属性赋值为坡面栅格，统计确定该子流域内沟道条数N_g，并对所有沟道长度求平均值得到沟道的平均概化长度len_g；

S24，对流域内每个河道和沟道栅格进行溯源遍历确定各栅格的汇流属性；

如果沟道栅格位于最低等高带，则给所有汇入该沟道栅格的栅格赋值最低等高带沟道汇流属性(F1)，同理，对所有汇入最低等高带河道栅格的栅格赋值最低等高带河道汇流属性(F2)，对所有汇入其他等高带沟道栅格的栅格赋值其他等高带沟道汇流属性(F3)，对所有汇入其他等高带河道栅格的栅格赋值其他等高带河道汇流属性(F4)；

S25,统计沟道和河道汇流栅格个数，并计算等高带坡面、沟道、河道径流分配系数；

最低等高带径流分配系数计算公式如下：

r_p＝0

r_g＝N_F1/(N_F1+N_F2)

r_r＝N_F2/(N_F1+N_F2)

其他等高带径流分配系数计算公式如下：

r_p＝(N_F1+N_F2)/(N_F1+N_F2+N_F3+N_F4)

r_g＝N_F3/(N_F1+N_F2+N_F3+N_F4)

r_r＝N_F4/(N_F1+N_F2+N_F3+N_F4)

其中，r_p，r_g，r_r表示坡面、沟道、河道的径流分配系数；N_F1为该等高带内最低等高带沟道汇流属性栅格个数；N_F2为该等高带内最低等高带河道汇流属性栅格个数；N_F3为该等高带内其他等高带沟道汇流属性栅格个数；N_F4为该等高带内其他等高带河道汇流属性栅格个数。

6.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S3中：

等高带内的梯田面积确定可有以下两种方法确定：

(1)收集地市统计年鉴中的梯田面积，并根据土地利用类型，按山区旱地斑块进行空间展布，并同等高带计算单元进行相交处理，得到等高带内梯田的面积；

(2)通过遥感手段，基于DEM、遥感影像资料，根据一定的标准，解译得出流域范围内的梯田面积；

当前等高带梯田对上一等高带坡面径流拦截比例由以下公式计算：

r_tj＝(2.2409-0.4075*A_sj/A_C1)*A_s/A_C0

其中，r_tj为当前等高带梯田对上一等高带坡面径流拦截比例，A_sj为当前等高带内梯田面积，A_Cj为当前等高带面积，A_C(j-1)为上一等高带面积，j表示等高带编号；如果当前等高带是最上面的等高带，则r_t1＝0。

7.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S4中：

统计淤地坝相关参数包括每个等高带内的淤地坝面积、最大洼地储留深、沟道径流拦截比例、淤地坝透水系数；

每个等高带内的淤地坝面积可根据DEM和遥感影像，基于一定的准则解译出子流域范围内的淤地坝范围，并同各等高带范围进行叠合处理，从而得出各等高带内的淤地坝面积；

淤地坝最大洼地储留深可根据参考文献设置；

淤地坝沟道径流拦截比例用于表征淤地坝对沟道径流的拦截程度，范围在0-1之间，是反映淤地坝位置和溢流信息的一个参数，如果实际情况下淤地坝大多位于沟道下游，则该比例系数较大，如果实际情况下淤地坝大多位于沟道上游，则该比例系数较小，如果淤地坝有溢流设施，则该比例相对较小；

淤地坝透水系数用于表示通过坝体渗透的水量占储水量的比例；

将沟道汇流过程均等分为N_g个，且根据实际情况确定子流域内是否有淤地坝；如果有淤地坝，则将流域内的所有淤地坝概化为N_g个大小一样的等效淤地坝，且假定位于各沟道出口，从而起到淤地坝拦蓄沟道径流的作用。

8.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S5中：

运动波方程模拟等高带的汇流过程计算公式为：

ql_i＝q_i/len_i；

S_fi＝S_0i；

其中,ql_i为汇流单宽流入量，q_i为参与汇流过程侧向径流量，len_i为汇流长度，A_i为汇流过水断面面积，Q_i为汇流断面出口流量，n_i为曼宁糙率系数，R_i为水力半径，S_0i为坡降，S_fi为摩擦坡降，t为时间，x为沿汇流方向的坐标，i表示坡面(＝pj,j为等高带编号)、沟道(＝g)、河道(＝r)不同汇流情况；

对每个等高带采用以下方法模拟梯田的产流过程：

H_j＝(Q_p(j-1)*r_tj*dt/A_sj+P_j-Inf_j)；

其中，H_j为当前等高带梯田内堆积的水深，H_max为梯田最大可接受的储留深，Q_p(j-1)为上一等高带流入当前等高带的流量(如果j为第一等高带，该值为0)，dt为模拟时间间隔，P_j为当前时段降水量，Inf_j为当前时段梯田部分入渗量，R_tj为当前等高带梯田部分产流深，其余参数意义同前面。

9.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S6中：

对于沟道汇流而言，只在子流域层面进行汇流演算，且认为子流域内有N_g条长度为len_g的沟道，对每条沟道分别进行汇流演算得到各沟道汇流量Q_g，再求和得到总的沟道汇流量Q_g*N_g；沟道的断面概化为倒等腰梯形，相关参数根据调查或文献获得；沟道汇流上游入流量＝0，参与侧向径流的q_g等于各等高带计算的坡面汇流量Q_pj*r_g之和，坡降S₀₂为外部给定的参数。

10.根据权利要求1所述的考虑梯田淤地坝影响的子流域汇流模拟方法，其特征在于，步骤S7中：

对河道采用运动波方程进行河道汇流模拟；对于河道而言，也只在子流域层面进行汇流演算，只有1个河段，长度为根据DEM提取的len_r,断面也概化为倒等腰梯形，相关参数根据河流上水文实测断面进行概化得出；河道上游入流量等于上游子流域出口径流量，参与侧向径流的q_r为各等高带计算的坡面汇流量Q_pj*r_r以及经过淤地坝调蓄后的沟道汇流量Q′_g*N_g总和。

Images