張海龍，張亞州，王 杰，冶鑫晨，王萬(wàn)瓊，李 嘉，張 萌，杜 旭

(1. 中國(guó)科學(xué)院新疆天文臺(tái)，新疆 烏魯木齊 830011；2. 中國(guó)科學(xué)院大學(xué)，北京 100049；3. 中國(guó)科學(xué)院射電天文重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210033；4. 國(guó)家天文科學(xué)數(shù)據(jù)中心，北京 100101)

射電天文學(xué)中，射頻干擾定義為任何可能影響天文觀測(cè)的有害信號(hào)。射電天文學(xué)發(fā)展過(guò)程中，抑制射頻干擾一直是天文學(xué)家研究的熱點(diǎn)。隨著信息技術(shù)飛速發(fā)展，人為造成的干擾信號(hào)導(dǎo)致天文臺(tái)站周圍電磁環(huán)境惡化，影響射電望遠(yuǎn)鏡正常觀測(cè)。

射頻干擾的來(lái)源多種多樣，包括外部來(lái)源和內(nèi)部來(lái)源。外部來(lái)源主要由天文臺(tái)址外的設(shè)備產(chǎn)生，如人造衛(wèi)星[1](北斗和全球定位系統(tǒng)導(dǎo)航衛(wèi)星等)、飛機(jī)、臺(tái)址附近的基站信號(hào)以及電視廣播信號(hào)等；內(nèi)部來(lái)源主要是天文臺(tái)址內(nèi)的電子設(shè)備，如計(jì)算機(jī)、視頻監(jiān)控、網(wǎng)絡(luò)交換設(shè)備、無(wú)線輸入設(shè)備等[2]。天文信號(hào)通常為寬頻帶且在時(shí)間尺度上平滑變化，而射頻干擾在時(shí)域和頻域的幅值強(qiáng)度高，絕大部分射頻干擾與天文信號(hào)相比有明顯差異[3]。射頻干擾與熱噪聲不同，它通常由通信系統(tǒng)、人造雷達(dá)或電子設(shè)備產(chǎn)生的干擾構(gòu)成，且具有復(fù)雜的時(shí)間和頻率結(jié)構(gòu)，常見通信信號(hào)的功率比天文信號(hào)高多個(gè)數(shù)量級(jí)并且隨時(shí)間變化，因此無(wú)法通過(guò)長(zhǎng)時(shí)間積分來(lái)降低強(qiáng)度，嚴(yán)重影響天文數(shù)據(jù)質(zhì)量。天文學(xué)家針對(duì)上述問(wèn)題已提出了多種射頻干擾抑制方法，但仍受到以下因素制約：

(1)隨著無(wú)線技術(shù)的飛速發(fā)展及人類活動(dòng)范圍的擴(kuò)大，地面及空間產(chǎn)生的人為干擾逐年遞增。

(2)隨著制造技術(shù)和信息技術(shù)的飛速發(fā)展，射電望遠(yuǎn)鏡向大口徑或大規(guī)模陣列方向發(fā)展，隨著多波束或相控陣饋源接收系統(tǒng)的投入使用，觀測(cè)設(shè)備的靈敏度及數(shù)據(jù)收集能力顯著提升，同時(shí)射頻干擾抑制問(wèn)題更加復(fù)雜。

在射電天文觀測(cè)以及數(shù)據(jù)處理過(guò)程中需要針對(duì)不同階段采取不同的射頻干擾應(yīng)對(duì)策略，本文將射頻干擾抑制策略分為主動(dòng)預(yù)防階段、預(yù)相關(guān)階段及后相關(guān)階段，不同階段的射頻干擾問(wèn)題可以采取不同的抑制方法。主動(dòng)預(yù)防階段若能采取有效防護(hù)手段，能夠杜絕通信基站、電視廣播等大部分射頻干擾進(jìn)入系統(tǒng)；預(yù)相關(guān)階段采取參考天線及空間濾波等方法可以應(yīng)對(duì)特定射頻干擾源；后相關(guān)階段可以采取閾值法處理幅值強(qiáng)度遠(yuǎn)大于天文信號(hào)的射頻干擾，或利用機(jī)器學(xué)習(xí)相關(guān)技術(shù)實(shí)現(xiàn)射頻干擾數(shù)據(jù)標(biāo)記。下面針對(duì)不同階段可采取的射頻干擾抑制策略進(jìn)行詳細(xì)分析。

1 射頻干擾主動(dòng)預(yù)防階段

1.1 電磁波寧?kù)o區(qū)

