行星際閃爍單站單頻數(shù)據(jù)處理

劉麗佳，彭 勃

(中國科學(xué)院國家天文臺，北京 100012)

行星際閃爍(Interplanetary Scintillation，簡稱IPS)是指來自遙遠(yuǎn)致密射電源的電磁波在穿過太陽系行星際空間時(shí)，被太陽風(fēng)等離子體散射形成的射電波強(qiáng)度和相位的一種隨機(jī)起伏現(xiàn)象。通過觀測這種閃爍現(xiàn)象，可以研究致密射電源的角徑，也可以測量太陽風(fēng)的速度和研究太陽風(fēng)等離子的不規(guī)則結(jié)構(gòu)[1-4]，是兼有理論研究和實(shí)際應(yīng)用價(jià)值的重要課題。

國際上英國、前蘇聯(lián)和日本都相繼開展了行星際閃爍測量用以研究太陽風(fēng)及瞬變現(xiàn)象[5-7]，印度上世紀(jì)70年代也開始了行星際閃爍觀測[8]。日本的三站觀測系統(tǒng)可直接測出投影太陽風(fēng)風(fēng)速，印度的Ooty射電望遠(yuǎn)鏡利用單站單頻模式觀測致密星體，通過譜擬合的方式得到的閃爍功率譜可測得太陽風(fēng)風(fēng)速及閃爍指數(shù)。

我國的行星際閃爍觀測研究開始于上世紀(jì)90年代末期[9]，北京天文臺利用密云綜合孔徑射電望遠(yuǎn)鏡的復(fù)合干涉儀對行星際閃爍觀測進(jìn)行了初步試驗(yàn)和研究，根據(jù)所得閃爍頻譜得到了太陽風(fēng)風(fēng)速及閃爍指數(shù)，并與同期印度Ooty觀測數(shù)據(jù)進(jìn)行了對比，符合較好。但由于采用單站單頻模式進(jìn)行觀測Fresnel膝不明顯，不好確定對應(yīng)的頻率，對估測速度造成困難。此種方式對模型選取依賴較大，存在一定誤差。

此后國際上發(fā)展了單站雙頻模式對行星際閃爍進(jìn)行觀測[10-11]。這種模式通過對歸一化互相關(guān)譜第一零點(diǎn)的測量來計(jì)算太陽風(fēng)風(fēng)速，同時(shí)也可以研究電子密度起伏的空間譜，比起單站單頻模式準(zhǔn)確性更高。

中國科學(xué)院國家天文臺正在密云基地50m射電望遠(yuǎn)鏡上建立一套單站雙頻行星際閃爍觀測系統(tǒng)。此系統(tǒng)服務(wù)于國家重大科學(xué)項(xiàng)目——子午工程。該系統(tǒng)建成投入使用后，將使國家天文臺密云基地成為國際上唯一的多組雙頻行星際閃爍觀測站。

由于執(zhí)行探月任務(wù)，目前密云基地只能在較高頻(S/X波段)進(jìn)行單站單頻的觀測。針對目前設(shè)備情況，作者研究編寫了開展行星際閃爍單站單頻觀測的數(shù)據(jù)處理軟件。本文的第一部分介紹行星際閃爍單站單頻觀測的相關(guān)理論，第二部分討論單站單頻數(shù)據(jù)處理，最后對數(shù)據(jù)處理軟件進(jìn)一步分析。

1 行星際閃爍單站單頻模式

