一种河网区一维模型计算单元自动划分与编码的方法与流程

本发明涉及水环境模拟技术领域，尤其涉及一种河网区一维模型计算单元自动划分与编码的方法。

背景技术：

在水环境模型领域，河流水质模型从20世纪30年代起步，目前发展成熟的一维水动力水质模型有qual系列、mike11等，可结合多种技术进行综合运用，分析多种变量。一维水动力水质模型的基本方程是一维平流-扩散物质输运和反应方程，将模拟通道划分为一系列恒定的非均匀流段，然后将每一流段划分为多个计算单元，功能简便强大，得到了广泛应用。

同时，随着水环境管理的现代化，河网编码技术得到进一步。河网编码能充分反映水系河网中各河流的流向、河流之间、河流与水体之间的拓扑关系，以及与相关属性数据之间的关系，致力于水的精细化管理。现有国内外河网编码技术可以归结为两大类：一类是以河流的水文分析为主，从便于分析各个集水区域的水流关系和河流空间位置管理的角度，建立河网水系的拓扑关系；另一类则以河网上的相关属性管理为主，在分析各个集水区域的水流关系和河网空间位置管理的基础上，进一步分析河流与相关属性的相互关系。这两类方法考虑的问题侧重点不同，编码思路和方法也大相径庭。

目前，常用的河网编码技术有拓扑属性表法、二叉树编码、多叉树编码、pfafstetter法、干支拓扑编码、nhd等，其中，拓扑属性表法采用自然数顺序编码方式，编码载体为数字，支持单一或多河段的汇合编码；二叉树编码是一种基于二叉树理论，将整个河网概化成以出口河段为根节点的二叉树结构，并以二元形式表示的河网编码方法，可由上游河段的二叉树编码直接计算出其紧邻下游河段的编码值，进而快速定位该河段，以进行信息传递；pfafstetter法是基于河网的拓扑关系以及河道的集水面积，对流域进行从大到小的逐级划分和编码的方法，对于给定的区域，可自动识别所有上游子流域和河段、或者下游河段；nhd采用动态分段技术对河段编码，是一种新的线性特征的动态分析、显示和绘图的技术，在传统gis数据模型的基础上利用线性参照系统和相应算法，动态地计算出属性数据的空间位置，能充分反映水系中各河流的流向，以及河流之间、河流与水体之间、甚至河段的拓扑关系。这些河网编码技术各有优缺点，侧重不一，需根据具体的研究内容和编码要求选择适宜的编码方法，由此保证编码的稳定性，提高计算机的信息管理效率，实现真实高效地分析提取所研究河网的拓扑关系。

然而，一维水质模型中的河段划分与编码仍需人工进行，首先识别上下游关系，区分干支流，打断河流进行子河段划分得到各个子河段，之后根据子河段编码原则依次对子河段进行编码，最后利用工具箱中的工具获取子河段属性。这一系列操作耗时耗力，但目前并没有一种简易的自动河段划分与编码技术能满足研究人员对模型河流的抓取、河段划分与编码的需求。

技术实现要素：

有鉴于此，本发明提供了一种河网地区的一维模型计算单元自动划分与编码方法，将一维水动力水质模型与河网编码技术结合，得到一种自动识别河网、区分主干流、划分河段并进行编码的方法，能满足一维水质模型的河段数据获取的系列流程，且操作过程简单明了，能快速获得一维综合水质模型的河段输入数据，为一维综合水质模型的应用提供了一定技术支持，从而推动了水环境管理领域的现代化和高效化。

本发明提供一种河网区一维模型计算单元自动划分与编码的方法，包括以下步骤：

s1、获取待研究流域的dem数据，通过填洼、流向判断、汇流累积计算判断流域源头及出口，进一步确定河网区河流上下游关系及水系走向，输出矢量化后的河网形文件；

s2、利用所述河网形文件，根据河网区每一河流所控制的上游集水面积之和判定河流的干支流关系，从而进行支流划分得到n级支流，并按照从上游至下游、从上级支流至下级支流的顺序依次对n级支流进行编码；