射電望遠(yuǎn)鏡觀測(cè)來(lái)自宇宙遙遠(yuǎn)天體的微弱信號(hào)，對(duì)臺(tái)址周邊電磁干擾有較高要求，選擇地理環(huán)境較好的臺(tái)址及建立無(wú)線電寧?kù)o保護(hù)區(qū)是抑制射頻干擾的重要一步，可以從根源上去除大部分干擾。 500 m口徑球面射電望遠(yuǎn)鏡(Five-hundred-meter Aperture Spherical radio Telescope, FAST)是世界最大單口徑、最靈敏的射電望遠(yuǎn)鏡，位于貴州省黔南布依族苗族自治州平塘縣克度鎮(zhèn)大窩凼，以臺(tái)址為中心，建立半徑為30 km的電磁波寧?kù)o區(qū)，分為不同要求的 3 個(gè)區(qū)域。其中以臺(tái)址為圓心，半徑5 km范圍為核心保護(hù)區(qū)，5～10 km的環(huán)帶為中間區(qū)，10～30 km的環(huán)帶為邊遠(yuǎn)區(qū)[4]，如圖1。電磁波寧?kù)o區(qū)能夠屏蔽絕大部分外來(lái)干擾源。

圖1 500 m口徑球面射電望遠(yuǎn)鏡無(wú)線電寧?kù)o保護(hù)區(qū) [5]Fig.1 FAST radio quiet zone[5]

新疆奇臺(tái)110 m全向可動(dòng)射電望遠(yuǎn)鏡(Qi Tai Radio Telescope, QTT)[6]位于新疆昌吉州奇臺(tái)縣，以臺(tái)址為中心，依據(jù)地形特征及當(dāng)?shù)貙?shí)際情況, 規(guī)劃了30 km的無(wú)線電寧?kù)o區(qū)。該寧?kù)o區(qū)共分為3個(gè)區(qū)域，分別為核心區(qū)、限制區(qū)和協(xié)調(diào)區(qū)，其中核心區(qū)為2.5 km × 4 km的矩形區(qū)域；限制區(qū)為10 km × 15 km的矩形區(qū)域[7-8]，如圖2。

圖2 新疆奇臺(tái)110 m全向可動(dòng)射電望遠(yuǎn)鏡

1.2 預(yù)留頻段

國(guó)際電信聯(lián)盟(International Telecommunication Union, ITU)通過(guò)劃分頻譜來(lái)協(xié)調(diào)8.3 MHz～300 GHz之間的無(wú)線電頻譜，將其分為幾個(gè)頻帶，僅允許指定服務(wù)運(yùn)行[9]。在分配的所有頻段中，約有70到80個(gè)頻段對(duì)應(yīng)于射電天文學(xué)(確切的數(shù)量取決于地區(qū)和當(dāng)?shù)胤?。針對(duì)射電天文頻段，一種是射電天文專用頻段(或與無(wú)源服務(wù)共享)，嚴(yán)格禁止無(wú)線電發(fā)射；另一種是與有源業(yè)務(wù)共享的頻段，只能部分施加保護(hù)。射電天文的波段選擇通常與科學(xué)目標(biāo)有關(guān)。例如，氫線位于1 420 MHz附近，因此1 400～1 427 MHz頻帶給射電天文學(xué)保留[10]。

即使是射電天文專用頻段，由于有源元件的無(wú)線電寧?kù)o保護(hù)區(qū)[7]

Fig.2 QTT radio quiet zone[7]

諧波和功率泄漏，仍然可能存在干擾。對(duì)于這些情況，國(guó)際電信聯(lián)盟在其 ITU-R RA.769建議報(bào)告中將最大可接受的干擾定義為在測(cè)量功率中引入不超過(guò)10%的誤差，這是天文學(xué)家普遍接受的閾值，該報(bào)告還為一些典型的天文波段提供了一份使用普通望遠(yuǎn)鏡和觀測(cè)參數(shù)的建議閾值清單[11]。

1.3 射頻干擾主動(dòng)預(yù)防策略的不足

射頻干擾主動(dòng)預(yù)防策略中最直接的方法是為天文臺(tái)選擇遠(yuǎn)離干擾源的偏遠(yuǎn)位置。自然形態(tài)(例如山脈)可以有效阻止外部干擾，但干擾源會(huì)因?yàn)槎鄰絺鞑ギa(chǎn)生反射和衍射，從而增加干擾抑制的難度。研究發(fā)現(xiàn)，在望遠(yuǎn)鏡臺(tái)址周圍種植松樹等針葉樹可以有效抑制射頻干擾，針葉中的水分可以吸收1 GHz以上的信號(hào)[12]。部分天文臺(tái)站建在高海拔地區(qū)以減少大氣對(duì)觀測(cè)的影響，遠(yuǎn)離人類活動(dòng)范圍，同時(shí)盡可能避免干擾源多徑傳播。面對(duì)內(nèi)部干擾，通常的辦法是在產(chǎn)生干擾的電子設(shè)備(計(jì)算機(jī)、微波和射頻元件)上覆蓋電磁屏蔽材料如導(dǎo)電箔，也可以將電子設(shè)備集中放入屏蔽室，將干擾控制在密閉空間，并不影響設(shè)備性能。

