基于影像数据的冯·卡门撞击坑形貌分析

郑晨，姚鸿泰

（1. 河南大学 数学与统计学院应用数学研究所，开封 475000；2. 中国科学院 国家天文台，北京 100012）

0 引 言

从古时的夜观天象到当今的“嫦娥”飞天，人类一直未停止过对浩瀚宇宙的观测步伐。而距离人类最近的天体——月球，也成为了人类探测的重要对象。从20世纪50年代美国和前苏联发射的探测器飞掠月球，到实现人类首次登月的“阿波罗号”，再到我国开展的“嫦娥系列”探月工程，人类对月球的认知在不断深入和细化。

2013年12月15日，我国“嫦娥3号”着陆器在虹湾地区顺利着陆，成为了世界上第3个成功实施月球软着陆的国家[1-2]，也实现了我国探月工程“绕、落、回”中“落”的目标。但是，目前成功软着陆于月球的人类探测器主要位于月球正面，还没有在月球背面成功着陆的先例。2014年以来，中国探月与航天工程中心经多方论证，确定了“嫦娥4号”着陆月球背面这一任务目标和相应的技术方案[3-4]，有望实现人类首次月球背面着陆及勘察。着陆区域则初步选定为月球背面南极–艾肯（South Pole-Aitken，SPA）盆地内的冯·卡门（Von Kármán）撞击坑，该区域属于高地地形[2]，为“嫦娥4号”着陆器的软着陆带来了困难。

着陆区地形地貌的分析是影响着陆设计的一个重要方面。为了辅助“嫦娥4号”探月工程的开展，本文根据已有冯·卡门撞击坑的影像观测数据，利用基于随机场的聚类分析方法，从聚类的角度对该区域的地形地貌进行了分析与表示，以期为“嫦娥4号”着陆区的选择提供有意义的参考，为探测任务的顺利开展提供更多的技术支持。

1 数据与方法

1.1 实验数据

主要采用中心点坐标为（44.8°S，175.9°E），40～48°S，172°E～180°～178°W范围内的LOLA（Lunar Orbiter Laser Altimeter）高程DEM数据[5]和LROC（Lunar Reconnaissance Orbiter Camera）的高分辨率影像数据[6]（图1），对冯·卡门撞击坑的内部形貌进行分析。

冯·卡门撞击坑属于SPA内的典型地貌类型，其内部的月球物质成分具有代表性[2,7-8]，对月球火山及月壳活动的研究具有重要价值。同时，它位于月球背面，是开展低频射电天文观测的理想地点。而且，从LOLA高程数据和LROC影像数据的直方图（见图2）可以发现，冯·卡门撞击坑的高程数据（图2（a））主

要分布在–6～–1 km之间，且撞击坑内部的起伏变化较小；而LROC的影像数据值呈近似正态分布（图2（b）），坑内的地形地貌整体也较为平坦。因此，该区域是“嫦娥4号”初步选定的着陆区之一[4]。但是，从图1（b）可以观测到，撞击坑内仍存有大量的局部中小型环状撞击坑区域，在直方图中也反映出异常的峰值点（图2红色虚线框标注）。由于LOLA高程数据和LROC影像数据的灰度值分布相对集中，造成局部形貌的数值变化较小。这种局部低对比度现象将直接影响人工判读的精度，对分析撞击坑内部的形貌造成潜在的影响。

图1 冯·卡门撞击坑的实验数据Fig. 1 Test data of the Von Kármán crater

1.2 实验方法

为解决观测数据中低对比度问题，提高观测数据的可读性，本文研究了基于马尔科夫随机场模型（Markov Random Field，MRF）[9]的聚类分析方法。根据已有研究，聚类分析是提高深空探测数据局部低对比度的一种有效方法[10]，而MRF模型是一类利用概率图模型进行聚类的方法。它通过在观测数据和聚类标记间建立概率模型，不仅可以凸显特征相近的不同地物形貌间的差异，而且还可以刻画相邻像素的空间关系，提高模型的纹理描述能力。

图2 冯·卡门撞击坑LOLA和LORC数据的直方图Fig. 2 Histogram of LOLA and LROC data for the Von Kármán crater

具体而言，假设大小为的观测影像数据定义在栅格位置集合上，其中表示像素处的观测数值。标记随机场也定义在位置集合上，其中每个表示像素处的聚类类别。此处，是一个随机变量，且，其中是聚类的类别总数。若是标记随机场的一个实现，那么，聚类问题在MRF模型中就被转化为最大后验概率的求解问题，即

