CN102034003B - 基于蓄水容量曲线和topmodel的流域水文模型的设计方法 - Google Patents

Abstract

本发明提供了一种基于蓄水容量曲线和TOPMODEL的流域水文模型的设计方法，其步骤为：地表产流的计算；地下产流的计算；蒸散发计算；水量平衡方程的应用；汇流机制；等流时线的生成；土壤湿度的动态可视化表达；NoahLSM的耦合。本发明其模拟确定性系数高，对日径流过程和次洪过程的模拟均优于半分布式模TOPMODEL和新安江模型以及分布式水文模型ESSI、HEC-HMS、SWAT，克服了传统流域水文模型往往不适于大尺度及干旱半干旱条件下水文过程精确模拟的缺点。本发明模型还成功地与NoahLSM耦合，提高了业务天气预报模式GRAPES-MESO对降水、土壤湿度和径流的预报能力。

Description

基于蓄水容量曲线和TOPMODEL的流域水文模型的设计方法

技术领域

本发明半分布式水文模型的设计方法，更具体的说是基于蓄水容量曲线和TOPMODEL的流域水文模型的设计方法。

背景技术

TOPMODEL是Beven和Kirkby于1979年提出的一个以地形为基础的半分布式水文模型，其主要特征是利用流域地形指数（lnα/tanβ）作为一个水文相似性的指数，所有具有相同地形指数值的点，对降雨具有相同的水文响应。作为该假定的结果，流域的地形可以用地形指数的空间分布来表示。具有较大地形指数值的点，将首先趋于饱和，因此暗示着潜在的地表和地下径流的产流面积。由流域湿润或干燥引起的产流面积的变化可由地形指数的空间分布来表示。

地形是降雨—径流动态物理过程中非常重要的影响因子。TOPMODEL以地形空间变化为主要结构，用地形指数（topography index,也称湿润指数，wetness index）的形式描述地下水的分布情况以及自由水排水沿坡向的运动，表达了变动面积产流的概念。根据这个概念，流域上产生的坡面流被认为是不均匀的。在一定程度上，变动产流面积可看作河道系统的延伸。

TOPMODEL参数比较少，共有五个，它们分别是：土壤孔隙尺度m，饱和导水率T0，根层最大缺水深SRmax，根层初始含水量SR0以及汇流速度ChVel。

新安江模型在国内湿润和半湿润地区得到广泛的应用。模型的核心是蓄水容量曲线,表示降雨分布均匀时的产流面积的变化情况。

TOPMODEL模型和新安江模型都是概念性模型，都是基于变动产源面积理论的模型，因此都是半分布式流域水文模型，而且两者都主要应用于湿润和半湿润地区。这两个模型结构简单，优选参数少(尤其是TOPMODEL)，而且与观测的物理水文过程有密切联系，模拟效果好，己被应用到水科学研究的各个方面，并不断发展、改进。由于这两个模型有很多相似之处，不少学者试图利用这两个模型的优点改进其中的一个模型或新建一个模型。例如，韩杰和张万昌（2005）将新安江模型三水源的结构应用到两水源TOPMODEL中，率先将其改造成三水源的模型，并在汉江子午河流域进行了应用，实现了TOPMODEL模型对三种径流成分（地表径流、壤中流及地下径流）的分割。文佩和陈喜（2006）也作了类似的改造工作，并用改进的TOPMODEL模型模拟了岳城等三个流域的日径流过程、次洪水过程和基流过程。结果表明，改进模型能较好地反映南方湿润地区产汇流过程及各种径流成份模拟。雍斌和张万昌(2007)用改进型SIMTOP（TOPMODEL的一个变种）的饱和产流方案与新安江水量平衡计算方法相结合，构建了一个具有地形指数尺度转换机制的陆面水文过程模型TOPX，并在汉江酉水街小流域进行了单点离线测试。结果表明TOPX模型在温润区小流域的日径流量模拟精度高，而且TOPX模型能够较好捕捉土壤含水量、饱和源面积等各种水文过程变量的时空动态变化过程。然而，到目前为止尚未见将新安江模型中的蓄水容量曲与地下水位波动关联，并与TOPMODEL地下产流方案相结合以构建新模型的报道。

陆面过程模式的发展历程和研究进展可以分为三个阶段：①“水箱”模式阶段；②SVAT模式阶段；③VIC和TOPMODEL大发展阶段。VIC模型是一个不断发展的陆面水文过程模型。该模型吸取了新安江模型的蓄水容量曲线的思想，考虑了网格系统内由于地形、土壤及植被的变化产生的入渗能力的变化，网格的入渗能力用一个抛物线型函数描述。VIC模型的产流思想又不断地被其它陆面参数化方案所引用，以改善以往对产流方案处理的不足。

Noah LSM由NCEP、NCAR、AFWA和UCLA共同发展，已在MM5V3.6、NCEP的Eta等中应用和验证。其特点为：多层土壤模式（通常用4层，可到2m深）计算土壤温度、湿度；显式考虑植被冠层的阻力；地表和次地表的径流过程、雪和冻土过程使用1-km的全球植被数据；1-km和5-min全球土壤类型图更新了土壤热传导系数、不协调的雪盖、雪密度、雪最大反照率数据库、积雪下土壤热通量的处理和雪的粗糙度；增加了地表感热、潜热、地表热通量、地表温度和地表径流等项；耦合了一个简单的城市下垫面模型。在Noah陆面模式中，地表径流参数化方案采用Schaake（1996）提出的简单水量平衡模型SWB。SWB模型是一个两层土壤的简单水量平衡模型，相当于两个水桶（bucket）模型。和bucket模型一样，未显式考虑植被对径流的影响，输入数据仅需要降水和潜在蒸发量。Beven（1979）曾指出，TOPMODEL疏于考虑的一个问题就是空间差异性。地形指数在进入模型计算的时候只是以统计的方式参与的，并没有考察其在空间分布上的差异。由此一来，就产生了一系列的问题：降水、下垫面情况、地形等空间上的差异就无法在水文响应中表现出来。而且，在将TOPMODEL与陆面过程模型LSM、CLM耦合时，是将地形指数拟合成三参数gamma分布函数或二参数指数分布函数的统计形式表达空间异质性，这与VIC模型非常像似。而且，TOPMODEL模型对地下水位的刻画上存在一个严重问题，即TOPMODEL模型刻画地下水动态变化时存在与实际情况偏差的问题（平均地下水位线高于地面），而且仅仅通过改进该模式很难解决这一问题。因此，非常有必要通过其他措施（比如模型与模型的“杂交”）解决这一问题。

新安江模型中只考虑土壤中含水量的变化，无法表达地下水位的动态变化。地下水位对大气模式和气候模式中水分和能量平衡有重要的影响。因此，解决新安江模型中地下水位动态表达的问题是非常必要的、紧迫的任务。

Noah LSM对水文过程特别是径流的描述还存在严重缺陷。其依据是它的径流参数化方案使用的是简单水量平衡模型（SWB）。SWB与TOPMODEL和新安江（或VIC）模型相比，其共同点是都能以统计分布函数刻画土壤湿度的空间不均匀性，不同点是，SWB不能有效地表达径流产源面积的变动情况，而且其产汇流参数化方案过于简单，不能精确地刻划陆面及浅地表水文过程。而TOPMODEL和新安江模型是基于变动产源面积理论的模型，而且对陆面及浅地表水文过程均有深入地刻划。因而用TOPMODEL和新安江模型或更先进的模型取代SWB模型，并能紧密地而不是松散地耦合到Noah LSM就显得迫切而重要。

发明内容

本发明的目的是提供基于蓄水容量曲线和TOPMODEL的流域水文模型的设计方法，紧密关联新安江模型中的蓄水容量分布曲线与TOPMODEL中关于地下水位线的表达，并在此基础上构建一个新的降雨—径流模型XXT；实现XXT模型与单点离线Noah LSM陆面过程模式的紧密耦合。

本发明的目的是通过以下研究方案来实现的：

基于蓄水容量曲线和TOPMODEL的流域水文模型（简称“XXT模型”）的设计方法，其步骤为：

A）地表产流的计算；

B）地下产流的计算；

C）蒸散发计算；

D）水量平衡方程的应用；

E）汇流机制；

F）等流时线的生成；

G）土壤湿度的动态可视化表达；

H）与Noah LSM的耦合。

在步骤A）中，将地下水水位线与新安江模型中的流域蓄水容量曲线关联形成新流域蓄水容量分布曲线。令流域上地下水水位线刚好位于地面0mm处的蓄水容量为

则新流域蓄水容量曲线定义为蓄水容量小于或等于

的流域面积。把全流域面积看作1，流域由无数个点组成，每个点都对应一个蓄水容量，其中最大的蓄水容量记为

α为饱和区面积与全流域面积之比，其计算公式为

