基于纹理分水岭的斜坡单元提取算法

程露 周波

摘 要：斜坡单元在以滑坡为主的地质灾害预防和评价中有着广泛的应用，其提取和划分是滑坡灾害风险评估的首要工作和重要基础。针对传统地理信息系统（GIS）方法提取的斜坡单元存在平行边界和误分割问题，提出了基于纹理分水岭的斜坡单元提取方法，通过分割地形图像划分斜坡单元。首先通过预处理地形数据得到数字高程模型（DEM）图像，利用灰度共生矩阵提取DEM纹理特征;然后计算融合灰度和纹理特征的梯度图像，对梯度图像进行基于标记的分水岭分割，使其能够准确获取山体和流域边界;最后，结合正负地形，对山体对象进行分水岭分割以实现斜坡单元的提取。实验结果表明，所提方法对不同地貌类型和分辨率的DEM图像都有良好的划分效果;相较于传统的GIS方法，该方法能够正确分割水平面和倾斜面，有效避免因洼地填平处理而产生的平行边界问题。

关键词：斜坡单元;数字高程模型;分水岭算法;正负地形;纹理梯度

中图分类号： TP391.41图像识别及其装置

文献标志码：A

Abstract： Slope units are widely used in the prevention and evaluation of landslide-based geological hazards， whose extraction and division are the primary target and important foundation for the risk assessment of landslide hazards. Considering the parallel boundaries and incorrect segmentation problems of the slope units extracted by traditional Geographic Information System （GIS） method， a slope unit extraction algorithm based on texture watershed was proposed， in which slope units were extracted by segmenting terrain images. Firstly， a Digital Elevation Model （DEM） image was obtained by the pretreatment of terrain data， and DEM texture features were extracted by gray level co-occurrence matrix. Then， the gradient image with gray level fused with texture features was calculated and segmented by marker-based watershed segmentation to accurately obtain mountain boundaries and watershed boundaries. Finally， combined with positive and negative terrains， the mountain objects were segmented by watershed segmentation to extract slope units. The experimental results show that the proposed method is pretty effective in segmentation for DEM images of different landform types and resolutions. Compared with traditional GIS method， horizontal planes and inclined planes can be segmented correctly， and the problem of parallel boundaries caused by filling of depressions can be effectively avoided through the proposed method.

Key words： slope unit; Digital Elevation Model （DEM）; watershed algorithm; positive and negative terrain; texture gradient

0 引言

數字高程模型（Digital Elevation Model， DEM）是空间数据基础设施的重要组成部分，蕴含着丰富的地形地貌信息，从DEM数据中自动提取其隐含的地形结构信息具有重要的研究意义[1]。斜坡单元是由山脊线和山谷线共同分割的图元，通常情况下，一个区域内相邻几个斜坡单元会有明显不同的地形特征，如地形线、流线、坡向、坡面等[2]。斜坡单元能够有效研究地形元素、分析单元间的联系、保留斜坡的整体性，贴近实际地貌，可以很好地体现研究区域中滑坡的实际发育状况。随着地质灾害风险评价研究的不断深入，斜坡单元在以滑坡为主的地质灾害评价中发挥着越来越重要的作用[3]，其提取和划分是开展大比例尺区域滑坡灾害易发性、危险性、风险评价的首要工作和重要基础[4]。

当前，斜坡单元的提取主要依赖于地理信息系统（Geographic Information System， GIS）的水系提取和子流域分割[5-7]，其基本原理是利用无洼地正负地形分别提取山脊线和山谷线（分别对应分水线和汇水线），把生成的集水流域与反向集水流域融合，再经后期人工修编不合理的单元，最终得到的由汇水线与分水线所划分的区域即为斜坡单元。其改进方法有：考虑栅格分辨率和集水面积阈值等因素进行改进的方法[8]、通过叠加曲率和反转曲率的流域边界进行改进的方法[2]、基于平均曲率进行改进的方法[9]等。

尽管基于水文分析的斜坡单元提取方法有着广泛的应用[1]，但依然存在着如下缺点：一是由于洼地填平对地形有一定的改变，极易形成大片平地，提取的斜坡单元往往会出现平行边界[10];二是由于一般不区分水平面和倾斜面，当地形中存在山间盆地或大型宽谷时，只能将水平面机械地并入到倾斜面中，不符合实际地貌[11]，从而导致提取的斜坡单元往往需要进行繁杂的人工修编工作，效率低下。