通常，主動(dòng)預(yù)防是處理干擾的最有效方法，可以作為射頻干擾的第1層防護(hù)。但仍然存在不足：

(1)只能為電磁頻譜的小部分提供保護(hù)；

(2)對(duì)于高靈敏度的設(shè)備，電磁屏蔽材料不能完全屏蔽干擾；

(3)衛(wèi)星通信以及飛機(jī)等產(chǎn)生的干擾，由于相對(duì)位置限制，無(wú)法通過(guò)主動(dòng)預(yù)防避免。

2 射頻干擾預(yù)相關(guān)階段

本文將觀測(cè)獲取的原始數(shù)據(jù)在未寫入磁盤形成可供科學(xué)研究處理之前的射頻干擾處理階段定義為預(yù)相關(guān)階段。對(duì)于單天線觀測(cè)可利用參考天線通過(guò)自適應(yīng)濾波方式去除特定干擾，對(duì)于陣列天線可采用空間濾波方法抑制射頻干擾。

2.1 自適應(yīng)濾波

自適應(yīng)濾波算法是由Barnbaum和Bradley在1998年首次引入射電天文領(lǐng)域，用來(lái)解決射頻干擾問(wèn)題[13]，在2005年，文[14]的現(xiàn)場(chǎng)試驗(yàn)證明，在射頻干擾存在的情況下，自適應(yīng)濾波方式能夠極大提高脈沖星觀測(cè)數(shù)據(jù)的質(zhì)量。

參考天線接收的信號(hào)可以表示為

R=Rnoise+Rrfi+Rastr，

(1)

其中，R為參考天線接收的信號(hào)；Rnoise為參考天線的系統(tǒng)噪聲；Rrfi為參考天線接收的射頻干擾；Rastr為參考天線接收的天文信號(hào)，天文信號(hào)通過(guò)參考天線的旁瓣接收，其強(qiáng)度微弱可忽略不計(jì)。同理，主望遠(yuǎn)鏡接收的信號(hào)可以表示為

M=Mnoise+Mrfi+Mastr，

(2)

主望遠(yuǎn)鏡和參考天線接收的系統(tǒng)噪聲不相關(guān)，而接收的射頻干擾具有相關(guān)性。為了消除主天線中的射頻干擾，使其對(duì)準(zhǔn)射頻干擾源，并通過(guò)不斷優(yōu)化參考天線的指向和極化方向提高參考天線中射頻干擾的信噪比。

自適應(yīng)濾波具體實(shí)現(xiàn)步驟：

(1)在自適應(yīng)相關(guān)環(huán)路中確定復(fù)增益系數(shù)g(該系數(shù)是不斷迭代優(yōu)化的)與參考天線接收的信號(hào)R相乘，逐步迭代gRrfi，使其大小接近主天線的射頻干擾信號(hào)，能最大程度地消除主天線某一特定方向的射頻干擾。

(2)當(dāng)g達(dá)到最優(yōu)值時(shí)，gRrfi基本接近Mrfi，天文信號(hào)與參考信號(hào)的差值即去除特定干擾后的天文信號(hào)。

(3)圖3為自適應(yīng)濾波器基本設(shè)計(jì)圖，將原始參考信號(hào)與濾波器輸出相關(guān)，當(dāng)相關(guān)系數(shù)為0時(shí)，達(dá)到最佳增益。

相關(guān)項(xiàng)計(jì)算公式為

C(g)=(M-g·R)·R，

(3)

當(dāng)C(g)=0時(shí)，可剔除主天線的射頻干擾，即

(4)

圖3 自適應(yīng)濾波器基本設(shè)計(jì)圖Fig.3 Adaptive filter basic design diagram

代入天線接收的各個(gè)分量得

(5)

因?yàn)閮蓚€(gè)天線的系統(tǒng)噪聲、天文信號(hào)、其他信號(hào)均與射頻干擾不相關(guān)，參考天線中的天文信號(hào)非常微弱，可以忽略不計(jì)，即Rastr=0，那么上式進(jìn)一步簡(jiǎn)化得

(6)

令Mrfi=ξRrfi，g可以寫為

(7)

(8)

在天文觀測(cè)中，主天線指向待觀測(cè)的天體，Mrfi在觀測(cè)期間可以認(rèn)為是恒定量，Srfi與干擾噪聲比成反比，即干擾噪聲比越大，Srfi越小，說(shuō)明提高射頻干擾在參考天線中的信噪比能提高主天線射頻干擾抑制能力。

影響自適應(yīng)濾波性能的幾個(gè)因素：

(1)混合射頻干擾源

