基于CLEAN算法的嫦娥四號低頻射電頻譜儀信號干擾抑制

劉晨迪，周建鋒，蘇 彥

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

嫦娥四號低頻射電頻譜儀放置在月球背面，由于月球屏蔽了大量來自地球的無線電干擾，因此具有得天獨厚的觀測條件，可以對太陽等宇宙天體的射電信號進(jìn)行觀測[1]。但是，由于著陸器平臺干擾信號比接收的低頻射電天文信號大很多，如果不分離平臺干擾信號，將無法得到有效的低頻射電天文信號數(shù)據(jù)。因此，本文提出一種有效分離平臺干擾信號的算法，實現(xiàn)對干擾的抑制。

在射電天文領(lǐng)域，傳統(tǒng)抑制干擾的方法有很多，包括對空白脈沖干擾的時間限制[2]，自適應(yīng)頻率選擇濾波[3]，參數(shù)信號相減[4]，自適應(yīng)波束成形[5]，時間自適應(yīng)濾波[6]，帶有參考天線的頻域后相關(guān)處理[7]，子空間投影的空間濾波[8]等等。截至目前，應(yīng)用比較廣泛的是單通道總功率變化檢測器。文[9]基于對功率變化的檢測，在RATAN600上實現(xiàn)了改進(jìn)的功率檢測器。文[10]提出在所有滯后的情形下使用量化相關(guān)性來測試干擾的存在。文[11]提出了一種基于小波分解的檢測器。這些都是單通道檢測器，它們不利于干擾的空間特性分析。文[12]提出了結(jié)合使用多個望遠(yuǎn)鏡以改善探測質(zhì)量的方法，此方法用于低頻干涉測量，文中還提出了一種基于交叉頻譜時間行為的魯棒數(shù)據(jù)檢測方法。文[6]考慮自適應(yīng)過濾技術(shù)，提出使用參考天線和最小均方(Least Mean Square, LMS)類型算法，消除綠岸望遠(yuǎn)鏡(Green Bank Telescope, GBT)的干擾[6]。文[13]提出一種針對韋斯特博克綜合射電望遠(yuǎn)鏡的射電干擾消除子系統(tǒng)。這些都是在射電干擾消除方面典型的例子。

射電干擾信號消除算法主要分為3類：(1)通過線性方式區(qū)分射電干擾信號的特點，比如奇異值分解(Singular Value Decomposition, SVD)[14]、主成分分析(Principal Component Analysis, PCA)[15]等等。如果射電干擾在時間和頻率上表現(xiàn)出重復(fù)模式，這些方法效果很好，但無法處理更多隨機信號，例如由衛(wèi)星引起的不規(guī)則信號。(2)基于閾值的算法，例如CUMSUM[8]和SUMTHRESHOLD[4]，這些算法簡單且可靠有效。(3)使用傳統(tǒng)的監(jiān)督機器學(xué)習(xí)技術(shù)(例如K最近鄰和高斯混合模型)對射電干擾信號進(jìn)行聚類[16]。

嫦娥四號低頻射電頻譜儀接收的平臺干擾信號很強，和太陽爆發(fā)峰值的量級10-15W/(m2·Hz)相當(dāng)，比普通的太陽爆發(fā)信號10-18W/(m2·Hz)大3個量級，比太陽寧靜時的信號10-22W/(m2·Hz)大7個量級。之前儀器設(shè)計時推薦使用的譜減法[17]，雖然可以將信號的信噪比提高約30 dB，但是仍然不足以探測有效的射電天文信號。為此，本文提出了基于CLEAN算法的平臺干擾信號分離方法，利用低頻射電頻譜儀3個天線同時觀測的信號有相關(guān)性的特點，并結(jié)合天線之間數(shù)據(jù)的互相關(guān)功率譜、傅里葉級數(shù)等工具，有效地分離平臺干擾信號，從而實現(xiàn)對干擾的抑制。

1 低頻射電頻譜儀的工作原理及接收信號的構(gòu)成

1.1 低頻射電頻譜儀簡介

宇宙中的來波信號是一個矢量信號。根據(jù)矢量場理論，來波信號可以分解為任意3個互相垂直的分量，只要將這3個分量接收并合成，就可以得到矢量的大小和方向。通過后期處理，我們可以獲得電場的頻譜和時變信息。

低頻射電頻譜儀的天線位置分布如圖1。3個分量有源天線接收宇宙的來波電場。3個天線分別標(biāo)注為天線A、天線B和天線C，兩兩正交，還有一個短天線D，主要用于接收平臺干擾信號,為天線數(shù)據(jù)處理做參考。

