digital

下载 SDSS 上的天文数据

Contents
  1. 1. 下载 CSV 文件
    1. 1.1. 使用搜索表格
    2. 1.2. 使用 SQL 查询
      1. 1.2.1. 登录
      2. 1.2.2. 编辑 SQL 语句
      3. 1.2.3. 提交并下载 CSV 文件
  2. 2. 下载数据
    1. 2.1. Wget 下载
    2. 2.2. Python 程序下载
  3. 3. 说明

学习天文相关的同学们应该都会怎么去下载天文数据，但对于门外汉来说，刚一接触的确有点懵，好在网络上也有一些只言片语，巡天网站的文档也写得足够详细，但简单的中文经验的分享似乎不多。故本文就介绍了如何下载斯隆数字巡天的天文数据。

摘一段维基百科的介绍：斯隆数字化巡天（英语：Sloan Digital Sky Survey，缩写为SDSS）是使用位于新墨西哥州阿帕契点天文台的 2.5 米口径望远镜进行的红移巡天项目。该项目开始于 2000 年，以阿尔弗雷德·斯隆的名字命名，计划观测 25% 的天空，获取超过一百万个天体的多色测光资料和光谱数据。斯隆数字化巡天的星系样本以红移0.1 为中值，对于红星系的红移值达到 0.4，对于类星体红移值则达到 5，并且希望探测到红移值大于6的类星体。

根据其官网，SDSS 最近发布了截至 2016 年 7 月所收集的数据，Data Release 14 (DR14)，用户可以在 Skyserver 上查看文档和使用有关工具。

下载 CSV 文件

这步要做的是，在网站上下载包含 plate, mjd, fiberid 的 CSV 文件。有两种不同方法，可根据需要选择其中一种。

使用搜索表格

SDSS 提供了一种方便，可视化的数据检索方法，Search Form，用户可以通过指导手册 来了解如何使用此表格。

按照条件选择好参数，勾选 “Output Format” 为 CSV，点击 “Submit Query to SkyServer” 即可下载文件。

使用 SQL 查询

另外一种可定制性更好的方法：使用 Casjob 服务器查询。

登录

首先访问 Casjob，点击右上角的 Login 进行登录。

用户名：jiangbin
密码：123456

点击 “Submit” 登录。希望大家合理使用，不要恶意修改密码等，谢谢。

之后便进入一个页面，页面大部分是编辑器文本框，我们就在文本框中输入查询命令。

观察到编辑器上方有 Context 、MyScratch、Table (optional)、Task Name 这几项，右侧有 Syntax、Plan、Quick、Submit 这四个按钮。

将上方 Context 的下拉菜单选择为 “DR14”。

编辑 SQL 语句

在大文本框中输入 SQL 查询命令，在此举几个例子：

例 1

1
2
3
4
-- 选取查询到的前 2000 个恒星数据，包括 plate、mjd、fiberid、和测光数据，要求改条数据有对应的光谱数据，要求 plate 值小于 3000
select top 2000 s.plate, s.mjd, s.fiberid, p.psfMag_u, psfMag_g, psfMag_r, psfMag_i, psfMag_z, modelMag_u, modelMag_g, modelMag_r, modelMag_i, modelMag_z, petroMag_r, petroMag_i, petroMag_z
from star p, specobj s
where p.objid=s.bestobjid and s.plate < 3000

之所以限制 plate 值小于 3000，是因为发现大于 3000，即 SEGUE-2 后的数据总是有 404 错误，这点不知为何。期待解答。

例 2

1
2
3
4
-- 选取查询到的前 2000 个星系数据，包括 plate、mjd、fiberid、和测光数据，要求改条数据有对应的光谱数据，要求 plate 值小于 3000
select top 2000 s.plate, s.mjd, s.fiberid, p.psfMag_u, psfMag_g, psfMag_r, psfMag_i, psfMag_z, modelMag_u, modelMag_g, modelMag_r, modelMag_i, modelMag_z, petroMag_r, petroMag_i, petroMag_z
from galaxy p, specobj s
where p.objid=s.bestobjid and s.plate < 3000