行星際閃爍單站單頻模式是指在一個(gè)觀測站使用單一頻率進(jìn)行IPS觀測，用這種模式測量太陽風(fēng)風(fēng)速有兩種方法，第一種稱為參數(shù)擬合法，第二種稱為特征頻率法。參數(shù)擬合法依賴于太陽風(fēng)的幾個(gè)主要參數(shù)：太陽風(fēng)電子密度功率譜的冪指數(shù)α，各向異性軸比AR，太陽風(fēng)風(fēng)速V，射電源角尺度θ0，選擇適當(dāng)?shù)奶栵L(fēng)模式，計(jì)算出模式譜，使其與實(shí)測譜吻合，從而確定太陽風(fēng)和射電源的上述參數(shù)。第二種方法是通過實(shí)測譜確定兩個(gè)特征頻率：Fresnel頻率fF和第一最小頻率fmin，根據(jù)兩個(gè)特征頻率中的任一個(gè)來計(jì)算太陽風(fēng)風(fēng)速，公式如下[5]：

(1)

(2)

其中λ是觀測頻率；Z是太陽與射電源視線方向之間的最短距離。圖1為一行星際閃爍單站單頻閃爍譜示例。

圖1 單站單頻閃爍譜示例Fig.1 IPS spectrum observed at a single station for single frequencies

2 單站單頻數(shù)據(jù)處理

2.1 數(shù)據(jù)處理介紹

圖2 數(shù)據(jù)處理流程圖Fig.2 Flow chart of data reduction

為了得到行星際閃爍譜，用Visual Fortran語言編寫了針對單站單頻觀測模式的數(shù)據(jù)處理軟件，具體流程如圖2。圖中所示的觀測數(shù)據(jù)采集和記錄系統(tǒng)是一個(gè)實(shí)時(shí)顯示系統(tǒng)。在觀測過程中數(shù)據(jù)可以被實(shí)時(shí)顯示在電腦顯示屏上，使觀測數(shù)據(jù)的質(zhì)量實(shí)時(shí)可見。在遇到強(qiáng)干擾時(shí)，可以通過調(diào)節(jié)望遠(yuǎn)鏡接收機(jī)終端上的步進(jìn)衰減或者輻射計(jì)的增益來控制數(shù)據(jù)的輸出數(shù)值范圍，將數(shù)據(jù)控制在線性范圍內(nèi)，以免數(shù)據(jù)因?yàn)楦蓴_而溢出產(chǎn)生失真。數(shù)據(jù)記錄軟件記錄范圍為0～65535，線性范圍為0～32000。為使數(shù)據(jù)不失真，當(dāng)天線不對準(zhǔn)目標(biāo)源(off-source)觀測時(shí)，調(diào)整步進(jìn)衰減和輻射計(jì)增益，使數(shù)據(jù)顯示在2000～5000的范圍內(nèi)；天線對準(zhǔn)目標(biāo)源(on-source)觀測時(shí)，用采集軟件觀察采集的數(shù)據(jù)，使數(shù)據(jù)不超過30000。

數(shù)據(jù)處理軟件濾波去干擾包括3個(gè)步驟：原始數(shù)據(jù)成圖、去掉緩變成分和去野點(diǎn)。圖3是一段密云50m射電望遠(yuǎn)鏡在S波段觀測3C48的原始數(shù)據(jù)圖示例。

選擇數(shù)據(jù)進(jìn)行去掉緩變成分處理方法是：根據(jù)采樣頻率、行星際閃爍現(xiàn)象的性質(zhì)和觀測數(shù)據(jù)情況，選取一定長度數(shù)據(jù)，用滑動平均方法計(jì)算出每段數(shù)據(jù)的平均值，然后用原始數(shù)據(jù)減去每段數(shù)據(jù)的平均值。在本圖中首先選取每段數(shù)據(jù)的長度為30s，進(jìn)行一次試去緩變成分處理。密云觀測時(shí)積分時(shí)間為0.1s，因此每段數(shù)據(jù)的長度為300個(gè)點(diǎn)。