近年来，基于纹理分析的地形特征研究开始应用于地表地貌类型的划分，将纹理分析与DEM的数字地形特征相结合，被认为是划分宏观地貌类型的突破点[12]。针对传统GIS方法的不足，本文尝试从图像分割的角度划分斜坡单元，将图像分割方法中的分水岭算法应用到地形分割中，并结合地形图像的纹理特征，提出了一种基于纹理分水岭[13]的斜坡单元提取方法，使得在提取斜坡单元的同时，能够克服传统GIS方法产生平行边界、不能有效区分水平面和倾斜面的缺点。

1 本文算法的原理及应用背景

分水岭算法是一种基于拓扑理论和模拟地形学的图像分割技术，因其运算简单、计算速度快，可以得到单像素宽、连通封闭的轮廓等特点，在图像分割领域有着广泛的应用[14]。当前最为普遍使用的一种分水岭算法实现方式是模拟浸没过程，其思想是将一幅灰度图像看作是测地学上的拓扑地貌[15]，图像中每个点的灰度值对应为地形表面的高程，每一个局部极小值及其影响区域为集水盆，而集水盆的边界则形成分水岭，即地形中的山脊线，如图1所示。在每个局部极小值处，刺穿一个小孔，让水均匀上升，从低到高直至淹没地形，在不同集水盆的水汇合处修建一条大坝，这些大坝就是分水岭。

从分水岭算法的原理可以看出，待分割图像中每个局部极小值对应分割结果的一个子区域，子区域的个数由局部极小值的数量决定。在DEM数据中，除真实地形（溶洞、湖泊、洼地等）对应的局部极小值外，在数据采集或数据内插生成DEM时也会产生大量的局部极小值，这些伪极小值会造成严重的过分割[7]。为解决过分割问题，本文在标记分水岭算法的基础上进行改进以提取斜坡单元。

DEM数据的纹理特征属于自然纹理的一种，平原、丘陵和山脉具有不同的地貌形态和高度起伏，反映在影像上则呈现具有方向和粗糙度差异的纹理特征。由于灰度梯度在表达纹理信息方面的不足[16]，考虑将DEM的纹理信息融合进灰度梯度中，构建一个混合的包含灰度差异和纹理特征的梯度图像，结合分水岭分割，使其能够更准确定位山体和流域边界。

分水岭分割后得到的子区域为地形中单个的山体对象或流域对象，需对这些子区域再次分割以提取斜坡單元。由于此时每一块子区域的极小值数量相比原图像已大大减少，对其作分水岭分割可提取地形中的山脊线，借鉴GIS中正负地形思想，考虑将子区域反转可提取山谷线，最后通过山脊线和山谷线的叠加划分出斜坡单元。

综上，本文考虑DEM图像的纹理特征，利用基于标记的分水岭算法解决过分割问题，提取山体边界线，获取山体对象，再结合正负地形，对山体对象的正负地形作分水岭变换，分别提取山脊线和山谷线，进而叠加获取斜坡单元。本文算法流程如图2所示，具体步骤如下：

1）地形数据预处理：对原始DEM数据进行灰度域转换生成DEM图像，随后对该图像作非线性滤波运算，以平滑图像、去除噪声。

2）结合纹理的分水岭分割：分别计算预处理后图像的纹理梯度和灰度梯度，构建一个混合梯度图像;利用前背景标记技术修正该混合梯度图像，使得仅标记过的位置出现极小值，随后对修正后的梯度图像作分水岭变换，提取山体边界线，获得山体区域。

3）结合正负地形的分水岭分割：对山体区域作形态学开闭重建，随后对重建后的山体图作分水岭变换，提取山脊线，将山体区域反转得负地形提取山谷线，最后叠加山脊线和山谷线划分出斜坡单元。

4）合并每一块山体区域的斜坡单元，可以得到输入DEM数据的斜坡单元划分。

2 本文算法的具体实现

2.1 地形数据预处理

2.1.1 灰度域转换

DEM为连续地表高程值的采样序列，可用一般函数形式Vi=（Xi， Yi， Zi）表达，其中，Xi、Yi为平面坐标，Zi为坐标（Xi， Yi）对应的高程[17]。和数字灰度图像不同，DEM数据的栅格属性值代表的是高程值，表达的地形空间表面点上的高程值范围较宽，尤其对于国家基础地理信息中心发布的5m、25m等高分辨率DEM数据，其栅格点高程值均为浮点类型。为提高数据的运算效率，在进行灰度共生矩阵的计算前，需先对DEM高程值进行相应的灰度域映射。DEM数据映射后，设其栅格对应的最大灰度量化级为Mz（最大取值255），则DEM的灰度量化集G = {0，1，…，Mz}，由此可以将DEM数据模拟的二维数字地形图像A表示为：A=Sx*Sy→G。

