劉宇軒，李春來，劉建軍

(1. 中國科學(xué)院國家天文臺，北京 100101；2. 中國科學(xué)院大學(xué)，北京 100049，3. 中國科學(xué)院月球與深空探測重點實驗室，北京 100101)

行星表面通常分布著大量的撞擊坑構(gòu)造，它們是行星表面最顯著的地貌特征，也是研究行星內(nèi)部物質(zhì)的窗口[1]。由于月球上沒有大氣，撞擊坑保存比較完好，因此研究月球形貌、相對地質(zhì)年齡、巖石構(gòu)造以及月球的起源有重要意義。撞擊坑識別的研究，從早年使用望遠(yuǎn)鏡觀察描繪撞擊坑的形貌，到利用遙感圖像數(shù)據(jù)目視識別撞擊坑，再到近年來通過計算機自動識別撞擊坑。隨著月球影像與地形數(shù)據(jù)的大量獲取，撞擊坑的自動識別提取成為一個新的研究熱點。近年來，國內(nèi)外眾多學(xué)者提出大量的研究方法，并建立多種尺度的撞擊坑數(shù)據(jù)庫[2-6]。通過對現(xiàn)有方法的分析總結(jié)，撞擊坑識別的方法大致可以歸納為4類：人工目視識別、基于形態(tài)擬合、基于機器學(xué)習(xí)以及基于地學(xué)信息分析的方法。人工目視識別法：利用影像和地形數(shù)據(jù)，建立識別撞擊坑的解譯標(biāo)識，在計算機環(huán)境下，采用人工交互的方式對撞擊坑逐個識別，這種方法是目前最準(zhǔn)確的撞擊坑識別方法，該方法的優(yōu)點是準(zhǔn)確性高，但費時費力，不適合于全月球和小尺度撞擊坑識別?；谛螒B(tài)擬合法：由于撞擊坑在圖像上有類圓形的幾何特點，在圖像分割的基礎(chǔ)上以霍夫變換和二次曲線為代表的形狀擬合算法用于模擬提取完整的撞擊坑邊界，其中文[7]利用組合霍夫變換算法，文[8]采用模糊霍夫變換算法對撞擊坑進行識別，文[9]提出一種基于弦中點霍夫變換算法，文[10]采用最小二乘法擬合邊緣橢圓的方法實現(xiàn)對撞擊坑的提取?；跈C器學(xué)習(xí)的方法：文[11]采用BP神經(jīng)網(wǎng)絡(luò)加訓(xùn)練樣本的算法，文[12]利用決策樹算法，文[13]利用面向?qū)ο蠓椒▽崿F(xiàn)對撞擊坑的識別。基于地學(xué)信息分析方法：文[14]利用坡度作為識別撞擊坑的要素，文[15-17]利用500 m分辨率的嫦娥一號數(shù)字高程模型數(shù)據(jù)和形態(tài)擬合的方法，進行撞擊坑的識別提取和研究，文[18]利用坡度數(shù)據(jù)、紋理信息和剖面曲率作為識別依據(jù)，采用水淹法對火星撞擊坑進行分析。

本文使用目前世界上覆蓋最全、精度最高的嫦娥二號數(shù)字地形數(shù)據(jù)，提出基于等高線地形分析的方法，針對小型撞擊坑進行識別提取研究，填補撞擊坑數(shù)據(jù)庫小型撞擊坑數(shù)量不足和準(zhǔn)確率不高的缺陷。

1 數(shù) 據(jù)