進(jìn)行去緩變成分處理之后是去野點(diǎn)處理。所謂野點(diǎn)是由于干擾等因素出現(xiàn)的高出正常值的數(shù)據(jù)。野點(diǎn)的存在會使閃爍譜結(jié)果不準(zhǔn)確。去野點(diǎn)的方法是：根據(jù)采樣頻率和實(shí)際觀測數(shù)據(jù)情況選擇一定長度的數(shù)據(jù)，計(jì)算這段數(shù)據(jù)的均方根(rms)，然后使這段數(shù)據(jù)中的每一個(gè)值同一個(gè)設(shè)定的閾值相比較(通常是選擇3rms或5rms)，絕對值大于這個(gè)閾值的點(diǎn)被其前后兩個(gè)點(diǎn)的平均值所取代。圖4是對原始數(shù)據(jù)去緩變并同3rms比較后所得。

圖3 密云50m射電望遠(yuǎn)鏡IPS觀測3C48原始數(shù)據(jù)Fig.3 Raw IPS data toward 3C48 obtained with the Miyun 50m radio telescope；x axis:recorder moments;y axis:ralative values of voltages

圖4 IPS數(shù)據(jù)去緩變、野點(diǎn)示例(3rms比較)Fig.4 Plot for IPS after the subtraction of slow changing component and the removal of excessive changing values (“wild point”)

去野點(diǎn)過程中選擇半分鐘的數(shù)據(jù)來計(jì)算rms值，閾值選定為3rms，比較后，共去掉638個(gè)點(diǎn)，占數(shù)據(jù)總量的1.47%。

軟件去干擾后是濾波再采樣、減去平均值和三角加權(quán)。設(shè)軟件去干擾后的數(shù)據(jù)為DATA。在減去平均值這個(gè)步驟中選擇一定長度的窗函數(shù)，計(jì)算這段數(shù)據(jù)的長度，然后減去平均值。在減去平均值時(shí)采用滑動平均的方法(同軟件去干擾中去掉緩變成分類似)?；瑒悠骄低瑯舆x擇一定長度的數(shù)據(jù)進(jìn)行平均。記減去平均值后的數(shù)據(jù)為DATA1。減去平均值后是三角加權(quán)，在這一步中選擇與上一步相同長度的三角函數(shù)，三角函數(shù)的兩端為0，中間數(shù)值最大，兩側(cè)對稱，使三角函數(shù)與DATA1相乘，所得數(shù)據(jù)記為DATA2。在下一步中DATA2將進(jìn)行傅里葉變換，數(shù)據(jù)從時(shí)域變換到頻域，頻域中的數(shù)據(jù)記錄為DATA3，DATA3再進(jìn)行自相關(guān)后即為閃爍譜。圖5是密云50m望遠(yuǎn)鏡觀測得到的一個(gè)閃爍譜。

圖5 密云50米射電望遠(yuǎn)鏡觀測到3C48的閃爍譜;2008年3月27日；S波段Fig.5 IPS Spectrum observed with the Miyun 50m radio telescope through the S band on March 27,2008

圖6 2005年9月Ooty望遠(yuǎn)鏡觀測3C279數(shù)據(jù)原始圖Fig.6 Raw data toward 3C279 observed by the Ooty radio telescope in September 2005

密云現(xiàn)有的輻射計(jì)積分時(shí)間過長且增益可調(diào)范圍小，且IPS現(xiàn)象一般在較低頻率觀測明顯，目前不適合進(jìn)行行星際閃爍觀測。在密云IPS專用輻射計(jì)研制期間，我們已開始申請烏魯木齊25m望遠(yuǎn)鏡開展試驗(yàn)觀測研究，取得了初步成功(在另文討論)。

2.2 數(shù)據(jù)處理驗(yàn)證

為了對自主編寫的數(shù)據(jù)處理軟件進(jìn)行驗(yàn)證，2009年3月，聯(lián)系了印度IPS組的Prof Monoharan P K,請他提供一些Ooty實(shí)測數(shù)據(jù)和處理結(jié)果，以檢驗(yàn)處理軟件。Prof Monoharan提供了2005～2007年間印度Ooty射電望遠(yuǎn)鏡觀測到的行星際閃爍原始數(shù)據(jù)及對應(yīng)的閃爍譜圖。