$\mathrm{\α}=1-{(1-\frac{W_m^{\′}}{W_{\mathrm{mm}}^{\′}})}^B$

式中B是表征蓄水容量分布不均匀性的参数，B愈大分布愈不均匀。据上式有：

$W_m^{\′}=W_{\mathrm{mm}}^{\′}\×[1-{(1-\mathrm{\α})}^{1/B}]$

令WM为流域平均蓄水容量，则WM等于未饱和层与饱和层的面积之和：

$\mathrm{WM}=\underset{0}{\overset{1}{\∫}}{W}_m^{\′}\mathrm{d\α}=\underset{0}{\overset{1}{\∫}}{W}_{\mathrm{mm}}^{\′}[1-{(1-\mathrm{\α})}^{1/B}]\mathrm{d\α}=W_{\mathrm{mm}}^{\′}/(1+B)$

新模型中土壤在垂直方向上分三层，即植被截留层、未饱和层及饱和层。植被（含根系）截留层有一最大储水容量SRmax，降雨在经过截留层时将被截留，只有该层蓄满后才会向下一层渗透。设植被截留层某一时刻的土壤水分亏缺深为SRZ（单位为m），则当P-SRZ≤0时，没有水从该层渗透到下面一层（未饱和层）；当P-SRZ>0时，有（P-SRZ）的水渗透到未饱和层。新模型土壤分层是有理论和实践依据的：即美国爱荷华州立大学（Iowa State University）的Iowa Stormwater Management Manual给出的土壤分层和水分渗透示意图。此图是在实践基础上的理论概括。XXT模型在垂直方向的土壤分层在本质上与该图保持一致。而且，XXT模型在垂直方向的土壤分层也与BTOP（TOPMODEL模型的一个变种）基本一致。

地表产流计算原理如下：

P为降水量，SRZ为植被截留层某一时刻的土壤水分亏缺深，P-SRZ是从植被截留层中渗透下来的降水量。α₁为（P-SRZ）的降水进入土壤未饱和层之前的饱和区面积与全流域面积之比，α₂为（P-SRZ）的降水进入该层之后的饱和区面积与全流域面积之比。

${\mathrm{\α}}_1=1-{(1-\frac{W_{m1}^{\′}}{W_{\mathrm{mm}}^{\′}})}^B$

${\mathrm{\α}}_2=1-{(1-\frac{(W_{m1}^{\′}+P-\mathrm{SRZ})}{W_{\mathrm{mm}}^{\′}})}^B$

令

代表流域平均土壤水分亏缺深，

为

对应的流域上某点的蓄水容量，

的计算公式推导过程如下：

$\stackrel{\‾}{Z}={\∫}_{{\mathrm{\α}}_1}^1{W}_m^{\′}\mathrm{d\α}-W_{m1}^{\′}\×(1-{\mathrm{\α}}_1)$

先解出积分式，再化简，得

$W_{m1}^{\′}=\mathrm{WM}\×(1+B)\×[1-{(\stackrel{\‾}{Z}/\mathrm{WM})}^{1/(1+B)}]$

地表产流的计算分三种情况：当P-SRZ≤0时，没有地表产流（Q_S），即：

Q_S=0

当 $0<P-\mathrm{SRZ}\&le;\mathrm{WM}\&times;(1+B)-W_{m1}^{\&prime;}$ 时：

$Q_S=(P-\mathrm{SRZ})\&times;{\mathrm{\&alpha;}}_2-[{\&Integral;}_{{\mathrm{\&alpha;}}_1}^{{\mathrm{\&alpha;}}_2}{W}_m^{\&prime;}\mathrm{d\&alpha;}-({\mathrm{\&alpha;}}_2-{\mathrm{\&alpha;}}_1)\&times;W_{m1}^{\&prime;}]$

化简后为

$Q_S=P-\mathrm{SRZ}-\stackrel{\&OverBar;}{Z}+\mathrm{WM}\&times;{\{1-(P-\mathrm{SRZ}+W_{m1}^{\&prime;})/[\mathrm{WM}\&times;(1+B)\left]\right\}}^{(1+B)}$

当 $P-\mathrm{SRZ}>\mathrm{WM}\&times;(1+B)-W_{m1}^{\&prime;}$ 时，易知地表产流为

$Q_S=P-\mathrm{SRZ}-\stackrel{\&OverBar;}{Z}$

由于

和SRZ是状态变量，故该产流方案中只有两个可调参数，WM和B。WM是整个流域很干燥时的缺水量，代表流域的干旱程度，是个气候因素。B代表蓄水容量空间分布的不均匀性，其值愈大蓄水容量空间分布愈不均匀，其值大小取决于地形地质条件。可见该产流方案显式地考虑了地形的空间异质性。

在步骤B）中，XXT模型中的地下产流引用Sivapalan等（1987）年提出的地下径流公式：

$Q_b=Q_0{e}^{-f\stackrel{\&OverBar;}{Z}}$

式中Q_b为地下径流量，Q₀表示暴雨前退水期地下水流量（稳定流流量），f是土壤水力传导率沿深度方向的衰减因子,

代表流域平均土壤水分亏缺深。XXT模型巧妙地通过变量

将蓄水容量曲线与地下产流有机关联，使其在描述地下水动态变化时具有明确的物理意义。此地下产流方案中有两个参数Q₀和f，可通过参数率定获得。

在步骤C）中，实际蒸散发（Q_E）的计算方法为先用下面的公式计算

Q_E＝E×(1-SRZ/SR_max)

上式中E为潜在蒸散发，其值可直接利用蒸发皿观测值，如果没有观测值可用气温、辐射和风等气象资料代入彭曼公式计算得到。当Q_E值计算之后还需作判断，如果Q_E＞SR_max-SRZ，则Q_E应重新赋值：

Q_E＝SR_max-SRZ

此蒸散发计算方案中只有一个参数需要率定，即SRmax。

在步骤D）中，地下水位线是XXT模型中的关键，其对应的变量

在计算下一个时间步长之前必须更新，依据下面的水量平衡方程：

${\stackrel{\&OverBar;}{Z}}_{t+1}-{\stackrel{\&OverBar;}{Z}}_t=Q_{\mathrm{Et}}+Q_{S_t}+Q_{b_t}-P_t$

式中

表示t+1时刻的流域平均土壤水分亏缺深，

Q_Et、

和P_t分别表示t时刻的流域平均土壤水分亏缺深、实际蒸散发、地表径流深、地下径流深和降水。

在步骤E）中，XXT模型的流域汇流程序与TOPMODEL类似：对子流域内部采用等流时线汇流，等流时线基于DEM提取得到。此部分只需要一个率定参数，即汇流速度(CHV)。河道汇流采用马斯京根法。

步骤F）中，假定径流在空间上相等，考虑流速由于受到坡面、粗燥度等的影响是变化的，把地表径流与基流合在一起计算，得出出水点的流量。根据到出水点的距离与一个地表水流速度参数ChVel计算即可得到等流时线。本发明提出了一个基于DEM提取等流时线的新方法，该方法可快速并可视化地提取等流时线。

基于DEM提取等流时线的主要步骤为：

（1）DEM填洼；（2）DEM流向计算，坡度计算；（3）计算每栅格到出水口的距离；（4）将距离矩阵按quantile分类方法进行分类，并在ARCGIS9.0中生成可视化图；（5）生成等流时线文件。

步骤G）中，XXT模型基于变动产流源面积和由DEM提取的地形指数可以实现对土壤湿度的动态可视化表达。其实现步骤主要包括：（1）DEM中“洼地”填充；（2）计算DEM的水流方向矩阵和汇流面积；（3）计算地形指数；（4）根据XXT模型模拟的变动产流面积和地形指数绘制土壤湿度空间分布图。

为了描述方便，将初始DEM记为RDEM，将对应的分类DEM矩阵记为CDEM。由于流域出口在DEM的边界上，因此应从分类值最小的一类开始填洼。对该类中CDEM的值大于其对应的RDEM的值的所有栅格作如下处理：如果某一栅格对应的RDEM的值大于或等于其在CDEM中相邻的8个栅格中的任意一个，则说明该栅格中的水可以流出去，该栅格不属于洼地，应当将该栅格在CDEM的值重新赋值为对应的RDEM的值；如果该栅格对应的RDEM的值比其在CDEM中相邻的8个栅格中的值都小，则将其在CDEM中的值重新赋值为其相邻的8个栅格中最小的值。重复上述过程，直到该类中不存在这样的栅格：即它对应的CDEM的值大于其对应的RDEM的值并且其对应的RDEM的值大于或等于其在CDEM中相邻的8个栅格中的任意一个。这种填洼算法的思想最早是由Moran和Ve′zina（1993）提出，Planchon和Darboux（2001）给出了此思想的一种实现算法。但此实现算法在每次对DEM进行扫描时都必须对所有栅格进行扫描。事实上，这是完全没有必要的，因为有大量栅格不属于洼，一旦确定它们不是洼地就不需要再处理了。而新算法启用两个栈（stack）来解决这一问题。