自適應(yīng)濾波方法消除單源射頻干擾效果較好，但是面對(duì)復(fù)雜的射頻干擾環(huán)境，其抑制射頻干擾的性能會(huì)大幅度下降。對(duì)于混合射頻干擾，當(dāng)兩個(gè)干擾信號(hào)處于相同頻域時(shí)，振幅大的干擾被消減。天文臺(tái)址附近射頻干擾來(lái)源復(fù)雜，且與天文信號(hào)相比具有較高的功率，在此類條件下自適應(yīng)濾波方法射頻干擾抑制效果一般[15]。

(2)多路徑傳播影響

多徑傳播是指由于信號(hào)到達(dá)接收天線的路徑不同，而存在同一信號(hào)多個(gè)副本的現(xiàn)象，通常是由大型陸地物體(山脈、 建筑物)、大型水體(湖泊、 瀉湖)和電離層反射產(chǎn)生的。

圖4為多路徑干擾信號(hào)傳播場(chǎng)景。因?yàn)榈竭_(dá)望遠(yuǎn)鏡和參考天線的多個(gè)副本有不同的時(shí)間延遲，如果延遲足夠大，那么添加的信號(hào)就可以視為完全不相關(guān)。在這種情況下，多路徑情況相當(dāng)于多個(gè)射頻干擾源，射頻干擾信號(hào)種類繁多，難以獲得濾波器多徑傳播的通用表達(dá)式。防止多徑傳播的最佳措施是將望遠(yuǎn)鏡和參考天線放置在盡可能近的位置，以減小傳播距離差異，并通過(guò)增加快速傅里葉變換點(diǎn)數(shù)來(lái)提高系統(tǒng)的頻率分辨率，或在快速傅里葉變換點(diǎn)數(shù)不變的條件下減小系統(tǒng)帶寬[15]。

圖4 射頻干擾多徑傳播情況Fig.4 RFI multipath propagation situation

(3)違反線性時(shí)不變(LTI)條件

假設(shè)干擾傳播路徑至少在濾波器收斂的時(shí)間內(nèi)滿足線性時(shí)不變條件，如果不滿足線性時(shí)不變條件，則表示射頻干擾路徑在某些地方呈現(xiàn)非線性或時(shí)變的傳播特性。然而，沒(méi)有任何自然介質(zhì)表現(xiàn)出足夠強(qiáng)的非線性，或具有快速變化的傳播特性而影響濾波器的性能，特別是考慮到最終算法收斂時(shí)間是毫秒級(jí)的。當(dāng)望遠(yuǎn)鏡跟蹤天體時(shí)，捕獲射頻干擾的旁瓣隨望遠(yuǎn)鏡一起移動(dòng)，被測(cè)射頻干擾的路徑特性隨時(shí)間逐漸變化。自適應(yīng)濾波方法需要具備相應(yīng)計(jì)算能力，不會(huì)因望遠(yuǎn)鏡緩慢移動(dòng)而受到明顯影響。

接收機(jī)前端的一些模擬元件是信號(hào)飽和的，如放大器和頻率混頻器。當(dāng)一個(gè)分量達(dá)到飽和時(shí)，輸入信號(hào)的幅值會(huì)被削去，從而產(chǎn)生諧波失真，影響系統(tǒng)的線性時(shí)不變條件。這種失真被視為附加的損壞信號(hào)，其頻率是原始輸入信號(hào)中心頻率的幾倍。此外，由于數(shù)字化后的混疊等原因，諧波信號(hào)往往被混入基帶，從而破壞天文數(shù)據(jù)。主信道射頻干擾諧波失真會(huì)嚴(yán)重破壞濾波器性能，附加的諧波信號(hào)不會(huì)出現(xiàn)在參考信道中，使得這類射頻干擾更難抑制[15]。

2.2 空間濾波

空間濾波方法主要針對(duì)多通道射電天文觀測(cè)中的射頻干擾抑制問(wèn)題[16-18]，該技術(shù)適用于干涉射電望遠(yuǎn)鏡陣列，如荷蘭的Westerbork綜合射電望遠(yuǎn)鏡(Westerbork Synthesis Radio Telescope, WSRT)、美國(guó)的甚大陣列(Very Large Array, VLA)，或配置了相控陣接收系統(tǒng)的單天線望遠(yuǎn)鏡。圖5為空間濾波技術(shù)在天線陣列系統(tǒng)中的應(yīng)用。

圖5 天線陣列系統(tǒng)在不同階段抑制射頻干擾