使用我們的軟件對多組印度Ooty望遠(yuǎn)鏡的觀測數(shù)據(jù)進(jìn)行了處理。處理結(jié)果與Ooty提供的處理結(jié)果符合很好，證明我們的數(shù)據(jù)處理達(dá)到了預(yù)期目標(biāo)。

圖6～8為處理的Ooty觀測數(shù)據(jù)對比圖。圖中觀測源為3C279，觀測時(shí)間是2005年9月2日，觀測頻率為327MHz。其中圖6為Ooty望遠(yuǎn)鏡所觀測的原始數(shù)據(jù)；圖7為該數(shù)據(jù)印度提供的線性坐標(biāo)閃爍譜圖與利用自主編寫的數(shù)據(jù)處理軟件得到的線性坐標(biāo)閃爍譜圖對比；圖8為利用自主編寫的數(shù)據(jù)處理軟件得到的對數(shù)坐標(biāo)閃爍譜圖及擬合圖，圖中實(shí)線為觀測所得閃爍譜圖，虛線為擬合圖。

(a) (b)

圖7 2005年9月Ooty望遠(yuǎn)鏡觀測3C279的IPS閃爍譜
(a)印度提供的閃爍譜；(b)自編寫數(shù)據(jù)處理軟件閃爍譜
Fig.7 IPS spectrum toward 3C279 observed by the Ooty radio telescope，(a) provided by the Indian group,(b) reduced by our own software

圖8 自編數(shù)據(jù)處理軟件的Ooty望遠(yuǎn)鏡數(shù)據(jù)閃爍譜及參數(shù)擬合結(jié)果(橫坐標(biāo)取對數(shù))觀測時(shí)間：2005年9月。Fig.8 IPS Spectrum from the data obtained with the Ooty radio telescope in September,2005 and reduced by our own software;x axis:frequenuy;y axis：logarithm of spectral power

由圖7可見，印度和我們的兩種數(shù)據(jù)處理方法得出的結(jié)果符合得很好。

對2005年Ooty望遠(yuǎn)鏡的觀測結(jié)果進(jìn)行了參數(shù)擬合，模型模擬結(jié)果如圖8。模型采用的參數(shù)為：AR=2.5，α=3.5，V=600km/s。由fF=1.1Hz和公式(1)得到V為672.4km/s,兩者在誤差允許范圍(～10%)內(nèi)吻合很好。 spaceweather.com公布的當(dāng)天太陽風(fēng)速為672.4km/s,與觀測的結(jié)果基本一致。

通過對Ooty射電望遠(yuǎn)鏡IPS觀測的數(shù)據(jù)處理，以及與spaceweather.com相應(yīng)時(shí)間觀測結(jié)果比較，可以得出以下結(jié)論：數(shù)據(jù)處理方法是正確的，數(shù)據(jù)處理結(jié)果也達(dá)到預(yù)期。

3 結(jié) 論

編寫的數(shù)據(jù)處理軟件適用于行星際閃爍單站單頻模式的數(shù)據(jù)處理。當(dāng)進(jìn)行單站雙頻模式觀測時(shí)，可以對此軟件進(jìn)行相應(yīng)的拓展：觀測得到的數(shù)據(jù)先分別進(jìn)行以上處理流程，在傅里葉變換后，得到的兩組數(shù)據(jù)不僅要做自相關(guān)(記為DATA11，DATA22)，同時(shí)還要進(jìn)行互相關(guān)運(yùn)算(記為DATA12，DATA21)。得到的4組數(shù)據(jù)再進(jìn)行減平均值的計(jì)算就可以得到歸一化互相關(guān)譜。通過計(jì)算歸一化互相關(guān)譜第一零點(diǎn)頻率可估算太陽風(fēng)風(fēng)速。