新填洼算法先建立两个栈（栈A和栈B）,把临时DEM中分位数最小的那些栅格的位置信息（i行j列）存储在栈A或栈B中，分位数最小的那一类中每个栅格（用c表示）进行循环处理，如果栅格单元属于洼地则令其等于其相邻栅格中最低的一个；若不是洼地则令其等于它对应的初始矩阵的值，并且将此栅格从栈中剔出，将没有被剔出的栅格存放到另一个栈中（栈B或栈A）。然后再对栈B或栈A进行类似的处理。重复此过程，直到栈中只剩下洼地为止，这些洼地的位置信息最终会保存在两个栈中的一个里。第一类洼地填充和平坦区处理的C++伪码如表1所示。

表1中符号δ(n)要么为0，要么代表微小高程。如果δ(n)取值为0,洼地被填成平地；如果δ(n)取一个微小高程，比如对于对角线方向取0.000141，非对角线方向取0.0001，则初始DEM被新方法处理后，没有平坦区的存在。

表1对分位数最小的一类进行填洼的算法伪码

第一类处理完后，在第二类处理之前还有一关键步骤：检查栅格中是否还存在数值等于该类分位数的栅格，如果存在，则将其赋值为下一类（即第二类）分位数，因为这些栅格虽然是洼地，但其在空间上与填洼后的第一类不相邻，需要进一步处理。之后进行第二类，第三类等填洼处理，直到所有类全部处理完，这时CDEM就是经过填洼和平坦区处理后的DEM，此时所有的洼地对应的位置信息存放在栈A或栈B两个中的一个之中，这非常有利用随后的水文模拟与分析。

本发明提出的新地形指数提取方法首先用分位数法将DEM分成8类，然后将DEM处理成无洼地及平坦区的DEM。然后，由于DEM中无洼地及平坦区，可直接根据D8算法指定水流方向。最后进行地形指数的计算，其具体计算方法参照Wolock and McCabe(2000)的方法。

步骤H）中，为了实现XXT与Noah LSM耦合，第一步需要用XXT模型中的地表径流方案取代Schaake的简单水量平衡模型SWB。XXT中的地表径流径流的公式表达式如下：

$Q_S=P-\mathrm{SRZ}-\stackrel{\&OverBar;}{Z}+\mathrm{WM}\&times;{\{1-(P-\mathrm{SRZ}+W_{m1}^{\&prime;})/[\mathrm{WM}\&times;(1+B)\left]\right\}}^{(1+B)}$

在Noah LSM中，Rsur=Q_S并且Pd=P-SRZ，代入得：

$R_{\mathrm{sur}}=P_d-\mathrm{DD}+\mathrm{WM}\&times;{\{1-(P_d+W_{m1}^{\&prime;})/[\mathrm{WM}\&times;(1+B)\left]\right\}}^{(1+B)}$

其中，

代表前一时刻在water table处的土壤蓄水容量值，可用下面的公式计算：

$W_{m1}^{\&prime;}=\mathrm{WM}\&times;(1+B)\&times;[1-{(\mathrm{DD}/\mathrm{WM})}^{1/(1+B)}]$

其中DD为四层土壤总的土壤水分亏缺深，在此方案中WM可由Noah LSM中的土壤属性计算得出，只剩下一个参数B是待定的。参数B是表示土壤蓄水容量空间分布不均匀程度的一个物理量，B愈大表示土壤蓄水容量空间分布愈不均匀。对于具体的某一个流域，B值可通过参数率定获得。对于较大的模拟区域（跨几个流域甚至更多），一种简单的办法即是将全模拟区域的B值取为一个值（常数），但这样做就导致新Noah LSM不具有空间分布式的特点。为解决这一问题，借鉴TOPMODEL模型中用DEM提取得到的地形参数而导出模型参数的方法，本发明利用DEM提取参数B。其提取原理是依据土壤蓄水容量空间分布公式，先对DEM进行相高对程的累积分布进行统计，然后用最小二乘法计算出参数B。

${\mathrm{\&alpha;}}_1=1-{(1-\frac{W_{m1}^{\&prime;}}{W_{\mathrm{mm}}^{\&prime;}})}^B$

由于Noah LSM是一维单点模式，并没有汇流方案。一种简单的办法就是将XXT模型中基于DEM的等流时线汇流程序直接和Noah LSM相联。但等流时线汇流方案不是基于格点的分布式汇流方案，为了解决这一问题，在原等流时线汇流方案的基础上，我们开发了一个新的基于格点的分布式汇流方案算法。新算法主要由以下步骤构成：

（1）将Noah LSM的每个格点的产流插值到100m分辨率的格点；

（2）将USGS 90m分辨率的DEM重采样到100m，并进行DEM洼地填充，进而生成流向和坡度矩阵；

（3）径流沿着由DEM计算的流向汇流，其速度为VEL=V0*(1+a*exp(S))，穿过第一个网格需要的时间是T=L/VEL，其中V0为初始速度；

（4）计算出每一个网格的径流到汇流至流域出口所需要时间；

（5）计算出每一时间步长所有格点汇流至流域出口的总径流量。

至此，XXT模型的产流方案及汇流方案已耦合到Noah LSM之中，并且均为空间分布式的，其参数并不需要用具体的流域资料进行率定，这使得新Noah LSM可移植性增强、适用性更广。

本发明与现有技术相比，其改进主要体现如下：

（1）基于新安江模型的蓄水容量分布曲线和TOPMODEL简单的框架结构构建了一个新型流域水文模型，命名为XXT。XXT模型的产流方程是在描述流域蓄水容量分布曲线的函数基础上嵌入了表征地下水位的物理量Z（土壤水分亏缺深），并经过严密地数学推导得到的。该产流方程通过变量Z和模型的地下产流方程相关联。XXT的地下产流方程和蒸发计算模块类似于TOPMODEL。XXT模型的汇流部分采用基于DEM的等流时线法。

XXT模型基于变动产流源面积和由DEM提取的地形指数可以实现对土壤湿度的动态可视化表达。从DEM提取地形指数首先遇到的一个难题是DEM中洼地和平坦区的处理方法。因为传统的方法速度太慢，不适用于模型的可视化表达要求以及XXT模型和气象模式耦合的需要。为此本章提出了一个新的高效的算法。传统方法相比，新方法显著提高了DEM洼地和平坦区的处理效率：平均而言，新方法要比著名的、目前几乎所有GIS软件所使用的Jenson和Domingue(1988)的方法快几千倍以上，比Planchon和Darboux(2001)的快速算法快7倍以上。而且，新方法可提供多种选择，洼地可以被填成“平地”，或被填成有微小高程差异的“坡地”。在此算法的基础上，提出了一个新的地形指数提取方法，使得地形指数的提取更简单、更高效。而且新地形指数提取方法完全不影响TOPMODEL模型对水文过程的模拟和描述。

（2）用XXT模型中的地表径流方案取代Schaake的简单水量平衡模型SWB，从而实现了XXT的产流方案与Noah LSM的耦合。新耦合方案中参数B是表示土壤蓄水容量空间分布不均匀程度的一个物理量，需要率定或能过其它方法获得。借鉴TOPMODEL模型中用DEM提取得到的地形参数而导出模型参数的方法，本文利用DEM提取参数B。其提取原理是依据土壤蓄水容量空间分布公式及最小二乘法。基于DEM获取参数B使得新耦合Noah LSM具有空间分布式特征。耦合分布式XXT模型的Noah LSM为空间分布式的，其参数并不需要用具体的流域资料进行率定，这使得新Noah LSM可移植性增强、适用性更广。

在原等流时线汇流方案的基础上，开发了一个新的基于格点DEM的分布式汇流方案算法。并将此新分布式汇流方案耦合到Noah LSM中。

有益效果：

（1）本发明建立一个新的降雨—径流模型

TOPMODEL和新安江模型都是已经被广泛的接受和应用的半分布式水文模型。TOPMODEL的主要特征是利用流域地形指数（lnα/tanβ）来反映水文现象、特别是径流运动的分布规律，模型的物理意义明确，模型框架结构简单，需要率定的参数较少。新安江模型的主要特点是采用流域蓄水容量分布曲线来表示流域单元中各点蓄水容量空间分布的不均匀性，模型效率高且适用性广。本发明基于新安江模型的蓄水容量分布曲线和TOPMODEL简单的框架结构构建了一个新型流域水文模型，命名为XXT。XXT模型的产流方程是在描述流域蓄水容量分布曲线的函数基础上嵌入了表征地下水位的物理量Z（土壤水分亏缺深），并经过严密地数学推导得到的。该产流方程通过变量Z和模型的地下产流方程相关联。XXT的地下产流方程和蒸发计算模块类似于TOPMODEL。XXT模型的汇流部分采用基于DEM的等流时线法。

