山脊线与山谷线提取实验

作为地形特征线的山脊线、山谷线对地形地貌具有一定的控制作用。它们与山顶点、谷底点以及鞍部等一起构成了地形起伏变化的骨架结构。同时由于山脊线具有分水性，山谷线具有合水特性，使得它们在地形分析中具有特殊的意义。

基于规则格网DEM是最主要的自动提取山脊线和山谷线的方法，从算法设计原理上来分，大致可以分为以下5种：基于图像处理技术的原理；基于地形表面几何形态分析的原理；基于地形表面流水物理模拟分析的原理；基于地形表面几何形态分析和流水物理模拟分析相结合的原理；平面曲率与坡形组合法。

本案例中以地形表面流水物理模拟分析方法来分别演示山脊线和山谷线的提取过程。使用的数据为上杭县域DEM数据（xianyudem）。

具体提取过程如下：

▷　步骤1：正负地形的提取。

首先，在ArcMap中加载上杭县域DEM数据xianyudem.grid。

然后，点击ArcToolbox中的“Spatial Analyst工具”—“邻域分析”—“焦点统计”，弹出“焦点统计”对话框（图3-21）。在“输入栅格”中输入xianyudem；在“输出栅格”中定义输出文件的路径（shiyan03文件夹下）和名称（tongji）；在“邻域分析”中选择“矩形”；“邻域设置”中高度为11，宽度为11，单位为像元；“统计类型”中选择“MEAN”。点击“确定”按钮，进行均值邻域统计。

其次，使用“栅格计算器”（图3-22），定义“地图代数表达式”为“xianyudem”-“tongji”；定义“输出栅格”的路径（shiyan03文件夹下）与文件名称（zhengfudixing）。点击“确定”按钮，执行栅格计算，得到规划研究区的正负地形的分布区域。

最后，使用“重分类”工具（图3-23）获取正负地形。定义“输入栅格”为zhengfudixing；点击“分类”按钮，弹出“分类”对话框，设置分类类别数为2，并定义中断值为0，点击“确定”按钮返回“重分类”对话框。重新定义“新值”，将大于0的区域赋值为1，小于0的区域赋值为0；在“输出栅格”中定义输出文件的路径（shiyan03文件夹下）和文件名称（zhengdixing）。同样方法，将大于0的区域赋值为0，小于0的区域赋值为1；在“输出栅格”中定义输出文件的路径（shiyan03文件夹下）和文件名称（fudixing）。

图3-21　“焦点统计”对话框

图3-22　“栅格计算器”对话框

图3-23　“重分类”对话框（正地形提取）

▷　步骤2：DEM数据的填洼、流向与流量分析。

首先，使用“Spatial Analyst工具”—“水文分析”—“填洼”工具进行原始DEM数据的洼地填充。在“填洼”对话框中（图3-24），定义“输入表面栅格数据”文件为xianyudem；并在“输出表面栅格”中定义输出文件的路径（shiyan03文件夹下）和名称（filldem）。“Z限制（可选）”用以定义要填充的凹陷点与其倾泻点之间的最大高程差。如果凹陷点与其倾泻点之间的Z值差大于Z限制，则不会填充此凹陷点。默认情况下将填充所有凹陷点（不考虑深度）。本例中采用默认值，即将规划研究区内的所有洼地都进行填充。点击“确定”按钮，执行填洼命令。

然后，使用“Spatial Analyst工具”—“水文分析”—“流向”工具进行无洼地DEM的水流方向计算。在“流向”对话框中（图3-25），定义“输入表面栅格数据”文件为filldem；并在“输出流向栅格数据”中定义输出文件的路径（shiyan03文件夹下）和名称（flowdir）。“输出下降率栅格数据（可选）”用于输出下降率栅格数据。下降率栅格用于显示从沿流向的各像元到像元中心间的路径长度的最大高程变化率，以百分比表示。本例采用默认设置（不选）。点击“确定”按钮，执行流向计算命令。