2.1.2 非线性滤波

在分水岭分割前，对原始DEM图像作滤波处理以降低图像噪声、平滑目标内部纹理，从而去除图像中部分伪局部极小值，缓解过分割。本文对原始DEM图像A先后作最小值滤波、最大值滤波运算（式（1）），得到滤波后的图像W，从而消除A中部分高斯噪声、模糊内部纹理，缓解过分割，而且边界依然清晰，基本不丢失A中的重要信息。

2.2 梯度图像的构建

考虑到DEM图像山体边界处的梯度值较大，山体区域内部梯度值较小，若用梯度图像代替原始DEM图像作分水岭变换，则可以有效提取出山体边界，正确区分水平面和倾斜面。此外，地形纹理作为表征地面固有自然特征，是人们描述和区分不同地形特征的重要依据之一。DEM的纹理特征反映了由地表高程变化引起的地形形态差异，故利用 DEM 的纹理信息研究地形形态特征，是当前提取地形纹理的重要方法之一[18]。本文综合考虑DEM的纹理特征和灰度梯度在表达山体对象各自的优势，研究一种混合梯度用于分水岭分割，以更好定位山体和流域对象的边界。

2.2.1 灰度梯度

在灰度梯度运算中，由于Sobel算子计算简单、运算速度快，可以产生较好的边缘检测效果[19]，所以常被用来计算图像的一阶梯度之近似值。考虑到Sobel算子比其他算子具有更强的抗噪性[19]，本文选用Sobel算子计算滤波后的图像W的梯度幅值（式（2）），得到灰度梯度图像Gra。

2.2.2 纹理梯度

为了构造纹理梯度，首先需要研究纹理的描述方法，目前纹理分析大致可分为统计型纹理分析方法和结构型纹理分析方法。研究[17]表明，统计型分析方法更适用于自然纹理。为此，本文选用统计型纹理分析方法中的灰度共生矩阵（Gray-Level Co-occurrence Matrix， GLCM）表达纹理图像。基于GLCM的纹理可以用14个特征参数表示，由于这些参数存在一定的冗余，可以用较少的参量表示纹理。因此，为了降低运算量且不失准确性，本文选取二阶角矩、对比度和熵3个纹理特征参数（如表1所示），并采用文献[13]方法计算提取到的纹理特征的纹理梯度TG。

2.2.3 混合梯度

由于纹理梯度可以很好检测纹理边缘，灰度梯度可以有效表达高程差异信息，为了提高DEM图像的分割效果，更好地定位山体和流域的边界，考虑将灰度梯度结合到纹理梯度中。本文采用文献[14]中的线性叠加方法计算混合梯度Gmix：

Gmix=ωTG（x，y）+（1-ω）Gra（x，y）

其中：TG（x，y）表示纹理梯度;ω为纹理梯度所占的权重;Gra（x，y）表示灰度梯度;1-ω为灰度梯度所占的权重。本文纹理梯度所占的权重ω取值为0.5。

2.3 标记分水岭分割

为进一步缓解过分割情况，对上述混合梯度图像进行基于标记的分水岭分割。本文先后对2.1节中滤波后的DEM图像W作形态学开闭重建运算和前背景标记处理，修正混合梯度Gmix，最后对修正过的梯度图像作分水岭变换，以提取山体和流域的边界线。

2.3.1 形态学开闭重建运算

数学形态学具有定位效果好、分割精度高、抗噪性能好等特点，将数学形态学方法引入分水岭变换中，可以有效减少极小值，很好缓解过分割现象[20]。本文通过对滤波后的DEM图像W构建形态学混合开闭重建运算（式（3），式（4）），得到图像F，以达到同时消除图像中明暗细节和噪声的目的。

2.3.2 前背景标记

使用前背景标记法标记DEM混合梯度图像，使得仅标记过的位置出现极小值。本文通过寻找重建后图像F中局部极大值（即地形中的山顶点）获取前景标记fgm，通过自动化阈值将图像F二值化后作距离变换运算，再对结果作分水岭变换，分水岭的边界线即为背景标记bgm，利用前景標记fgm和背景标记bgm修正混合梯度图像Gmix，把bgm中无关的极小值屏蔽掉（式（5）），得到修正后的梯度图像Gmark。

