黄萌萌,李 灿,郭凤娇
(1.重庆市勘察规划设计有限公司,重庆 401121;2.重庆欣荣土地房屋勘测技术研究所有限责任公司,重庆 401121;3.重庆市巴南区规划和自然资源局,重庆 401000)
地形中沿山脊走向布设的路线即为山脊线。山脊线与山谷线作为地形骨架线,共同决定着地形的高低起伏。在地形表达中,山脊线和山谷线对其他地理要素起着控制性的作用。因此,山脊线常被作为专题进行研究,山脊线的提取和分级方法则是地形综合、制图表达等领域的重要研究内容。目前,根据基础数据,山脊线的提取方法可分为基于等高线[1]的方法和基于DEM分析的方法。基于DEM分析的方法主要包括地表径流模型(反地形)[2-3]和基于流域分析[4]的提取方法,前者已有较成熟的算法,并集成于GIS软件中,可快速自动化提取山脊线,但存在无法精确提取平坦地区山脊线,且出现平行效应的问题;后者是由郭万钦[4]等提出的一种利用流域边界和坡向差提取山脊线的方法,能较精准地提取山脊线,有效解决平坦地区的平行效应,但该方法需要设定汇流累积量、坡向差、坡向距离等6个阈值参数,人工干预性较大、自动化程度较低。针对上述问题,本文提出了一种基于流域边界与河网空间关系分析的山脊线提取方法,可有效避免平坦地区的平行效应问题,解决了阈值参数过多的问题,提高了山脊线提取的自动化效率。
1.1.1 流域及其组成
流域是指流经其中的水流或其他物质从不同的方向流过公共的出水口,这些水流或其他物质流经的由分水线所包围而形成的闭合集水区域[5-6]。从尺度上看,不同的水系等级都有自己的流域,流域本身又可按照水系级别细分为多个子流域。相邻流域之间公用的边界为分水岭,也是集水区域的边界线。水流或其他物质流经出水口的路径形成河网,从水流方向定向的特点来看,河网的形态为树结构,树的根节点是出水口,即为该流域内的最低点[7-9]。流域的组成如图1所示。
图1 流域的组成示意图
1.1.2 流域边界搜索
流域分割是以出水点(该集水区域内的最低点)为基础,结合水流方向矩阵,逆流而上搜索该出水点上游所有流经该点的栅格单元,直到该集水区域内所有栅格单元都确定了位置,即搜索到指定出水口点的流域边界[10]。为了确保山脊线网络与河网的匹配性,可采用河流链接数据(Stream Link)作为出水口数据,其中隐含着每条河网弧段的联结信息,包括河段的起点和终点,河段的终点即为该汇水区域出水口所在的位置[11-12]。同时,河网提取采用的阈值也决定了该河段所处汇水区域的详略程度,即提取河网的阈值越大,子流域面积也越大。
如图2所示,通过D8算法计算无洼地DEM每个栅格单元的水流方向,并计算流入每个栅格单元的汇流总量,获得汇流累积量矩阵;图2d为将栅格单元(4,4)指定为出水口点,逆向搜索所有流经该点的栅格单元(图中的蓝色单元),所有蓝色单元形成了出水口(4,4)的流域边界。
1.2.1 流域边界到山脊线的转换
流域边界除了分水岭界线外,还有一部分是位于负地形的边界线。因此,可利用流域边界线与山脊线之间的包含关系,剔除位于平坦河谷地带等负地形的流域边界,保留的则是需要提取的山脊线网络。从流域边界与河网的空间位置关系来看,位于负地形的流域边界与河网相交,可通过以下思路粗略提取山脊线:首先将栅格的流域边界转换为矢量数据,并由面图层转换为线图层;然后合并所有的线要素,并在节点(线与线的交点)处打断;最后通过河网空间位置查询,选取位于负地形的流域边界线要素,删除该部分要素,形成粗略的山脊线网络。
1.2.2 山脊线精细化处理
1)闭合环检查。受DEM分辨率、流域提取算法或其他不明原因的影响,删除负地形的流域边界后,山脊线网络内部尚存在闭合环。为提高山脊线的使用效率,确保山脊线网络的拓扑正确性,应对网络内部的闭合环进行破环处理。对删除位于负地形的边界后剩下的流域边界(简称图层L)进行要素转面操作,可将所有的闭合图形转为面要素,未形成闭合图形的线要素则不能转面成功,形成面图层P。面图层P中的每个面要素则是待处理的环,对每个面要素的RID进行编号,RID是每个环的唯一标识码。
2)闭合环分析。由于面图层P与图层L中的部分要素存在共线的关系,将面图层P的唯一编号RID通过空间链接赋值到图层L,空间位置关系为:SHARE_A_LINE_SEGMENT_WITH。若链接要素中具有与目标要素共线的要素,则匹配这些要素,赋值到图层L的RID字段,并将空间链接成功要素的BZ字段赋值为环要素,其余未链接成功的赋值为非环要素,RID为0。图层L的主要字段如表1所示。
表1 图层L的主要字段
根据环要素与非环要素的空间位置关系,可将环分为两种类型:①多要素环,由多条要素组成,每个节点由环与非环要素相连接,如图3a所示,多边形(环)ABCDE由边AB、BC、CD、DE、EA组成,每个 节点均是环与非环要素的交点;②单要素环,由一条边形成,单环与非环要素只有一个交点,如图3c所示,多边形与非环要素相交于V0点,V0既是边界的起点,也是终点。
图3 多要素环和单要素环破环前后示意图
3)破环方法。通过对图层L的RID字段进行频数分析,可识别位于同一环的各线要素,RID相同的要素为同一环的边界,RID的频数可区分环的类型,即
为保持破环后,山脊线网络整体的连通性、完整性和拓扑正确性,可分别对单要素环和多要素环进行以下处理:①针对多要素环由多条线要素构成的几何特征,可对位于边界上的各线要素进行对比分析,选取长度最长的要素(DE),并删除;②针对单要素环由一条线要素组成的几何特征,线要素的起点和终点重合,可从起点开始,计算每个节点至起点的欧氏距离,标记距离最长的节点P,删除该节点的下一个节点至终点的所有节点,点P为终点。
实验区位于108.051°~109.493 E、35.635°~ 37.206 N之间,海拔高度为780~1 800 m。实验数据来源于ASTER GDEM 30 m分辨率的DEM数据。本文提出的数据处理方法基于ArcGIS DeskTop10.2提供的工具箱以及ArcGIS Engine10.2,采用C#语言在Microsoft Visual Studio 2012开发平台下实现。
2.2.1 河网的提取
首先计算实验区的水流方向,利用ArcGIS Toolbox里的水文分析工具,双击Flow Direction工具,选择实验区域的DEM数据,生成水流方向图层FlowDir;然后采用Sink工具对FlowDir进行洼地检测,若存在洼地,则利用Fill工具对原始DEM进行填洼处理,若无洼地,则采用Flow Accumulation计算汇流累积量数据,生成FlowAcc数据。本文通过不断试验和利用现有地形图等其他数据辅助检验的方法来确定能满足研究需要且符合研究区域地形地貌条件的合适的阈值。为满足不同尺度的河网和山脊线数据要求,选取了3组阈值(阈值T分别设置为500、200、100)进行实验,结果如图4所示。通过Raster Calculator工具,将汇流累计量大于等于阈值的属性值设定为1,这些栅格就是河网的潜在位置;将小于阈值的属性值设定为No Data,形成StreamNet栅格数据。采用Stream to Feature完成对河网的栅矢转换,形成Stream矢量图层。
图4 不同阈值的河网提取结果
2.2.2 流域边界的提取
首先双击Hydrology工具集中的Stream Link工具,在Input Stream Raster文本框中选择StreamNet,在Input Flow Direction Raster文本框中选择fdirfill,输出数据为Stream Link;然后双击Watershed工具,打开集水区域(贡献区域)计算的对话框,分别在水流方向数据和出水口数据的输入文本框中选择fdirfill和Stream Link数据,输出数据为Watershed,并利用Raster to Polygon工具将Watershed转换为矢量图层。由图5可知,流域的边界能很好地概括山脊线的走向,并与相同阈值的河网数据形成良好的衔接关系。
图5 基于不同阈值河网提取的流域边界
2.2.3 流域边界与河网空间分析
首先选择阈值相匹配的河网S t r e a m图层和Watershed图层;然后通过Feature to Polyline将Watershed转换为线图层;最后将目标图层设置为Watershed,源要素设置为Stream,空间位置关系选择与源要素相交,通过空间位置查询选取Watershed中位于负地形的流域边界,并将其删除,剩下的流域边界则是粗略的山脊线Ridge_R。受DEM分辨率、流域边界提取算法等因素的影响,提取的山脊线内部存在部分小环,如图6所示。
图6 山脊线内部存在的闭合环
为了获取科学合理、实用性强的山脊线,需对粗略提取的山脊线网络进行精细化处理。首先采用Feature to Polygon工具将Ridge_R图层转为面图层,确定Ridge_R内部的闭合环,并对所有环要素的RID字段编码,可直接赋值为要素的FID+1;然后利用空间链接将RID链接回图层Ridge_R,通过频数分析确定闭合环的类型,并利用前文叙述的环处理思路,通过代码对山脊线内部环进行破环处理。由图7可知,与传统的基于反地形D8算法的提取方法相比,本文方法的实验结果能有效地与地形相吻合,且在山顶平坦位置或洼地位置避免了因填洼出现的平行效应;同时,在参数设置和自动化提取方面,较参考文献[4]也有一定的简化。
图7 D8算法与本文方法提取山脊线的对比图
本文利用基于流域边界分析的山脊线提取方法,采用3组阈值分析了提取的山脊线与相同密度的山谷线的衔接情况,如图8所示,可以看出,山脊线与相同密度的山谷线的衔接较好,能构建出不同层次的地形网络骨架,为有效划分地形层次结构提供了实践基础。
图8 本文算法提取的山脊线与相同阈值的山谷线衔接情况
针对山脊线提取存在的平行效应、自动化程度较低的问题,本文提出了基于流域边界分析的山脊线提取方法。该方法能够得到连续的山脊线,且与DEM可见的山脊线走向保持一致,避免了脊线段断续和平坦地段出现平行效应等问题;同阈值山脊线与河网能相互有效衔接,形成可概括地形的骨架线网络;不同汇流阈值的河网可获得不同面积大小的流域边界,从而形成详略程度不一的山脊线网络,为地形综合研究提供了新的思路。目前还没有科学有效的山脊线提取评价机制,山脊线的提取还应配合行之有效的分级方法,才能更好地为地形综合和水文分析提供理论基础。今后笔者将针对山脊线的分级做进一步研究。