[发明专利]一种基于实时雨情的动态水系提取方法

权利要求书

1.一种基于实时雨情的动态水系提取方法，其特征在于，包括以下步骤：

S1、基于栅格数字高程模型DEM数据将流域划分为若干个正交的栅格单元，并且计算每一个栅格单元中的坡度，得到坡度栅格；

S2、估算流域的水系阈值，提取出流域中的水系，并提取出流域中的河源汇流单元；

S3、计算该流域的多年汛期平均逐小时降雨量，并计算在该小时降雨量下每一个栅格单元中的产流，得到汛期平均径流深栅格；

S4、结合实时逐小时降雨资料，在每一个栅格单元中进行产流计算，得到实时径流深栅格；

S5、基于水系阈值、汛期平均径流深栅格、坡度栅格以及实时径流深栅格估算逐小时变化的动态水系阈值，并根据该动态水系阈值提取逐小时变化的动态水系；

所述步骤S1中，计算每一个栅格单元中的坡度，得到坡度栅格采用的方法是：

以栅格单元Cell为中心，通过周围栅格单元的高程值与该栅格单元的高程值的对比，找出与其相比最低的栅格单元Cell_D，并计算Cell和Cell_D之间的高程差DH_max和水平投影距离Dis，结合DH_max和Dis计算栅格单元Cell的坡度S：

S＝DH_max/Dis

然后按照以上方法遍历流域中每一个栅格单元，从而得到坡度栅格Raster_Slope；

所述步骤S2中，估算流域的水系阈值，提取出流域中的水系，并提取出流域中的河源汇流单元所采用的方法包括以下步骤：

S2.1、基于步骤S1得到的坡度栅格Raster_Slope获得流向栅格Raster_Dir和汇流累计栅格Raster_Acc：

将Cell作为出流栅格单元，Cell_D作为入流栅格单元，入流栅格单元汇流累计值加1，逐栅格循环，计算出每一个栅格单元中的汇流累计值Acc，同时依据Cell和Cell_D之间的相对位置关系，采用D8流向法确定Cell中的流向，按照以上方法遍历流域中每一个栅格单元，从而得到流向栅格Raster_Dir和汇流累计栅格Raster_Acc；

S2.2、估算流域的水系阈值，首先初步假定水系阈值TA，将Raster_Acc中Acc高于TA的栅格单元判定为河道栅格单元Cell_River，并将提取得到的河道栅格单元的位置与GoogleEarth上的卫星影像中的实际河道相对比，反复调整TA，使得提取出来的河道栅格单元的位置与卫星影像中的实际河道的位置重合，从而去确定汇流阈值TA，同时提取出流域中由所有河道栅格单元Cell_River组成的水系栅格Raster_River；

S2.3、基于水系栅格Raster_River，提取出所有一级河流的终点Points，以Points中的一个点Point为出口点，结合Raster_Dir，提取出所有汇流到Point的栅格单元，并将这些栅格单元所在的区域作为Point所对应的河源汇流单元Subaisn，遍历Points中的点，从而提取出流域中所有的河源汇流单元；

所述步骤S3中，计算该流域的多年汛期平均逐小时降雨量，并计算在该小时降雨量下每一个栅格单元中的产流，得到汛期平均径流深栅格的步骤包括：

S3.1、计算流域的多年汛期平均逐小时降雨量；

式中：P汛_h为多年汛期平均逐小时降雨量；YE为流域降雨数据的总年数，y为流域降雨数据以年时间为单位的编号，从1到YE；XD为流域汛期的总天数，d为一个汛期内流域降雨数据以天为时间单位的编号，从1到XD；P_y，d为第y年的从汛期起算开始第d天的降雨量；

S3.2、以多年汛期平均逐小时降雨量为输入，采用蓄满产流计算方法计算每一个栅格单元中的径流深，得到汛期平均径流深栅格Raster_Rve；

所述步骤S4中，结合实时逐小时降雨资料，在每一个栅格单元中进行产流计算，得到实时径流深栅格的步骤包括：

S4.1、基于站点观测的逐时段降雨数据，通过反距离权重插值得到流域中降雨的空间分布栅格径流深栅格Raster_Ph；

式中：Ph_k是编号为k的栅格单元内经过插值得到的降雨量；Ph_r为编号为r的站点的观测降雨，r为站点的编号从1到R，R为站点的总数目；d_r，k为编号为k的栅格单元到编号为r的站点之间的距离，x_k为编号为k的栅格单元的横坐标，x_r为编号为r的站点的横坐标，y_k为编号为k的栅格单元的纵坐标，y_r为编号为r的站点的纵坐标；

S4.2、以实时的逐小时降雨数据为输入，采用蓄满产流计算方法计算每一个栅格单元中的径流深，得到逐小时变化的实时径流深栅格Raster_Rtime；

所述步骤S5中，基于水系阈值、汛期平均径流深栅格、坡度栅格以及实时径流深栅格估算逐小时变化的动态水系阈值，并根据该动态水系阈值提取逐小时变化的动态水系的步骤包括：

S5.1、初步假定动态水系阈值TH，将Raster_Acc中Acc高于TH的栅格单元判定为河道栅格单元Cell_River，提取出流域中由所有河道栅格单元Cell_River组成的水系栅格Raster_River，基于水系栅格Raster_River，提取出所有一级河流的终点Points，以Points中的一个点Point为出口点，结合Raster_Dir，提取出所有汇流到Point的栅格单元，并将这些栅格单元所在的区域作为Point所对应的沟道汇流单元G_Subaisn；

S5.2、在沟道汇流单元中构建动态平衡方程：

式中：NS是流域中河源汇流单元的总数目，b是流域中河源汇流单元的编号，从1到NS；TA是水系阈值，即河源汇流单元中栅格单元的总数目，i是流域河源汇流单元中栅格单元的编号，从1到TA；S_i，b是编号为b的河源汇流单元中第i个栅格单元的坡度值；R_汛h为以多年汛期平均逐小时降雨量为输入，计算得到的流域中每一个栅格单元内的径流深；GS是流域中沟道汇流单元的总数目，g是流域中沟道汇流单元的编号，从1到GS；TH是动态水系阈值，即沟道汇流单元中栅格单元的总数目，j是流域沟道汇流单元中栅格单元的编号，从1到TH；S_j，g是编号为g的沟道汇流单元中第j个栅格单元的坡度值；R_时h，j，g为以实时的逐小时降雨数据为输入，计算得到的流域中编号为g的沟道汇流单元中第j个栅格单元的径流深；

反复调整TH的值，取得相应的GS的值，使得沟道汇流单元中的动态平衡方程成立，此时的TH即为该时段下的动态水系阈值；

S5.3、将Raster_Acc中Acc高于TH的栅格单元判定为河道栅格单元Cell_River，同时提取出流域中由所有河道栅格单元Cell_River组成的水系栅格Raster_River，并且以实时降雨为输入，逐时段反复执行步骤S2到步骤S5之间的操作，从而得到流域中随实时降雨数据逐时段变化的动态水系。