Gmark=imimposemin（Gmix，bgm|fgm）

（5）

2.3.3 分水岭变换

对修正后的梯度图像Gmark进行分水岭变换（式（6）），得到单像素宽、封闭的分割线Edge1，即为地形中的山体边界线，从而可以获得单个山体对象图像I。值得注意的是，从分水岭算法的原理中可知，分水岭为地形上的山脊线，但此处对于DEM图像而言，分割线Edge1对应图像中的对象轮廓线，即为山体边界线。

Edge1=watershed（Gmark）

（6）

2.4 结合正负地形的分水岭分割

上述基于标记的分水岭分割为DEM图像的第一次分割，此时图中的每一块子区域为地形的山体或流域对象，并非斜坡单元，需对子区域再次分割以提取斜坡单元。对子区域而言，其伪极小值和噪声相对原图已大大减少，可以通过对其正负地形作分水岭分割，来获取山脊线和山谷线，再叠加提取出斜坡单元。本文通过对山体对象的灰度图像I作形态学开闭重建运算得I′（式（3）、式（4）），以去除图像噪声，缓解过分割。对重建后的图像I′直接作分水岭变换得到Edge2（式（7）），提取图中的山脊线;反转山体对象的DEM图像获取负地形，提取山谷线。最后，叠加山脊线和山谷线划分斜坡单元。

3 实验与结果分析

本文算法在软件Matlab R2017a上实现，为了验证纹理梯度的正确性及算法的有效性，本文对不同地貌类型和分辨率的DEM数据进行提取斜坡单元实验，并与传统GIS方法的提取结果进行比较。

3.1 灰度梯度与纹理梯度的比较

为验证传统灰度梯度在表达地形纹理特征的不足，本文选择某山区1∶50000比例尺、25m分辨率的DEM数据进行实验。使用本文预处理的方法，先后对原始DEM数据作灰度域映射和非线性滤波得到的局部尺度图如图3（a）所示，分别计算该图像的灰度梯度和纹理梯度，实验结果的局部放大图分别如图3（b）、（c）所示。从实验结果可以看出，与传统灰度梯度相比，纹理梯度进一步加强了纹理区域的边缘信息，模糊了山体内部区域的纹理，更有效表达了DEM的纹理特征。在纹理梯度的基础上，结合分水岭分割，就能够更加准确地定位山体和流域对象的边界。

图片

图3 灰度梯度与纹理梯度的对比

Fig. 3 Comparison of gray level gradient and texture gradient

3.2 本文算法的分割结果

图4（a）和图5（a）分别为1∶50000比例尺、25m分辨率的某山区和某平原地区DEM图，图6（a）为1∶10000比例尺、5m分辨率的某山区DEM图，对应的山体阴影图（可看作分割的参考图像）分别如图4（d）、图5（d）和图6（d）所示。为验证算法的有效性，利用本文算法分别对上述实验样区的DEM数据进行斜坡单元的提取，实验结果分别如图4（c）、图5（c）和图6（c）所示。从实验结果可以看出，本文算法有效可行，可以实现斜坡单元的提取，且坡面提取较为完整，分割边界定位较为准确，符合实际地形。同时，对于不同分辨率和地貌类型的DEM数据，本文算法都可以得到单像素宽且封闭的分割区域，分割边界清晰、定位准确。由此，本文算法分割性能较好，具有良好的鲁棒性。

3.3 与传统GIS方法比较

当前斜坡单元的提取算法主要集中在基于GIS的水文分析上，故将本文算法与传统GIS方法作对比，进而验证本文方法的有效性和优越性。为此，分别使用传统GIS方法和本文算法对图6（a）所示的山区DEM图像进行斜坡单元提取。利用ArcGIS软件的水文分析工具，通过无洼地正负地形子流域边界的划分，实现该山区的斜坡单元提取，实验所取正地形流量阈值为800，负地形流量阈值为500，实验结果局部放大區域如图7（a）所示，用本文算法提取斜坡单元的实验结果局部放大区域如图7（b）所示。