2008年2月，國務(wù)院批準(zhǔn)探月工程二期立項，二期工程主要目標(biāo)是實現(xiàn)月面軟著陸，開展月面就位探測與自動巡視勘察。嫦娥二號作為二期工程的技術(shù)先導(dǎo)星，首要的任務(wù)是對全月球和嫦娥三號的月面著陸區(qū)進行詳查，精細(xì)測繪全月球的地形地貌。嫦娥二號衛(wèi)星于北京時間2010年10月1日18時59分57秒在西昌衛(wèi)星發(fā)射中心成功發(fā)射。10月24日16時49分，嫦娥二號衛(wèi)星上的兩線陣CCD立體相機首次開機工作，并成功獲取了前視與后視兩個視角的影像數(shù)據(jù)。截至2011年5月20日，嫦娥二號衛(wèi)星CCD立體相機完成全部任務(wù)，共獲取607軌圖像數(shù)據(jù)。全月球影像圖制作過程中選取的100 km軌道高度成像的384軌數(shù)據(jù)(其中南北緯70°范圍344軌，極區(qū)40軌)，用于全月球數(shù)字高程模型數(shù)據(jù)產(chǎn)品生產(chǎn)。數(shù)據(jù)處理流程包括：數(shù)據(jù)準(zhǔn)備、測區(qū)數(shù)據(jù)組織、平差點準(zhǔn)備、全球平差、測區(qū)初級產(chǎn)品生產(chǎn)、標(biāo)準(zhǔn)產(chǎn)品生產(chǎn)、全球數(shù)據(jù)鑲嵌等。結(jié)合嫦娥二號衛(wèi)星軌道與姿態(tài)測量數(shù)據(jù)、5個月面絕對控制點進行全球平差，經(jīng)過4年的數(shù)據(jù)處理，成功制作了7 m、20 m、50 m分辨率的數(shù)字高程模型數(shù)據(jù)。嫦娥二號數(shù)字高程模型數(shù)據(jù)在空間分辨率、數(shù)據(jù)一致性和完整性鑲嵌精度等方面優(yōu)于國際同類全月球數(shù)字產(chǎn)品。本文選擇7 m分辨率的數(shù)字高程模型數(shù)據(jù)進行小型撞擊坑的識別提取。

2 方 法

2.1 撞擊坑模型分析

根據(jù)基本形態(tài)，撞擊坑可分為簡單撞擊坑(碗型小型撞擊坑)、復(fù)雜撞擊坑(距中央峰的較大型撞擊坑)、多環(huán)撞擊坑(撞擊盆地)3類[19]。如圖1，一個理想狀態(tài)下的復(fù)雜撞擊坑模型主要由坑唇邊緣(rh)、坑壁(w)、坑底(f)、撞擊培巖(m)和中央峰或中央環(huán)(p)5種形態(tài)要素組成。撞擊濺射堆積物由坑唇邊緣(rh)、丘狀沉積物(rr)和輻射線沉積物(rc)構(gòu)成。撞擊溶巖(m)覆蓋在其他形態(tài)要素之上。簡單碗型撞擊坑(c)一般分布在丘狀沉積物(rr)之上。撞擊坑識別是一種特征識別，對影像或地形數(shù)據(jù)進行分析處理，如圖2，識別撞擊坑的坑唇邊緣(rh)，并提取撞擊坑的直徑(D)、中心位置和高程等信息。

圖1 理想的的復(fù)雜撞擊坑模型(Don wilhelms, 1987)
Fig.1 Ideal complex impact crater model

等高線是地形圖上高程相等的相鄰各點連成的閉合曲線，是最為常用的表達地形形態(tài)特征和起伏情況的形式，如圖3。等高線從數(shù)字高程模型中提取出來，科學(xué)地反映月表的高程 、坡度、坡向、高地走勢等基本月表地形地貌。嫦娥二號數(shù)字高程模型數(shù)據(jù)的空間分辨率很高，月表高程在二維地理空間上得到了連續(xù)表達，利用該數(shù)字高程模型數(shù)據(jù)較為方便、精細(xì)地提取月球撞擊坑地形等高線，表達月球撞擊坑形態(tài)，并對其進行分析和處理。

從等高線的形態(tài)分析，如圖4撞擊坑具有3個典型的空間特點：(1)撞擊坑坑唇邊緣等高線通常為近似圓形或橢圓形；(2)撞擊坑坑唇邊緣等高線內(nèi)側(cè)通常存在多條嵌套的相似等高線，多為同心圓或橢圓，且間距較小，外側(cè)等高線形狀不定，且間距較遠(yuǎn)，與坑唇邊緣有明顯區(qū)別；(3)撞擊坑坑底高程小，坑邊緣高程大，其剖面線可以用下凹的拋物線形態(tài)擬合。因此，通過撞擊坑地形等高線的特征以及與其他等高線之間的相互關(guān)系，可以分析并識別月表的撞擊坑。因此，只要選取合適的等高距，撞擊坑坑唇邊緣的閉合近圓形等高線就可以逼近真實邊緣，可以準(zhǔn)確識別并提取撞擊坑的地理位置、大小和深度等信息。