射電望遠(yuǎn)鏡接收的信號(hào)由3部分組成，射電源信號(hào)(擬觀測(cè)目標(biāo))、系統(tǒng)噪聲和射頻干擾。系統(tǒng)噪聲由宇宙背景噪聲、大氣噪聲、接收機(jī)噪聲等噪聲分量組成。根據(jù)中心極限定理，系統(tǒng)噪聲是時(shí)間獨(dú)立且高斯分布的信號(hào)，在短時(shí)間間隔上是短時(shí)平穩(wěn)的，且所有的射電源假定在統(tǒng)計(jì)學(xué)上彼此獨(dú)立[12]。

射電天文學(xué)中的射頻干擾特征是各種各樣的。本文選擇中心窄帶高斯平穩(wěn)、二階非周期和循環(huán)平穩(wěn)信號(hào)3種不同的信號(hào)特性對(duì)射頻干擾建模。

天線陣列共有M個(gè)天線，其中x(t)=[x1(t)…xM(t)]T表示M× 1個(gè)天線陣列的輸出，其中，第k個(gè)天線在t時(shí)刻的輸出為xk(t)，那么x(t)可以寫為

x(t)=Ac·c(t)+Ar·r(t)+n(t)，

(9)

其中，c(t)=[c1(t)…cNc(t)]T表示在t時(shí)刻射電源的特征向量，大小為Nc× 1；Ac=[c1(t,θc1,φc1)…cNc(t,θcNc,φcNc)]是射電源的空間特征向量矩陣，大小為M×Nc，其中cn(t,θcn,φcn)=[cn,1(t,θcn,φcn)…cn,M(t,θcn,φcn)]T是第n個(gè)射電源的空間特征向量。同理，Ar和r(t)分別代表射頻干擾的空間特征向量矩陣以及在t時(shí)刻射頻干擾的信號(hào)特征向量，n(t)=[n1(t)…nM(t)]T是系統(tǒng)噪聲在t時(shí)刻的向量，大小為M× 1。

通過(guò)協(xié)方差矩陣描述射電源在時(shí)間維度上的相關(guān)性，即

(10)

其中，R(t,τ)為在t時(shí)刻τ范圍內(nèi)射電源信號(hào)的相關(guān)性；xH代表x的共軛轉(zhuǎn)置。假設(shè)射頻干擾、射電源和系統(tǒng)噪聲三者相互獨(dú)立，即互不相關(guān)，那么

(11)

從線性代數(shù)的觀點(diǎn)來(lái)看，相控天線陣列射電望遠(yuǎn)鏡數(shù)據(jù)協(xié)方差矩陣可以看作是線性變換矩陣。這個(gè)矩陣生成一個(gè)數(shù)據(jù)向量空間，該向量空間由射頻干擾子空間、射電源子空間和系統(tǒng)噪聲子空間組成。圖6展示了一個(gè)無(wú)噪聲場(chǎng)景下的二維數(shù)據(jù)向量空間的示例。紅色和黑色矢量分別代表與射電源和射頻干擾源相關(guān)的向量。

使用正交投影技術(shù)減少干擾，將數(shù)據(jù)向量空間投影到與射頻干擾子空間正交的子空間上，投影的射頻干擾子空間完全為0，由此產(chǎn)生的數(shù)據(jù)只包含射電源。

從圖6可以清楚地看出，在投影后，恢復(fù)的射電源能量取決于兩個(gè)向量的夾角大小，即

(12)

其中，P2為射電源的功率；θac,arn為射頻干擾向量與射電源向量的夾角大小，在0～90°范圍內(nèi)，θac,arn越大，最終投影得到的射電源功率越大。但是射電源子空間和干擾子空間之間的夾角難以測(cè)量，需要借助投影矩陣間接計(jì)算。

圖6 二維數(shù)據(jù)向量空間的正交投影

設(shè)投影矩陣為P，它的特征值是0或1。特征值為0對(duì)應(yīng)的特征向量生成投影核子空間，特征值為1對(duì)應(yīng)的特征向量生成值域子空間。因此投影矩陣的秩等于它的值域子空間的維數(shù)。如果H是投影范圍子空間基，K是投影核子空間基，則有

P·H=H，P·K=0.

(13)

由基矩陣H生成的投影矩陣定義為

P=H(HTH)-1HT，

(14)

其值域是正交于H的子空間的正交投影為P⊥，定義為

P⊥=I-H(HTH)-1HT.

(15)

在多干擾的情況下，干擾子空間是多維的，這個(gè)子空間是由存儲(chǔ)在矩陣Ar中的獨(dú)立射頻干擾向量集合生成。因此，投影矩陣可以表示為

(16)

在預(yù)相關(guān)階段，可用數(shù)據(jù)為天線陣列輸出數(shù)據(jù)向量x(t)，將其投影到射頻干擾子空間的正交空間，修正后的數(shù)據(jù)向量xclean(t)為

xclean(t)=Pr⊥x(t).

(17)

