Bo Zhang's Homepage
..The universe is unfolding as it should..

2016-11-3

ALFALFA巡天数据处理系统学习手记

归档于: 天文空间科学, 天文软件 @ 4:44 pm

ALFALFA巡天是阿雷西博望远镜从2004秋季一直进行到这个十年之初的河外HI谱线“盲巡”项目,旨在对邻近宇宙中的中性氢作系统的普查。FAST的早期巡天无论从科学目标还是工作模式上都应该与ALFALFA很相似,所以入职以来就在集中精力学习ALFALFA的数据处理流程。虽然这套程序是专门为ALFALFA的特点而量身定做,其他望远镜不便直接照搬,但个中的原理还是相通的。

ALFALFA的这套程序是在IDL环境下开发的,最好在Linux系统下运行。程序的配置不算很麻烦,首先将代码包解压(具体下载地不明,因为这些东西都是组内互相流通的……),然后将这个文件拷贝到家目录下,并将其中第7行的“/home/dorado3/galaxy/”替换为解压目录。然后打开解压目录中的idl_alfa/子目录,找到“alfinit”文件,将解压目录添入其中的“.comp fsc_droplist.pro”一行并保存:

.comp /untarred directory/eggidl/coyote/fsc_droplist.pro

再将解压目录中的aoidl/gen/aodefdir.pro打开,将“return,’/home/dorado3/galaxy/aoidl/”改为:

return, ‘/untarred directory/aoidl/’

最后打开家目录下的.cshrc文件,加入如下一行就大功告成:

setenv IDL_STARTUP ~/.idlstartup

然后在csh环境下进入IDL,并输入如下两行命令,调出相关程序包即可。

@wasinit
@alfinit

整个数据处理流程大致分为4步,分别是在观测当天完成的格式转化、第一阶段的流量校准/带通改正/基线校准/射频干扰去除、第二阶段的DataCube生成与信号提取,还有第三阶段的星表编纂/科学研究/虚拟天文台对接/论文发表等等。对于射电谱线数据处理来说,这样的步骤应该是普适的,只是各软件所用的具体算法可能有所不同。由于手头数据所限,现在本人只能亲身体验第二阶段生成DataCube之后的步骤,不过前面的流程还是在此大致概述一下。

ALFALFA巡天采用的是射电观测中很常见的“漂移”模式,望远镜固定不动,而天体在东升西落,期间接收机每秒记录一次数据。这样设置的好处是机械结构简单,毋须考虑跟踪精度问题。由于天空每秒转过15角秒,而ALFA接收机的波束宽度足有3个角分还多,所以每次数据记录期间,天空的位移是可以忽略不计的。在正常情况下,巡天的一次漂移持续10分钟,然后发射一次为时1秒、流量已知的人工噪声作为校准的基准,再开始下一次漂移。这最初的数据格式转化,除了要将原始的FITS文件转化为IDL的.sav格式以便进行中间处理之外,更重要的目的就是截取校准数据并为其编目。

校准数据是以3个一组的.sav文件的形式呈现的,包括持续1秒的人工噪声,以及作为参考的前一次漂移最后1秒和后一次漂移第1秒的观测数据。所以每次漂移最后都会生成4个文件,包括3个校准文件以及1个观测数据文件。这个数据文件有多大呢?ALFA接收机有7个波束外加1个空闲波束,每个波束对应2个偏振分量,每个偏振分量又分别引出4096个频谱通道,所以一次漂移就要记下(7+1)×2×4096×600个浮点数据,连带与之相伴的头文件。

文件格式转化的工作通常在天文台即已完成,随后生成的.sav文件被送往康奈尔大学,进行后续的处理。由于第一阶段的校准工作事关重大,通常是由ALFALFA巡天的领导者之一、康奈尔大学的教授Martha Haynes亲自上阵的,所以其他人能够拿到的数据大抵是已经进行过初期处理的加工品了。