圖2 簡單撞擊坑和復(fù)雜撞擊坑剖面圖
Fig.2 Impact crater structure

圖3 等高線示意圖
Fig.3 Contour diagram

圖4 撞擊坑等高線示意圖
Fig.4 Crater contour diagram

2.2 指標(biāo)體系和參數(shù)的確定

撞擊坑等高線往往表現(xiàn)為相互嵌套的類圓環(huán)狀。因此，初步判別撞擊坑特征的思路是根據(jù)等高線的圓度、矩形度、體態(tài)比等指標(biāo)參數(shù)保留滿足條件的等高線，再分析等高線之間的關(guān)系，聚合嵌套的等高線并標(biāo)示外層等高線；然后通過剖面線的凹凸剔除穹隆等正向地形區(qū)域，保留負(fù)向地形的撞擊坑區(qū)域；最后通過霍夫變換擬合撞擊坑的邊緣，得到撞擊坑的中心位置和大小。因此，指標(biāo)體系和參數(shù)包括圓度、矩形度、體態(tài)比、剖面線凹凸度。

(1)圓形度是圖像處理中重要的概念，用于特征提取與描述，其值反映被測量邊界的復(fù)雜程度，公式為

e=4πA0/P2,

(1)

其中，A0為面積；P為周長。這個特征對圓形面域取最大值1，坑唇邊緣越接近1的等高線越符合條件。

(2)矩形度是常用的目標(biāo)物體矩形擬合因子，其公式為

e1=A0/AR，

(2)

其中，A0為面積；AR為最小外接矩形的面積；e1表示某物體對其最小外接矩形的充滿程度。對于矩形物體，e1取最大值1；對于圓形物體，e1取0.79；對于其他細(xì)長物體，e1取值變小。矩形擬合因子的值限定在0～1之間。

(3)體態(tài)比定義為

e2=W/L，

(3)

其中，W為最小外接矩形的寬；L為最小外接矩形的長。圓形面域的體態(tài)比取值為1，這個參數(shù)可以把細(xì)長不規(guī)則的等高線剔除。

(4)剖面線凹凸性：因為月表還存在其他隆起的地形，如穹隆也會產(chǎn)生近圓形的閉合等高線。潛在撞擊坑還存在，對此情況的處理方法是判斷潛在撞擊坑所在的數(shù)字高程模型剖面線是否為下凹形態(tài)，如果剖面線為下凹形態(tài)，則是撞擊坑；若剖面線為上凸形態(tài)，則不是撞擊坑。根據(jù)內(nèi)等高線與等高線的數(shù)值大小關(guān)系確定是否為真實撞擊坑。

2.3 撞擊坑邊界的擬合

對滿足條件的等高線，利用霍夫變化的方法擬合其邊界，最終得到中心點和半徑大小。

霍夫變換算法是圖像處理中識別幾何圖形以及其特征的一種常用方法?；舴蜃儞Q的思想是：對二值圖像I上任意一個邊緣點I(x,y)，如果它在以(a,b)為圓心，r為半徑的圓上，則滿足：

(x-a)2+(y-b)2=r2.

(4)

從(4)式可以看出，每個邊緣點I(x,y)映射成空間 (a,b,r)上一個二次曲面。將(4)式寫成參數(shù) 形式：

a=x-rcosθ,b=y-rsinθ，

(5)

其中，θ∈[0, 2π)。設(shè)待檢測圓的半徑r∈[R1,R2]，則按照霍夫變換，對任意一個邊緣點I(x,y)，θ，r分別以Δθ，Δr的步長遍歷[0, 2π)和[R1,R2]，得到參數(shù)空間(a,b,r)的子空間，并對相應(yīng)的累加器單元AC(a,b,r)加1。對圖像上所有的點遍歷之后，得到累加器數(shù)組AC，其中任意一點AC(a,b,r)表示在圓(a,b,r)上的邊緣點的數(shù)目。相同數(shù)目的邊緣點在不同半徑的圓上被認(rèn)為是圓的可能性不同，半徑越小，可能性越大。因此定義圓度的概念：

p=n/r，

(6)

