本发明属于水土保持和国土资源监测技术，具体涉及一种基于高分辨率遥感影像和dem的黄土浅沟自动化提取方法。
背景技术
浅沟是黄土高原丘陵沟壑区一种较为常见的微地貌现象。在降雨发生时，由于地形起伏引起坡面径流的汇集而产生股流，随着径流继续汇集，其侵蚀能力加大，沟道侵蚀加重，沟深和沟宽都逐渐扩展发育，形成了细沟、浅沟、切沟、冲沟、干沟及河沟等沟壑纵横的独特地貌形态。黄土浅沟作为细沟发展为切沟的过渡类型，在沟道侵蚀系统中具有特殊意义，其潜在危害极大。它是造成黄土高原水土流失严重的主要原因之一，其分布面积可占到沟间地的35％左右，侵蚀量占坡面侵蚀量的35％～75％。
野外调查和遥感影像目视解译是最传统和最常用的浅沟提取方法。由于黄土地形地貌的复杂性，实地测量获得大范围的浅沟比较困难。用遥感影像目视解译浅沟的特征信息，数据的采集操作复杂且低效，采集的数据精度受人为因素影响大。因此，亟需提出一种准确、快捷、自动的获取浅沟提取方法。
再加上由于浅沟较细长，宽度在0.5米到1米之间，在高分辨率的影像中也只占1-2个像元，属于遥感影像中的微弱信息。因此，准确、快捷、自动化的提取浅沟在技术方面也一直是个难题。
有鉴于以上问题，有必要提出一种黄土浅沟的自动化提取方法。为水土保持监测和国土资源调查等有关部门提供基础数据和决策支持。
技术实现要素：
发明目的：本发明的目的在于解决现有技术中存在的不足，提供一种基于高分辨率遥感影像和dem的黄土浅沟自动化提取方法，为水土保持监测和国土资源调查等有关部门提供基础数据和决策支持。
技术方案：本发明的一种基于高分辨率遥感影像和dem的黄土浅沟自动化提取方法，依次包括以下步骤：
(1)获取研究区域的高分辨率遥感影像数据，并利用已有的底图资料裁剪影像，得到浅沟所在的坡面区域图；然后对其进行边缘检测，边缘检测后得到一个像元值由0和1组成的边缘图像，其中1代表检测到的边缘；
(2)对浅沟坡面的数字高程模型dem进行流向判断得到流向栅格，并对流向进行重新编码使其像元值在1-8之间，分别代表右、右下、下、左下、左、左上、上和右上8个流向方向；其中右是指东方、右下是指东南方向、下是指南方、左下是指西南方向、左是指西方、左上是指西北方向、上是指北方和右上是指东北方向；
(3)对步骤(1)中得到的边缘图像进行定向检测，遍历整个边缘图像矩阵，当目标像元值为0时，跳过该像元；当目标像元值为1时，根据步骤(2)流向栅格中的方向，选择对应的定向检测方式生成候选浅沟栅格；
(4)沿相应定向检测方向对步骤(3)中得到的相邻的候选浅沟像元进行连接和唯一值标识，并记录每一条连接栅格的长度，长度的计算公式为：
式中，l为连接栅格的长度，n为该线段中水平相连或竖直相连的像元的个数；m为该线段中对角相连的像元的个数；d为像元大小；
(5)根据黄土浅沟是沿流向方向较长的线性要素的特征，以连接栅格长度的平均值为阈值，过滤掉检测中残存的伪浅沟信息(即低于阈值长度的浅沟均视为伪浅沟)，仅保留长度超过该阈值的连接栅格；另外，此处的分割阈值也可以使用实地调查的一般浅沟长度，或者从遥感影像中人工解译得出的一般浅沟长度；
(6)栅格转线，将步骤(5)中得到的连接栅格转换成矢量线，即得到黄土浅沟线。
进一步的，所述步骤(1)中用于选边缘检测的边缘检测算子为canny算子。
进一步的，所述步骤(2)采用d8流向算法来判断浅沟所在坡面的流向。
进一步的，所述步骤(3)中对边缘图像中值为1的目标像元进行定向检测的具体方法为：
若该点的流向为1(向栅格的像元值为1时)，依次检测与该目标像元邻接的右方、右上方、右下方像元；若右方像元值为1，则该右方像元为候选浅沟像元，并以该像元为新的目标像元开始新的定向检测；若右方像元为0，则检测右上方像元；若右上方像元值为1，则该右上方像元为候选浅沟像元，并以该像元为新的目标像元开始新的定向检测；若右上方像元为0，则检测右下方像元；若右下方像元值为1，则该右下方像元为候选浅沟像元，并以该像元为新的目标像元开始新的定向检测；若右下方像元为0，则停止该检测；
若流向为其他值时，检测方法与流向为1时相同，仅优先判断顺序不同：值为2时优先判断顺序为右下方、右方、下方；值为3时优先判断顺序为下方、右下方、左下方；值为4时优先判断顺序为左下方、下方、左方；值为5时优先判断顺序为左方、左下方、左上方；值为6时优先判断顺序为左上方、左方、上方；值为7时优先判断顺序为上方、左上方、右上方；值为8时优先判断顺序为右上方、上方、右方。
有益效果：首先，通过对浅沟坡面的高分辨率影像进行边缘检测得到浅沟边缘图像；其次，根据dem对浅沟坡面的流向就行判断，并进行重新编码，生成流向栅格；再次，根据流向对边缘图像进行定向检测，从边缘图像的每个边缘像元出发，沿流向检测下一个像元是否为边缘像元，若是，则连接并重复上一步的检查，若不是，则检测终止；最后，对上一步结果使用适当的阈值剔除结果中的伪浅沟信息和栅格转线，即得到提取区域的黄土浅沟。本发明提供了自动化的黄土浅沟提取方法，节省了野外调查的人力物力，能为水土保持监测和国土资源调查等提供更丰富、便利的基础数据。
附图说明
图1是本发明提供的黄土浅沟提取方法的流程图；
图2是实施例样区高分辨率遥感影像；
图3是实施例样区边缘检测结果图像；
图4是实施例样区dem；
图5是本发明提供的流向重编码示意图；
图6是实施例样区的流向栅格；
图7是本发明提供的定向检测算法示意图；
图8是本发明提供的黄土浅沟提取结果图。
具体实施方式
下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。
如图1所示，本发明的一种基于高分辨率遥感影像和dem的黄土浅沟自动化提取方法，依次包括以下步骤：
(1)获取研究区高分辨率遥感影像数据，并利用已有的底图资料裁剪影像，得到浅沟所在的坡面区域图；然后对其进行边缘检测，边缘检测后得到一个像元值由0和1组成的边缘图像，其中1代表检测到的边缘；
(2)对浅沟坡面的dem进行流向判断得到流向栅格，并对流向进行重新编码使其像元值在1-8之间，分别代表右、右下、下、左下、左、左上、上和右上8个流向方向；其中右是指东方、右下是指东南方向、下是指南方、左下是指西南方向、左是指西方、左上是指西北方向、上是指北方和右上是指东北方向；
(3)对步骤(1)中得到的边缘图像进行定向检测，遍历整个边缘图像矩阵，当目标像元值为0时，跳过该像元；当目标像元值为1时，根据步骤(2)流向栅格中的方向，选择对应的定向检测方式生成候选浅沟栅格；
(4)沿定向检测方向对步骤(3)中得到的相邻的候选浅沟像元进行连接和唯一值标识，并记录每一条连接栅格的长度，长度的计算公式为：
式中，l为连接栅格的长度，n为该线段中水平相连或竖直相连的像元的个数；m为该线段中对角相连的像元的个数；d为像元大小；
(5)根据黄土浅沟是沿流向方向较长的线性要素的特征，以连接栅格长度的平均值为阈值，过滤掉检测中残存的伪浅沟信息，仅保留长度超过该阈值的连接栅格；另外，此处的分割阈值也可以使用实地调查的一般浅沟长度，或者从遥感影像中人工解译得出的一般浅沟长度；
(6)栅格转线，将步骤(5)中得到的连接栅格转换成矢量线，即得到黄土浅沟线。
实施例：
针对如图2所示的样区高分辨率遥感影像，采用以下步骤：
步骤一，首先获取高分辨率影像，并利用已有的底图资料进行裁剪，得到浅沟所在坡面，然后对其进行边缘检测获取浅沟的边缘图像；本例中使用的高分辨率遥感影像为0.5m分辨率的googleearth影像，裁剪后的浅沟坡面影像如图2所示；使用canny算子进行边缘检测，得到边缘图像如图3所示。
步骤二，使用已有的d8流向算法对浅沟所在坡面的dem进行流向判断，并对流向进行重新编码使其像元值在1-8之间的，分别代表右(东)，右下(东南)，下(南)，左下(西南)，左(西)，左上(西北)，上(北)，右上(东北)等8个流向方向；样区dem如图4所示，新的流向编码如图5所示，生成的流向栅格如图6所示。
步骤三，对步骤一中得到的边缘图像进行定向检测，遍历整个边缘图像矩阵，当目标像元值为0时，跳过该像元；当目标像元值为1时，根据流向栅格中的方向，选择对应的检测方式生成候选浅沟栅格；
如图7所示，具体检测算法为：
若该点的流向为1，依次检测与该目标像元邻接的右方、右上方、右下方像元，若右方像元值为1，则该右方像元为候选浅沟像元，并以该像元为新的目标像元开始新的定向检测；若右方像元为0，则检测右上方像元。若右上方像元值为1，则该右上方像元为候选浅沟像元，并以该像元为新的目标像元开始新的定向检测；若右上方像元为0，则检测右下方像元。若右下方像元值为1，则该右下方像元为候选浅沟像元，并以该像元为新的目标像元开始新的定向检测；若右下方像元为0，则停止该检测。
若流向为其他值时，检测方法与上述基本相同，仅优先判断顺序不同：值为2时优先判断顺序为右下方、右方、下方；值为3时优先判断顺序为下方、右下方、左下方；值为4时优先判断顺序为左下方、下方、左方；值为5时优先判断顺序为左方、左下方、左上方；值为6时优先判断顺序为左上方、左方、上方；值为7时优先判断顺序为上方、左上方、右上方；值为8时优先判断顺序为右上方、上方、右方；
步骤四，沿定向检测的方向对步骤三中得到的相邻的候选浅沟像元进行连接，并对每一条连接栅格进行唯一值标识和记录其的长度；其中长度计算公式为：
式中，l为连接栅格的长度，n为该线段中水平相连或竖直相连的像元的个数；m为该线段中对角相连的像元的个数；d为像元大小，此处的像元大小为0.5m。
步骤五，根据黄土浅沟是沿流向方向较长的线性要素的特征，此例以连接栅格长度的平均值13.4为阈值，仅保留长度大于该阈值13.4的连接栅格；
步骤六，栅格转线，将步骤(5)中得到的连接栅格转换成矢量线，即得到黄土浅沟线。本例中使用的栅格转线工具是arcgis10.2中提供的栅格转线工具；最终结果如图8所示，通过与原始影像叠加显示，可以发现其匹配精度较好，几乎所有的黄土浅沟均被成功的提取。
对于最终的黄土浅沟提取结果，采用完整性和冗余度进行定量评价，公式如下：
完整性＝正确提取的浅沟总长度/参考浅沟总长度
冗余度＝错认作浅沟的总长度/参考浅沟总长度
其中参考浅沟由人工目视解译得到。
定量评价结果如表1所示，本发明所提取的浅沟，完整性为93.28％，冗余度为6.50％；即本发明具有较高的精度和实用性。
表1实施例中浅沟提取结果精度评价