圖1 著陸器低頻射電頻譜儀的天線位置分布

低頻射電頻譜儀主要由4個接收天線、前置放大器、電子學(xué)單元、電纜組件等部分組成。有源天線安裝于衛(wèi)星本體外，其余組件安裝于衛(wèi)星本體內(nèi)。電子學(xué)單元由控制器、配電器、基準(zhǔn)時鐘模塊、四通道接收機、內(nèi)定標(biāo)模塊、存儲單元、通訊接口、數(shù)傳接口等組成。低頻射電頻譜儀系統(tǒng)組成及內(nèi)部各單元與數(shù)管分系統(tǒng)、供配電分系統(tǒng)之間的接口關(guān)系如圖2。

圖2 低頻射電頻譜儀系統(tǒng)組成及接口關(guān)系示意圖

低頻射電頻譜儀有6種工作模式，分別為(1)待機模式；(2)內(nèi)定標(biāo)模式；(3)時頻對比模式；(4)頻譜巡查模式；(5)三天線模式；(6)時域模式。其中應(yīng)用最廣泛的是三天線工作模式，在此模式下，低頻射電頻譜儀分別采集4個天線4個通道的科學(xué)數(shù)據(jù)。本文對三天線模式下的數(shù)據(jù)進(jìn)行處理。

低頻射電頻譜儀的科學(xué)數(shù)據(jù)分為0A級、0B級、01級、2A級、2B級和2C級共6級。預(yù)處理數(shù)據(jù)經(jīng)過分路解幀和分包，生成0A級數(shù)據(jù)產(chǎn)品；0A級數(shù)據(jù)經(jīng)過排序、去重復(fù)、兩站優(yōu)化拼接、去源包包頭，根據(jù)采樣周期形成低頻射電頻譜儀數(shù)據(jù)塊，生成0B級數(shù)據(jù)產(chǎn)品；多個0B級數(shù)據(jù)文件拼合成一個探測周期，生成01級數(shù)據(jù)產(chǎn)品；對01級數(shù)據(jù)進(jìn)行傅里葉變化、定標(biāo)等操作，生成2A級數(shù)據(jù)產(chǎn)品；2A級數(shù)據(jù)加入太陽位置等信息，生成2B級數(shù)據(jù)產(chǎn)品；對2B級數(shù)據(jù)進(jìn)行消噪聲處理，生成2C級數(shù)據(jù)產(chǎn)品。其中1級數(shù)據(jù)為時域數(shù)據(jù)，2級數(shù)據(jù)為頻域數(shù)據(jù)。本文所用的原始數(shù)據(jù)均為1級數(shù)據(jù)。

1.2 接收信號的構(gòu)成

低頻射電頻譜儀接收的信號主要由3部分構(gòu)成：(1)平臺干擾信號；(2)接收機噪聲；(3)低頻射電天文信號。

平臺干擾信號的特點是幅度相對較強，兩天線之間的相關(guān)性高。如果干擾源位置確定，各天線接收信號之間的強度比例不隨時間變化。來自電子設(shè)備的干擾信號，一部分是固定頻率的發(fā)射信號，如通訊發(fā)射機、時鐘信號等，另一部分是寬帶隨機噪聲信號，這是由電磁兼容(Electro Magnetic Compatibility, EMC)問題帶來的。最終的平臺干擾信號既有確定變化的成分，也有寬帶隨機噪聲信號。

低頻射電天文信號具有一定的相關(guān)性，除此之外，還有一些天體源特有的性質(zhì)。比如，低頻太陽射電信號有瞬時和頻譜漂移等特點，具有很強的頻譜可識別特征[18]，如圖3。再比如，快速射電暴(Fast Radio Burst, FRB)FRB200428的持續(xù)時間約0.61 ms，能量通量(Fluence)達(dá)到1.5 MJyms[19]，輻射強度足以被低頻射電頻譜儀探測到。另外，由于儀器指向隨月球的自轉(zhuǎn)而改變，低頻射電天文信號在低頻射電頻譜儀各天線之間的強度比例隨時間呈現(xiàn)有規(guī)律的變化。所有這些特性有助于在強平臺干擾的背景下識別有用的低頻射電天文信號。

圖3 低頻太陽射電信號頻譜

接收機噪聲的特點是天線兩兩之間無相關(guān)性。我們可以利用該特性，通過計算天線之間數(shù)據(jù)的相關(guān)性，評估平臺干擾信號的分離程度。