從以上內(nèi)容可知，空間濾波需要預(yù)先知道干擾源的來(lái)源，以便確定投影矩陣，對(duì)于復(fù)雜的干擾環(huán)境，測(cè)定干擾源的具體位置是非常困難的。

3 射頻干擾后相關(guān)階段

由于射頻干擾信號(hào)強(qiáng)度通常比天文信號(hào)大幾個(gè)數(shù)量級(jí)，閾值法是一種有效抑制射頻干擾的方法。閾值水平通常根據(jù)統(tǒng)計(jì)量決定，可以用一段數(shù)據(jù)的均值、均方根等作為閾值參考值，超出這個(gè)范圍的值標(biāo)記為射頻干擾。閾值法不足之處是會(huì)誤將天文信號(hào)標(biāo)記為射頻干擾，在后相關(guān)階段降低射頻干擾檢測(cè)的錯(cuò)誤率是研究射頻干擾抑制問(wèn)題的一個(gè)重要方向。

3.1 VarThreshold算法

在頻域和時(shí)域上，射頻干擾往往影響相鄰的數(shù)據(jù)，因?yàn)槠鋸?qiáng)度比較低，不會(huì)被標(biāo)記為射頻干擾，增加了檢測(cè)射頻干擾的錯(cuò)誤率。VarThreshold是一種組合閾值算法，其思想是當(dāng)組合樣本的某一統(tǒng)計(jì)屬性超過(guò)某個(gè)限值時(shí)，將標(biāo)記樣本的組合為射頻干擾[19-20]。 假設(shè)A和B是相鄰的樣本采樣點(diǎn)，常規(guī)閾值法是分別查看樣本A和B是否超過(guò)設(shè)定的閾值水平線，如果超過(guò)，則標(biāo)記為射頻干擾。但對(duì)于組合閾值來(lái)說(shuō)，只有樣本A和B都超過(guò)閾值，才將這一組合的樣本都標(biāo)記為射頻干擾。如果樣本A和B沒(méi)有被標(biāo)記為射頻干擾，它們將和下一個(gè)相鄰的采樣點(diǎn)C進(jìn)行組合，繼續(xù)與閾值進(jìn)行對(duì)比，迭代多次完成射頻干擾標(biāo)記。算法公式為

flagυM(υ,t)=?i∈{0,…，M-1}:?j∈{0,…，M-1}:|R(υ+aΔυ,t)|>χM，

(18)

(18)式表達(dá)了從t時(shí)刻起，M個(gè)數(shù)的范圍內(nèi)，如果該范圍內(nèi)所有采樣點(diǎn)的絕對(duì)值大于χM，那么將這M個(gè)數(shù)都標(biāo)記為射頻干擾。

VarThreshold算法用一組嚴(yán)格遞減的閾值序列χi(i=1,…,N)判定某個(gè)采樣點(diǎn)是否應(yīng)被標(biāo)記為射頻干擾，公式為

(19)

其中，χ1是單采樣點(diǎn)閾值，當(dāng)ρ=1.5時(shí)效果最好。χ1的取值由數(shù)據(jù)的統(tǒng)計(jì)水平量決定，即選取的閾值使射頻干擾的假陽(yáng)率最小。

3.2 SumThreshold算法

SumThreshold算法也是一種組合閾值算法[19-20]，與VarThreshold算法求解閾值序列類似，只是ρ值不同，當(dāng)ρ=1.2時(shí)效果最好。求出一組閾值T={T1,T2,T4, …}，通過(guò)迭代的方式與樣本進(jìn)行比較，即T1與單樣本比較，T2與兩個(gè)樣本的均值比較，Tm與m個(gè)樣本的均值比較，當(dāng)樣本均值大于閾值時(shí)，標(biāo)記為射頻干擾，且不參與后續(xù)的迭代計(jì)算過(guò)程。以一維序列為例，設(shè)序列為{1, 7, 4, 3, 2, 1, 1, 6, 1, 1, 3, 2, 3, 2}，閾值χ為{χ1=6,χ2=3,χ4=2}。

圖7、圖8和圖9展示了 SumThreshold算法的過(guò)程。從算法角度考慮，該算法的時(shí)間復(fù)雜度為

圖7 第1輪，閾值為χ1=6，只有7被標(biāo)記

圖8 第2輪，閾值為χ2=3Fig.8 In the second round, the threshold value is χ2=3

圖9 第3輪，閾值為χ4=2Fig.9 In the third round, the threshold value is χ4=2

O(n2)。為了減小計(jì)算速度，當(dāng)閾值M取1, 2, 4, 8, 16, 32, 64, 128, 256 … 時(shí)，時(shí)間復(fù)雜度降至O(nlogn)。

3.3 奇異值分解(Singular Value Decomposition, SVD)

奇異值分解是線性代數(shù)中一種重要的矩陣分解，能夠提取信息，簡(jiǎn)化數(shù)據(jù)，去除噪聲等，可將數(shù)據(jù)映射到低維空間[21]。奇異值可以分解為

