2. 地理学国家级实验教学示范中心(陕西师范大学), 陕西 西安 710119;
3. 华北水利水电大学艺术与设计学院, 河南 郑州 450046
2. National Demonstration Center for Experimental Geography Education(Shaanxi Normal University), Xi'an 710119, China;
3. School of Art and Design, North China University of Water Resources and Electric Power, Zhengzhou 450046, China
撞击坑是小行星、彗星、陨石等小天体撞击行星表面而形成的一种凹状地形单元[1]。因月球无大气、无水、无风等环境因素，月表撞击坑得以较好的完整保存[2]。作为月球表面的一种主要地貌形态，撞击坑的精准识别与统计分析对于人类研究月球地质运动和演化历史等具有重要的意义[3]。由于月表撞击坑的数目繁多、形态不一，目视解译耗时高，效率低，因此撞击坑的自动识别成为当前月球科学领域的研究基础和关键。
撞击坑的识别经历了人工观测手绘、遥感影像目视判读和高分辨率数据自动提取3个阶段[4]。从数据源来看，可将撞击坑的识别方法分为3种类型：①基于月表遥感影像数据；②基于月表DEM数据；③基于影像与DEM融合数据。其中，第一种方法主要是利用撞击坑在遥感影像上的灰度差异，在经过边缘检测[5-6]、纹理分析[7-8]等数据预处理后，运用面向对象[9-11]、Hough变换及其改进[12-13]、模板匹配[14-16]、机器学习[17-18]等方法自动提取撞击坑[19]。这种方法对于灰度差异比较明显的撞击坑都能很好地提取出来。然而由于遥感影像数据固有的光学误差及月表实际撞击坑复杂的地形特征等影响，导致撞击坑在遥感影像上出现几何学畸变，增加了影像识别的不确定性。第二种方法主要是利用撞击坑实际地形信息，使用Hough变换[1, 20-21]、傅里叶变换[22]、机器学习[23-24]等方法提取撞击坑[25]。相比于影像数据，DEM数据能够再现撞击坑的三维地形结构[26]且不受光照影响，基于DEM数据的撞击坑识别方法能根据月表真实的地形形态正确提取撞击坑。然而，由于早期月表DEM数据分辨率低及与之相关的技术水平的局限，许多撞击坑不能准确有效地提取出来。近年来，随着月表DEM数据分辨率的提高及地理信息系统技术的日趋成熟，诸多学者开始采用融合影像和DEM数据的方法提取撞击坑[27-30]，这种方法虽然在一定程度上提高了撞击坑识别的准确性，但是由于影像数据和DEM数据不同的数据格式及不同分辨率的影响，二者不能有效地融合，并且通常适用于影像识别的方法不适合DEM数据的提取。
综上所述，撞击坑作为月表的一种典型月貌单元，高程是定量表达其形态特征最为基本的几何量。因此，月表DEM数据所蕴含的撞击坑结构和地形特征较影像数据更为丰富，DEM数据更适用于撞击坑的自动提取[31]。本文基于美国LRO飞行器获取的100 m分辨率月表DEM数据，采用线状窗口邻域分析方法准确提取撞击坑边界。
1 研究样区与研究方法 1.1 研究样区月球表面主要由两种地貌单元组成：月陆和月海[32]。基于地貌类型的完整性和算法的普适性原则，本文分别选择月陆和月海两个研究样区(如图 1所示)。其中，月陆样区位于LQ-6洪堡海和LQ-7达朗贝尔[33]交界的莫斯科海和洪堡海之间地区，与韦尔斯、达贝尔、达朗贝尔等大型撞击坑相邻，经纬度范围为35°N—50°N、117°E—140°E，面积大约为314 400 km2。此研究样区位于月球背面，撞击事件频繁，以大型撞击坑为主。月海样区位于LQ-4虹湾和LQ-11哥白尼交界的雨海地区，与风暴洋、澄海、虹湾等地区相邻，经纬度范围为29°N—44°N、10°W—33°W，面积大约为317 300 km2。此研究样区位于月球正面，撞击事件较背面少，撞击坑多是由月海玄武岩经岩溶填充而成，以中小型撞击坑为主。
|
| 图 1 研究样区概况 |
本文利用线状窗口邻域分析法(如图 2所示)[34]提取撞击坑。地形特征线微分单元走向大致可分为4个方向：东西方向、南北方向、东北-西南方向、西北-东南方向。在DEM中，分析窗口可设定为水平方向上的1×n和垂直方向上的n×1的线状邻域分析窗口(n=3，5，7，…)。在基于DEM进行线状窗口分析时，线状分析窗口将遍历整幅DEM，导致水平方向的分析窗口对大致南北方向的地形特征线具有过滤作用，而垂直方向的分析窗口对大致东西走向的地形特征线具有过滤作用。最后，对分别使用1×n和n×1的线状邻域分析窗口过滤后的结果求和即可得到潜在撞击坑的坑唇线。
|
| 图 2 线状窗口邻域分析原理 |
撞击坑自动提取流程如图 3所示。首先，对原始样区DEM数据进行水平和垂直方向上线状窗口的邻域均值分析。给定一个阈值，将原始样区DEM分别与水平和垂直方向均值分析后DEM进行栅格减运算并使其大于给定阈值；其次，将水平与垂直方向上的结果进行栅格加运算，将和值进行重分类二值化得到初步结果；最后，将二值化结果进行数学形态学处理，得到最终的月表撞击坑提取结果。
|
| 图 3 算法的技术路线 |
针对两个研究样区，首先分别使用1×n和n×1的水平和垂直方向上的线状分析窗口对样区DEM进行邻域均值分析，将原始样区DEM与均值分析后的DEM作栅格减运算。然后给定一个阈值，对水平与垂直方向上的相减运算结果分别取相同阈值大小进行过滤，将二者相加求和后二值化，得到该阈值下1×n和n×1分析窗口的结果。最后依次改变分析窗口大小及提取阈值试验，直至得到坑唇线栅格突出、连续，区域杂质较少的撞击坑最佳提取效果。试验结果表明，随着分析窗口的扩大，撞击坑的提取效果依次增强，但是杂质噪音也随之增多。需要指出的是，不同大小的分析窗口下都有一个最佳的提取阈值，阈值的选取标准与原始DEM减去均值分析后的DEM高程差值相关。最终将各个窗口的最佳提取效果进行叠加比对，从中选择最佳的线状分析窗口，以及该窗口下的最佳提取阈值。本文经过撞击坑的坑唇线提取效果与杂质噪音二者综合取舍试验，确定本文的月海样区最佳提取窗口为1×7和7×1，该窗口下的最佳提取阈值为4.5；月陆样区最佳提取窗口为1×9和9×1，该窗口下的最佳提取阈值为6(如图 4所示)。
|
| 图 4 线状窗口邻域分析最佳结果 |
研究样区存在一定的高低起伏，导致提取的撞击坑存在一定伪坑，这些伪坑主要是月表山脉、月溪、月谷、断裂带等月貌单元[33]。由图 4可知，本文算法对撞击坑进行正确提取的同时，存在局部撞击坑边界不连续，没有形成闭合曲线等问题。为此需要对提取结果进行数学形态学处理[35-36]。
本文利用Matlab软件选择形态学开运算，对提取结果进行细化成线、去除孤立对象、膨胀连接等操作，得到坑唇线连续、闭合的撞击坑提取结果(如图 5所示)。
|
| 图 5 数学形态学处理结果 |
由图 5可知，数学形态学处理后的效果仍然存在少量杂质。为得到最佳的提取效果，必须将上述形态学处理后的栅格坑唇线转换为矢量坑唇线，以原始DEM数据晕渲图为映衬，对矢量坑唇线辅以人工手动编辑处理，去除剩余的杂质，连接少量未闭合的坑唇线，得到正确完整提取、正确部分提取、错误提取和未提取4种结果(如图 6所示)。为便于目视识别，本文只保留直径大于1 km的撞击坑。
|
| 图 6 两试验样区撞击坑提取结果 |
本文以基于DEM坡度变率为基准、以晕渲图为背景，辅以坡度、坡向等地形因子，综合判断识别直径大于1 km的撞击坑为标准结果，对撞击坑提取结果从以下3方面进行评价：①提取精度；②提取状态；③与前人提取结果比对。
3.1 提取精度评价本文采用Shufelt[37]提出的算法评价因子来评价算法的提取精度(见表 1)。其因子指标有：构建检测百分比D(building detection percentage)=100TP/(TP+FN)，分支系数B(branching factor)=FP/TP，质量百分比Q(quality percentage)=100TP/(TP+FP+FN)。其中，TP为识别到的撞击坑且为真实撞击坑的个数；FN为未能识别到的真实撞击坑个数；FP为识别到的撞击坑但为非真实的撞击坑个数。构建检测百分比D可以作为撞击坑识别精度的度量，分支系数B可以作为算法的正确性度量，质量百分比Q可以作整个算法质量的度量。由于存在未能完全提取的撞击坑，本文以提取部分占所在撞击坑的1/2及以上作为真实提取的撞击坑，反之为未能提取的撞击坑。
| 研究样区 | TP | FN | FP | D /(%) | B | Q /(%) |
| 月海 | 85 | 24 | 11 | 77.98 | 0.13 | 70.83 |
| 月陆 | 164 | 47 | 21 | 77.73 | 0.13 | 70.69 |
评价结果显示，本文月海和月陆的撞击坑识别精度评价因子D达到77%以上，算法的分支系数B为0.13，整个算法的提取质量Q在70%以上。由于月海的撞击坑多以简单撞击坑为主而月陆的撞击坑多以复杂撞击坑为主，但是月海撞击坑在直径上却远远小于月陆撞击坑，导致本文算法在提取精度上对于月海和月陆样区的评价因子数值比较接近，表明本文算法对月海撞击坑和月陆撞击坑在提取精度上虽然都不够高但是彼此差别不大，因此可将本文算法应用于全月球撞击坑的初步自动提取。
3.2 提取状态评价由提取的结果图显示，撞击坑提取状态有正确提取TE(true extraction)、未提取NE(not extraction)和错误提取FE(false extraction)3种情况，其中正确提取和错误提取又可分为完全提取CE(complete extraction)、部分提取PE(partial extraction)两种情况。因此将正确提取分为正确完整提取(TCE)和正确部分提取(TPE)。为了更加精确地评价算法，在正确部分提取基础上又细分为提取部分占整个撞击坑比例的1/2及以上和1/2以下两种情况。分别计算4种提取状态下撞击坑个数占手动勾绘标准撞击坑个数的百分比(见表 2)，从而进一步评价算法的提取效果，展示提取的撞击坑的完整情况。
| 研究样区 | TCE | TPE | NE | FE | 标准 | TCEZ占比/(%) | TPE占比/(%) | NE占比/(%) | FE占比/(%) | ||
| ≥1/2 | < 1/2 | ≥1/2 | < 1/2 | ||||||||
| 月海 | 64 | 21 | 9 | 15 | 11 | 109 | 58.72 | 19.27 | 8.26 | 13.76 | 10.09 |
| 月陆 | 120 | 44 | 16 | 31 | 21 | 211 | 56.87 | 20.85 | 7.58 | 14.69 | 9.95 |
评价结果显示，通过与目视识别的标准结果进行比对计算，月海和月陆样区的撞击坑正确完整提取率达到56%以上；正确部分提取率在27%以上，其中提取部分占整个撞击坑比例的1/2及以上的达到20%左右；未提取率在15%以下；错误提取率在10%左右，表明两个样区的各项状态指标在数值上比较接近。如果将正确部分提取中提取部分占整个撞击坑比例的1/2及以上算作理想的提取结果，加上实际的正确完整提取的结果，两个样区的撞击坑有效提取率将分别达到77.99%和77.72%，既验证了前面部分对撞击坑提取精度的评价，也说明了本文算法的提取结果满足对全月表撞击坑的初步自动提取。
3.3 与前人提取结果比对评价通常不同算法的提取结果评价应当建立在相同分辨率的数据类型和研究样区的基础之上，这样的结果才具有可比性。虽然目前撞击坑的提取主要是基于遥感影像和DEM数据两种不同数据类型而进行的，但遗憾的是目前无法找到基于相同数据类型且相同分辨率的研究，而且研究样区也都不一致，因而无法进行完全公平的检测结果比对。为了比较客观地评价算法的提取结果，本文选择与Zuo W等[38]基于DEM数据等高线的撞击坑自动识别的检测方法的提取结果进行评价因子比对(见表 3)。选择原因有二：①在Zuo W等的研究中，研究数据采用嫦娥一号500 m分辨率和嫦娥二号30 m分辨率DEM数据，与本文100 m分辨率DEM的数据类型相同且有很好的过渡性；②Zuo W等的研究也将DEM数据应用于两个研究样区，且与本文研究样区的撞击坑特征相似，二者具有极高的可比性。
| 研究样区 | 空间分辨率/m | TP | FN | FP | D /(%) | B | Q /(%) |
| area#1 | 500 | 96 | 16 | 18 | 85.71 | 0.19 | 73.85 |
| area#2 | 30 | 145 | 34 | 39 | 81.01 | 0.27 | 66.51 |
| 月海 | 100 | 85 | 24 | 11 | 77.98 | 0.13 | 70.83 |
| 月陆 | 100 | 164 | 47 | 21 | 77.73 | 0.13 | 70.69 |
比对结果显示，Zuo W等两个样区的提取结果中，检测百分比D都达到80%以上，高于本文算法提取结果，但分支系数B高于本文结果。说明其结果虽然在提取精度上高于本文结果，但是其提取结果的错误率也要高于本文结果。从质量百分比Q来看，本文100 m分辨率提取结果对于其500 m分辨率和30 m分辨率提取结果，两个算法的整体提取质量相当。
4 结论与讨论本文利用100 m分辨率DEM数据，根据撞击坑坑唇相对凸起的地形特征，提出一种线状窗口邻域均值分析方法自动提取月表撞击坑，通过在月海与月陆两种不同样区的研究，在综合分析提取方法和提取结果的基础之上，主要得出以下结论：
(1) 本文提出的基于DEM线状窗口邻域分析的方法，能够根据月表撞击坑实际地形形态对其进行自动检测。由提取结果的评价因子显示，该方法具有一定的普适性。
(2) 从试验的过程来看，不同大小的线状分析窗口对撞击坑的提取结果不同，且在相同大小的线状分析窗口下，不同的提取阈值也会产生不同的提取结果。因此，不同大小的分析窗口下都有一个最佳的提取阈值。通过不同大小分析窗口及相同大小分析窗口下不同阈值的试验比对，本文的月海样区的最佳提取窗口为1×7和7×1，该窗口下的最佳提取阈值为4.5；月陆样区的最佳提取窗口为1×9和9×1，该窗口下的最佳提取阈值为6。
(3) 从提取结果的状态来看，本文提取结果与目视识别结果比对显示，两个样区的完整提取、部分提取、未提取和错误提取占比分别约为58%、28%、14%和10%，其中部分提取中提取部分占整个撞击坑比例的1/2及以上的约占20%，如果将其算作理想提取结果，则本文撞击坑有效提取将达到78%以上。
(4) 从本文结果与Zuo W等基于DEM数据等高线的撞击坑自动识别的检测方法结果比对情况来看，虽然Zuo W等检测结果在提取精度上高于本文结果，但是其提取错误率也要高于本文结果。两个算法的提取质量上各有优势，表明本文方法在月表撞击坑的自动提取方面有一定的应用价值。
然而，由于月表实际地形及撞击坑几何形态组合的复杂多样性，该方法未能完全实现提取过程的自动化，特别是未提取和错误提取达到了24%左右，这对月表撞击坑的自动提取来说显然不是最佳的。未提取的撞击坑主要有两类：一是坑唇退化的撞击坑，其退化原因主要是月球自身地质活动演化和后期撞击坑的随机撞击结果所致；二是地势较低的撞击坑，由于其高程差值不大，与地势较高的撞击坑组合在一起时，算法难以确定一个能完全囊括不同地势高低起伏的提取阈值。错误提取的撞击坑实为伪撞击坑，包括由于多个撞击坑边缘彼此相交连接而形成的伪坑及其他月貌单元，如山脉、月岭、月溪、月谷等。因此，在后期的研究中，因充分结合坡度、坡向、曲率等地形因子综合判断，设计不同的线状分析窗口大小，实现对撞击坑的正确有效提取。
| [1] | LUO L, WANG X, GUO H, et al. Automatic Detection of Lunar Craters Based on CE-1 DEM Data[J]. Journal of Remote Sensing, 2014(1): 105–116. |
| [2] | 欧阳自远. 月球科学概论[M]. 北京: 中国宇航出版社, 2005: 58-64. |
| [3] | HUANG Q, PING J, SU X, et al. New Features of the Moon Revealed and Identified by CLTM-s01[J]. Science in China:Physics, Mechanics and Astronomy, 2009, 52(12): 1815–1823. DOI:10.1007/s11433-009-0284-x |
| [4] | 王娇, 程维明, 周成虎. 全月球撞击坑识别、分类及空间分布[J]. 地理科学进展, 2015, 34(3): 330–339. DOI:10.11820/dlkxjz.2015.03.008 |
| [5] | LEROY B, MEDIONI G, JOHNSON E, et al. Crater Detection for Autonomous Landing on Asteroids[J]. Image & Vision Computing, 2001, 19(11): 787–792. |
| [6] | 李超, 王心源, 骆磊, 等. 基于Apollo图像的月表撞击坑自动提取[J]. 国土资源遥感, 2012, 24(4): 71–75. DOI:10.6046/gtzyyg.2012.04.12 |
| [7] | BARATA T, ALVES E I, SARAIVA J, et al. Automatic Recognition of Impact Craters on the Surface of Mars[J]. Lecture Notes in Computer Science, 2004, 3212(3212): 489–496. |
| [8] | BANDEIRA L, WEI D, STEPINSKI T F. Detection of Sub-kilometer Craters in High Resolution Planetary Images Using Shape and Texture Features[J]. Advances in Space Research, 2012, 49(1): 64–74. DOI:10.1016/j.asr.2011.08.021 |
| [9] | 岳宗玉, 刘建忠, 吴淦国. 应用面向对象分类方法对月球撞击坑进行自动识别[J]. 科学通报, 2008, 53(22): 2809–2813. DOI:10.3321/j.issn:0023-074X.2008.22.020 |
| [10] | 陈伟涛, 闫柏琨, 张志. 基于嫦娥一号CCD数据空间特征的特定目标识别[J]. 国土资源遥感, 2009, 21(4): 40–44. DOI:10.6046/gtzyyg.2009.04.08 |
| [11] | BARNES A H. An Object-oriented Classification of Impact Craters Using Lunar Reconnaissance Orbiter Data[D]. San Francisco, California: San Francisco State University, 2016. |
| [12] | DING M, CAO Y, WU Q. Method of Passive Image Based Crater Autonomous Detection[J]. Chinese Journal of Aeronautics, 2009, 22(3): 301–306. DOI:10.1016/S1000-9361(08)60103-X |
| [13] | HONDA R, AZUMA R. Crater Extraction and Classification System for Lunar Images[R]. Kochi: Kochi University, 2000: 13-22. |
| [14] | 魏士俨, 张建利, 彭松, 等. 虹湾地区月面撞击坑自动提取[J]. 计算机仿真, 2013, 30(8): 74–77. |
| [15] | BANDEIRA L, SARAIVA J, PINA P. Impact Crater Recognition on Mars Based on a Probability Volume Created by Template Matching[J]. IEEE Transactions on Geoscience & Remote Sensing, 2007, 45(12): 4008–4015. |
| [16] | BANDEIRA L P C, SARAIVA J, PINA P. Development of a Methodology for Automated Crater Detection on Planetary Images[J]. Springer Berlin Heidelberg, 2007, 4447(4): 193–200. |
| [17] | WANG Y, YANG G, GUO L. A Novel Sparse Boosting Method for Crater Detection in the High Resolution Planetary Image[J]. Advances Space Research, 2015, 56(5): 982–991. DOI:10.1016/j.asr.2015.05.014 |
| [18] | XIN X, DI K C, WANG Y X, et al. Automated Detection of New Impact Sites on Martian Surface from HiRISE Images[J]. Advancesm Space Research, 2017, 60(7): 1557–1569. DOI:10.1016/j.asr.2017.06.044 |
| [19] | 丁萌, 李海波, 曹云峰, 等. 基于光学图像的撞击坑识别研究综述[J]. 深空探测学报, 2015, 2(3): 195–202. |
| [20] | MICHAEL G G. Coordinate Registration by Automated Crater Recognition[J]. Planetary & Space Science, 2003, 51(9/10): 563–568. |
| [21] | BUE B D, STEPINSKI T F. Machine Detection of Martian Impact Craters from Digital Topography Data[J]. IEEE Transactions on Geoscience & Remote Sensing, 2006, 45(1): 265–274. |
| [22] | HARADA N, HAYASHI T, HIRATA N, et al. Recognition Algorithm for Topographic Features[C]//17th IEEE International Conference on Computer and Information Technology(CIT 2007). [S. l. ]: CIT, 2007: 685-689. |
| [23] | DI K, LI W, YUE Z, et al. A Machine Learning Approach to Crater Detection from Topographic Data[J]. Advances in Space Research, 2014, 54(11): 2419–2429. DOI:10.1016/j.asr.2014.08.018 |
| [24] | COHEN J P. Automated Crater Detection Using Machine Learning[D]. Boston: University of Massachusetts Boston, 2016. |
| [25] | 刘宇轩, 刘建军, 牟伶俐, 等. 撞击坑识别方法综述[J]. 天文研究与技术, 2012, 9(2): 203–212. |
| [26] | 张锋, 邹永廖, 郑永春, 等. 月表撞击坑自动识别与提取的新方法及其应用[J]. 地学前缘, 2012, 19(6): 118–127. |
| [27] | 罗中飞, 康志忠, 刘心怡. 融合嫦娥一号CCD影像与DEM数据的月球撞击坑自动提取和识别[J]. 测绘学报, 2014, 43(9): 924–930. |
| [28] | KANG Z, LUO Z, HU T. Automatic Extraction and Identification of Lunar Impact Craters Based on Optical Data and DEMs Acquired by the Chang'E Satellites[J]. IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2016, 8(10): 4751–4761. |
| [29] | SALAMUNIC'CAR G, LONČARIC' S, PINA P. Integrated Method for Crater Detection from Topography and Optical Images and the New PH9224GT Catalogue of Phobos Impact Craters[J]. Advances in Space Research, 2014, 53(12): 1798–1809. DOI:10.1016/j.asr.2013.11.006 |
| [30] | 周增坡, 程维明, 周成虎, 等. 基于"嫦娥一号"的月表形貌特征分析与自动提取[J]. 科学通报, 2011, 58(1): 18–26. |
| [31] | 贺力. 基于"DEM"的月球撞击坑判识方法研究[D]. 南京: 南京师范大学, 2012. |
| [32] | 王琛智, 汤国安, 袁赛, 等. 基于DEM纹理特征的月貌自动识别方法探究[J]. 地球信息科学学报, 2015, 17(1): 45–53. |
| [33] | 王梁, 丁孝忠, 韩坤英, 等. 月球数字地质图的编制与研究[J]. 中国地质, 2015(1): 331–341. |
| [34] | 周毅, 汤国安, 张婷, 等. 基于格网DEM线状分析窗口的地形特征线快速提取方法[J]. 测绘通报, 2007(10): 67–69. DOI:10.3969/j.issn.0494-0911.2007.10.021 |
| [35] | 陶兴龙. 数学形态学和分形理论在图像去噪和边缘提取中的应用[D]. 南京: 南京信息工程大学, 2013. |
| [36] | 袁悦峰, 朱培民, 赵娜, 等. 基于数学形态学的月海圆形撞击坑自动识别方法[J]. 中国科学(物理学力学天文学), 2013, 43(3): 324–332. |
| [37] | SHUFELT J. Performance Evaluation and Analysis of Monocular Building Extraction from Aerial Imagery[J]. IEEE Transactions on Pattern Analysis & Machine Intelligence, 1999, 21(4): 311–326. |
| [38] | ZUO W, ZHANG Z, LI C, et al. Contour-based Automatic Crater Recognition Using Digital Elevation Models from Chang'E Missions[J]. Computers & Geosciences, 2016, 97(C): 79–88. |