校准工作最重要的目的是将接收机收到的电信号强度转化为有意义的物理单位。对于ALFALFA的第一阶段处理来说,这个单位选为天线温度,不过可以在第二阶段进一步转化为流量密度。这一步首先要整理校准文件的列表,以每个session(包括数十次漂移)为单位进行。期间只要观测无中断,校准文件就要三个一组按照时间顺序排好。

随后大致的步骤是用calib1生成校准用数据结构(内含校准用人工噪声信号,以及对噪声发射前后1秒的参考数据的平均)-> calib2给出接收机各个波束和偏振方向的流量转化系数 -> bpdgui校准接收机带通,扣除基线 -> flagbb人工检查射频干扰,并为数据设置权重 -> strip_pv平滑光谱。当然,这些东西说着轻松,真的要实际操作起来,需要留意的细节还是颇多的。然后就迎来了整个流程中的重头戏:第二阶段之始的DataCube生成。

因为ALFALFA巡天是全程无跟踪的,在进行谱线信号提取或生成天空图之前,自然要将彼此错位的数据对齐并重新拼合在一起。在对齐之前,首先要对原始数据进行反卷积运算,这是因为我们接收到的影像实际上是天体的真实情况同仪器的波束图之间卷积的结果。一般常用的是以高斯、sinc或两种函数的乘积来模拟单个波束的形态,按相应的分布在一定像素范围内取权重,计算出理论上的真实天空图。然后先前在flagbb中设置的权重在这里发挥了用场,被标为干扰的区域不再参与反卷积,而是被当成坏点忽略掉。此外,这一阶段还要进行第二次的单位转化——修正馈源增益相对望远镜天顶角的起伏,并将天线温度转化为流量密度。

数据对齐的方法示意

对齐过程需要的是多个session的数据输入,在耗时数小时的计算后,输出的是3维数组DataCube,也就是在2维空间坐标的基础上加上1维频谱信息。出于减轻计算负担的考虑,在默认情况下,ALFALFA生成DataCube的进程grid_prep会输出4个一组的文件,每组覆盖2.4度见方的天区,其在赤经与赤纬方向的步长皆为1角分。各个文件名对应相应天区的中心点坐标,分称grid_XXXX+XXa/b/c/d.sav,分别在频谱方向覆盖每秒-2000到3300、2500到7900、7200到12800,以及12100到17900千米的速度范围,各涵盖1024个频谱通道。

然后就可以正式进行信号提取了。首先导入所选文件,以“2156+07”为例:

restore, ‘grid_2156+07X.sav’

以上命令中的“X”请用a、b、c或d自行替换。然后使用grid_base进一步扣除数据对齐带来的额外基线,也就是衬在谱线之下的连续谱背景:

grid_base, grid, grid2

这里grid表示从grid_2156+07X.sav文件中提取的数据,grid2是扣除基线后的结果。这时屏幕上显示的是各个赤纬截面上每个偏振的情况。图中横轴为光谱通道,纵轴为赤经,其下方还配有相应的权重图,而最下方和右方的曲线分别是相应方向上的平均流量分布,所以最下方给出的就是平均频谱:

同时,终端中还会显示一系列的选项:

这其中最常用的是“a”,也就是通过点击相应区域的左下与右上角,来设置不参与拟合的部分。基线拟合是不需要将银河系本身发出的HI谱线包括在内的,而且相应的射频干扰也要一并去除,比如在阿雷西博数据的低频部分往往存在San Juan机场阴魂不散的空管雷达信号……

考虑从平均频谱图来看,这里不存在明显弯曲的背景连续谱,所以选择一阶多项式拟合基线,在终端中输入“f 1”即可。拟合完毕,一般在数据不太离谱的情况下,后面选择一个H加速,询问是否有需要跳过的频谱通道时再次选择银河系HI和射频干扰所在的范围,让程序自动处理剩下的数据就是了,否则就要逐个审阅144个赤纬截面×2个偏振分量上的288张图片……

不要忘记将扣除基线后的数据重新保存,文件名一般设为“gridb……”:

grid = grid2
save, grid, file=’gridb_2156+07X.sav’