A=UΣVT.

(20)

如(20)式，奇異值分解將原始矩陣A分解為3個(gè)矩陣：U,Σ和VT，其中U和V是正交矩陣，Σ是對(duì)角矩陣，對(duì)角線上的元素為奇異值。絕大部分射頻干擾信號(hào)在時(shí)域或頻域上的振幅較大。用A代表接收的天文信號(hào)，將其進(jìn)行奇異值分解，即A=UΣVT，將Σ中最大特征值置為0，逆

向恢復(fù)后，得到新的矩陣A＾，那么原矩陣A中幅值較大的信號(hào)(射頻干擾)會(huì)被抑制。

當(dāng)射頻干擾信號(hào)足夠強(qiáng)，且σRFIσastro=0時(shí)，算法效果最好。該方法僅適用于寬帶射頻干擾信號(hào)，對(duì)于遵從高斯分布的射頻干擾信號(hào)并不適用[19-20]。

4 基于機(jī)器學(xué)習(xí)的射頻干擾抑制方法

隨著機(jī)器學(xué)習(xí)和深度學(xué)習(xí)等技術(shù)快速發(fā)展以及天文數(shù)據(jù)的急速增長(zhǎng)，兩者的結(jié)合已成必然趨勢(shì)。通過(guò)對(duì)大量數(shù)據(jù)進(jìn)行學(xué)習(xí)和訓(xùn)練后，能有效地標(biāo)記天文信號(hào)中的射頻干擾，并進(jìn)行特征提取。

4.1 主成分分析法(Principal Component Analysis, PCA)

主成分分析是一種用于數(shù)據(jù)特征降維的方法，常用于減少數(shù)據(jù)集的維數(shù)，同時(shí)保持?jǐn)?shù)據(jù)集對(duì)方差貢獻(xiàn)最大的特征。這樣低階成分往往能夠保留數(shù)據(jù)的重要信息，方便特征提取和分類[22-23]。

文[24]對(duì)生活中常見的9種瞬時(shí)射頻干擾來(lái)源進(jìn)行分類測(cè)試實(shí)驗(yàn)，在時(shí)域范圍內(nèi)分別使用標(biāo)準(zhǔn)主成分分析和核主成分分析進(jìn)行測(cè)試。之后用聚類分離的方法比較了兩種方法的聚類精度，結(jié)果如圖10和圖11，發(fā)現(xiàn)核主成分分析在區(qū)分來(lái)源方面比標(biāo)準(zhǔn)主成分分析更好，能有效區(qū)分瞬時(shí)射頻干擾來(lái)源。

圖10 標(biāo)準(zhǔn)主成分分析[24] Fig.10 Standard PCA[24]

圖11 核主成分分析[24]Fig.11 Kernel PCA[24]

4.2 支持向量機(jī)(Support Vector Machine, SVM)

在機(jī)器學(xué)習(xí)中，支持向量機(jī)是在分類與回歸分析中分析數(shù)據(jù)的監(jiān)督式學(xué)習(xí)模型與相關(guān)的學(xué)習(xí)算法[25]。

文[26]選取時(shí)間長(zhǎng)度為300 μs的非射頻干擾和射頻干擾信號(hào)訓(xùn)練，提取均方根、均值、方差、均值與方差之比、偏度、峰度、最大振幅、最小振幅和峰值個(gè)數(shù)等信息，導(dǎo)入支持向量機(jī)訓(xùn)練模型，之后對(duì)測(cè)試數(shù)據(jù)進(jìn)行射頻干擾和非射頻干擾分類，流程如圖12。結(jié)果表明，即使干擾噪聲比低和射頻干擾占空比小的情況下，該方法仍能準(zhǔn)確地檢出射頻干擾，效果良好。

圖12 基于支持向量機(jī)的射頻干擾抑制過(guò)程Fig.12 RFI suppression process based on SVM

4.3 神經(jīng)網(wǎng)絡(luò)

神經(jīng)網(wǎng)絡(luò)的基本結(jié)構(gòu)包含輸入層、隱藏層和輸出層，由大量節(jié)點(diǎn)相互連接構(gòu)成，是一種非線性統(tǒng)計(jì)數(shù)據(jù)建模工具，常用來(lái)對(duì)輸入和輸出間復(fù)雜的關(guān)系進(jìn)行建模，或用來(lái)探索數(shù)據(jù)的模式[27]。