例 3

1
2
3
4
-- 选取查询到的前 2000 个类星体数据，包括 plate、mjd、fiberid、和测光数据，要求改条数据有对应的光谱数据，要求 plate 值小于 3000
select top 2000 s.plate, s.mjd, s.fiberid, p.psfMag_u, psfMag_g, psfMag_r, psfMag_i, psfMag_z, modelMag_u, modelMag_g, modelMag_r, modelMag_i, modelMag_z, petroMag_r, petroMag_i, petroMag_z
from photoobj p, specobj s
where p.objid=s.bestobjid and s.plate < 3000 and (s.class = 'QSO')

例 4

1
2
3
4
-- 选取查询到的前 2000 个类星体数据，包括 plate、mjd、fiberid、和测光数据，要求改条数据有对应的光谱数据，要求 plate 值小于 3000，要求信噪比均值大于 5，红移值大于 0
select top 2000 s.plate, s.mjd, s.fiberid, p.psfMag_u, psfMag_g, psfMag_r, psfMag_i, psfMag_z, modelMag_u, modelMag_g, modelMag_r, modelMag_i, modelMag_z, petroMag_r, petroMag_i, petroMag_z
from photoobj p, specobj s
where p.objid=s.bestobjid and s.plate < 3000 and (s.class = 'QSO') and s.snmedian > 5 and s.z > 0

提交并下载 CSV 文件

输入好语句后，点击我们提到过的，文本框右上方的 Quick 按钮。提示 “Query complete!” 后发现页面下方出现了查询到的表格预览，最下方的 “Save As” 右侧有个下拉菜单，选择 “CSV”，然后点击 “Save As” 按钮，即可开始下载 CSV 文件。

下载数据

如果我们用单纯的文本编辑器打开下载好的 CSV 文件，我们会发现是类似这样的：

1
2
3
4
5
6
plate,mjd,fiberid,psfMag_u,psfMag_g,psfMag_r,psfMag_i,psfMag_z,modelMag_u,modelMag_g,modelMag_r,modelMag_i,modelMag_z,petroMag_r,petroMag_i,petroMag_z
324,51666,81,20.1458,18.91237,18.43639,18.22727,18.20093,20.14658,18.91115,18.42017,18.22268,18.19052,18.42331,18.26658,18.19617
324,51666,98,23.16555,19.7103,18.21049,17.35237,16.80854,23.2576,19.72584,18.22088,17.33814,16.81207,18.30366,17.35972,16.86364
324,51666,44,21.64006,19.30029,17.79194,16.48456,15.72073,21.64766,19.30349,17.79343,16.46885,15.7269,17.82307,16.48844,15.73929
324,51666,46,23.19713,19.87924,18.46723,17.16783,16.52462,23.20551,19.86484,18.43119,17.13482,16.4612,18.47289,17.20836,16.46233
324,51666,88,19.13607,18.05983,17.60652,17.40258,17.30837,19.12605,18.03544,17.57149,17.37517,17.27461,17.60952,17.42687,17.31413

前三列可以作为定位的信息，下载到相应的光谱数据，后面的列都是自己在 SQL 语句里指定要查询的测光数据。

Wget 下载

Linux 下一般都预装 wget 了，Windows 需要下载可执行文件，并在命令提示符下调用。下面以 Linux 为例，介绍下载方法。

首先使用 LiberOffice 打开 CSV 文件，只留下前三列的数据，另存为新的 CSV 文件。

访问 DAS，首先在表格里选择 All Data，然后点击文本框上方的 “Choose File” 按钮上传编辑好的 CSV 文件。然后点击 “Submit Request” 按钮提交。