其中，n為圓的邊緣點的數(shù)目 ；r為圓的半徑；p為圓度。將累加器數(shù)組AC中的元素按(6)式計算圓度，得到三維圓度數(shù)組C。圓形目標(biāo)對應(yīng)于圓度數(shù)組的局部峰值，即對應(yīng)的C(a,b,r)為局部最大值。

3 結(jié) 果

3.1 識別流程

在嫦娥二號7 m全月數(shù)字高程模型數(shù)據(jù)中隨機選擇G227區(qū)、G228區(qū)進行實驗，基于ENVI和IDL，按照章節(jié)2中的原理和方法，編寫程序進行識別提取，如圖5。首先，設(shè)置合理的等高距在嫦娥二號數(shù)字高程模型數(shù)據(jù)中提取等高線；然后，通過圓形度、矩形度、體態(tài)比、剖面凹凸性篩選出初選撞擊坑區(qū)域；最后，通過霍夫變換提取撞擊坑邊界信息。將提取的信息導(dǎo)入ArcGIS，確定撞擊坑的地理位置和半徑大小，并通過高程差計算得到撞擊坑的深度信息，生成具有位置、大小、深度的小型撞擊坑數(shù)據(jù)庫。

圖5 IDL程序計算演示
Fig.5 IDL program diagram

3.2 識別結(jié)果

基于數(shù)字高程模型數(shù)據(jù)等高線自動識別結(jié)果如圖6，具體統(tǒng)計結(jié)果見表1。

本文的方法在G227區(qū)、G228區(qū)分別識別撞擊坑728個和989個。為了評價本文方法的精度，使用Shufelt[20]提出的算法評價因子：檢測百分比為(Detection percentage)D=100TP/(TP+FN)，分支系數(shù)為B=FP/TP，質(zhì)量百分比為Q=100TP/(TP+FP+FN)。TP表示正確識別的撞擊坑數(shù)量；FN表示遺漏識別的撞擊坑數(shù)量；FP表示錯誤識別的撞擊坑數(shù)量。目前最好的評價方法是采用人工識別的結(jié)果作為對比參照。

圖6 撞擊坑識別結(jié)果示意圖
Fig.6 The crater detection result diagram

表1撞擊坑識別結(jié)果統(tǒng)計表
Table1Craterdetectionresultstatistics

TP/個FN/個FP/個D/%B/%Q/%G227697533192.94.489.2G228918817191.97.781

從表1可知，本文方法在兩個實驗區(qū)的評價因子D分別為92.9%、91.9%，均值為92.4%，B分別為4.4%、7.7%，均值為6.1%，Q分別為89.2%、81%，均值為85.1%。從評價因子可以看出，本文的方法對于小型撞擊坑有非常好的檢測結(jié)果，正確識別率高，且錯分率低，存在小部分遺漏識別的情況。

3.3 討 論

本文提出一種基于高精度數(shù)字高程模型數(shù)據(jù)的派生數(shù)據(jù)等高線，通過計算機自動識別提取小型撞擊坑的方法。相對于傳統(tǒng)的逐點分析坡度指數(shù)，提取圖像邊緣轉(zhuǎn)換成對地形線的分析，可以有效避免地形皺褶和碎石群的干擾，識別流程更為簡單高效，識別速度較快。設(shè)置不同的等高距，分批次分區(qū)間識別不同尺度的撞擊坑，其正確識別率高，遺漏識別率低，說明該方法是可行有效的。不過需要注意的是，對于部分退化嚴(yán)重的撞擊坑，坑緣破壞嚴(yán)重，坑壁平緩，等高線不能明顯表達和準(zhǔn)確勾畫出撞擊坑的形貌輪廓，無法被準(zhǔn)確識別，導(dǎo)致漏分情況的發(fā)生。但通過對比同位置的影像數(shù)據(jù)，通過人工識別分析是可以分辨的。因此，今后的研究方向和重點在于如何識別退化嚴(yán)重的撞擊坑，這也是所有撞擊坑識別方法中的薄弱環(huán)節(jié)和有待解決的問題。在實際應(yīng)用中，可用該方法的識別結(jié)果補充現(xiàn)有月球撞擊坑數(shù)據(jù)庫中的小尺度撞擊坑，豐富現(xiàn)有的月球撞擊坑數(shù)據(jù)庫。