上式中是实现的集合，即根据贝叶斯公式，式（1）可表示为

其中第2个等号成立的原因是由于是已知的观测数据，不影响最终的聚类结果。所以，在MRF模型的聚类过程中，需要首先确定似然函数和标记场联合概率分布的具体函数形式。

对于联合概率MRF模型假设标记随机场具有空间马氏性，即其中表示与像素空间相邻的像素点组成的集合。由Hammersley–Clifford定理可知[9]，此时服从Gibbs分布

上式中是归一化常数，是实现的能量函数。其中，是定义在势团上的势函数。本文模型中，为了在保持聚类精度的前提下降低计算量，我们只考虑具有二阶势团的多层逻辑模型[11-12]，此时势函数为

对于似然函数假设在给定标记随机场实现的条件下，各像素点的观测特征是相互独立的，即根据图2的直方图，尤其是LROC的影像数据分布，本文模型采用正态分布来描述各像素点的似然函数分布，即在的条件下其中n是观测数据的波段数，参数分别是类别的正态均值和方差。

在确定了似然函数和标记场的联合概率后，通过不动点的思路来逐像素迭代更新（2）式的后验概率，实现最终的聚类分析，即

具体的算法请见表1。其中，关于式（7）中的正态分布参数，可以根据EM算法[13]估计如下

式中，表示类别为的像素点的数目。

表1 MRF算法流程Table 1 Algorithm of the MRF model

③根据公式（9），以及标记实现更新(7)式的参数及似然函数

④根据公式（5）和标记实现更新（6）式的标记联合概率

⑤根据公式（8）和最大后验准则，更新像素处的标记

⑥更新标记实现则输出作为聚类结果；否则，令t=t+ 1，返回步骤③。

2 实验结果

在上述的MRF模型，定义邻域集合为8邻域，分别对LROC影像数据和LOLA高程数据进行了聚类分析，具体结果如下。

2.1 LROC实验结果

图3 LROC影像数据K = 3时的聚类结果和对应直方图Fig. 3 Clustering results of LROC data with K = 3 and its corresponding histogram.

对于LROC影像数据，首先在K= 3的类别数目下进行了聚类分析，结果如图3（a）所示。此时， MRF模型的结果主要分为了撞击坑内较平坦的区域、向光区域和阴影区域这3种类别，它们对应的影像灰度值分别位于图2（b）中的峰值和两边的尾部区域，如图3（b）的直方图和表2的定量数值所示。由于冯·卡门撞击坑内的中小型撞击坑具有明显的向光和阴影区域，因此根据该结果，可以更直观快速地判读出这些局部撞击坑。然而，不同撞击坑间的差异却难以得到更精细的表示。

表2 K = 3时影像各类别的定量指标Table 2 Quantitative indexes of images with K = 3

为了进一步凸显不同撞击坑间的差异，考察了K=7时的聚类结果，如图4所示。此时聚类结果中每种类别对应的地形地貌被进一步细化，其中坑内的平坦区域被细分为两类（橙色和粉色区域）、向光区域被细分为3种不同程度的类别（青色、黄色和绿色区域）、阴影区域则被细分为蓝色和紫色两种类别。同时，相比于K= 3，此时每种类别内的标准差会相应降低，各类别的灰度均值则分别位于图2（b）中的一些局部异常峰值点附近，如表3所示。随着聚类的细化，图4的聚类结果可以描述出不同中小型撞击坑间的差异。

图4 LROC影像数据K = 7时的聚类结果Fig. 4 Clustering results of LROC data with K = 7