从以下两个角度对比实验结果：一是分割结果区域的单元形态，圆形、正方形和正三角形更贴合实际地貌，长条形则需要繁杂的人工修编;二是平均提取准确率，即提取到的完整坡面（区域面积达到坡面的70%以上）/分割区域（斜坡单元）数，平均提取准确率越高说明分割结果越接近实际地形的边界线。图7（a）中共生成516个分割区域，其中提取到较为完整坡面有23个，平均提取准确率为4.46%。图7（b）中共生成81个分割区域，其中提取到较为完整的坡面有59个，平均提取准确率为72.84%。从图7（a）可以看出，传统的GIS方法由于洼地填平对地形的改变，造成了大片平地，提取的斜坡单元出现了很多长条形的单元形态和平行的边界线;此外，由于不能有效区分水平面和倾斜面，山区的山体和流域部分都出现了明显的错误分割，坡面提取不完整，不符合实际地形。从图7（b）可以看出，分割结果区域的单元形态少有长条形，更贴合实际地貌，表明本文算法克服了传统方法的缺点，能够有效区分水平面和倾斜面，正确分割出了山体对象和流域对象，且避免了因洼地填平处理产生的平行边界，免去了繁杂的修编工作。相比于传统的GIS方法，本文提取的斜坡单元的边界线与实际地形的边界线更为相符，更加贴近真实地貌，且本文算法的平均提取准确率相比于传统GIS方法大大提高。

4 结语

本文结合DEM纹理信息和正负地形，提出了一种改进的基于标记分水岭的斜坡单元提取方法，其基本过程为：利用灰度共生矩阵求取纹理梯度，结合灰度梯度构建混合梯度图像，对此进行标记分水岭分割以准确进行山体和流域边界的定位，得到山体对象。随后结合正负地形，对山体对象进行分水岭分割获取斜坡单元。实验结果表明，本文算法能够实现斜坡单元的提取，对于不同地貌类型、不同分辨率的DEM数据都具有良好的有效性和适应性，提取结果的边界清晰、定位准确，比较贴合实际地形地貌。此外，本文算法能够克服传统GIS方法容易产生平行边界、不能正确区分水平面和倾斜面等缺点，避免后期繁杂的人工修编工作，大大提高提取效率。

由于形态学结构元素形状和大小、灰度梯度算子的选取等关键因素对划分结果的影响难以定量考量，实验过程中参数的选取缺乏定量计算分析，考虑这些参数对划分结果影响，建立较完整的斜坡单元定量划分算法将是下一步的研究工作。

参考文献 （References）

[1] 洪莹，王继周，李昂.地形特征提取的一种简易算法[J].测绘科学，2009，34（6）：127-129.（HONG Y， WANG J Z， LI A. A simple terrain feature extraction algorithm [J]. Science of Surveying and Mapping， 2009， 34（6）： 127-129.）

[2] 谢谟文，江崎哲郎，周国云.基于边坡单元的三维滑坡灾害评价的GIS方法[J].岩石力学与工程学报，2003，22（6）：969-976.（XIE M W， TETSURO E， ZHOU G Y. GIS method for slope-unit-based 3d landslide hazard evaluation [J]. Chinese Journal of Rock Mechanics and Engineering， 2003， 22（6）： 969-976.）

[3] 颜阁，梁收运，赵红亮.基于GIS的斜坡单元划分方法改进与实现[J].地理科学，2017，37（11）：1764-1770.（YAN G， LIANG S Y， ZHAO H L. An approach to improving slope unit division using GIS technique [J]. Scientia Geographica Sinica， 2017， 37（11）： 1764-1770.）

[4] 张曦，陈丽霞，徐勇，等.两种斜坡单元划分方法对滑坡灾害易发性评价的对比研究[J].安全与环境工程，2018，25（1）：12-17.（ZHANG X， CHEN L X， XU Y， et al. Comparison of two methods for slope unit division in landslide susceptibility evaluation [J]. Safety and Environmental Engineering， 2018， 25（1）： 12-17.）

[5] 谢谟文，江崎哲郎，邱骋，等.空间三维滑坡敏感性分区工具及其应用[J].地学前缘，2007，14（6）：73-84.（XIE M W， TETSURO E， QIU C， et al. A spatial three-dimensional landslide susceptibility mapping tool and its applications [J]. Earth Science Frontiers， 2007， 14（6）： 73-84.）

[6] ALVIOLI M， MARCHESINI I， REICHENBACH P， et al. Automatic delineation of geomorphological slope units with r.slopeunits v1.0 and their optimization for landslide susceptibility modeling [J]. Geoscientific Model Development， 2016， 9（11）： 3975-3991.