XXT模型基于变动产流源面积和由DEM提取的地形指数可以实现对土壤湿度的动态可视化表达。从DEM提取地形指数首先遇到的一个难题是DEM中洼地和平坦区的处理方法。因为传统的方法速度太慢，不适用于模型的可视化表达要求以及XXT模型和气象模式耦合的需要。为此本章提出了一个新的高效的算法。传统方法相比，新方法显著提高了DEM洼地和平坦区的处理效率：平均而言，新方法要比著名的、目前几乎所有GIS软件所使用的Jenson和Domingue的方法快几千倍以上，比Planchon和Darboux的快速算法快7倍以上。而且，新方法可提供多种选择，洼地可以被填成“平地”，或被填成有微小高程差异的“坡地”。在此算法的基础上，本章以提出了一个新的地形指数提取方法，使得地形指数的提取更简单、更高效。而且新地形指数提取方法完全不影响TOPMODEL模型对水文过程的模拟和描述。

新模型测试选择位于温暖湿润气候区的汉江上游的酉水街流域。结果表明，因1981-1983年为丰水年，径流系数较大，三个模型在率定期模拟效率都很高，均在83%以上，其中XXT模拟效率最高；而1984-1985年为枯水年，径流系数小，三个模型在验证期的模拟效率都较低，其中XXT模拟效率明显高于另外两个模型。TOPMODEL模拟的流域平均土壤水分亏缺深出现了负值，显然不符合实际情况。而XXT模拟的土壤水分亏缺深则未出现负值，从此角度上看XXT在模拟地下水位深时比TOPMODEL更合理。XXT模型结构简单、参数少，效率高，便于与大尺度气候模式及数值天气预报业务模式耦合。

（2）XXT适用于多尺度多生态条件下的日径流和洪水模拟

温暖湿润气候区的汉江上游的汉中流域、褒河流域（马道控制站和江口控制站），半湿润半干旱区的淮河流域的子流域—沂沭河流域（沂临控制站）以及位于干旱、半干旱区的黑河流域（鹰落峡控制站）为研究区，对新构建的半分布式水文模型XXT进行包括洪水预报和日径流模拟的较全面的测试和评估，并与半分布式模型TOPMODEL和新安江以及分布式水文模型ESSI、SWAT、HEC-HMS进行对比验证；最后将XXT模型应用于基于TIGGE超级集合的洪水预报试验，以评价新构建的模型对不同流域尺度和不同气候和生态条件下的适用性。研究结果表明：

a）在温润区的江口流域使用了半分布式模型XXT、新安江和TOPMODEL以及分布式水文模型ESSI、HEC-HMS进行3场次洪过程模拟。从洪水过程模拟结果的对比分析来看，总体而言，XXT模型表现最突出，相比其他几个模型而言，XXT模型在江口流域次洪过程模拟中的确定性系数最高，模拟的洪峰流量和峰现时间较准确，模拟的洪水退水过程与观测值最为接近。对于江口流域次洪模拟，半分布式模型优于分布式模型。所有参与比较的模型对峰现时间的预报较好，确定性系数均在0.7以上，均达到洪水预报标准乙等以上，能够正确描述洪峰过程的变化趋势，但均低估了洪峰流量。

b）在温暖湿润气候区的汉江上游选择了不同尺度的四个子流域：汉中流域（集水面积9329km2）、褒河流域（马道控制站，集水面积3415km2）、褒河流域（江口控制站，集水面积2501km2）和酉水街流域（集水面积908km2，见第二章）进行XXT模型对日径流模拟能力的测试，并与新安江和TOPMODEL进行了比较。结果表明在所有这个流域XXT模型的确定性系数均最高，特别是在大尺度的汉中流域，验证期XXT的确定性系数显著高于TOPMODEL和新安江，分别比这两个模型高21.5%和7.5%。在江口流域，模拟结果还与分布式水文模型ESSI、SWAT模拟的日径流进行了对比，结果表明半分布式水文模型TOPMODEL、XXT和新安江在江口流域的日径流模拟能力高于分布式水文模型ESSI和SWAT，其中XXT表现尤其突出。在褒河流域（马道控制站）将XXT与TOPMODEL对地下水位的模拟能力作了比对，结果显示TOPMODEL模拟的流域平均地下水位竟然超出了地面，这显然与实际情况不符，而XXT模拟的土壤水分亏缺深未出现负值，表明XXT在模拟地下水位深时克服了TOPMODEL这方面的缺陷。