s3、首先对干流的起点和终点、支流的起点进行编码，然后确定支流入汇点并将其作为河段划分节点，根据空间位置关系，按照从上游至下游、从干流至支流、从上级支流至下级支流的顺序依次对所述河段划分节点进行编码；

s4、依次检索编码后的河段划分节点，根据所述河段划分节点自动对河网区河流进行划分形成河段，并按照河段编码原则进行编码；

s5、根据河段编码顺序，按照子河段节点确定原则，确定子河段划分节点；

s6、根据子河段划分节点对河流进行划分形成子河段，再根据子河段编码原则对形成的子河段进行编码，得到河网区的完整编码；

s7、根据河网区一维模型的模拟需求，确定计算单元的大小并进行划分，再根据子河段空间位置关系，从上游至下游按照子河段的编码顺序以递增形式完成计算单元的编码；

s8、建立子河段属性表，根据河网形文件检索河流空间信息，得到计算单元的相关属性。

进一步地，所述步骤s1的具体过程为：

s101、获取流域的dem数据，进行填洼和削峰；

s102、根据d8算法确定dem数据中每个栅格单元的流向，生成水流方向栅格数据；

s103、根据水流方向栅格数据，对任一栅格单元，确定流入所述任一栅格单元的所有累积上游单元格数目，生成流域汇流能力栅格数据，进一步根据设置的nip阈值得到河流网格；

s104、对流域汇流能力栅格数据进行单值化处理得到单值化后的栅格数据；

s105、根据单值化后的栅格数据以及水流方向栅格数据，得到矢量化后的河网形文件，完成河网提取。

进一步地，步骤s2的具体过程为：

s201、根据流域汇流能力栅格数据得到河网中的最大积流量网格，从而确定干流出口，并由此向上游检索，利用各网格的集水面积之和判定其与相邻网格的干支流关系，其中，集水面积总和最大的网格为干流所在网格，由此划分得到干流，并将其编码为m；

s202、对干流上游输入口的网格进行判断，并将其作为1级支流出口，据此向上游检索，根据其上游各个网格的集水面积之和判定其与相邻网格的干支流关系，集水面积总和最大的网格为1级支流所在网格，以此划分得到1级支流；同样地，对所述1级支流上游输入口的网格进行判断，并将其作为2级支流出口，向上检索，根据其上游各个网格的集水面积之和判定其与相邻网格的干支流关系，集水面积总和最大的网格为2级支流所在网格，以此划分得到2级支流；以此类推，得到1至n级支流；

s203、按照从上游至下游、从上级支流至下机支流的原则进行编码，其中，对于同一级别的支流，按照上游至下游的顺序递增编码。

进一步地，所述河段编码原则为：按照从上游至下游、从干流至支流的原则对河段依次编码；从干流起点开始，若检索到河段划分节点，则对所述河段划分节点与干流起点之间的河段进行编码；然后由所述河段划分节点向对应的支流上游进行检索，从支流上游起点开始对河段进行编码，若检索到河段划分节点，则对之间的河段采用递增形式进行编码，若未检索到河段划分节点，则直接对该河段进行编码；重复上述过程，直至无支流存在。

进一步地，所述子河段节点确定原则将河流起始点、行政区划边界与河段交点、大型水库、闸泵站、水质监测断面、重要取水口、控制单元以及河流拐弯处作为子河段划分节点。

进一步地，所述子河段节点确定原则包括：将河段交点作为1级节点，将闸泵站作为2级节点，将行政区划或控制单元与河段的交汇处、以及控制断面作为3级节点，将水质监测断面、重要取水口作为4级节点，将河道拐点作为5级节点；

且根据所述子河段节点确定原则得到的子河段划分节点满足：子河段划分节点之间的距离适当，以确保子河段长度合理；河流上的大型水库自成一段，其入口、出口处为子河段划分节点。