相比于LROC影像数据的直接判读，结合图3和图4的聚类结果，冯·卡门撞击坑内的地形地貌，尤其是撞击坑等局部低对比度区域，可以得到更合理的表示，如图5所示。具体而言，对于该区域内较大型的撞击坑，它的向光与阴影区域是渐变的，因此在K= 3的聚类结果中会从左至右地表现出向光区域与阴影区域的交替（一般都为两次，分别对应坑的边缘和内部光照变化）；而在K= 7的聚类结果中，大型撞击坑会含有多种类别的聚类结果，如在图5（a）中示例的撞击坑就由两种阴影区域类别和3种向光区域类别共同构成。对于坑内的中型撞击坑，由于坑的边缘不是很明显，在K= 3的聚类结果中通常只表现出坑内部的一次光照变化情况，如图5（c）所示，而K= 7的聚类结果中出现的类别数目也会少于大型撞击坑。对于小型和微型撞击坑，由于其直径很小，导致在K= 3的聚类结果中只表现出小块的阴影类别，其在K= 7的聚类结果中也一般只含有两类表示阴影区域的类别，如图5（b）所示。除了上述3种不同大小的环形撞击坑，冯·卡门撞击坑内还具有一些长条型的凹陷区域，本文称之为长型撞击坑。这种类型的撞击坑在光照变化方面与中型撞击坑类似，一般只含有一次的阴影与光照区域交替，所以类别特征的分布是相似的，其差别主要体现在坑的形状上。根据上述特征，冯·卡门撞击坑内典型形貌的统计结果总结如表4所示（注：因小型撞击坑数目过多，因此采用程序自动地统计了聚类结果，其数目显示为3 994，因程序统计可能存在一定误差，所以下表记为约4 000个）。

表3 K = 7时影像各类别的定量指标Table 3 Quantitative indexes of images with K = 7

图5 不同类型撞击坑在K = 3和K = 7时的特征Fig. 5 Characteristics of different craters with K = 3 and K = 7

表4 冯·卡门撞击坑内各典型形貌的统计结果Table 4 Statistical results of various craters in the Von Kármán crater

2.2 LOLA实验结果

对于LOLA高程数据，由于冯·卡门撞击坑内数据差异较小，导致整个区域的可视化特征非常接近，如图6（b）所示。这种低对比度的影像数据不能有效地凸显撞击坑内不同地物的差异，而一些在LROC数据中可以观测到的典型形貌也难以在高程数据中进行判读。为了描述不同地物形貌在LOLA高程数据中的差异，我们仍采用MRF模型对其进行聚类分析，具体结果请见图6。

从图6可以发现，当K值设定较小，如图6（c）中K= 3时，除了坑中部小型凸峰的峰值区域外，整个撞击坑区域被划分成一种类别，此时聚类结果主要表现出数据的宏观特征，即坑内区域整体较平坦。当进一步增大K值时，如图6（d）中K= 5时，坑内的典型地物，如中部的凸峰就可以得到较完整的显示。当K= 10时，聚类结果中可以显示出该区域内整体的高程变化和一些较深的中小型撞击坑。而当K= 22时，一些从图6（b）中难以目视判读的小型撞击坑，也可以在结果中得到表现；同时，整个冯·卡门撞击坑内的高程变化也可以细分为左上部分和右下部分两种类别，这和图4中橙色和粉色的平坦区域划分是相对应的。最后，当K= 28时，坑内右下部分中一块略有凹陷的区域可以得到表现。

纵观上述结果，通过MRF模型的聚类分析，LOLA高程数据中的低对比度区域可以得到更合理的表示，一些典型地物形貌也能在聚类结果中得以彰显。

因此，相比于目视判读，基于聚类的表示方法可以有效地提高低对比度区域的可视化水平，能辅助相应区域的地形分析。同时，本文使用的马尔科夫随机场模型因其有效的空间描述能力，可以进一步优化聚类的结果。为了突显本文模型的特点，在图7中对比了本文模型和经典Kmeans算法[10]在LOLA高程数据的聚类结果，可以看出，在K= 5的条件下，经典Kmeans算法因为仅考虑了各像素的高程数据，没有考虑不同像素间的空间关系，所以它的聚类结果（图7（a））没有完整地识别出中央的凸峰。本文采用的马尔科夫随机场模型通过似然函数和联合概率分别考虑了影像数据中的像素高程数据和各像素间的空间关系，更充分地使用了影像信息，所以聚类结果图7（b）中较完整地识别出了中央的凸峰。

图6 LOLA高程数据聚类结果Fig. 6 Clustering result of LOLA data

图7 本文方法和Kmeans方法的聚类结果对比Fig. 7 Comparison of clustering results between our method and the Kmeans

3 结 论