[7] SCHLGEL R， MARCHESINI I， ALVIOLI M， et al. Optimizing landslide susceptibility zonation： effects of DEM spatial resolution and slope unit delineation on logistic regression models [J]. Geomorphology， 2018， 301： 10-20.

[8] 黄启乐，陈伟，唐绪波，等.区域地质灾害评價中斜坡单元划分方法研究[J].自然灾害学报，2017，26（5）：157-164.（HUANG Q L， CHEN W， TANG X B， et al. Study on the method of slope unit zoning in regional geo-hazards risk assessment [J]. Journal of Natural Disasters， 2017， 26 （5）： 157-164.）

[9] ROMSTAD B， ETZELMüLLER B. Mean-curvature watersheds： a simple method for segmentation of a digital elevation model into terrain units [J]. Geomorphology， 2012， 139/140： 293-302.

[10] 叶爱中，夏军，王纲胜，等.基于数字高程模型的河网提取及子流域生成[J].水利学报，2005，36（5）：531-537.（YE A Z， XIA J， WANG G S， et al. Drainage network extraction and subcatchment delineation based on digital elevation model [J]. Journal of Hydraulic Engineering， 2005， 36（5）： 531-537.）

[11] GUZZETTI F， CARRARA A， CARDINALI M， et al. Landslide hazard evaluation： a review of current techniques and their application in a multi-scale study， central Italy [J]. Geomorphology， 1999， 31（1/2/3/4）： 181-216.

[12] 刘凯，汤国安，陶旸，等.基于灰度共生矩阵的DEM地形纹理特征量化研究[J].地球信息科学学报，2012，14（6）：751-760.（LIU K， TANG G A， TAO Y， et al. GLCM based quantitative analysis of terrain texture from DEMs [J]. Journal of Geo-Information Science， 2012， 14（6）： 751-760.）

[13] HILL P R， CANAGARAJAH C N， BULL D R. Image segmentation using a texture gradient based watershed transform [J]. IEEE Transactions on Image Processing， 2003， 12（12）： 1618-1633.

[14] 杨家红，刘杰，钟坚成，等.结合分水岭与自动种子区域生长的彩色图像分割算法[J].中国图象图形学报，2010，15（1）：63-68.（YANG J H， LIU J， ZHONG J C， et al. A color image segmentation algorithm by integrating watershed with automatic seeded region growing [J]. Journal of Image and Graphics， 2010， 15（1）： 63-68.）

[15] 李小红，武敬飞，张国富，等.结合分水岭和区域合并的彩色图像分割[J].电子测量与仪器学报，2013，27（3）：247-252.（LI X H， WU J F， ZHANG G F， et al. New color image segmentation based on watershed and region merging [J]. Journal of Electronic Measurement and Instrument， 2013， 27（3）： 247-252.）

[16] 沈晶，杨学志.一种新的基于纹理分水岭的纺织品缺陷检测方法[J].中国图象图形学报，2009，14（10）：1997-2003.（SHEN J， YANG X Z. A new method for the fabric defect defection based on texture watershed [J]. Journal of Image and Graphics， 2009， 14（10）： 1997-2003.）

[17] 李志林，朱庆.数字高程模型[M].2版.武汉：武汉大学出版社，2003：19-226.（LI Z L， ZHU Q. Digital Elevation Model [M]. 2nd ed. Wuhan： Wuhan University Press， 2003： 19-226.）

[18] 王琛智，汤国安，袁赛，等.基于DEM纹理特征的月貌自动识别方法探究[J].地球信息科学学报，2015，17（1）：45-53.（WANG C Z， TANG G A， YUAN S， et al. A method for identifying the lunar morphology based on texture from DEMs [J]. Journal of Geo-Information Science， 2015， 17（1）： 45-53.）

[19] 袁春兰，熊宗龙，周雪花，等.基于Sobel算子的图像边缘检测研究[J].激光与红外，2009，39（1）：85-87.（YUAN C L， XIONG Z L， ZHOU X H， et al. Study of infrared image edge detection based on Sobel operator [J]. Laser and Infrared， 2009， 39（1）： 85-87.）

[20] 任玉刚，张建，李淼，等.基于分水岭算法的作物病害叶片图像分割方法[J].计算机应用，2012，32（3）：752-755.（REN Y G， ZHANG J， LI M， et al. Segmentation method for crop disease leaf images based on watershed algorithm [J]. Journal of Computer Applications， 2012， 32（3）： 752-755.）