进一步地，所述子河段编码原则为：按照河段编码从小到大、从上游到下游，依次对子河段划分节点之间的子河段采用递增形式进行编码。

进一步地，所述计算单元为河网区一维模型模拟水文水质的最小模拟单元。

进一步地，所述子河段属性表包括河长、河段编码、子河段起始点坐标、子河段起始点高程、以及是否为源头子河段。

本发明提供的技术方案带来的有益效果是：基于arcgis软件python的二次开发功能，获得流域河流拓扑关系，自动划分计算单元并进行编码，操作简洁，可拓展应用到全国各大河网水系的一维模型构建过程中，确保模型构建的高效率与准确率，有助于拓展河网中的水质水生态环境研究的模型体系，为我国水环境质量精细化管理提供支撑。

附图说明

图1是本发明实施例提供的河网区一维模型计算单元自动划分与编码的方法流程图；

图2是本发明实施例提供的支流编码示意图；

图3是本发明实施例提供的河段划分节点的编码示意图；

图4是本发明实施例提供的河段编码的流程图；

图5是本发明实施例提供的河段编码示意图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地描述。

请参考图1，本发明的实施例提供了一种河网区一维模型计算单元自动划分与编码的方法，包括以下步骤：

s1、河网识别：获取流域的dem数据，得到河网区河流的水文特征，通过填洼、流向判断、汇流累积计算等手段判断流域源头及出口，进一步确定河网区河流上下游关系及水系走向，输出矢量化后的河网形文件，以反映河流的上下游信息、汇流关系、干支流信息以及河网要素信息等流域拓扑关系。

步骤s1的具体过程为：

s101、获取流域的dem数据，利用hydrology工具包中的fill工具进行填洼和削峰。

s102、根据d8算法确定dem数据中每个栅格单元的流向，生成水流方向栅格数据。

具体地，计算每一栅格单元与其相邻的8个单元之间的坡度，然后按照最陡坡度原则设定所述栅格单元的水流流向，所述最陡坡度原则是指，任一栅格单元的坡度的最佳代表值为根据相邻8个单元计算得到的坡度中的最大数值，所述栅格单元的水流方向即为最陡坡度对应的方向，利用数字进行表示。其中，两个相邻栅格单元i、j之间的坡度为θij＝arctan|(hi-hj)/d|，hi、hj分别表示两个栅格单元i、j的高程值，d表示两个栅格单元中心之间的距离。

为唯一确定栅格单元的水流流向，对每个方向赋以代码表示，例如，分别采用数字1、2、4、8、16、32、64、128表示东、东南、南、西南、西、西北、北、东北8个方向，若栅格单元的水流流向与其相邻的8个单元中的某个单元，则将水流方向设为该单元表示的方向代码，从而生成水流方向栅格数据。

s103、将水流方向栅格数据作为输入，利用hydrology下的flowaccumulation()函数，对任一栅格单元，确定流入该单元的所有累积上游单元格数目(nip)，从而生成流域汇流能力栅格数据，进一步根据设置的nip阈值确定河流网格，其中，大于nip阈值的单元格为沟谷线上的点，连接各个沟谷线上的点，形成河流网格。

s104、利用singleoutputmapalgebra命令,以flowaccumulation生成的流域汇流能力栅格数据作为输入，输出单值化后的栅格数据。具体地，在mapalgebraexpression中填入condition语句，语句格式为con(<condition>，<trueexpression>，<falseexpression>)，其中<falseexpression>为可选项，可选项为空时,若<condition>语句返回值为假，则单元格将被默认赋值为null。

s105、利用streamtofeature命令，以singleoutputmapalgebra生成的栅格数据作为“inputstreamraster”的输入，以水流方向栅格数据作为“inputflowdirectionraster”的输入，输出矢量化后的河网形文件(即shape文件)，完成河网提取。

s2、干支流划分及编码：利用所述河网形文件，根据河网区每一河流所控制的上游集水面积之和判定河流的干支流关系，从而进行支流划分得到n级支流，并按照从上游至下游、从上级支流至下级支流的顺序依次对n级支流进行编码。其中，集水面积总和最大的河流作为干流。