然后是通过grid_flatfield扣除平场。这个“平场”似乎与光学观测中的平场定义不太一样,可以简单粗暴地理解为空间方向上的基线,固定赤纬,用多项式拟合赤经方向上的起伏。这样的起伏很大程度上是由数据对齐过程中的权重不均匀引发的,而且扣除方法也是做减法而非除法。

grid_flatfield, grid, grid2

输入这个命令,稍候片刻,屏幕上就会弹出这样的画面:

这其中,左上的窗口显示的是一条平均的频谱。如果自己没有理解错的话,应该就是对DataCube覆盖天区内各频谱通道图上每一点的数据取平均而得的结果。这里要用户人工选择需要进行自动拟合平场的通道范围,最主要的建议就是避开明显的谱线,如河内HI线附近。选好之后,这些区域内的图像一般就毋须操心了。

剩下的通道图就要人工逐一检视并确定拟合方案了,显示上都是一式两份,左侧是拟合前的情况,右侧是拟合后的。如果存在需要进一步扣除的干扰信号之类,也可以在此阶段去掉。在选择非拟合区域时,记得要在右侧选择,然后相应的选区会出现在左侧的方框中。

扣除完平场的文件一般命名为“gridbf……”并保存:

grid = grid2
save, grid, file=’gridbf_2156+07X.sav’

对a、b、c、d四个文件都完成扣除基线和平场的操作之后,可以运行ex3dh_cube,用模板匹配滤波法提取可能的河外HI信号:

ex3dh_cube, ‘2156+07′, SN=4.3

这一步毋须太多的人工干预,不过相对比较耗时。由此生成的是“cdh……”文件,内里记录下了找到的各源的位置与置信度。之后使用ex3d_d人工检查搜索结果并生成对应天区内的星表文件,又一次逐个检视赤纬切面:

restore,’gridbf_2156+07X.sav’
restore,’cdh_2156+07X.sav’
ex3d_d, grid, sources, res, axis=’dec’, step=2, catalog=’ex3d_2156+07X.cat’

这期间,所有由ex3dh_cube找到的HI源均按置信度标示在图上,颜色越冷则置信度越高,如下图中的26号源就是一个相对可靠的例子。一般来说,一个河外星系在多个频谱通道上都会有辐射,而且存在强度的变化;另外很多源在空间上也是具有延展形态的。所以如果频谱分布过窄或空间延展不足,可以将这个源视作不可信,对其进行删除(实为降低置信度)操作;而若认为程序有所疏漏,也可以人工增添新源。一条不错的建议是,在不确定的情况下不轻举妄动,宁可将可疑的源留到后面作进一步的分析。另外,河内HI辐射(下图中的0号源)以及HI线附近的疑似高速云(如下图中的11和12号源)也不要删掉,它们都有着各自的研究价值。

然后要写下相关日志:

restore,’gridbf_2156+07X.sav’
grid_q, grid, q

对4个文件都做完同样的操作之后,就可以调用gridview进行数据可视化了:

gridview, grid, cat=’ex3d_2156+07X.cat’

这里的ex3d_2156+07X.cat就是由ex3d_d生成的星表文件。切记,gridview这一步骤一定要使用大屏幕显示器!否则程序界面很可能会显示不完全。上图中,左上是相应天区在各个频谱通道上的天空图,右上是数据拼合过程中产生的权重图,右下是经由斯隆等星表调出的射电源光学影像,中下是天空中某点的辐射谱。由此还可以调用galflux测量星系的谱线流量之类,嗯,本人还没有试到这一步,日后酌情补充吧。

至于整个流程中的具体算法,最近可能还要再细细啃一啃,其他的射电数据处理软件也可能会尝试一下。为什么射电频段的数据处理各自为政现象如此严重,连阿雷西博内部都没能统一呢?还是光学波段的一个IRAF或是高能的一套HEASOFT一统天下省心啊……

No Comments

No comments yet.

RSS feed for comments on this post.

Sorry, the comment form is closed at this time.

首页 | 天文 | 科学 | 摄影 | 模型 | CV | 版权声明 | 联系站长
京ICP备05002854号-2 Powered by WordPress Version 2.0.6
Licensed under Creative Commons Licenses

porno izle