圖4展示了低頻射電頻譜儀B天線觀測的一個頻譜實例。由圖4可知，頻譜在1 MHz附近的流量密度大約為10-15W/(m2·Hz)，與太陽爆發(fā)峰值流量密度相當(dāng)，比太陽爆發(fā)流量密度大2～3個量級，說明著陸器平臺干擾比我們希望探測的太陽爆發(fā)信號大許多。如果不通過一定的方法分離平臺干擾信號，將無法有效提取低頻射電天文信號。

圖4 天線B采集的高頻波段的一個頻譜實例，對應(yīng)第23月晝第2道數(shù)據(jù)

2 低頻射電頻譜儀平臺干擾信號的分離方法

低頻射電源與平臺干擾源在低頻射電頻譜儀的3個天線上激發(fā)出互相關(guān)聯(lián)的接收信號。本文的目的是從兩組信號的相關(guān)性出發(fā)，借助CLEAN算法以及互相關(guān)功率譜、傅里葉級數(shù)等工具，把低頻射電頻譜儀觀測信號中的強相關(guān)成分與部分相關(guān)成分分離，以便進(jìn)行后續(xù)的科學(xué)分析。分離之后，原始的低頻射電頻譜儀時域信號被切分為兩部分，即強相關(guān)的CLEAN模型信號和部分相關(guān)的殘余信號。在這個過程中，原始信號的信息沒有丟失或損傷。

CLEAN模型信號主要由平臺干擾信號組成。當(dāng)?shù)皖l強射電爆發(fā)信號到來時，進(jìn)入CLEAN模型信號。由于數(shù)據(jù)傳輸?shù)南拗?，低頻射電頻譜儀大約以1 s為間隔采集一道時域數(shù)據(jù)，采集間隔并不嚴(yán)格相等。每道數(shù)據(jù)由4 096個數(shù)據(jù)組成，采樣頻率為100 MHz，采樣持續(xù)時間大約40 μs。因此，不同道時域數(shù)據(jù)之間有明顯的差別。

為了尋找有效的低頻強射電爆發(fā)信號，我們需要提取每道時域數(shù)據(jù)的CLEAN模型信號，并對其中的確定成分和寬帶隨機成分進(jìn)行處理。例如，對于CLEAN模型信號中的確定時域變化結(jié)構(gòu)或頻譜成分，我們把它們當(dāng)作系統(tǒng)誤差，借助CLEAN算法扣除，剩下的CLEAN模型信號，由寬帶隨機信號組成，通過平均多道數(shù)據(jù)的方式降低漲落水平，提高探測低頻強射電天文信號的靈敏度。

殘余信號主要由不相關(guān)的接收機噪聲、未分離的平臺干擾信號以及常規(guī)低頻射電天文信號組成。殘余信號的特征類似于隨機白噪聲，沒有確定的信號成分，可以借助數(shù)理統(tǒng)計工具進(jìn)行分析。例如，借鑒微波輻射計的思想，通過積分平均的方式進(jìn)一步提高殘余信號的探測靈敏度。

我們的最終目標(biāo)是從嫦娥四號低頻射電頻譜儀的觀測數(shù)據(jù)中搜索低頻射電天文信號，并進(jìn)行深入的科學(xué)研究。此時，無論是CLEAN模型信號還是殘余信號都缺一不可。

時標(biāo)短、幅度強的低頻強射電爆發(fā)往往隱藏在CLEAN模型信號中，對這類射電爆發(fā)信號的提取離不開寧靜狀態(tài)的CLEAN模型信號的分類整理。特別強的低頻射電天文爆發(fā)，在單道數(shù)據(jù)中有明顯異常的表現(xiàn)。此時，可以借助交叉相關(guān)直接扣除該道數(shù)據(jù)對應(yīng)的寧靜狀態(tài)的CLEAN模型信號，獲得射電爆發(fā)信號的時域、頻譜特征。如果射電爆發(fā)強度與平臺干擾中的寬帶隨機成分接近，先扣除CLEAN模型信號中的確定成分，再通過多道數(shù)據(jù)積分平均的方式提高探測靈敏度，最終獲得低頻射電天文爆發(fā)信號的時域、頻譜特征。