步骤s2的具体过程为：

s201、根据河网中的最大积流量网格，确定干流出口，并由此向上游检索，利用各网格的集水面积之和判定其与相邻网格的干支流关系，其中，集水面积总和最大的网格为干流所在网格，由此划分得到干流，并将其编码为m。

s202、划分支流级别l，对干流上游输入口的网格进行判断，并将其作为1级支流出口，据此向上游检索，根据其上游各个网格的集水面积之和判定其与相邻网格的干支流关系，集水面积总和最大的网格为1级支流所在网格，赋值l＝1，以此划分得到1级支流；同样地，对所述1级支流上游输入口的网格进行判断，并将其作为2级支流出口，向上检索，根据其上游各个网格的集水面积之和判定其与相邻网格的干支流关系，集水面积总和最大的网格为2级支流所在网格，以此划分得到2级支流，赋值l＝2；以此类推，以此得到1至n级支流。

s203、按照从上游至下游、从上级支流至下机支流的原则进行编码，其中，根据支流级别l，对于同一级别的支流，按照上游至下游的顺序递增编码。

请参考图2，本实施例中1级支流的编码前缀为t1，若共有k条1级支流，则按照从上游至下游的顺序，分别命名为t11、t12、……t1k；2级支流的编码前缀为t2，若共有j条2级支流，从t11支流开始编码，从上游至下游，第一条汇入t11的2级支流命名为t21，第二条汇入t11的2级支流命名为t22，直至无2级支流汇入t11，接着继续顺着1级支流的编码t12、……t1k向下游检索，直至最后一条汇入1级支流的2级支流编码为t2j；以此类推，n级支流的编码前缀为tn，依次得到1-n级支流的编码。

s3、河段节点划分及编码：首先对干流的起点和终点、支流的起点进行编码，然后确定支流入汇点并将其作为河段划分节点，根据空间位置关系，按照从上游至下游、从干流至支流、从上级支流至下级支流的顺序依次对所述河段划分节点进行编码。

请参考图3，将干流起点编码为ms，终点编码为me；接着对支流的起点进行编码，在支流编码后加上后缀s，例如，将1级支流t11的起点编码为t11s，将2级支流t21的起点编码为t21s，以此类推，完成对支流起点的编码；然后确定干流与支流的交点、上下级支流之间的交点，也即支流入汇点，将其作为河段划分节点，并保存为点状要素矢量文件；利用所述点状要素矢量文件，从干流开始检索，按照从上游至下游、从干流至1级支流、2级支流、……n级支流的顺序依次对河段划分节点进行编码。

干流与支流的交点的编码前缀为mj，将干流与1级支流的交点从上游至下游按照1级支流的顺序依次编码为mj1、mj2、……mjk；1级支流与2级支流的交点的编码前缀为t1j，从第一条1级支流t11开始编码，将第一条汇入t11的2级支流t21与t11的交点编码为t1j1，将第二条汇入t11的2级支流t22与t11的交点编码为t1j2，直至没有2级支流汇入t11，接着继续顺着1级支流的编码t12、……t1k向下游检索，直至最后一条汇入1级支流的2级支流t2j与所述1级支流的交点编码为t1jj；n-1级支流与n级支流的交点的编码前缀为tn-1j，以此类推，完成对所有河段划分节点的编码。

s4、河段划分及编码：依次检索编码后的河段划分节点，根据河段划分节点自动对河网区河流进行划分形成河段，并按照河段编码原则进行编码。

优选地，调用arctoolbox工具箱中的datamanagementtools→features→splitlineatpoint，在inputfeatures中输入要处理的河流文件，在pointfeatures中输入河段划分节点文件，设置输出路径，进行处理，自动在河段划分节点处打断河流，完成河段划分。