本文基于LOLA高程数据和LROC高空间分辨率影像数据，通过马尔科夫随机场模型对月球背面SPA内的冯·卡门撞击坑进行了地形地貌的聚类表示和分析。本文的贡献主要体现在以下两点：①模型的似然函数刻画了影像数据的概率分布，而标记随机场则描述了地形地貌间的空间关系，结合这两点的本文模型通过聚类可以有效地突显出影像数据中低对比度区域内地形地貌的差异；②通过分析聚类结果可知，冯·卡门撞击坑的地形从整体来看，除了在中上部有一个凸峰外，其余区域是相对平坦的，高低起伏也较缓慢；从局部来看，该区域内仍存在大量的环状撞击坑。而这些撞击坑又可以根据聚类结果分为大型、中型、小型和长型撞击坑，其中大型撞击坑不仅深度较大，而且坑的边缘还有明显的地势起伏，但数量相对较少；中型和长型撞击坑的边缘相对平缓，但是撞击坑自身仍具有一定的深度，且分布于冯·卡门撞击坑内部各处；而小型撞击坑的高程变化较小，但数量却很多，广布于拟着陆区域。在工程可控的前提下，着陆区域可考虑选择凸峰西侧的下方。因为此处不仅相对平坦，而且还接近各类型撞击坑，便于展开后续的科学探测工作。

致谢

本文实验LROC数据和LOLA数据分别从Java Missionplanning and Analysis for Remote Sensing（JMARS）https://jmars.mars.asu.edu/和网站http://imbrium.mit.edu/DATA/LOLA_GDR/CYLINDRICAL/FLOAT_IMG/下载，在此表示感谢！

参考文献

[1]SUN Z Z，JIA Y，ZHANG H. Technological advancements and promotion roles of Chang’e-3 lunar probe mission[J]. Science China（Technological Sciences），2013，43（11）：2702-2708.

[2]李飞，张熇，吴学英，等. 月球背面地形对软着陆探测的影响分析[J].深空探测学报，2017，4（2）：143-149.LI F，ZHANG H，WU X Y，et al. Influence analysis of terrain of the farside of the Moon on soft-landing[J]. Journal of Deep Space Exploration，2017，4（2）：143-149.

[3]YE P J，SUN Z Z，ZHANG H，et al.An overview of the mission and technical characteristics of Change’4 Lunar Probe [J]. Science China（Technological Sciences），2017，60（5）：658-667.

[4]吴伟仁，王琼，唐玉华，等. “嫦娥4号”月球背面软着陆任务设计[J].深空探测学报，2017，4（2）：111-117.WU W R，WANG Q，TANG Y H，et al. Design of Chang’e-4 lunar farside soft-landing mission [J]. Journal of Deep Space Exploration，2017，4（2）：111-117.

[5]SMITH D E，ZUBER M T，JACKSON G B，et al. The Lunar orbiter laser altimeter investigation on the lunar reconnaissance orbiter mission[J]. Space Science Reviews，2010，150 （1-4）：209-241.

[6]ROBINSON M S，BRYLOW S M，TSCHIMMEL M，et al. Lunar reconnaissance orbiter Camera（LROC）instrument overview [J]. Space Science Reviews，2010，150 （1-4）：81-124.

[7]PETERSON C A，HAWKE B R，BLEWETT D，et al. Geochemical units on the Moon：the role of South Pole-Aitken basin[C]//Lunar and PlanetaryScience XXXIII. USA：[s.n.]，2002.

[8]Peterson C A，Hawke B R，Lucey P G，et al. Anorthosite on the lunar farside and its relationship to South Pole-Aitken Basin[C]//Lunar and Planetary Science XXXI. USA：[s.n.]，2000.

[9]LI S Z. Markov random field modeling in computer vision. 3ed ed.[M]. New York：Springer-Verlag，2009.

[10]ZHENG C，PING J S，WANG M Y. Hierarchical classification for the topography analysis of Asteroid（4179）Toutatis from the Chang’E-2 images [J]. ICARUS，2016，278：119-127.

[11]DERINH，ELLIOTT H. Modeling and segmentation of noisy and textured images using gibbs random fields [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence，1987，9（1）：39-55.

[12]DERINH，COLE W S. Segmentation of textured images using Gibbs random fields [J]. Computer Vision，Graphics and Image Processing，1986，35：72-98.

[13]DEMPSTER A P，LAIRD N M，RUBIN D B. Maximum likelihood from incomplete data via the EM algorithm [J]. Journal of the Royal Statistical Society（series b），1977，39（1）：1-38.

[14]JAIN A K，MURTY M N，FLYNN PJ. Data clustering：A review [J].ACM Computing Surveys，1999，31（3）：264-323.