行星際閃爍的單站單頻模式，無論從觀測、數(shù)據(jù)采集設(shè)備，還是數(shù)據(jù)流程來講都較簡單，易于操作，被國際上多數(shù)觀測站采用。但是，這種觀測模式要求有較高信噪比(大于25dB)和靈敏度。在觀測中幾乎所有觀測到的閃爍都可以找到fF，這是Fresnel過濾因子與電子密度冪函數(shù)譜相互作用的結(jié)果。進(jìn)行模式擬合時(shí)，由fF可知選擇合適的太陽風(fēng)速度值。另外，由于路徑積分效應(yīng)，F(xiàn)resnel振蕩現(xiàn)象已不明顯，極小值退變成平臺，只能觀測到第一極小值fmin，而其它的極小值總是被淹沒。當(dāng)各向異性軸比AR增大時(shí)，fF就會變得不明顯，從而對估測速度造成困難。同時(shí)fmin也容易受到噪聲、AR和ε的影響，因此參數(shù)擬合的準(zhǔn)確性容易受到所選參數(shù)的影響而存在一定的不確定性[12]。

相對于單站單頻模式，單站雙頻模式具有以下幾個(gè)優(yōu)點(diǎn)：更高的靈敏度，更精確的測量太陽風(fēng)速，對于各種太陽風(fēng)參數(shù)具有更高的穩(wěn)定性[6]。但是這種模式對觀測系統(tǒng)和接收系統(tǒng)有較高要求，同時(shí)數(shù)據(jù)處理流程也較單站單頻模式復(fù)雜。

目前國家天文臺密云基地50m射電望遠(yuǎn)鏡開展的單站雙頻模式的接收機(jī)系統(tǒng)已處于研制過程中。

致謝：感謝印度Prof Monoharan P K 提供Ooty望遠(yuǎn)鏡觀測數(shù)據(jù)，特別感謝張喜鎮(zhèn)老師的有益討論和指導(dǎo)。

[1] Hewish A,Symonds M D.Radio investigation of the solar plasma.Planetary and Space Science[J].Nature,1969,17:313.

[2] Armstrong J W,Coles W A.Analysis of three-station interplanetary scintillation[J].Geophysics Res,1972,77:4602-4610.

[3] Watanabe T,Shibasaki K,Kakinuma T.Latitudinal distribution of solar wind velocity and its relation to solar EUV corona[J].Geophys Res,1974,79:3841-3843.

[4] Hewish A,Scott P F,Wills D.Radio investigation of the solar plasma[J].Nature,1964,203:1214.

[5] Purvis A,Tappin S J,Rees W G,et al.The Cambridge IPS survey at 81.5 MHz[J]. Royal Astronomical Society,Monthly Notices,1987,229:589-619.

[6] Vitkevich V V，Glushaev A A,Iliasov Iu P,et al.Antenna equipment of the Lebedev Institute BSA radio telescope facility[J].Radiofizika，1976，19:1594-1606.

[7] Kojima M，Asai，Kozuka Y，et al.Velocity observations at high latitude and the acceleration phenomena[J].Advances in Space Research，1995,16(9):101-110.

[8] Swarup G，Sarma N V G，Joshi M N,et al.Large Steerable Radio Telescope at Ootacamund,India[J].Nature Physical Science,1971,230:185-188.

[9] Wu J H,Zhang X Z,Zheng Y J.IPS Observations at Miyun Station,Bao[J].ApSS,2001，278:189-192.

[10] Scott L,Rickett B J,Armstrong J W.The velocity and the density spectrum of the solar wind from simultaneous three-frequency IPS observations[J].Astronomy and Astrophysics,1983,123:191-206.

[11] Zhang X Z.A Study on the Technique of Observing Interplanetary Scintillation with Simultaneous Dual-frequency Measurement[J].Chin J Astron & Astrophys,2007, 7(5):712-720.

[12] 葉品中，邱育海.行星際閃爍(IPS)單站觀測對太陽風(fēng)速度的診斷[J].天體物理學(xué)報(bào)1996,16(4):389-394.