所述河段编码原则为：按照从上游至下游、从干流至支流的原则对河段依次编码，具体地，从干流起点开始，若检索到河段划分节点，则对所述河段划分节点与干流起点之间的河段进行编码；然后由所述河段划分节点向对应的支流上游进行检索，从支流上游起点开始对河段进行编码，若检索到河段划分节点，则对之间的河段采用递增形式进行编码，若未检索到河段划分节点，则直接对该河段进行编码；重复上述过程，直至无支流存在。

请参考图4，对河段进行编码的具体过程为：

(1)从干流起点ms开始，沿着干流的水流方向进行检索，找到干流上第一个河段划分节点mj1，将两点之间的河段编码为r1，并设定节点mj1为该河段的终止节点t11e；

(2)以与干流相交的交点所在的1级支流(t11)为分段对象，从支流t11的起点t11s开始，沿着该支流的水流方向进行检索，若该支流上存在河段划分节点，找到第一个2级支流(t21)与所述支流(t11)的交点t1j1，将t11s至t1j1之间的河段编码为r2，并设定节点t1j1为2级支流(t21)的终点t21e；

(3)以与1级支流(t11)相交的第一个交点所在的2级支流(t21)为分段对象，从2级支流(t21)的起点t21s向下游检索，并判断是否存在河段划分节点，若存在，找到第一个3级支流(t31)与所述2级支流(t21)的交点t2j1，将2级支流(t21)的起点t21s与2、3级支流交点t2j1之间的河段编码为r3，并设定节点t2j1为3级支流(t31)的终点t31e；

(4)在3级支流(t31)上，从河段划分节点t2j1(即t31e)向3级支流上游进行检索，若不存在河段划分节点，直接将该支流编码为r4；若存在河段划分节点，将起点t31s与河段划分节点之间的河段编码为r4，并将该河段划分节点设置为t41e，继续向上游检索，重复上述过程，直至对最小级别(n级)支流(tn1)检索完成；

(5)返回上一级支流，按照从上游至下游、上级河流至下级河流的原则，对tn1e与tn-1j1之间的河段依次进行递增编码；

(6)在n-1级支流tn-1j上，从河段划分节点tn-1j1向下游检索，判断是否还存在河段划分节点，若存在则对两个河段划分节点之间的河段进行编码，再由新出现的河段划分节点向上游检索，重复步骤(4)；若不存在则重复步骤(5)，直至由mj1汇入干流上游的所有河段均编码完毕；

(7)由mj1向干流下游检索，找到河段划分节点mj2，确定两点之间的河段编码，设定mj2为2级支流的终止节点t12e，之后重复步骤(2)-(6)，由此将所有的干支流交汇节点对应的入汇河段编码完毕，直到出口端点me，结果如图5所示。

s5、子河段划分节点识别：根据河段编码顺序，按照子河段节点确定原则，确定子河段划分节点。所述子河段节点确定原则将河流起始点(n)、行政区划边界与河段交点(a)、大型水库(r)、闸泵站(p)、水质监测断面(m)、重要取水口(t)、控制单元(c)以及河流拐弯处(w)等位置作为子河段划分节点。

本实施例中，所述子河段节点确定原则为：河段交点为1级节点，闸泵站为2级节点，行政区划或控制单元与河段的交汇处、控制断面为3级节点，水质监测断面、重要取水口为4级节点，河道拐点为5级节点。

按照所述子河段节点确定原则，根据河段编码顺序，确定子河段划分节点，其中，所述子河段划分节点满足：子河段划分节点之间距离适当，以确保子河段长度合理；河流上的大型水库自成一段，其入口、出口处为子河段划分节点。

本实施例中，对于行政区划边界与河段交点(a)或控制单元与河段交点(c)处，首先利用python调用arctoolbox→datamanagermenttools→features→polygontoline，将行政区划边界面要素或控制单元面要素转化为线要素，再调用arctoolbox工具箱→analysistool→overlay→intersect，提取二者交点输出为点要素文件，即可得到行政区划边界与河段交点或控制单元与河段交点，作为子河段划分节点；