一些時標(biāo)長(超過月量級)、幅度弱的低頻射電天文信號隱藏在殘余信號中。為了探測這類弱信號，首先借助星載定標(biāo)裝置或其他外部穩(wěn)定的定標(biāo)源，確認(rèn)接收機工作穩(wěn)定。之后，可以直接從殘余信號出發(fā)，同樣借助微波輻射計的思想，通過積分平均提高殘余信號的探測靈敏度。理論上，低頻射電頻譜儀3個天線殘余信號的積分值隨月球自轉(zhuǎn)發(fā)生調(diào)制變化,這樣的調(diào)制信息甚至可以用來對全天進(jìn)行粗略成像。

2.1 CLEAN算法分離平臺干擾信號的基本思路

設(shè)低頻射電頻譜儀中一個天線接收的信號為S(t)，它包含平臺干擾信號I(t)、宇宙射電信號C(t)以及接收機噪聲N(t)。這些信號之間的關(guān)系為

(1)

其中，SA(t)，SB(t)和SC(t)分別為天線A，B和C接收的信號；NA(t)，NB(t)和NC(t)分別為天線A，B和C對應(yīng)接收機的本地噪聲，它們之間相互獨立；αA(t)，αB(t)和αC(t)為平臺干擾信號在天線A，B和C的投射系數(shù)；βA(t)，βB(t)和βC(t)為宇宙射電信號在天線A，B和C的投射系數(shù)。

平臺干擾信號比低頻射電天文信號大3～10個量級，因此，分離平臺干擾是進(jìn)行有效天文觀測的基本前提。

經(jīng)典的CLEAN算法在臟圖或殘余圖中搜索最大值，此處包含部分真實目標(biāo)，減去點擴散函數(shù)(Point Spread Function, PSF)提取這部分真實目標(biāo)，不斷迭代，直到殘余圖中只剩噪聲結(jié)構(gòu)[20]。CLEAN算法同樣可以應(yīng)用到時域信號，比如利用CLEAN算法對時間序列進(jìn)行頻譜分析[21]，識別聲源[22]，以及提取一些特定結(jié)構(gòu)的信號(周期信號、準(zhǔn)周期信號)等[23]。

我們利用CLEAN算法提取低頻射電頻譜儀的平臺干擾信號。算法的輸入數(shù)據(jù)為一組天線的觀測數(shù)據(jù)。首先，計算每個天線的頻譜，以及每對天線的互相關(guān)功率譜，并將所有的互相關(guān)功率譜疊加。由于天線的接收機噪聲N(t)相互獨立，因此互相關(guān)功率譜中的噪聲水平降低，而疊加后互相關(guān)功率譜的信噪比進(jìn)一步增強，有利于提取平臺干擾信號。

疊加的互相關(guān)功率譜包含平臺干擾信號對應(yīng)的頻譜成分，同時包括振幅和相位信息。根據(jù)這些信息，我們借助CLEAN算法和傅里葉級數(shù)展開公式，可以逐一將某個頻率成分對應(yīng)的三角函數(shù)信號從時域數(shù)據(jù)中扣除。考慮到每個天線及接收機的增益不一致，在扣除時每個天線的CLEAN因子是不一樣的。

離散傅里葉變換對應(yīng)的時域及頻域信號都是周期信號。如果信號首尾的值不一致，會在傅里葉變換時出現(xiàn)頻譜泄漏，導(dǎo)致頻率估計出現(xiàn)偏差。通常，實際采集的信號不滿足首尾連續(xù)的條件，一般可以采用加窗函數(shù)避免頻譜泄漏。但是加窗函數(shù)會導(dǎo)致信號扭曲，幅度值發(fā)生變化，假如我們直接采用傅里葉逆變換的方式進(jìn)行信號分離，則導(dǎo)致較大的偏差。為此，我們引入CLEAN算法，通過加窗函數(shù)的方式獲得準(zhǔn)確的頻譜估計；然后在信號分離時，逐步扣除某一個諧波分量，避免信號扭曲，從而把首尾不連續(xù)的干擾信號精確分離。

提取大部分平臺干擾信號后，殘余信號由接收機噪聲主導(dǎo)，此時可以停止CLEAN迭代算法。來自宇宙的射電信號幅度往往小于噪聲水平，加上小部分平臺干擾信號，低頻射電天文信號仍然隱藏在接收機噪聲中。之后，可以通過對時域信號積分平均，提高殘余信號的信噪比，并利用信號隨時間的變化特性進(jìn)一步辨別、分離平臺干擾信號和低頻射電天文信號，從而實現(xiàn)對干擾的抑制。

2.2 CLEAN算法的流程

