孫浩民,鄧 輝,梅 盈,衛(wèi)守林,戴 偉,王 鋒
(1. 廣州大學(xué)天體物理中心,廣東 廣州 510006;2. 昆明理工大學(xué),云南 昆明 650051)
近十年來我國在射電天文領(lǐng)域得到了飛速發(fā)展,21 cm望遠鏡陣列(The 21 CentiMeter Array, 21CMA[1])、新疆的天籟陣列[2]、內(nèi)蒙古的明安圖太陽頻譜日像儀(Mingantu Ultrawide Spectral Radioheliograph, MUSER)[3]以及中國甚長基線干涉測量(Very Long Baseline Interferometry, VLBI)網(wǎng),這些射電望遠鏡的建成在我國天文科學(xué)研究、數(shù)據(jù)處理、經(jīng)驗積累、人才培養(yǎng)等方面都做出了重要貢獻。平方千米陣(Square Kilometre Array, SKA)是人類正在建設(shè)的最大規(guī)模射電干涉望遠鏡,由13個成員國共同建設(shè),中國是其中7個發(fā)起國之一。
觀測數(shù)據(jù)的存儲和交換是天文觀測和研究的基本要求。FITS(Flexible Image Transport System, FITS)一直是天文數(shù)據(jù)存儲交換的標準格式,針對射電觀測數(shù)據(jù),在FITS基礎(chǔ)上發(fā)展出了UVFITS和FITSIDI等格式。近年來,測量集應(yīng)用越來越廣泛,逐漸成為射電天文數(shù)據(jù)處理分析的標準格式,被CASA、WSCLEAN[4]等射電天文數(shù)據(jù)處理軟件廣泛支持。測量集在國內(nèi)射電天文領(lǐng)域應(yīng)用相對較少,國內(nèi)射電望遠鏡往往根據(jù)各自接收機的特點,自行定義相應(yīng)的原始數(shù)據(jù)存儲格式。如明安圖太陽頻譜日像儀采用祼二進制的方式保存觀測文件,以大幅度降低存儲空間,需要進行數(shù)據(jù)交換時,通過格式轉(zhuǎn)換軟件轉(zhuǎn)換為UVFITS或FITSIDI。
射電天文模擬校準成像庫是目前正在研制中的平方千米陣的算法參考庫(Algorithm Reference Library, ARL)。為了實現(xiàn)射電天文模擬校準成像庫與主流射電天文數(shù)據(jù)處理軟件(如CASA)的數(shù)據(jù)對接,需要解決測量集格式輸出問題。本文結(jié)合實際需求,系統(tǒng)討論并實現(xiàn)了利用Python和Python-casacore生成測量集格式文件,并將Python代碼集成到射電天文模擬校準成像庫,為平方千米陣科學(xué)數(shù)據(jù)處理研究提供支撐,對其它射電望遠鏡數(shù)據(jù)處理工作也有較高的參考價值。
測量集是一個遵從射電干涉測量方程(Radio Interferometer Measurement Equation, RIME[5])的文件格式標準,在AIPS ++ Note 191中被正式定義,用于規(guī)范校準前的射電天文觀測數(shù)據(jù)的存儲。測量集標準發(fā)布以后,CASA團隊和歐洲VLBI網(wǎng)團隊的多個天文軟件開發(fā)小組進行了代碼實現(xiàn)。由于CASA采用測量方程作為其基本校準方案,測量集很自然地成為CASA觀測數(shù)據(jù)的存儲標準。隨著CASA成為阿塔卡瑪大型毫米波天線陣和甚大陣(Very Large Array, VLA)的指定數(shù)據(jù)處理分析軟件包[6],測量集也成為它們數(shù)據(jù)分析中的缺省數(shù)據(jù)格式。但是,阿塔卡瑪大型毫米波天線陣和甚大陣的原始數(shù)據(jù)存儲格式分別為阿爾馬科學(xué)數(shù)據(jù)模型(ALMA Science Data Model, ASDM)和科學(xué)數(shù)據(jù)模型(Science Data Model, SDM),因此也都開發(fā)了相應(yīng)接口,實現(xiàn)ASDM/SDM與測量集的轉(zhuǎn)換。
測量集適用于目前射電天文學(xué)中的所有用例,包括單碟、少數(shù)天線構(gòu)成的簡單干涉儀以及成千上萬個天線構(gòu)成的大型射電干涉陣。測量集借鑒了關(guān)系型數(shù)據(jù)庫的建模方法來降低數(shù)據(jù)的冗余,采用主表和子表、主鍵和外鍵,把干涉得到的可見度函數(shù)或單天線總功率測量值及其時間戳保存在主表(建主鍵),反復(fù)使用的元數(shù)據(jù)保存在子表(建外鍵),通過外鍵實現(xiàn)對元數(shù)據(jù)的引用。主表中有多個數(shù)據(jù)列以及關(guān)聯(lián)到子表的鍵值,主表中必須有DATA列,用于存放干涉陣的可見度數(shù)據(jù),或者FLOAT_DATA列,用于存放單天線功率值。
CASA遵循了測量集第2版本(MSv2)[7],為了確保與CASA的數(shù)據(jù)兼容,本文在測量集第2版本基礎(chǔ)上開展研究工作。
子表中存儲了測量集的關(guān)鍵字,表1列出了CASA采用的測量集第2版的所有子表,其中括號里的子表為可選子表。在實際應(yīng)用中,必須生成非可選子表,即必須有的子表,子表內(nèi)容可以為空。每個測量集文件一定要有一個主表(MAIN),表中包含數(shù)據(jù)列和各子表的鍵。
表1 測量集的子表[7]
顯然,與FITS相比,測量集格式要復(fù)雜很多。在實際應(yīng)用中,根據(jù)射電望遠鏡的不同觀測數(shù)據(jù)需求,生成相應(yīng)的子表并存入相關(guān)的數(shù)據(jù),最終構(gòu)成一個完整的測量集格式文件目錄樹結(jié)構(gòu),原則上必須生成所有非可選表。
與其它天文文件存儲格式不同,測量集格式采用多級目錄多文件保存的方法,各個表都以CASA表格存儲。這意味著一個表格包含了多個文件,整個測量集格式文件也不是單個文件,而是一個由多級目錄構(gòu)成的目錄樹。一般來說,主表位于第1級目錄,各個子表位于第2級目錄。每一級目錄均包含實際數(shù)據(jù)存放位置信息的table.info, table.f0, table.f1等文件。
每個測量集文件必須有一個主表(MAIN TABLE),結(jié)構(gòu)見表2。受限于篇幅,表2僅列出了主表的部分字段,完整的主表字段請參考測量集技術(shù)規(guī)范[7]。
表2 測量集的主表結(jié)構(gòu)[7]
注:Nc為相關(guān)器數(shù);Nf為頻率通道數(shù);Ncat為分類數(shù);*表示不同的天線類型有不同的取值,對于本文討論的射電干涉陣來說,都是相關(guān)參數(shù)的數(shù)量。
由表2可見,在數(shù)據(jù)存儲方面,測量集格式與FITS格式基本類似,需要保存的數(shù)據(jù)類型也包括整型(int)、浮點(float)、雙精度(double)、字符串(string)等。與FITS文件的頭定義相比,測量集格式文件字段設(shè)計更為復(fù)雜,有3種類型的字段:
(1)關(guān)鍵字(Keywords):MS_VERSION,用來標識所保存的測量集格式文件遵從哪一個版本的規(guī)范。
(2)鍵(Key):相當于關(guān)系型數(shù)據(jù)庫中的主鍵,用來和子表進行關(guān)聯(lián)。如TIME給出了觀測時刻,ANTENNA1,ANTENNA2指定兩個相關(guān)天線。
(3)非鍵屬性(Non-key attributes):根據(jù)實際需要定義的重要參數(shù)或?qū)傩浴?/p>
其它表結(jié)構(gòu)與主表類似,也各自規(guī)定了相應(yīng)的保留字、鍵值和參數(shù)。生成測量集格式文件,必須首先明確各字段的內(nèi)容、格式、單位等。
從AIPS++發(fā)展到CASA后,CASA軟件的開發(fā)采用了多計算機語言混合編程的方法,其核心代碼來自于C/C++,用戶接口與調(diào)用部分基本上采用Python實現(xiàn),C/C++部分構(gòu)成了CASA核心庫,即所謂的Casacore。Casacore是目前最完整的射電天文數(shù)據(jù)處理軟件包,也是唯一實現(xiàn)了測量集文件的讀寫操作的軟件包。Python-casacore(1)https://github.com/casacore/casacore實現(xiàn)了Python對Casacore的調(diào)用。本文深入分析了Python-casacore的基本函數(shù),研究其數(shù)據(jù)表的讀寫方法,完成了一個測量集格式輸出對象,并集成到射電天文模擬校準成像庫中,實現(xiàn)射電天文模擬校準成像庫的測量集格式輸出功能。
Python通過Python-casacore訪問Casacore,表3列出了 Python-casacore的重要測量集表寫入函數(shù)。
表3 Python-casacore的測量集表操作函數(shù)
調(diào)用表3中相應(yīng)的函數(shù)或者類,可以生成一個完整的數(shù)據(jù)子表,主要代碼見表4:
(1)生成矩陣列col1,定義 ‘UVW’ 字段,浮點類型,單位為m,該列存儲一維數(shù)組。
(2)生成矩陣列col2,定義 ‘FLAG’ 字段,布爾類型,該列存儲二維數(shù)組。
(3)生成表定義desc,該表有4列,分別為col1, col2, col6, col7。然后按照表定義desc創(chuàng)建空表tb。
(4)把數(shù)據(jù)寫入表tb,循環(huán)次數(shù)等于波段數(shù)量。
本文利用Python-casacore,針對測量集第2版開發(fā)實現(xiàn)了測量集文件輸出模塊,封裝成WriteMS類,結(jié)構(gòu)見圖1。
表4 子表生成主要代碼
圖1 WriteMs類圖
Fig.1 The class of WriteMs
目前,開發(fā)的程序已經(jīng)集成到平方千米陣的算法參考庫RASCIL(https://gitlab.com/timcornwell/rascil/tree/master/rascil/processing_components/visibility),同時,在射電天文模擬校準成像庫中給出了一個利用WriteMS輸出測量集格式的例子,供使用者參考,見test_export_ms_rascil.py。圖2是該應(yīng)用實例的流程圖,包括了模擬觀測到測量集輸出的完整過程,主要步驟解釋如下:
(1)Input parameters。輸入?yún)?shù)主要包括觀測位置、觀測時間、觀測頻率、頻帶寬度、相位中心、輸出文件的名稱以及干涉陣列天線配置。
(2)Input the image of M31。實例中使用了阿塔卡瑪大型毫米波天線陣和其它模擬程序通常采用的M31星系圖作為模擬觀測圖像,見圖3(a)。
(3)Create block of visibilities。根據(jù)步驟(1)所輸入的參數(shù)和天線配置生成初始可見度數(shù)組,預(yù)測寬視場成像參數(shù),得到合適的柵格尺寸大小。
(4)Sample M31。對M31圖像進行模擬采樣,得到觀測可見度數(shù)組。
(5)Export visibilities to MS files。把重采樣得到的可見度函數(shù)按測量集格式輸出。
為了檢驗本文方法的正確性,使用CASA軟件對所生成的測量集文件成圖。即把圖2流程所生成的測量集格式的可見度數(shù)據(jù)導(dǎo)入CASA,利用CASA成圖(圖3(b)),可以看出所成的臟圖與原始圖像(圖3(a))輪廓基本一致。因為CASA軟件在讀取測量集格式數(shù)據(jù)的過程中會進行較多的校驗操作,因此,如果所生成的測量集格式數(shù)據(jù)通過了CASA的處理,得到了正確的結(jié)果,說明所生成的測量集各個子表和字段是合理的(雖然部分字段的取值可能與真實情況有區(qū)別,但這不影響成像的處理)。
圖2 測量集文件生成應(yīng)用實例流程圖
Fig.2 The flow chart of exporting MS format file
圖3 M31的原始圖像與模擬觀測圖像的比較。(a) 原始圖像;(b) 模擬觀測臟圖
Fig.3 The original M31 image and the observed M31 image. (a) The original image; (b)The simulated image
雖然測量集文件規(guī)范制定較早,但在我國的射電領(lǐng)域應(yīng)用較少。一方面是因為測量集格式文件占用較多的空間,另一方面是因為生成測量集格式文件一直依賴于Casacore這一底層軟件包,開發(fā)比較困難。本文結(jié)合平方千米陣工程建設(shè)的需要,研究了測量集格式文件的定義、結(jié)構(gòu)、字段以及利用Python-casacore寫入數(shù)據(jù)的方法,并基于Python-casacore開發(fā)實現(xiàn)了測量集格式輸出軟件包,集成到射電天文模擬校準成像庫。整體來看,本文的工作在平方千米陣算法參考庫的研制過程中起到了重要的支撐作用,為后續(xù)數(shù)據(jù)模擬與文件存儲提供了保障,也對其它射電天文數(shù)據(jù)的測量集格式文件生成有較好的參考作用。
致謝:感謝國家天文臺-阿里云天文大數(shù)據(jù)聯(lián)合研究中心對本項工作的支持。