对于河流拐弯处(w)，综合考虑拐点角度、拐点连接的前向河流弧段和后向河流弧段的长度来判断拐点的显著性，确保河流拐点弯度较合理且前后拐点间的河流程度适宜，进行人工设置；

对于其他子河段划分节点，比如大型水库(r)、闸泵站(p)、水质监测断面(m)、重要取水口(t)等位置，直接作为子河段划分节点。

s6、子河段划分及编码：根据子河段划分节点自动对河流进行划分形成子河段，再根据子河段编码原则对形成的子河段进行编码，得到河网区的完整编码。

优选地，本实施例首先利用python调用arctoolbox工具箱中的datamanagementtools→features→splitlineatpoint，在inputfeatures中输入要处理的河流文件，在pointfeatures中输入子河段划分节点文件，设置输出路径进行处理，自动在子河段划分节点处打断河流，完成子河段的划分。

所述子河段编码原则为：按照河段编码从小到大、从上游到下游，依次对子河段划分节点之间的子河段采用递增形式进行编码。

对子河段进行编码的具体过程为：从干流河段r1的起始节点开始，向下游检索找到第一个子河段划分节点，将该子河段划分节点与起始节点之间的子河段编码为sr1；接着继续向下游检索，按照编码递增的形式将找到前后两个子河段划分节点之间的子河段依次编码为sr2、sr3、……，直至河段r1的终止节点，完成河段r1上的子河段编码；转向河段r2，从河段r2的起始节点开始，向下游检索，继续按照递增形式对找到的前后两个子河段划分节点之间的子河段进行编码，直至河段r2的终止节点，完成河段r2上的子河段编码，转向河段r3……以此类推，直至最后一个河段的终止节点，完成所有子河段的编码，若共有n个子河段，则最后一个子河段编码为srn。

s7、计算单元划分及编码：根据河网区一维模型的模拟需求，确定计算单元的大小并进行划分，再根据子河段空间位置关系，从上游至下游按照子河段的编码顺序以递增形式完成计算单元的编码。其中，所述计算单元是一维模型模拟水文水质的最小模拟单元。

具体地，从子河段sr1的起始节点开始，向下游检索找到第一个计算单元，将其编码为e1；继续向下游检索，按照编码递增的形式将依次找到的计算单元编码为e2、e3、……，直至子河段sr1的终止节点，完成子河段sr1上的计算单元编码；转向子河段sr2，从子河段sr2的起始节点开始，向下游检索，继续按照编码递增的形式对找到的计算单元进行编码，直至子河段sr2的终止节点，完成子河段sr2上的计算单元编码，转向子河段sr3……以此类推，直至最后一个子河段srn的终止节点，若共有n’个计算单元，则最后一个计算单元编码为en’，完成对计算单元的编码。

s8、确定计算单元属性：建立子河段属性表，包括是否为源头子河段、河长、河段编码、子河段起始点坐标(经纬度)、子河段起始点高程等属性，通过arcgis对每个属性进行处理。

具体地，采用arcgis的地理几何计算工具对计算单元的长度进行计算，从而对河流长度进行赋值，单位为km；采用featureverticestopoints工具，选择start或者end，分别得到计算起点和终点，并通过地理几何计算工具，得到每个计算单元的起点、终点坐标，从而对子河段起始点坐标进行赋值，其中，计算单元的起点、终点坐标代表其上下游坐标值；基于dem文件，将dem的高程值扩展赋值给计算单元节点文件，选择的工具为spatialanalysttools的extractvaluestopoints，转换得到的图层中，即得到计算单元起点、终点的高程值，其中，计算单元起点、终点的高程值代表了其上下游的高程值。

上面结合附图对本发明的实施例进行了描述，但是本发明并不局限于上述的具体实施方式，上述的具体实施方式仅仅是示意性的，而不是限制性的，本领域的普通技术人员在本发明的启示下，在不脱离本发明宗旨和权利要求所保护的范围情况下，还可做出很多形式，这些均属于本发明的保护之内。