按照2.1節(jié)的思路，用CLEAN算法分離平臺干擾信號的流程如圖5。為了更簡明清晰地展示算法，我們先考慮兩個天線(分別記為A和B)的情形，并且只考慮連續(xù)的時域信號。在某一道觀測數(shù)據(jù)中，天線A的時域數(shù)據(jù)記為fA(t)，天線B的時域數(shù)據(jù)記為fB(t)。

圖5 CLEAN算法流程圖

相關(guān)的公式為

(2)

(3)

(4)

我們假定天線A和B的時域信號fA(t)和fB(t)同步。若由于線路補償沒有做好，兩個時域信號之間存在時間延遲，那么需要通過交叉相關(guān)的方法計算時間延遲，并對時域信號做延遲修正，然后再用CLEAN算法提取平臺干擾信號。

通過CLEAN算法分離平臺干擾信號之后，得到的殘余信號主要成分是接收機噪聲，以及隱藏在其中的宇宙射電信號，還有少部分沒有扣除的平臺干擾信號。此時兩個天線殘余信號之間的相關(guān)系數(shù)比較低。通過計算信噪比發(fā)現(xiàn)，殘余信號中的平臺干擾信號和接收機噪聲基本在同一量級。這樣可以通過積分平均的方式提高殘余信號的信噪比，這與傳統(tǒng)輻射計的工作原理一致。至于區(qū)分平臺干擾信號和宇宙射電信號，則要借助3個天線信號隨頻譜、時間和來波方向變化的規(guī)律。

3 低頻射電頻譜儀數(shù)據(jù)處理實例

3.1 低頻射電頻譜儀數(shù)據(jù)處理流程

低頻射電頻譜儀原始數(shù)據(jù)由嫦娥四號地面應(yīng)用系統(tǒng)數(shù)據(jù)發(fā)布中心提供。數(shù)據(jù)存儲在嫦娥四號低頻射電頻譜儀載荷相關(guān)文件夾下，按照0A，0B，01，2A，2B和2C進(jìn)行分級，每一級數(shù)據(jù)按照每個月晝一個文件夾進(jìn)行排列。每個月晝通常有14個文件，其中7個是數(shù)據(jù)文件，7個是數(shù)據(jù)采集參數(shù)文件，7個數(shù)據(jù)文件對應(yīng)一個月晝中7次開機時間段所采集的科學(xué)數(shù)據(jù)。我們在研究不同時間段的數(shù)據(jù)時，分別讀取這7個文件，并結(jié)合相應(yīng)的7個數(shù)據(jù)采集參數(shù)文件，確定每個時段的采集參數(shù)。

以第23月晝高頻1級時域數(shù)據(jù)的第1道數(shù)據(jù)為例，對應(yīng)的數(shù)據(jù)文件名為CE4_GRAS_LFRS-TR_SCI_P_20201014020000_20201014140000_0135_B.01，其中GRAS代表地面應(yīng)用系統(tǒng)；SCI代表科學(xué)數(shù)據(jù)；0135代表這個文件是第135個數(shù)據(jù)文件；數(shù)據(jù)結(jié)尾的01擴展名代表數(shù)據(jù)是1級數(shù)據(jù)，為時域數(shù)據(jù)。從文件名可以看出，數(shù)據(jù)采集的時間是從2020年10月14日2點0分0秒到2020年10月14日14點0分0秒，這個時間是世界時，換算為北京時間還需要在小時數(shù)上加8。

在數(shù)據(jù)文件中，數(shù)據(jù)按照道的順序依次排列，每道數(shù)據(jù)由數(shù)據(jù)頭(79字節(jié))和科學(xué)數(shù)據(jù)(8 192字節(jié))組成。其中數(shù)據(jù)頭包含儀器的相關(guān)參數(shù)，科學(xué)數(shù)據(jù)共有4 096個數(shù)據(jù)點，每個數(shù)據(jù)2個字節(jié)，由16位無符號整數(shù)組成，高字節(jié)在前，低字節(jié)在后。讀取后還需要將整數(shù)電壓轉(zhuǎn)化為浮點數(shù)電壓，具體數(shù)據(jù)形式參考內(nèi)部技術(shù)文檔 《探月工程嫦娥四號任務(wù)著陸器低頻射電頻譜儀數(shù)據(jù)預(yù)處理方法設(shè)計》。第1道數(shù)據(jù)一共采集了4 096個點，高頻數(shù)據(jù)的采樣間隔為0.01 μs，低頻數(shù)據(jù)的采樣間隔為0.2 μs。