图3-24　“填洼”对话框

图3-25　“流向”对话框

最后，使用“Spatial Analyst工具”—“水文分析”—“流量”工具进行汇流累积量的计算。在“流量”对话框中（图3-26），定义“输入流向栅格数据”文件为flowdir；并在“输出蓄积栅格数据”中定义输出文件的路径（shiyan03文件夹下）和名称（flowacc）。“输入权重栅格数据（可选）”用于对每一像元应用权重的可选输入栅格。配权数据一般是表示降水、土壤以及植被等造成径流分布不平衡的因子。如果未指定权重栅格，则将默认的权重值1应用于每个像元，计算出来的汇流累积量的数值就代表着该栅格位置流入的栅格数的多少。“输出数据类型（可选）”中提供了FLOAT浮点型（默认设置）或INTEGER整型。本例均采用默认设置。点击“确定”按钮，执行流量计算命令。(www.xing528.com)

图3-26　“流量”对话框

▷　步骤3：山脊线的提取。

首先，使用“栅格计算器”提取汇流累积量为0值的区域，通过定义地图代数表达式为"flowacc"＝＝0来获取。定义输出栅格的路径（shiyan03文件夹下）与文件名称（flowacc0）。

在ArcMap中打开flowacc0文件，可以发现很多地方不是山脊线，用户还需要对此数据进行进一步的处理。

然后，使用“Spatial Analyst工具”—“邻域分析”—“焦点统计”工具，对flowacc0数据进行3×3的邻域分析，求均值，使数据变得光滑，处理后的数据文件名称为linyuflowacc0。

其次，使用“Spatial Analyst工具”—“表面分析”—“等值线”工具和“山体阴影”工具，分别生成原始DEM数据的等值线图（dengzhixian）和山体阴影晕渲图（hillshade）。

再次，在“图层属性”对话框中，对linyuflowacc0进行重新分级，将数据分为两级，这时需要不断的调整分级临界点（阈值），并以等值线图（dengzhixian）和山体阴影晕渲图（hillshade）作为辅助判断，最终确定的分界阈值为0.55。

然后，按照该阈值将linyuflowacc0进行重分类，大于0.55的重新赋值为1，其余的赋值为0，并定义文件名称为reclassflow。

最后，使用“栅格计算器”将文件reclassflow和正地形数据zhengdixing相乘，得到shanjixian数据文件，这样就消除了那些存在于负地形区域中的错误山脊线。然后，再将shanjixian数据进行重分类，将所有属性值不为1的栅格重新赋值为NoData，得到规划研究区的山脊线reclassshjx（图3-27、图3-28）。

图3-27　“重分类”对话框

图3-28　计算得出的山脊线（黑色区域）

▷　步骤4：山谷线的提取。

首先，使用“栅格计算器”获取规划研究区的反地形。

在“栅格计算器”对话框中，定义地图代数表达式为Abs（1769-“xianyudem”）（1769为研究区最高海拔值）；定义输出栅格的路径（shiyan03文件夹下）与文件名称（fandixing）。

点击“确定”按钮，执行栅格计算命令，得到反地形DEM数据。

其次，按照前面山脊线提取的步骤，直到得到最终利用重分类方法将重新分级的邻域分析后的结果二值化为止（通过不断调整，最终分级阈值为0.65）。这里不需要对反地形进行DEM的填洼分析。分别定义得到的文件名称为fanflowdir，fanflowacc，fanflowacc0，linyufanacc0，reclassfanf1。

然后，使用“栅格计算器”将文件reclassfanf1和负地形数据fudixing相乘，得到shanguxian数据文件，这样就消除了那些存在于正地形区域中的错误山谷线。

最后，再将shanguxian数据进行重分类，将所有属性值不为1的栅格重新赋值为NoData，得到规划研究区的山谷线reclassshgx（图3-29）。

图3-29　计算得出的山谷线（黑色区域）