接下来跳转到另一个页面，点击 “Download selection” 下的 “this form” 链接。在后续的页面中，确认勾选了 “spSpec”，下方是 “wget” 后，点击 “request” 按钮。

会下载到一个 .lis 格式的文件，移动该文件到想要下载包含光谱数据的 .fit 文件的文件夹中，在此工作目录下打开终端，运行 wget 进行下载：

1
wget -i sdss-wget-XXXXXX.lis

Python 程序下载

采用 wget 方法下载有诸多不便，比如不能照顾到测光数据与下载到的光谱数据之间的联系，对应起来是比较麻烦的。故利用一些现成的模块库，我编写了一个下载小程序。注意运行程序需要安装相应库：

1
pip3 install --user numpy pandas astroML

程序如下：

download.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
import numpy as np
from pandas import DataFrame, read_csv
from astroML.datasets import fetch_sdss_spectrum
import urllib
import pickle

# Read CSV to Dataframe

print("Loading csv files...\n")

str_df = read_csv('./stars.csv', sep=',', header=0)
glx_df = read_csv('./galaxies.csv', sep=',', header=0)
qso_df = read_csv('./quasars.csv', sep=',', header=0)

str_loc = str_df.values
glx_loc = glx_df.values
qso_loc = qso_df.values
all_loc = np.vstack((str_loc, glx_loc, qso_loc))
all_unq = np.unique(all_loc, axis=0)                   # remove the possible duplicate rows 

print("Completed!\n")

# Download and preprocess the data

print("Downloading spectrum...\n")

all_label = []
all_data = []
spec_i = 0
warn_d = 0

for all_ins in all_loc:
    try:
        plate, mjd, fiberid, psfu, psfg, psfr, psfi, psfz, modelu, modelg, modelr, modeli, modelz, petror, petroi, petroz = all_ins
        if plate < 3000:
            spec = fetch_sdss_spectrum(round(plate), round(mjd), round(fiberid))
            data_list = [spec.hdulist.filename(), round(psfu,5), round(psfg,5), round(psfr,5), round(psfi,5), round(psfz,5), round(modelu,5), round(modelg,5), round(modelr,5), round(modeli,5), round(modelz,5), round(petror,5), round(petroi,5), round(petroz,5)]
            all_data.append(data_list)
            spec_i = spec_i + 1
            print("Downloaded " + str(spec_i) + " spectra(s) ! ")
    except urllib.error.HTTPError or urllib.error.URLError:
        warn_d = warn_d + 1
        spec_i = spec_i + 1
        print("\nPassed a error at: " + str(spec_i) + "\n")

print("\nTried to download: " + str(spec_i) + " spectra(s).")
if warn_d > 0:
    warn_d = str(warn_d)
    print("Warning: " + warn_d + " spectra(s) failed to download !\n")
else:
    print("And there's no error.\n")
print("Completed!\n")

# Store the data to pickle file

with open('./data/listfile.data', 'wb') as filehandle: 
    # store the data as binary data stream
    pickle.dump(all_data, filehandle)

事实上，fetch_sdss_spectrum 不仅可以下载 fit 文件，也会读取该文件。有一点问题就是，如果首次运行程序，all_data 的 filename 那列有些会变成 <class> 这类的对象描述，而不是 fit 文件的路径，都下载好再运行一次（由于下载好了，有缓存，运行会很快）方可正常，不知道是为什么，期待解答。

这样的话 listfile.data 里就存了一个 list，其中每个项目都有光谱文件路径和对应的测光数据信息，使用会非常方便：

read.py
1
2
3
4
5
6
7
8
9
10
11
import pickle
import numpy as np

# Read Files
print("Reading fits files...\n")

with open('/home/chaos/Projects/SVM/python/final/data/listfile_01.data', 'rb') as filehandle: 
    # read the data as binary data stream
    list_data = pickle.load(filehandle)

print("Completed!\n")

说明

本文到此就结束了。文字原创，转载请注明出处。

Mastodon