我們通過Python編程按照道的順序依次讀取數(shù)據(jù)，將每道4 096個點排成一行，生成一個科學(xué)數(shù)據(jù)矩陣。程序提取每行(每道)數(shù)據(jù)，分別用CLEAN算法進(jìn)行處理。

3.2 低頻射電頻譜儀部分處理結(jié)果展示

CLEAN算法處理前的時域數(shù)據(jù)見圖6，由圖6可以看出，天線A和天線B的信號范圍基本接近，均為-0.5～0.5 V，天線C數(shù)據(jù)比天線A和天線B略大。3個天線的信號隨時間變化的結(jié)構(gòu)具有較高的相似性，A和B之間的相關(guān)系數(shù)高達(dá)0.965，B和C之間的相關(guān)系數(shù)為0.823，A和C之間的相關(guān)系數(shù)為0.915。同時，數(shù)據(jù)中存在的震蕩結(jié)構(gòu)主要成分不是接收機噪聲，因為接收機噪聲之間無相關(guān)性。天線A，B和C信號的均值分別為-0.014 9 V，-0.014 3 V和-0.003 2 V，均接近于0。它們的標(biāo)準(zhǔn)差σ分別為0.189 V，0.138 V和0.409 V，存在明顯的差別。這種差別可能是由于3個天線的增益不同，或者平臺干擾在3個天線上的投射系數(shù)不一致，具體原因有待后續(xù)細(xì)致的數(shù)據(jù)分析。從CLEAN算法處理前第2道時域數(shù)據(jù)(圖7)中可以看到，第1道數(shù)據(jù)與第2道數(shù)據(jù)時域信號有較為明顯的差別。

圖6 天線A，B，C和D CLEAN算法處理前的第23月晝第1道原始時域數(shù)據(jù)。(a)天線A時域數(shù)據(jù)；(b)天線B時域數(shù)據(jù)；(c)天線C時域數(shù)據(jù)；(d)天線D時域數(shù)據(jù)

圖7 天線A，B，C和D CLEAN算法處理前的第23月晝第2道原始時域數(shù)據(jù)。(a)天線A時域數(shù)據(jù)；(b)天線B時域數(shù)據(jù)；(c)天線C時域數(shù)據(jù)；(d)天線D時域數(shù)據(jù)

從4個天線的數(shù)據(jù)來看，天線D的響應(yīng)值比其余3個天線小很多，這是由于天線D長度只有20 cm，比天線A，B和C的長度(5 m)小很多造成的。同時，天線D信號與其余3個天線信號在形狀上也有明顯的差別，說明天線D的頻譜響應(yīng)明顯異于其余3個天線，并且這種差異是非線性的。從天線之間的相關(guān)系數(shù)也可以看出，天線A，B，C與天線D的數(shù)據(jù)相關(guān)系數(shù)分別為0.440，0.397和0.298，顯著低于天線A，B和C之間的相關(guān)系數(shù)。因此，我們放棄使用天線D數(shù)據(jù)進(jìn)行平臺干擾信號的分離。

從天線A和天線BCLEAN算法處理前的頻譜(圖8)可以看到，頻譜的幅度譜強度分布在10-1～10-5V之間，低頻成分強度較大，并且頻譜中呈現(xiàn)復(fù)雜的結(jié)構(gòu)。將此時的電壓頻譜轉(zhuǎn)化為流量密度頻譜(圖9)后，兩天線較低頻段的流量密度峰值超過10-16W/(m2·Hz)，比太陽爆發(fā)信號的流量密度10-18W/(m2·Hz)大2個量級，說明直接讀取的時域信號如果不經(jīng)過去除平臺干擾，我們很難發(fā)現(xiàn)宇宙射電爆發(fā)信號。

圖8 天線A和天線B CLEAN算法處理前第23月晝第1道數(shù)據(jù)的頻譜。(a)天線A時域數(shù)據(jù)的頻譜；(b)天線B時域數(shù)據(jù)的頻譜

圖9 天線A和天線B第23月晝第1道數(shù)據(jù)CLEAN算法處理前的流量密度。(a)天線A的流量密度；(b)天線B的流量密度

兩天線在CLEAN算法處理前的互相關(guān)功率譜(圖10)顯示，其中存在一些較為復(fù)雜的結(jié)構(gòu)。此外，與單天線頻譜相比，兩天線互相關(guān)功率譜的信噪比更高。對比圖10和圖8可以發(fā)現(xiàn)，單天線頻譜中，低頻的峰值比高頻的平均值大約高1個數(shù)量級，而在兩天線互相關(guān)功率譜中，低頻的峰值比高頻的平均值大約高2個數(shù)量級。因此，從互相關(guān)功率譜出發(fā)，更有利于提取平臺干擾信號。