c）在黑河山区流域(莺落峡控制站，集水面积10009km²）选择了四次洪水过程对XXT进行洪水过程模拟能力评价，并与新安江、TOPMODEL和ESSI进行比较。所有模型模拟的确定性系数均在0.8以上，达到了业务洪水预报标准，预报水平均在乙等以上，与温润区的中尺度江口流域相比表现更优。与温润区的中尺度江口流域相比，分布式水文模型ESSI对洪水的预报能力显著增强。总体而言，XXT和ESSI表现较优，新安江次之，TOPMODEL表现最差。相比其他几个模型而言，XXT模型在黑河流域次洪过程模拟中的确定性系数较高，模拟的洪峰流量出现一致性的偏低，没有一次高估洪峰流量，而其他几个模型或高或低地估计了洪峰流量。

在该流域用1990-2000年共11年的降水、蒸发和径流数据对XXT模型日径流模拟能力进行评估，并与新安江、TOPMODEL、ESSI和SWAT模拟结果进行比较。从总体而言，所有参与比较的模型在黑河山区流域的模拟效果均不如江口流域，模拟效果也远不如该流域对次洪过程的模拟。其中XXT模型表现最优，新安江和分布式水文模型ESSI和SWAT表现较好，TOPMODEL表现最差。

d）在沂沭河流域(临沂控制站，集水面积10315km2）选择了四次洪水过程对XXT进行洪水过程模拟能力评价，并与新安江和TOPMODEL进行比较。所有模型模拟的确定性系数均在0.7以上，达到了业务洪水预报标准，预报水平均在乙等以上。与温润区的中尺度江口流域相比表现相差不大，与黑河山区流域（鹰落峡控制站）相比总体表现较差。总体而言，XXT表现最优，新安江次之，TOPMODEL表现较差，这可能与TOPMODEL主要适用于中小尺度流域有关。相比其他几个模型而言，XXT模型在林沂流域次洪过程模拟中的确定性系数较高，模拟的洪峰流量出现“稳定”的偏低。

在该流域用2001-2003年（参数率定期）和2006-2007年（模型验证期）的降水、蒸发和径流数据对XXT模型日径流模拟能力进行评估，并与新安江和TOPMODEL模拟结果进行比较。XXT模型表现较优，新安江次之，TOPMODEL表现较差。XXT模型模拟值较其他模型更接近观测值，更能捕捉几次主要洪峰，而新安江和TOPMODEL模型过高估计了主要洪峰。从总体而言，参与比较的这三个模型在临沂流域的模拟效果均不如江口流域和汉中流域，但比黑河流域要好。

总的来看，XXT模型在温暖湿润气候区的汉江上游的汉中流域、褒河流域（马道控制站和江口控制站）和酉水街流域，半湿润半干旱区的淮河流域的子流域—沂沭河流域（沂临控制站）以及位于干旱、半干旱区的黑河流域（鹰落峡控制站）的水文过程模拟都取得了成功，说明在新蓄水容量曲线和简单的类似TOPMODEL模型的简单框架结构的机制下，XXT模型不仅可以完成不同气候条件下流域的水文过程模拟，并具有令人满意的模拟精度，而且可以完成从短期次洪过程到长期日径流过程模拟，实现流域的洪水预报和长时段的径流过程模拟和预测，其简单的结构，较少的参数个数，较高的模拟效率势有利于耦合到陆面过程模式中，进而与数值天气预报系统或气候模式耦合。

e）TIGGE是一项计划，提供了多个集合预报系统的预报结果。TIGGE资料给早期洪水预报预警提供了一个很好的新途径。本发明基于TIGGE和水文模型XXT对洪水集合预报进行了有益的尝试。以沂沭河流域(临沂控制站）为研究区，所用资料为全球三大数据归档中心（CMA、ECMWF和NCEP）的TIGGE资料。对TIGGE资料进行插值后直接驱动XXT进行洪水预报试验。结果表明，TIGGE数据集和XXT水文模型可以成功预报洪水过程，能较准确地预报出峰现时间和洪峰流量所带的风险概率。与确定性洪水预报相比，基于TIGGE和XXT的洪水集合预报能给用户（包括政府部门）提供更多的有价值的信息和风险概率，是一条颇具前景的新方法。

综上所述，XXT模型在多个流域尺度和多气候生态条件下表现均很突出，模拟的确定性系数高，总体而言，对日径流过程和次洪过程的模拟均优于半分布式模型新安江和TOPMODEL以及分布式水文模型ESSI、HEC-HMS、SWAT。另外，相比其他几个模型而言，XXT模型在江口流域、黑河流域和临沂流域次洪过程模拟的洪峰流量出现“一致性”的偏低，这一点在集合预报中具有重要应用价值，该模型能给出集合预报的一个“下限”，而且XXT模拟的这个下限与观测值比较接近，其他几个模型都出现或高或低地估计洪峰流量，因而都不具备这一优点。下一步工作中的一个重要研究课题就是通过对多个水文模型的产汇流方案和模型框架结构进行研究，以期构建一个能稳定地略微偏高估计洪峰流量的水文模型，以为集合预报提供一个“上限”成员。

（3）XXT与Noah LSM耦合成功

首先将XXT的产流方案和汇流方案分别耦合到Noah LSM中，并进行off-line测试，然后将耦合XXT的Noah LSM应用到数值天气预报模式GRAPES-MESO中并进行on-line测试，探索新耦合方案对GRAPES-MESO降水预报的影响，再结合实测土壤湿度资料和由AMSR-E遥感影像反演的土壤湿度对GRAPES-MESO土壤湿度预报效果进行分析和检验，最后考察新耦合方案对GRAPES-MESO径流预报的影响。得出以下结论：

a）用XXT模型中的地表径流方案取代Schaake（1996）的简单水量平衡模型SWB，从而实现了XXT的产流方案与Noah LSM的耦合。新耦合方案中参数B是表示土壤蓄水容量空间分布不均匀程度的一个物理量，需要率定或能过其它方法获得。借鉴TOPMODEL模型中用DEM提取得到的地形参数而导出模型参数的方法，本文利用DEM提取参数B。其提取原理是依据土壤蓄水容量空间分布公式及最小二乘法。基于DEM获取参数B使得新耦合NoahLSM具有空间分布式特征。耦合分布式XXT模型的NoahLSM为空间分布式的，其参数并不需要用具体的流域资料进行率定，这使得新Noah LSM可移植性增强、适用性更广。

b）在原等流时线汇流方案的基础上，开发了一个新的基于格点DEM的分布式汇流方案算法。并将此新分布式汇流方案耦合到Noah LSM中。

c）为了评价耦合XXT后的新Noah LSM的模拟性能，首先对其进行离线测试。研究区位于美国Illinois州的Champaign，测试及验证数据为此区域的地面通量观测站(40.01°N,88.37°W)的1998年通量观测数据及气象数据。结果表明：新Noah LSM模式显著提高了地表潜热通量和表潜热通量的模拟能力。新Noah LSM提高了地表径流量。结合已有研究报道可知新Noah LSM表面优于原模式。为了进一步证明此结论，对新旧Noah LSM模拟的水量平衡残差进行了分析，水量平衡残差越大表明模式对水文过程描述能力越差。而新Noah LSM模拟的水量平衡残差远小于原Noah LSM模式，这进一步证实了新耦合模式比原Noah LSM在刻画水文过程时具有明显的优势。

d）新Noah LSM被耦合到中国气象局自主开发的业务数值天气预报模式GRAPES-MESO中进行在线测试。为了验证耦合了XXT的Noah LSM在线模拟效果，本文选沂沭河流域上游的临沂子流域（东经117.35°-118.85°，北纬34.9°-36.40°）为研究区；个例选择为山东临沂地区2007年8月17日-19日的一次暴雨过程，起报时间为2007年8月17日00时（世界时，本文中如未特别声明即为北京时），作72小时预报。所使用GRAPES-MESO为中国气象局气象中心业务运行版本,其空间分辨率为0.15°×0.15°。设计了3个试验方案：（1）陆面过程为热扩散方案的GRAPES-MESO模拟；（2）陆面过程为原Noah LSM的GRAPES-MESO模拟；（3）陆面过程为耦合XXT的新NoahLSM的GRAPES-MESO模拟。结果表明：带Noah LSM的GRAPES-MESO预报的降水比带热扩散方案的GRAPES-MESO预报的降水更接近观测值，不仅表现在降水数量级上，还表现在空间分布上。新Noah LSM方案的GRAPES-MESO预报的降水显著优于比带原Noah LSM方案的GRAPES-MESO。

将GRAPES-MESO预报的第一层土壤湿度与气象台站测定的10cm土壤湿度进行比对的结果显示：带原Noah LSM的GRAPES-MESO预报土壤湿度偏小，而带耦合了XXT的Noah LSM的GRAPES-MESO预报土壤湿度比带原Noah LSM的GRAPES-MESO预报土壤湿度在空间上分布更合理，在量值上更接近气象台站实测值，尤其在山东省和辽东半岛效果更明显。而带有耦合了分布式XXT的GRAPES-MESO预报土壤湿度比耦合半分布式XXT（参数B在全国取常数，即取在临沂流域率定的值）的GRAPES-MESO预报土壤湿度更接近观测值。

利用AMSR-E遥感反演的土壤湿度对GRAPES-MESO模拟的土壤湿度进行了进一步评估，结果表明：AMSR-E遥感反演土壤湿度值小于带有新Noah LSM方案的GRAPES-MESO模拟的土壤湿度，但在土壤湿度的空间分布形态上两者比较相似。带新Noah LSM方案的GRAPES-MESO模拟的土壤湿度在空间分布比带原Noah LSM方案的GRAPES-MESO更接近AMSR-E遥感反演土壤湿度的分布型，这进一步说明了耦合了XXT的Noah LSM能提高GRAPES-MESO对土壤湿度的预报能力。

带耦合了XXT的Noah LSM方案的GRAPES-MESO模拟的径流深在临沂有较大值，与该区域降水呈较好的正相关，而带原Noah LSM方案的GRAPES-MESO未能模拟出该区域的径流深，这说明不论是带半分布式XXT的Noah LSM方案的GRAPES-MESO还是带分布式XXT的Noah LSM方案的GRAPES-MESO对径流深的预报都要优于带原Noah LSM方案的GRAPES-MESO。尽管GRAPES-MESO本身对这次暴雨过程的降水预报不够准确，但经过分布式汇流之后，模拟的洪峰流量与实测值比较接近，基本捕捉到最大洪峰，这说明耦合了新分布式XXT的GRAPES-MESO能够成功地作洪水预报。

目前中国气象局国家气象中心使用的业务GRAPES-MESO模式中的陆面过程方案为SLAB方案，该方案没有土壤湿度的预报，没有地表径参数化方案，不能进行径流和洪水预报。而耦合分布式XXT的GRAPES-MESO模式不仅改善了降水预报，提高了土壤湿度预报能力，还能成功地作洪水径流预报。

说明书附图

图1为本发明的技术路线图；

图2为新降雨-径流模型XXT框架结构图；

图3为DEM洼地及平坦区处理新算法流程图。

具体实施方式

以下通过实施例对本发明作进一步描述：

实施例

小流域酉水街的集水面积约为908km²，其地理位置为陕南的中部地区，处在陕西省安康市和汉中市接壤地带。控制站点位于汉江支流酉水河流域的酉水街水文站（33°17′N，107°46′E）。流域内地貌条件复杂，气候变化多样，垂直差异明显，为典型的亚热带湿润气候区，平均气温变化在7.7-15.7℃之间，常年降水变化在730-1230mm之间。土壤类型比较单一，主要有黄褐土、黄棕壤等，土壤质地以颗粒较大的山地土壤为主。研究区内海拔范围425-2960m，植物类型较多，植物分布以温带成分为主，具有暖温带和亚热带植被相互过度的自然特征，属于暖温带落叶林地带和北亚热带常绿落叶混交林地带。植被类型可分为森林、野生植物和传统经济栽培作物，森林主要树种为华山松、油松和栎等，主要土地覆被为林地、灌木、草地和农田等。

单点离线测试的主要目的是检验水文模型本身对水文过程描述是否合理、模型结构是否存在大的问题、是否能反映陆面水文过程内在机理的主要特征等。因此为了更好地分析和检测XXT模型对重要水文过程机制的模拟能力，我们接下来将详细分析XXT的径流量、蒸散发等模型输出量，以及各层土壤湿度、饱和面积时空动态变化等重要水文过程中间变量的模拟效果。

表2给出了XXT、TOPMODEL和新安江三个模型1981-1983年在酉水街流域率定的参数及NE值。表3、表4和表5分别给出了这三个模型的参数说明。由表2可知，XXT模型有7个率定参数：f、Q₀、gama、SRmax、WM和B的物理意义前面已说明，CHV是汇流速度。XXT模型比TOPMODEL多1个，多出的两个参数（WM和B）来自于蓄水容量曲线。新模型比新安江少4个参数，因此其参数率定比新安江模型容易得多。在参数率定期XXT的模型效率高达0.906，高于TOPMODEL和新安江模型。

表2：三个模型在酉水街流域率定的参数（率定期为1981-1983年）

表3：TOPMODEL模型参数

表4：新安江模型主要可调参数

表5：XXT模型参数

Saulnier和Datin（1997）研究指出TOPMODEL模拟的流域平均土壤水分亏缺深常常出现负值的情况，这与理论假定及实际田间观测不符。TOPMODEL在酉水街流域的模拟结果也出现了这种问题。但XXT模拟的土壤水分亏缺深大于或等于零，有效地解决了TOPMODEL这一问题。

小流域（酉水街流域）的测试结果表明，因1981-1983年为丰水年，径流系数较大，三个模型在率定期模拟效率都很高，均在83%以上，其中XXT模拟效率最高；而1984-1985年为枯水年，径流系数小，三个模型在验证期的模拟效率都较低，其中XXT模拟效率明显高于另外两个模型，分别比TOPMODEL和新安江模型高5.2%和3.2%，这在湿润区来说是很不容易的（已有大量报道说湿润区水文模型模型效率差别不大）。TOPMODEL和XXT模型的地表产流及地下产流方案均使用了变量Z，亦即土壤水分亏缺深，并通过该变量将地表产流和地下产流关联，这一点使该模型在表达地下水位时具有明确的物理意义，而新安江模型不具备此优点。而且，TOPMODEL模拟的流域平均土壤水分亏缺深出现了负值，说明TOPMODEL模拟的流域平均地下水位高于地面，显然不符合实际情况，前人也指出TOPMODEL存在此缺陷。而XXT模拟的土壤水分亏缺深则未出现负值，从此角度上看XXT在模拟地下水位深时比TOPMODEL更合理。

新Noah LSM的Off-line测试：

为了评价耦合XXT后的新Noah LSM的模拟性能，首先对其进行离线测试。研究区位于美国Illinois州的Champaign，测试及验证数据为此区域的地面通量观测站(40.01°N,88.37°W)的1998年通量观测数据及气象数据。该通量观测站数据由NOAA/ARL的Tilden Meyers维护，提供了1998年的间隔30分钟一次的大气强迫数据、地表能量能量数据及土壤温度数据。本发明所用Noah LSM为耦合了XXT的Noah 2.7.1版，系一维土柱的陆面过程模式。

耦合了XXT的新Noah LSM与原Noah LSM模式模拟的地表感热通量对比可以看出观测的感热通量在初夏相对较低，峰值为150-180W/m2，其原因在于这个季节该区域正是植被生长旺盛时期，蒸散量大。用新Noah LSM模式模拟的地表感热通量比较接近观测值，而原Noah LSM模拟值低估了观测值。

模拟地表潜热难度大一些，因为植被的生理过程和蒸散发的相互作用机理非常复杂。与原Noah LSM相比，新Noah LSM模式显著提高了地表潜热通量的模拟能力，观测的最大潜热通量值发生1998年5月28日，其值约为W/m2。尽管新旧Noah LSM都高估了潜热通量，但新方案明显更接近观测值。

显示了新旧Noah LSM模式模拟的径流深对比。总体而言，新Noah LSM提高了地表径流量。已有不少研究报道陆面过程模式（如LSM、Noah LSM、CLM等）中的地表产流方案常常低估了观测值。因此，对地表产流而言，新Noah LSM表面优于原模式。为了进一步证明此结论，对新旧Noah LSM模拟的水量平衡残差进行了分析。下式是水量平衡残差计算公式。

WRES=(Precip+Dew-(Runoff+Evap))*DT

式中WRES代表水量平衡残差，Precip为未经植被截留的直接降水，Dew为植被截留中进入地表的水，Runoff为总径流量（地表径流和地下径流），DT为时间步长。

水量平衡残差越大表明模式对水文过程描述能力越差。从中可以看出新Noah LSM模拟的水量平衡残差远小于原Noah LSM模式，这进一步证实了新耦合模式比原Noah LSM在刻画水文过程时具有明显的优势。

新Noah LSM的on-line测试：

由于Noah LSM已广泛应用于众多数值天气预报模式(ETA、MM5、WRF和GRAPES等)中，仅对Noah LSM做离线测试不足以说明其优势。故有必要将耦合XXT的新Noah LSM放入数值天气预报模式进行检验。结合国家气象行业重大专项″面向TIGGE的集合预报关键应用技术研究″（GYHY200706001）项目的需要，我们选择了将新Noah LSM放入中国气象局自主开发的业务数值天气预报模式GRAPES-MESO中进行测试。

由于数值预报在气象业务服务中的基础性作用以及数值模式在大气与地球科学研究中的重要地位，自主发展我国的数值天气预报系统成为世纪之交我国大气科学界与管理部门共同关注的一个热点。2001年中国气象局把数值预报系统创新技术的发展列为气象科技创新的最重要项目之一，组织了数值预报技术创新的专门队伍与机构，同时，国家科技部在“十五”国家重点科技攻关计划中启动了“中国气象数值预报系统技术创新研究”项目，开始自主研究开发我国新一代全球/区域同化预报系统（GRAPES，Global/Regional Assimilation PrEdiction System）。经过5年的联合攻关，完全依靠中国科学家自己的力量，研究建立了一套新的数值预报系统，它以变分同化和全可压非静力平衡模式的理论为基础，并将全球/区域的资料同化、预报模式“一体化（unified）”的思想贯穿于系统设计和工程实现之中。目前GRAPES系统已在全国业务、科研、教学等许多单位得到了多方面的应用。

GRAPES是Global/Regional Assimilation and Prediction System的简称，它从理论方案设计、程序编写，到系统运行都是按“一体化（unified model）”的思想进行构造的，它可以通过简单的参数置换选用区域GRAPES MESO和全球GRAPES Global两个版本。GRAPES Global是在GRAPES-MESO基础上发展起来的全球中期数值预报系统，相对来说，它增加了包括近极地半拉格朗日处理方法、极区滤波、极区Arakawa-C网格设计及由此带来的Helmholtz稀疏矩阵的计算方法等一些特殊的处理技术。

GRAPES系统的核心部分包括模式动力框架、经过优化选取和改进的物理过程参数化方案、三维变分资料同化系统（下一步将更新为四维变分同化系统）、模式标准初始化系统。

●GRAPES系统构成

GRAPES系统的核心包含3个部分，分别是GRAPES三维变分资料同化系统、模式标准初始化系统以及预报模式系统。为了在业务环境下实时运行，除了同化和预报系统，还需要解决其它一些相关的环节：如实时观测资料从要素库的检索、观测资料的质量控制、同化所需背景场和模式积分所需边界条件的选取和接入，分析系统与模式的标准化接口、模式预报结果的后处理、模式产品的入库、模式产品图形显示和模式结果检验以及整个系统的作业管理等。另外，预报模式部分的计算是并行实现的，系统运行结合了计算机平台的并行环境建设。

●GRAPES物理过程

GRAPES的物理模式主要包括积云对流、微物理、辐射、行星边界层、陆面面过程等物理过程参数化。参数化物理过程之间有着直接相互作用：行星边界层影响地表要素风速V，温度T，水汽QV；地表面的感热、潜热输送影响行星边界层；大气向下长波、短波辐射影响地表的热量平衡；地表的向上反射、长波辐射影响大气的辐射过程；积云覆盖量和云效应影响大气的辐射过程；积云的卷出量影响微物理过程。另外，这些物理过程以及扩散过程可以改变大气环境场的分布，从而产生间接的相互作用。

●GRAPES中的陆面过程方案

GRAPES中的陆面过程主要有两种方案：（1）SLAB方案将土壤分为5层，每层土壤均考虑向上、向下的热通量，并通过热平衡方程对每一层土壤的温度进行预报。其初始条件通过地表温度（由大尺度背景场，如T213，NCEP等提供）与土壤深层温度（NCAR提供的气候平均值）经过简单的线性插值得到。但该陆面过程模式没有土壤湿度的预报。（2）Noah LSM方案它相对比较复杂，它包括了一个4层土壤的模块和一层植被冠层的植被模块。不仅能够提供土壤温度的预报，还可以预报土壤湿度、地表径流等。其初始温度、湿度场均由大尺度场提供的信息经过插值得到。目前中国气象局国家气象中心使用的业务GRAPES模式中的陆面过程方案为SLAB方案，该方案没有土壤湿度的预报，没有地表径参数化方案，不能进行径流预报。

●研究区及数据

为了验证耦合了XXT的Noah LSM在线模拟效果，本文选沂沭河流域上游的临沂子流域（东经117.35°-118.85°，北纬34.9°-36.40°）为研究区；个例选择为山东临沂地区2007年8月17日-19日的一次暴雨过程，起报时间为2007年8月17日00时（世界时，本文中如未特别声明即为北京时），作72小时预报。所使用GRAPES-MESO为中国气象局气象中心业务运行版本,其空间分辨率为0.15°×0.15°。本研究GRAPES-MESO的模拟范围有两个：一是沂沭河区域（东经113°-122°18’；北纬31°-44°18’，62×62格点）；二是中国区域（东经75°-145°18’，北纬15°-64°30’，共330×502格点）。GRAPES-MESO的初始场和侧边界条件由T213提供。Noah LSM的土壤湿度初始场数据也来自于T213。T213是第四代全球中期数值天气预报系统，于2002年9月1日起在中国气象局投入业务使用。其预报模式的水平分辨率约为60km左右，垂直分辨率为31层（σ）。每天可制作一次10天的中期数值天气预报，侧重于对全球大气环流演变和时间空间尺度较大的天气系统发展变化的预报，或者说较早地（如提前一周）对转折性、灾害天气过程预报，提供中期天气过程的预报参考，同时也为全国各有限区域短期数值天气预报模式提供初估场和侧边界条件。对降水预报效果的检验是用临沂及周边气象台站实测逐小时降水数据，经反距离加权法插值成面雨量。

●GRAPES-MESO/XXT在沂沭河模拟区域的预报结果与分析

由于原中国气象局气象中心业务运行的GRAPES-MESO版本中陆面过程方案并不是Noah LSM，而是热扩散方案（Thermal diffusion scheme），因此不仅要对耦合XXT的新Noah LSM与原Noah LSM进行比较，还应当将Noah LSM与热扩散方案作对比。鉴于此，我们设计了3个试验方案：（1）陆面过程为热扩散方案的GRAPES-MESO模拟；（2）陆面过程为原Noah LSM的GRAPES-MESO模拟；（3）陆面过程为耦合XXT的新Noah LSM的GRAPES-MESO模拟。在PC机上的linux操作系统下运行GRAPES-MESO。此GRAPES-MESO未带有数据同化系统。

●降水预报结果与分析

为了进一步测试新Noah LSM对GRAPES-MESO预报能力的影响，在中国气象局气象中心的IBM工作站上进行GRAPES-MESO模拟试验，此GRAPES-MESO带有数据同化系统，和气象中心所用GRAPES-MESO版本完全相同。仍然设计了3个试验方案：（1）陆面过程为热扩散方案的GRAPES-MESO模拟；（2）陆面过程为原Noah LSM的GRAPES-MESO模拟；（3）陆面过程为耦合XXT的新Noah LSM的GRAPES-MESO模拟。

带Noah LSM的GRAPES-MESO预报的降水比带热扩散方案的GRAPES-MESO预报的降水更接近观测值，不仅表现在降水数量级上，还表现在空间分布上。对临沂的这次强降水过程，在起报6小时后带热扩散方案的GRAPES-MESO几乎没有预报出来（空报），新Noah LSM方案的GRAPES-MESO预报的降水显著优于比带原Noah LSM方案的GRAPES-MESO，降水量值更接近观测值，范围和落区也与观测值更吻合。

●土壤湿度预报结果与分析

Noah LSM中将土壤分为4层，厚度分别为0.1m,0.3m,0.6m,1m，总厚度为2m。气象台站测定的土壤湿度有10cm、20cm、50cm等几个土壤层次，其中10cm数据在全国范围内较为齐全，也刚好同Noah LSM第一层土壤厚度相同。因此本试验将GRAPES-MESO预报的第一层土壤湿度与气象台站测定的10cm土壤湿度进行比对，以判别模式的优劣。

带有不同Noah LSM的GRAPES-MESO预报土壤湿度空间分布（2007年8月17日）。可见带原Noah LSM的GRAPES-MESO预报土壤湿度偏小，而带耦合了XXT的Noah LSM的GRAPES-MESO预报土壤湿度比带原Noah LSM的GRAPES-MESO预报土壤湿度在空间上分布更合理，在量值上更接近气象台站实测值，尤其在山东省和辽东半岛效果更明显，如图中红色圆圈所圈的部分所示。而带有耦合了分布式XXT的GRAPES-MESO预报土壤湿度比耦合半分布式XXT（参数B在全国取常数，即取在临沂流域率定的值）的GRAPES-MESO预报土壤湿度更接近观测值。本次个例中，耦合了XXT分布式方案的Noah LSM达到了最好的模拟效果，这说明XXT耦合了XXT分布式方案的Noah LSM不仅能够改善GRAPES-MESO降水预报效果，还显著改善了GRAPES-MESO对土壤湿度的模拟性能。

前面用气象台站实测土壤湿度数据对GRAPES-MESO土壤湿度的预报能力进行了评估。尽管测土壤湿度数据比较可靠，但由于在全国范围内我们能获得的壤湿度数据的站点个数不够多，且分布很不均匀，更重要的是每旬只测定一次，时效性太差。为了解决这一问题，我们用AMSR-E遥感反演的土壤湿度对GRAPES-MESO模拟的土壤湿度进行评估。之所以采用AMSR-E遥感数据，是因为它具有空间分布式特点且时效性好（1-3天即有一次观测值）。

AMSR-E遥感反演土壤湿度值小于带有新Noah LSM方案的GRAPES-MESO模拟的土壤湿度。原因之一是AMSR-E遥感反演所测定的土壤湿度是地表顶部很薄的一层土壤的土壤湿度，而GRAPES-MESO模拟的土壤湿度是距表层10cm处的土壤湿度，故而在数量上两者肯定存在差别。但在土壤湿度的空间分布形态上两者比较相似。带新Noah LSM方案的GRAPES-MESO模拟的土壤湿度在空间分布比带原Noah LSM方案的GRAPES-MESO更接近AMSR-E遥感反演土壤湿度的分布型，这进一步说明了耦合了XXT的Noah LSM能提高GRAPES-MESO对土壤湿度的预报能力。

●径流预报结果与分析

Noah LSM中水文过程描述的优劣直接反映在径流深的模拟上。新旧Noah LSM的GRAPES-MESO模拟的地表径流深。尽管没有实测的中国区域径流深资料，但由于本次降水过程之前已有一次较大的降水，前期研究区土壤基本趋于饱和，地表产流理应与降水有较大的正相关性。带耦合了XXT的Noah LSM方案的GRAPES-MESO模拟的径流深在临沂有较大值，与该区域降水呈较好的正相关，而带原Noah LSM方案的GRAPES-MESO未能模拟出该区域的径流深，这说明不论是带半分布式XXT的Noah LSM方案的GRAPES-MESO还是带分布式XXT的Noah LSM方案的GRAPES-MESO对径流深的预报都要优于带原Noah LSM方案的GRAPES-MESO

对于这次发生在临沂研究区的暴雨过程，GRAPES-MESO预报的降水与观测值有较大出入。GRAPES-MESO预报的降水出现多个峰值，但未能抓住主峰，预报的降水强度偏小。这说明GRAPES-MESO模式本身对定量降水的预报能力还有待提高。

用新汇流方案对GRAPES-MESO/XXT预报的径流深进行汇流，汇流到流域出口（临沂市水文站）的洪水流量与观测值的比较。尽管GRAPES-MESO本身对这次暴雨过程的降水预报不够准确，但经过分布式汇流之后，模拟的洪峰流量与实测值比较接近，基本捕捉到了最大洪峰，这说明耦合了新分布式XXT的GRAPES-MESO能够成功地作洪水预报。

之所以这次洪水预报不是很准确，还有一个重要原因，就是在研究区内人工水库众多，对洪水产流及汇流都会产生较大影响。临沂站以上流域有大大小小水库几十个，水库的存在可以有效地调节洪水，降低洪水的风险，但同时也给洪水预报工作带来了挑战。本次个例研究中对出口流量的模拟值偏大，其中重要原因之一是未能在模式中考虑水库的蓄水作用。

Claims (1)

1.一种基于蓄水容量曲线和TOPMODEL的流域水文模型的设计方法，其步骤为：

A）地表产流的计算；

将地下水水位线与新安江模型中的流域蓄水容量曲线关联形成新流域蓄水容量分布曲线；令流域上地下水水位线刚好位于地面0mm处的蓄水容量为

则新流域蓄水容量曲线定义为蓄水容量小于或等于

的流域面积；把全流域面积看作1，流域由无数个点组成，每个点都对应一个蓄水容量，其中最大的蓄水容量记为

α为饱和区面积与全流域面积之比，其计算公式为

$\mathrm{\&alpha;}=1-{(1-\frac{W_m^{\&prime;}}{W_{\mathrm{mm}}^{\&prime;}})}^B$

式中B是表征蓄水容量分布不均匀性的参数，B愈大分布愈不均匀，据上式有：

$W_m^{\&prime;}=W_{\mathrm{mm}}^{\&prime;}\&times;[1-{(1-\mathrm{\&alpha;})}^{1/B}]$

令WM为流域平均蓄水容量，则WM等于未饱和层与饱和层的面积之和：

$\mathrm{WM}=\underset{0}{\overset{1}{\&Integral;}}{W}_m^{\&prime;}\mathrm{d\&alpha;}=\underset{0}{\overset{1}{\&Integral;}}{W}_{\mathrm{mm}}^{\&prime;}[1-{(1-\mathrm{\&alpha;})}^{1/B}]\mathrm{d\&alpha;}=W_{\mathrm{mm}}^{\&prime;}/(1+B)$

地表产流计算如下：

P为降水量，SRZ为植被截留层某一时刻的土壤水分亏缺深，P-SRZ是从植被截留层中渗透下来的降水量；α₁为P-SRZ的降水进入土壤未饱和层之前的饱和区面积与全流域面积之比，α₂为P-SRZ的降水进入该层之后的饱和区面积与全流域面积之比；

${\mathrm{\&alpha;}}_1=1-{(1-\frac{W_{m1}^{\&prime;}}{W_{\mathrm{mm}}^{\&prime;}})}^B$

${\mathrm{\&alpha;}}_2=1-{(1-\frac{(W_{m1}^{\&prime;}+P-\mathrm{SRZ})}{W_{\mathrm{mm}}^{\&prime;}})}^B$

令

代表流域平均土壤水分亏缺深，

为

对应的流域上某点的蓄水容量，的计算公式推导过程如下：

$\stackrel{\&OverBar;}{Z}={\&Integral;}_{{\mathrm{\&alpha;}}_1}^1{W}_m^{\&prime;}\mathrm{d\&alpha;}-W_{m1}^{\&prime;}\&times;(1-{\mathrm{\&alpha;}}_1)$

先解出积分式，再化简，得

$W_{m1}^{\&prime;}=\mathrm{WM}\&times;(1+B)\&times;[1-{(\stackrel{\&OverBar;}{Z}/\mathrm{WM})}^{1/(1+B)}]$

地表产流的计算分三种情况：当P-SRZ≤0时，没有地表产流Q_S，即：

Q_S=0

当 $0<P-\mathrm{SRZ}\&le;\mathrm{WM}\&times;(1+B)-W_{m1}^{\&prime;}$ 时：

$Q_S=(P-\mathrm{SRZ})\&times;{\mathrm{\&alpha;}}_2-[{\&Integral;}_{{\mathrm{\&alpha;}}_1}^{{\mathrm{\&alpha;}}_2}{W}_m^{\&prime;}\mathrm{d\&alpha;}-({\mathrm{\&alpha;}}_2-{\mathrm{\&alpha;}}_1)\&times;W_{m1}^{\&prime;}]$

化简后为

$Q_S=P-\mathrm{SRZ}-\stackrel{\&OverBar;}{Z}+\mathrm{WM}\&times;{\{1-(P-\mathrm{SRZ}+W_{m1}^{\&prime;})/[\mathrm{WM}\&times;(1+B)\left]\right\}}^{(1+B)}$

当 $P-\mathrm{SRZ}>\mathrm{WM}\&times;(1+B)-W_{m1}^{\&prime;}$ 时，易知地表产流为

$Q_S=P-\mathrm{SRZ}-\stackrel{\&OverBar;}{Z}$

由于

和SRZ是状态变量，故该产流方案中只有两个可调参数，WM和B；WM是整个流域很干燥时的缺水量，代表流域的干旱程度，是个气候因素，B代表蓄水容量空间分布的不均匀性，其值愈大蓄水容量空间分布愈不均匀；

B）地下产流的计算；

地下产流引用地下径流公式：

$Q_b=Q_0{e}^{-f\stackrel{\&OverBar;}{Z}}$

式中Q_b为地下径流量，Q₀表示暴雨前退水期地下水流量，f是土壤水力传导率沿深度方向的衰减因子,

代表流域平均土壤水分亏缺深；

C）蒸散发计算；

实际蒸散发Q_E的计算方法为先用下面的公式计算

Q_E＝E×(1-SRZ/SR_max)

上式中E为潜在蒸散发，其值可直接利用蒸发皿观测值，如果没有观测值可用气温、辐射和风等气象资料代入彭曼公式计算得到；当Q_E值计算之后还需作判断，如果Q_E＞SR_max-SRZ，则Q_E应重新赋值：

Q_E＝SR_max-SRZ

此蒸散发计算方案中只有一个参数需要率定，即SRmax；

D）水量平衡方程的应用；

依据下面的水量平衡方程：

${\stackrel{\&OverBar;}{Z}}_{t+1}-{\stackrel{\&OverBar;}{Z}}_t=Q_{\mathrm{Et}}+Q_{S_t}+Q_{b_t}-P_t$

式中表示t+1时刻的流域平均土壤水分亏缺深，Q_Et、

和P_t分别表示t时刻的流域平均土壤水分亏缺深、实际蒸散发、地表径流深、地下径流深和降水；

E）汇流机制；

对子流域内部采用等流时线汇流，等流时线基于DEM提取得到；此部分只需要一个率定参数，即汇流速度CHV；河道汇流采用马斯京根法；

F）等流时线的生成；

假定径流在空间上相等，考虑流速由于受到坡面、粗燥度等的影响是变化的，把地表径流与基流合在一起计算，得出出水点的流量；根据到出水点的距离与一个地表水流速度参数ChVel计算即可得到等流时线，

基于DEM提取等流时线的步骤为：

1）DEM填洼；2）DEM流向计算，坡度计算；3）计算每栅格到出水口的距离；4）将距离矩阵按quantile分类方法进行分类，并在ARCGIS9.0中生成可视化图；5）生成等流时线文件；

G）土壤湿度的动态可视化表达；

XXT模型基于变动产流源面积和由DEM提取的地形指数可以实现对土壤湿度的动态可视化表达；其实现步骤主要包括：1）DEM中“洼地”填充；2）计算DEM的水流方向矩阵和汇流面积；3）计算地形指数；4）根据XXT模型模拟的变动产流面积和地形指数绘制土壤湿度空间分布图；

H）与NoahLSM的耦合；

将XXT模型与Noah LSM耦合，第一步需要用XXT模型中的地表径流方案取代Schaake的简单水量平衡模型SWB；XXT中的地表径流径流的公式表达式如下：

$Q_S=P-\mathrm{SRZ}-\stackrel{\&OverBar;}{Z}+\mathrm{WM}\&times;{\{1-(P-\mathrm{SRZ}+W_{m1}^{\&prime;})/[\mathrm{WM}\&times;(1+B)\left]\right\}}^{(1+B)}$

在Noah LSM中，Rsur=Q_S并且Pd=P-SRZ，代入得：

$R_{\mathrm{sur}}=P_d-\mathrm{DD}+\mathrm{WM}\&times;{\{1-(P_d+W_{m1}^{\&prime;})/[\mathrm{WM}\&times;(1+B)\left]\right\}}^{(1+B)}$

其中，

代表前一时刻在water table处的土壤蓄水容量值，可用下面的公式计算：

$W_{m1}^{\&prime;}=\mathrm{WM}\&times;(1+B)\&times;[1-{(\mathrm{DD}/\mathrm{WM})}^{1/(1+B)}]$

其中DD为四层土壤总的土壤水分亏缺深，在此方案中WM可由Noah LSM中的土壤属性计算得出，只剩下一个参数B是待定的；参数B是表示土壤蓄水容量空间分布不均匀程度的一个物理量，B愈大表示土壤蓄水容量空间分布愈不均匀；对于具体的某一个流域，B值可通过参数率定获得；对于较大的模拟区域，利用DEM提取参数B，依据土壤蓄水容量空间分布公式，先对DEM进行相高对程的累积分布进行统计，然后用最小二乘法计算出参数B；

${\mathrm{\&alpha;}}_1=1-{(1-\frac{W_{m1}^{\&prime;}}{W_{\mathrm{mm}}^{\&prime;}})}^B$

基于格点的分布式汇流方案算法；新算法主要由以下步骤构成：

（1）将Noah LSM的每个格点的产流插值到100m分辨率的格点；

（2）将USGS 90m分辨率的DEM重采样到100m，并进行DEM洼地填充，进而生成流向和坡度矩阵；

（3）径流沿着由DEM计算的流向汇流，其速度为VEL=V0*(1+a*exp(S))，穿过第一个网格需要的时间是T=L/VEL，其中V0为初始速度；

（4）计算出每一个网格的径流到汇流至流域出口所需要时间；

（5）计算出每一时间步长所有格点汇流至流域出口的总径流量；

产流方案及汇流方案已耦合到Noah LSM之中。

Images