劉 海,岳云鵬,韓 峰,孟 旭,周 斌,方廣有
(1. 廣州大學土木工程學院, 廣東廣州 510006; 2. 廈門大學電子科學與技術學院, 福建廈門 361005;3. 中國科學院空天信息創(chuàng)新研究院, 北京 100094)
月球探測是人類進行太空探測的開端,太空探測可加深人類對地球、月球和太陽的認識,帶動一系列基礎科學的創(chuàng)新[1-2]。在月球探測技術中,探月雷達是其中一種重要方法,也是探地雷達的一種。由于月球表面沒有液態(tài)水,介質的電阻率很高,電磁波穿透能力強,雷達可以探測到月面以下數百米深的次表層[3]。雷達不僅能夠實現(xiàn)月表的巡視探測,而且能在地球基地和繞月衛(wèi)星上對月球表面進行大范圍遙感測繪[4]。因此,在月球、火星等其他外星球探測中,雷達已經成為一種用于淺表層和次表層地質結構探測的有效工具,在現(xiàn)階段的月球探測中發(fā)揮著重要作用[1]。
我國發(fā)射的嫦娥五號探測器已于2020年12月31日于月球表面著陸,嫦娥五號探測器負責嫦娥三期工程“采樣返回”任務,為我國實施載人登月做技術儲備[5-6]。嫦娥五號月壤結構雷達探測儀采用非屏蔽的高頻天線,其工作環(huán)境相比嫦娥三號更加復雜,天線耦合、著陸器地板和其他金屬構件的強反射以及底板與地面之間的多次反射信號等都會對月面以下目標的反射信號形成強烈干擾,甚至淹沒目標信號[7]。根據嫦娥五號的任務安排,在返回雷達數據后,需要快速地對復雜探月雷達信號進行數據處理并實現(xiàn)高分辨率成像,及時為月壤鉆取任務探測鉆頭下方區(qū)域2 m深度范圍內的月壤結構和月巖分布。因此,應進一步研究數據處理方法,抑制雜波和多次波的干擾,才能保證數據解譯的準確性;并進一步研究適合強電磁干擾環(huán)境下的探月雷達快速高精度成像算法,對算法保證計算效率的同時,兼顧成像精度。
基于此,本文開展月壤分層結構的建模、電磁響應計算與分析以及嫦娥五號探月雷達數據的快速高精度成像研究,并通過地面試驗進行驗證,為月壤、月巖物理參數的準確分析提供基礎,為嫦娥五號探測器的月壤取樣提供科學依據,也為進一步的月球以及深空探測提供方法和技術。
嫦娥五號著陸器上搭載一套由12個Vivaldi天線組成的多發(fā)多收探月雷達系統(tǒng),如圖1所示。該雷達系統(tǒng)由中國科學院空天信息研究院團隊研制,其中心頻率為2 GHz。嫦娥五號探月雷達將在靜止的狀態(tài)下對著陸點下2 m深度范圍內的月壤結構進行精細探測,為后續(xù)鉆取采樣過程的順利實施提供月壤分層結構及可能存在的月巖等信息,并通過分析現(xiàn)場探測數據與返回樣品實驗室分析數據之間的聯(lián)系,為研究月壤形成的原因和演化機理提供數據支撐。
圖1 嫦娥五號月壤結構探測儀立體結構示意圖(藍色區(qū)域表示雷達探測儀天線和月壤鉆頭的布局)
由于嫦娥五號探測器所搭載的載荷眾多、系統(tǒng)復雜,內部電子器件和金屬構件對探月雷達數據干擾嚴重。此外,由于著陸器整體優(yōu)化的需要,配置給雷達系統(tǒng)的空間有限,因此天線設計采用尺寸較小的空氣耦合Vivaldi天線,天線布局方式如圖2所示。由于雷達天線距月球表面上方僅有90 cm(見圖3),月表以下的目標反射信號受到天線耦合信號、著陸器金屬底板和其他構件反射信號,以及地面多次反射信號等的嚴重干擾,信噪比低。
(a) 立面圖
圖3 嫦娥五號雷達探測儀的復雜干擾環(huán)境示意圖
考慮到嫦娥五號著陸器底板和其他金屬構件等對月壤結構反射信號的嚴重干擾,有必要對復雜環(huán)境下探月雷達的數據處理方法進行優(yōu)化,研究適合嫦娥五號月壤結構雷達探測儀的數據處理流程,在復雜干擾環(huán)境中提取有效的弱目標信號,為月表實測雷達數據的正確解譯提供依據。此外,根據嫦娥五號的任務安排,從探月雷達數據傳回地球到完成月壤采樣只有短短3個小時,給探月雷達數據的快速處理與成像提出了很大挑戰(zhàn)。因此,研究探月雷達數據處理和快速成像方法,是成功實現(xiàn)嫦娥五號月壤結構雷達探測儀科學探測任務的關鍵問題之一。
由于嫦娥五號月壤結構探測儀的數據成像是通過發(fā)射天線輻射波場和接收天線的逆時波場的互相關進行成像。在實際仿真模擬過程中,如果將天線甚至整個著陸器模型納入到仿真模型中,將大大增大計算量,計算時間無法承受。實際計算過程需要在仿真模型中將天線近似成點源,這需要對系統(tǒng)進行精確校準滿足近似條件。這些系統(tǒng)校準包括子波測試、點掃描組件及電纜延時校正、天線相位中心校準和背景干擾去除。
子波是指由安裝在著陸器上的雷達天線實際輻射到月壤介質中的電磁波波形,并用于模擬發(fā)射天線的輻射波場。本文采用如圖4所示的實驗方法采集各個發(fā)射天線的子波,將整個著陸器在實驗室內吊高至離地面約3 m高的位置,在地表放置吸波材料,然后在著陸器鉆頭的正下方架設一個接收天線,采集不同發(fā)射天線向下輻射的電磁波。探月雷達系統(tǒng)的天線陣列單元采用的Vivaldi天線,Vivaldi是行波天線。傳輸的激勵信號從饋點沿著天線元件傳播時會被延遲。在仿真計算中,應將天線簡化為理想源點,該源點可看作球面輻射場的起源,通常該點被稱作天線視在相位中心。
圖4 天線子波測試實驗示意圖
月壤結構探測儀的天線陣列搭載在著陸器的底板上,而且著陸器大部分組件都是由金屬構成,著陸器底板和其他金屬組件將對天線輻射的電磁波產生很強的干擾,若不對這些背景干擾進行處理,來自月壤內部的目標回波信號將被完全淹沒在這些背景干擾中,很可能造成成像失敗。因此必須預先采集這些背景干擾,在數據處理過程中進行去除。在實驗室內采用的背景干擾去除實驗方法是將著陸器放置在探測區(qū)域的上方,采集包含著陸器的背景干擾、環(huán)境干擾和目標回波信號的數據,然后在地表放置一層20 cm厚的吸波材料,再次采集背景和環(huán)境干擾信號。兩次數據相減,能夠抑制背景和環(huán)境干擾信號,從而增強目標回波信號,提高信號干擾比。
在實際月壤表面環(huán)境中,通過實驗室方式在月表鋪設吸波材料測量背景干擾的方法是不現(xiàn)實的,因此將搭載月壤結構探測儀的著陸器放置在航天五院大型微波暗室以及在繞月飛行階段開機采集背景干擾數據,儲存后進行實測數據的背景去除處理。
系統(tǒng)校準后,通過速度譜分析法對嫦娥五號探月雷達探測器所采集的多偏移距數據進行分析,獲取月壤分層結構不同月壤層的介電常數和厚度[8-9]。速度譜分析方法的基本原理如圖5所示。圖5(a)為地下單個水平反射界面進行多偏移距測量得到探月雷達剖面示意圖,通過給定的試速度生成一系列的雙曲線,然后沿著這些雙曲線疊加不同偏移距通道的信號能量,可得到試速度與零偏移距雙程走時的速度譜如圖5(b)所示。由圖5(b)可知:如果選擇的速度非常接近真實的速度,則疊加的能量將在速度譜上有一個最大的峰值,可以選擇相應的速度和零偏移距雙程走時,并認為其為地下介質的電磁波波速和電磁波到達水平反射界面的雙程走時,從而計算出水平反射界面的深度。
(a) 多偏移距探月雷達剖面圖
如果地下介質中含有多個水平反射界面,在共中心點剖面上會出現(xiàn)多個反射波,則共中心點速度譜上出現(xiàn)多個能量峰值,根據速度譜上能量峰值提取的速度是相應的疊加速度。當收發(fā)天線距離很小時,疊加速度近似等于均方根速度,為了計算地下各層介質的層速度,可應用Dix公式將均方根速度轉化為層速度。Dix公式(1)中第i層層間速度vint,i的計算表達式為
(1)
式中,vst,i+1為第i+1層的均方根速度,ti+1為從地表到i+1層的雙程走時,vst,i為第i層的均方根速度,ti為地表到i層的雙程走時。層速度是與地層巖性密切相關的速度,在各種速度中占有十分重要的地位,計算了各層的層速度之后,結合各層的雙程到達時間,則各層的厚度都可計算出來,第i層介質的厚度為
di=vint,i(ti+1-ti)
(2)
由于嫦娥五號的月壤結構探測儀固定在著陸器上,采用靜止探測的工作方式,僅由其所裝配的12個Vivaldi天線采集132道雷達數據,地下目標回波信息量很少,采用常規(guī)的基于路徑追蹤原理的常規(guī)雷達成像算法很難滿足高分辨率成像需求。逆時偏移技術是近十年發(fā)展起來的一種最新的基于波動方程的偏移成像技術,已經在地震勘探領域取得了巨大的成功[10-11]。將逆時偏移技術應用于雷達成像,是基于以下的電磁波波動方程:
(3)
相對其他偏移方法,逆時偏移算法采用雙程電磁波動方程對波場進行延拓避免了對波動方程的近似。它沒有傾角限制,并可以利用轉換波、棱鏡波或多次反射波進行成像,獲得更精確的振幅等動力學信息,實現(xiàn)保幅成像,還可以更好地對復雜速度場進行更細化更精確的估計。逆時偏移的主要流程包括激勵源波場的正演模擬、接收波場的逆時外推和成像條件的應用,其基本流程如圖6所示。
圖6 逆時偏移算法的基本流程
接收波場的逆時外推是波的逆時傳播過程,沿時間軸的負方向遞推,每完成一個時間步長的外推,則根據一定的成像條件,將該時刻的波場與預先存儲的正演模擬中對應時刻的波場進行一次成像,得到的成像值累加到成像剖面,直到波場外推過程完成就可以得到最終的偏移成像剖面。常用的成像條件為疊前逆時偏移互相關成像條件,其數學表達式為
(4)
式中,S為激勵源波場,R為接收波場,t和Txi分別表示時間和炮數。這種成像條件實際是利用發(fā)射波場和接收波場的互相關。
由于嫦娥五號工作任務要求在3個小時內完成鉆取采樣任務,希望探月雷達數據處理與成像的時間控制在分鐘量級。為提高計算效率,課題組提出了一種分層介質的格林函數的頻域逆時偏移算法[12],選取了零時滯互相關通常用作成像條件:
(5)
通過計算探月雷達信號的因果關系,公式(5)可以按照文獻[13]的卷積形式進行重寫:
(6)
由于時域中的卷積與頻域中頻譜的相乘相對應,方程(6)可以寫成:
(7)
(8)
(9)
分層介質格林函數是成像域中每個網格的每個頻率的二階張量。探月雷達通常在垂直寬邊模式下工作,并使用線性極化天線進行發(fā)送和接收。這樣,探月雷達僅記錄電場的一個極化分量,則
(10)
圖7 分層介質格林函數的頻域逆時偏移成像算法流程
為驗證速度譜分析法可使嫦娥五號能夠準確獲取月壤層狀介質信息,通過正演GPRMax2.0正演模擬軟件建立了如圖8(a)所示的仿真模型,采用時域有限差分方法(FDTD)模擬電磁場的傳播。該分層模型包含三層月壤介質,相對介電常數從上至下分別為2,2.5和3,電導率均為0.000 01 S/m。天線簡化為點源進行模擬,共12個,其空間布局和高度都設置成與嫦娥五號探測器一致。3個方向的網格的尺寸均為7.5 mm,饋電脈沖采用中心頻率為2 GHz的Ricker子波。
將正演模擬后得到的探月雷達數據仿真結果進行天線延時矯正、增益、道間均衡等預處理,預處理后得到的結果如圖8(b)所示,圖中展示了 1-11號天線作為發(fā)射天線時,所采集的110道信號,可清楚地分辨月壤表面、第一層月壤界面和第二層月壤界面的雷達反射信號,隨著收發(fā)天線偏移距的增加,反射信號的雙程走勢沿雙曲線增長。
由該正演數據計算得到的速度譜如圖8(c)所示,可以提取月壤表面、第一層月壤界面和第二層月壤界面的雷達反射信號的疊加速度及雙程零偏移距走時,根據公式(1)、(2),可計算得到月表上方空氣層、第一層月壤介質和第二層月壤介質的電磁波速度及厚度,如表1所示。與其真實值比較接近,最大誤差僅為11.2%。
圖8 月壤分層結構的正演模擬及速度反演
表1 月壤分層結構的介電常數及厚度估算結果
為驗證筆者課題組所提出的頻域逆時偏移算法準確性,本文使用多輸出多輸入探月雷達系統(tǒng)進行二維數值仿真實驗進行分析,并將計算結果與基于時域有限差分的時域逆時偏移算法進行比較。建立月壤中孤立月巖模型,如圖9(a)所示,地下層(風化石)和圓柱形巖石物體的相對介電常數分別設定為2.5和5,對應電導率分別設定為10-5S/m和10-4S/m。為與嫦娥五號月球探測雷達系統(tǒng)中使用的天線陣列相類似,本模型同樣輪流使用12個天線作為發(fā)射天線,剩余的11個天線作為接收天線,從而可得到132道合成探月雷達數據。使用中心頻率為2 GHz的雷克子波作為源激勵,且將接收天線的接收信號中的直達波信號濾除。
二維成像區(qū)域的尺寸為1.6 m×1.8 m。仿真圖像區(qū)域中的場點的空間采樣間隔在x和z方向均為4 mm。對于頻域逆時偏移成像,源頻譜和接收頻譜均以20 MHz的頻率步進從20 MHz到4 GHz進行采樣。
圖9(b)、(c)分別為時域逆時偏移和頻域逆時偏移使用半空間初始模型重建的二維圖像結果。除源位置處的存在相位差之外,兩個重建圖像之間沒有區(qū)別,導致源位置處的相位差的原因是時域有限差分方法不能精確地模擬非常接近源點的區(qū)域中的電場。被掩埋的巖石和地表都可以清晰地成像,巖石的頂部和底部界面可被準確識別,巖石的頂部位置大致在深0.3 m處與實際的位置相當吻合。相比時域逆時偏移算法,頻域逆時偏移算法消耗的內存不到五十分之一,計算時間成本只有百分之一。由于可以提前計算分層介質格林函數,因此可以在快速地根據嫦娥五號采集的探月雷達數據生成二維地表圖像,滿足嫦娥五號的任務要求。
圖9 時域逆時偏移和頻域逆時偏移效果對比
本文通過建立火山坑地面模型試驗對提出的頻域逆時偏移成像算法進行驗證,如圖10所示,實驗室采用7 m×3 m×2.5 m的火山灰填充坑模擬月壤,火山灰的相對介電常數約為2.5,其介電性能與月壤比較接近,電導率可忽略不計。選取深1 m處的大理石板模擬月巖結構,放置的大理石板的尺寸為60 cm×60 cm×3 cm,相對介電常數大約為8,在灰坑的底部(2.5 m深處),放置整塊金屬板。
(a) 嫦娥五號月壤結構探測著陸器的實際尺寸模型
裝備實際尺寸模型的嫦娥五號著陸器放置在火山坑表面以檢測地下目標,并對實測數據進行三維逆時偏移成像處理。首先根據估計的表面反射信號到達時間,從接收記錄的132道探月雷達數據中去除直達波。頻域逆時偏移和時域逆時偏移算法均采用半空間速度模型。由于本文只對比地下區(qū)域的成像結果,因此頻域逆時偏移算法不考慮計算地面上的分層介質格林函數。相反,由于時域有限差分方法必須從源點外推波場,時域逆時偏移方法須考慮地上部分。地下區(qū)域的空間采樣步長均為4 mm,共有440×220×340個網格。源頻譜和接收頻譜采用20 MHz的頻率步進從20 MHz至4 GHz采樣。
由兩種逆時偏移算法獲得的三維地下成像結果如圖11所示,由于兩種算法成像結果的數值僅代表對比度的不同,因此對地下區(qū)域的數值進行了歸一化處理,通過不同角度剖面所得圖像對比可知,兩種算法重建的三維圖像基本沒有區(qū)別,且兩種算法的成像結果的位置尺寸與實際情況相當吻合。
(a) 時域逆時偏移
從圖11可以看出,其中y方向的成像分辨率比x方向的成像分辨率稍差,這是由于在y方向上部署的只有兩個Vivaldi天線所產生的孔徑有限。表2分別給出了頻域逆時偏移和時域逆時偏移的計算成本。三維時域有限差分模擬計算成本相對較高,其消耗了較大的存儲空間和近三天的計算時間。由于分層介質格林函數可以預先計算并保存在磁盤中,因此實際頻域逆時偏移算法的三維成像所需的計算時間約為1 h,可以滿足嫦娥五號任務的要求。
表2 時域逆時偏移與頻域逆時偏移計算成本比較
本文提出了一種適合嫦娥五號探月雷達的數據處理方法,并實現(xiàn)月壤結構的三維高精度反演成像。通過對嫦娥五號月壤結構探測儀系統(tǒng)進行子波測試、天線相位中心校準和背景干擾去除等校正,有效地去除了嫦娥五號月壤結構雷達探測儀在探測過程的干擾雜波;利用速度譜分析方法獲取了月壤分層結構模型的介電常數和厚度,計算值與真實值最大計算誤差僅為11.2%,為逆時偏移成像提供了初始模型。地面驗證試驗結果表明:同等精度的頻域逆時偏移算法在對層狀結構進行成像時,計算時間僅為時域逆時偏移的1.5%,占用內存僅為其22.2%,計算效率有顯著的優(yōu)越性。綜上所述,本文所提的頻域逆時偏移算法可滿足嫦娥五號搭載的探月雷達系統(tǒng)在分鐘量級時間內對著陸器下方的月球表面可能的巖石進行高精度圖像的要求,為嫦娥五號鉆取采樣任務提供了重要的信息支持,并可為未來的月球探測任務提供技術方法支撐。