圖10 天線A和天線B CLEAN算法處理前的互相關(guān)功率譜

我們進(jìn)行CLEAN算法處理時選取的δ增益因子為0.2，迭代次數(shù)設(shè)置為10萬次。經(jīng)過CLEAN算法處理后，天線A和天線B的殘余信號顯示在圖11中。對比圖8可以看出，CLEAN算法處理前天線時域信號中存在的一些強峰結(jié)構(gòu)已經(jīng)消失。CLEAN算法處理后天線A數(shù)據(jù)的標(biāo)準(zhǔn)差為6.58 × 10-5V，天線B數(shù)據(jù)的標(biāo)準(zhǔn)差為5.36 × 10-5V，可見天線數(shù)據(jù)的波動性已經(jīng)減小。天線信號在不同時間段的幅度近似相等，這與接收機白噪聲本底的特點相符合。

圖11 天線A和天線B CLEAN算法處理后的殘余信號。(a)天線A的殘余信號；(b)天線B的殘余信號

從兩天線CLEAN算法處理后得到的模型信號分布(圖12)可以看出，CLEAN算法處理后的模型信號與初始信號很接近，說明原始信號中大部分是平臺干擾信號。同時，兩天線模型信號的相關(guān)系數(shù)為0.966，接近1，也說明模型信號中大部分是相關(guān)的平臺干擾信號。

圖12 天線A和天線B CLEAN算法處理后的模型信號。(a)天線A的模型信號；(b)天線B的模型信號

CLEAN算法處理后兩天線的頻譜(圖13)幅度峰值由原來約10-1V降低到約10-5V，降低4個數(shù)量級。殘余信號的頻譜也不再有明顯的頻譜結(jié)構(gòu)，在各個頻率分量的強度大致相等，可以看作白噪聲譜，說明天線信號中的強平臺干擾成分得到了非常有效的去除。

圖13 天線A和天線B CLEAN算法處理后頻譜。(a)天線A的頻譜；(b)天線B的頻譜

CLEAN算法處理后的互相關(guān)功率譜(圖14)同樣印證了算法的有效性。由圖14可以看出，兩天線的互相關(guān)功率譜中已經(jīng)沒有明顯的峰值結(jié)構(gòu)，各個頻段的互相關(guān)功率譜峰值接近10-9V2，對比圖10的結(jié)果，降低了約9個數(shù)量級。

圖14 天線A和天線B CLEAN算法處理后的互相關(guān)功率譜

圖15給出了兩天線信號之間的相關(guān)系數(shù)隨迭代次數(shù)變化的曲線。CLEAN算法處理前初始的相關(guān)系數(shù)為0.965，由強平臺干擾信號主導(dǎo)；CLEAN算法處理10萬次后的相關(guān)系數(shù)接近0.2，根據(jù)此相關(guān)系數(shù)計算的信噪比，可以反映此時殘余信號中的平臺干擾信號和接收機噪聲基本在同一水平，此時信號由接收機噪聲主導(dǎo)。當(dāng)然，CLEAN算法處理結(jié)束后，兩天線信號的相關(guān)系數(shù)并沒有顯著接近0，說明殘余信號中可能包含來自宇宙的低頻射電信號，以及沒有徹底清除的平臺干擾信號。

圖15 天線A和天線B CLEAN算法處理時相關(guān)系數(shù)隨迭代次數(shù)分布曲線

3.3 CLEAN算法效果分析

(5)

當(dāng)相關(guān)系數(shù)從最初的0.965降低到最后接近0.2時，由(5)式計算得到的信噪比從5.3降低到0.5，說明在末態(tài)殘余信號中，相干的信號成分和不相干的噪聲成分幾乎在同一水平。

通過將末態(tài)信號轉(zhuǎn)化為流量密度，得到天線B末態(tài)信號的流量密度隨頻率的分布(圖16)。從圖16可以看出，末態(tài)流量密度基本降低到10-23W/(m2·Hz)以下，與太陽寧靜時的信號水平相當(dāng)。后續(xù)可以在此基礎(chǔ)上對強射電爆發(fā)進(jìn)行研究。

圖16 天線B CLEAN算法處理后殘余信號的流量密度

4 討 論