文[28]提出了全卷積神經(jīng)網(wǎng)絡(luò)(Deep Fully Convolutional Neural Network, DFCN)用于圖像的分割，解決的核心問(wèn)題是圖像像素級(jí)別的分類。文[29]利用全卷積神經(jīng)網(wǎng)絡(luò)，如圖13，H和W分別對(duì)應(yīng)輸入的時(shí)間和頻率可見性維數(shù)，F(xiàn)為濾波器層數(shù)，L對(duì)應(yīng)輸入和全卷積層之間的層數(shù)總和，以振幅和相位為特征，瀑布圖進(jìn)行特征提取，并標(biāo)出射頻干擾。

圖13 全卷積神經(jīng)網(wǎng)絡(luò)架構(gòu)[28]Fig.13 DFCN architecture[28]

文[30]提出的U-Net神經(jīng)網(wǎng)絡(luò)，是在卷積神經(jīng)網(wǎng)絡(luò)的基礎(chǔ)上增加卷積的擴(kuò)展路徑，可以對(duì)小樣本的數(shù)據(jù)集進(jìn)行較快、有效的圖像分割，常用于醫(yī)學(xué)細(xì)胞圖像分割。文[3]成功應(yīng)用了U-Net深度神經(jīng)網(wǎng)絡(luò)模型來(lái)識(shí)別單天線射電望遠(yuǎn)鏡數(shù)據(jù)中的射頻干擾，其結(jié)構(gòu)如圖14。

圖14 U-Net架構(gòu) [3]Fig.14 U-Net architecture [3]

文[5]發(fā)現(xiàn)使用卷積神經(jīng)網(wǎng)絡(luò)對(duì)500 m口徑球面射電望遠(yuǎn)鏡數(shù)據(jù)進(jìn)行標(biāo)記射頻干擾時(shí)，常常標(biāo)記錯(cuò)誤，需要進(jìn)一步人工檢查，帶來(lái)了顯著的額外工作量。為了克服這一問(wèn)題，文[5]提出了一種稱為RFI-Net的神經(jīng)網(wǎng)絡(luò)模型，如圖15，結(jié)果表明，RFI-Net模型優(yōu)于U-Net, KNN(K-nearest Neighbour)和SumThreshold算法。

圖15 RFI-Net架構(gòu) [5]Fig.15 RFI-Net architecture [5]

循環(huán)神經(jīng)網(wǎng)絡(luò)(Recurrent Neural Network, RNN)是一類用于處理序列數(shù)據(jù)的神經(jīng)網(wǎng)絡(luò)，對(duì)處理上下文關(guān)聯(lián)的語(yǔ)音[31]數(shù)據(jù)有很好效果。天文信號(hào)也是一種時(shí)序數(shù)據(jù)，非常適合使用循環(huán)神經(jīng)網(wǎng)絡(luò)對(duì)其進(jìn)行處理，文[32]利用循環(huán)神經(jīng)網(wǎng)絡(luò)對(duì)射電望遠(yuǎn)鏡干涉陣列數(shù)據(jù)進(jìn)行射頻干擾檢測(cè)，根據(jù)巨型米波射電望遠(yuǎn)鏡(Giant Metrewave Radio Telescope, GMRT)610 MHz的射頻干擾振幅信息，區(qū)分射頻干擾和非射頻干擾數(shù)據(jù)。

5 總 結(jié)

本文系統(tǒng)闡述了射電天文中抑制和標(biāo)記射頻干擾算法和技術(shù)，分析了在主動(dòng)預(yù)防階段、預(yù)相關(guān)階段、后相關(guān)階段射頻干擾抑制方法的優(yōu)勢(shì)和不足。采用有效的屏蔽手段可以在主動(dòng)預(yù)防階段阻止大部分射頻干擾進(jìn)入系統(tǒng)；預(yù)相關(guān)階段基于各個(gè)天線中的天文信號(hào)具有相關(guān)特性抑制射頻干擾；在后相關(guān)階段大部分算法基于閾值標(biāo)記射頻干擾。通過(guò)機(jī)器學(xué)習(xí)相關(guān)技術(shù)對(duì)射頻干擾進(jìn)行標(biāo)記和識(shí)別是目前的研究熱點(diǎn)，基于海量天文數(shù)據(jù)進(jìn)行訓(xùn)練能夠極大提高標(biāo)記射頻干擾的正確率，但不足是前期手動(dòng)標(biāo)記訓(xùn)練樣本中射頻干擾和天文信號(hào)十分耗費(fèi)時(shí)間，目前可用的訓(xùn)練樣本有限。

射頻干擾抑制問(wèn)題需要在不同階段采用多種方法協(xié)同解決，從射電望遠(yuǎn)鏡臺(tái)址周邊的無(wú)線電環(huán)境保護(hù)到最終天文數(shù)據(jù)處理過(guò)程各個(gè)環(huán)節(jié)都需要考慮射頻干擾問(wèn)題。各天文臺(tái)站需要結(jié)合實(shí)際需求和所處環(huán)境的電磁干擾情況，合理選擇射頻干擾抑制方法。