經(jīng)過CLEAN算法處理后，嫦娥四號低頻射電頻譜儀各天線的探測靈敏度(未積分)從10-15W/(m2·Hz)提高到10-23W/(m2·Hz)，在此靈敏度下可以對太陽爆發(fā)進(jìn)行監(jiān)測。太陽爆發(fā)的流量密度一般在10-18W/(m2·Hz)，極端爆發(fā)可以達(dá)到10-15W/(m2·Hz)。除此之外，這樣的靈敏度也足以探測部分強的快速射電暴，其輻射流量密度可以達(dá)到幾十萬央斯基(Jy)，1Jy=10-26W/(m2·Hz)。如果再對時域數(shù)據(jù)進(jìn)行積分平均，可以進(jìn)一步提高探測靈敏度。例如，如果對1道4 096個數(shù)據(jù)進(jìn)行積分平均，靈敏度可以提高64倍，大大拓展低頻射電頻譜儀的科學(xué)探測目標(biāo)。

CLEAN算法分離的平臺干擾信號也很有用。首先可以用來分析平臺干擾信號的頻譜、空間和時間變化特征，并確定其來源，以便在今后開發(fā)類似設(shè)備時，改進(jìn)電磁兼容設(shè)計，大幅度降低平臺干擾水平。其次，在全面了解平臺干擾信號特征的基礎(chǔ)上，充分利用這些知識，發(fā)現(xiàn)特征異常的信號，這些異常信號往往來自宇宙。

分離并抑制干擾信號的最終目的是進(jìn)行有效的天文觀測。如何在強烈的平臺干擾背景下尋找并確認(rèn)宇宙射電信號呢？主要思路是尋找異常信號。平臺干擾往往有固定的一組頻譜特征，一旦有例外的頻譜信號出現(xiàn)，可能意味著記錄了來自宇宙的射電爆發(fā)信號。如果這些頻譜還具備色散等特征，那么幾乎可以肯定它來自宇宙。再例如，天線A，B和C信號的強度比例系數(shù)，若平臺干擾信號產(chǎn)生的空間分布是固定的，那么不管其在時域如何變化，信號之間的強度比例系數(shù)也是固定的。倘若這些比例系數(shù)有明顯的改變，則往往意味著天線接收到來自宇宙的強射電爆發(fā)信號。若射電爆發(fā)信號的時標(biāo)長達(dá)幾天量級，這些比例系數(shù)還會隨著月球自轉(zhuǎn)發(fā)生規(guī)律性變化，有助于我們更加可靠地確定爆發(fā)信號的來源。當(dāng)然，具體如何搜尋異常信號，還得依賴后續(xù)細(xì)致的數(shù)據(jù)分析工作。

本文提出的算法也可以應(yīng)用到其他射電觀測設(shè)備，尤其是地基綜合射電干涉陣。目前地基干涉陣受到越來越多的來自地面設(shè)備以及空間衛(wèi)星的無線電干擾。借助CLEAN算法，我們利用兩天線之間的互相關(guān)功率譜，在信號層面提取干擾信號，再對殘余信號數(shù)據(jù)相關(guān)獲得可見度數(shù)據(jù)，對可見度數(shù)據(jù)作成圖等處理工作。同樣，我們也可以對提取的干擾信號進(jìn)行類似的操作，獲得詳細(xì)的頻譜、空間分布和時域變化等信息。

5 結(jié) 論

本文提出一種CLEAN算法，從兩組信號的相關(guān)性出發(fā)，并借助互相關(guān)功率譜、傅里葉級數(shù)等工具，把低頻射電頻譜儀天線A，B和C的時域觀測數(shù)據(jù)分為強相關(guān)的CLEAN模型信號和部分相關(guān)的殘余信號。其中CLEAN模型信號主要由平臺干擾信號和可能的低頻強射電爆發(fā)組成；殘余信號由接收機噪聲、未扣除的平臺干擾信號和常規(guī)低頻射電信號組成。將該算法應(yīng)用到實際數(shù)據(jù)中，結(jié)果表明，嫦娥四號低頻射電頻譜儀的未積分靈敏度可以提高約8個數(shù)量級，達(dá)到10-23W/(m2·Hz)。

在此基礎(chǔ)上，基于對平臺干擾信號中確定成分和寬帶隨機成分的分類處理，再借助低頻射電爆發(fā)信號和平臺干擾信號在功率譜上不同的表現(xiàn)，以及常規(guī)低頻射電信號受月球自轉(zhuǎn)調(diào)制等信息，將來工作的重點是進(jìn)一步分析CLEAN模型信號和殘余信號，以期搜索低頻強射電爆發(fā)信號，乃至對全天區(qū)進(jìn)行粗略成像。