下载宽场无缝光谱 (WFSS) 数据#
本笔记本使用 Python 的 astroquery.mast Observations 类，通过 MAST API 查询特定程序的特定数据产品。我们正在寻找 NGDEEP 程序（ID 2079）的 NIRISS 成像和 WFSS 文件。观测使用三种 NIRISS 滤光片：F115W、F150W 和 F200W，采用 GR150R 和 GR150C 光栅。一个 WFSS 观测序列 通常由一幅直接图像和同一阻挡滤光片下的光栅观测组成，以帮助识别场中的源。在程序 2079 中，曝光序列遵循以下模式：直接图像 -> GR150R -> 直接图像 -> GR150C -> 直接图像。
用例：使用 MAST 下载数据产品。
数据：来自程序 2079 观测 004 的 JWST/NIRISS 图像和光谱。
工具：astropy、astroquery、glob、matplotlib、numpy、os、pandas、(yaml)
跨仪器：所有
内容
作者：Camilla Pacifici (cpacifici@stsci.edu) & Rachel Plesha (rplesha@stsci.edu) & Jo Taylor (jotaylor@stsci.edu)
首次发布：2024 年 5 月
本笔记本受到 关于 MAST 的 JWebbinar 会议 的启发。
导入#
from astropy.io import fits # 导入FITS文件处理库
from astroquery.mast import Observations # 导入用于查询MAST（微波天文科学数据中心）的库
from matplotlib import pyplot as plt # 导入绘图库用于数据可视化
import numpy as np # 导入NumPy库用于数值计算
import os # 导入操作系统接口库
import glob # 导入用于文件路径匹配的库
import pandas as pd # 导入Pandas库用于数据处理和分析
查询观测数据#
在astroquery.mast中的观察类用于下载JWST数据。使用元数据（metadata）函数查看可用的搜索选项及其描述。
请注意，对于JWST，仪器名称具有特定格式。有关更多信息，请访问：https://outerspace.stsci.edu/display/MASTDOCS/JWST+Instrument+Names
# 获取观测数据的元数据
Observations.get_metadata("observations") # 调用Observations类的get_metadata方法，获取指定观测的元数据
下载特定数据集的两种最常见方法是使用 提案 ID (proposal ID) 或使用 观测 ID (observation ID)。
使用提案ID搜索#
# 选择提案ID、仪器和一些有用的关键词（在这种情况下是滤光片）。
obs_table = Observations.query_criteria(obs_collection=["JWST"], # 查询JWST观测集合
instrument_name=["NIRISS/IMAGE", "NIRISS/WFSS"], # 查询NIRISS图像和WFSS仪器
provenance_name=["CALJWST"], # 查询已执行的观测
filters=['F115W', 'F150W', 'F200W'], # 查询特定的滤光片
proposal_id=[2079], # 查询特定的提案ID
)
print(len(obs_table), 'files found') # 打印找到的文件数量
# 查看查询结果中感兴趣的几个列名的内容
obs_table[['obs_collection', 'instrument_name', 'filters', 'target_name', 'obs_id', 's_ra', 's_dec', 't_exptime', 'proposal_id']] # 显示特定列的数据
使用观测ID搜索#
观测ID（obs_id）允许通过提案ID和观测ID灵活搜索，因为JWST文件名的结构使然。有关JWST文件命名约定的更多信息，请访问：https://jwst-pipeline.readthedocs.io/en/latest/jwst/data_products/file_naming.html。在本笔记本系列中，我们将仅使用程序2079中的两个观测之一（004）。
此外，在搜索条件中使用通配符也提供了灵活性。例如，我们可以指定“NIRISS*”，而不是同时指定“NIRISS/IMAGE”和“NIRISS/WFSS”，这样就可以同时匹配这两种文件模式。通配符在obs_id中也有效，因此我们不必列出所有不同的ID。
# 从特定的观测ID列表中获取下载文件的列表
obs_id_table = Observations.query_criteria(instrument_name=["NIRISS*"], # 查询使用NIRISS仪器的观测
provenance_name=["CALJWST"], # 查询已执行的观测数据
obs_id=['jw02079-o004*'], # 搜索PID 2079的观测004
)
# 这个数字会随着JWST管道和参考文件的更新而变化
print(len(obs_id_table), 'files found') # 打印找到的文件数量，预计约为613个文件
过滤和下载产品#
如果要下载的文件过多，API 将会超时。相反，最好将观测数据分批下载，一次下载一个。
batch_size = 5 # 每次最多下载5个文件，以最大化下载速度。
# 将文件列表 ``obs_id_table`` 按照批次大小进行分组。
obs_batches = [obs_id_table[i:i+batch_size] for i in range(0, len(obs_id_table), batch_size)]
print("有多少个批次?", len(obs_batches)) # 输出批次数量
single_group = obs_batches[0] # 获取第一个批次，方便检查所获得的文件
print("检查第一个批次以确保它符合你想要下载的预期:")
# 显示所需的列以检查数据
single_group['obs_collection', 'instrument_name', 'filters', 'target_name', 'obs_id', 's_ra', 's_dec', 't_exptime', 'proposal_id']
选择所需的产品类型。不同的级别如下：
未校准文件 (uncalibrated files)
productType=[“SCIENCE”]
productSubGroupDescription=[‘UNCAL’]
calib_level=[1]
速率图像 (rate images)
productType=[“SCIENCE”]
productSubGroupDescription=[‘RATE’]
calib_level=[2]
光谱和成像的二级关联 (level 2 associations for both spectroscopy and imaging)
productType=[“INFO”]
productSubGroupDescription=[‘ASN’]
calib_level=[2]
成像的三级关联 (level 3 associations for imaging)
productType=[“INFO”]
productSubGroupDescription=[‘ASN’]
dataproduct_type=[“image”]
calib_level=[3]
下载前的数据过滤#
# 创建一个字典，用于存储上述信息，以便检查过滤函数
file_dict = {
'uncal': { # 未校准数据
'product_type': 'SCIENCE', # 产品类型为科学数据
'productSubGroupDescription': 'UNCAL', # 产品子组描述为未校准
'calib_level': [1] # 校准级别为1
},
'rate': { # 速率数据
'product_type': 'SCIENCE', # 产品类型为科学数据
'productSubGroupDescription': 'RATE', # 产品子组描述为速率
'calib_level': [2] # 校准级别为2
},
'level2_association': { # 级别2关联数据
'product_type': 'INFO', # 产品类型为信息
'productSubGroupDescription': 'ASN', # 产品子组描述为关联
'calib_level': [2] # 校准级别为2
},
'level3_association': { # 级别3关联数据
'product_type': 'INFO', # 产品类型为信息
'productSubGroupDescription': 'ASN', # 产品子组描述为关联
'calib_level': [3] # 校准级别为3
},
}
# 查看每个不同级别的现有文件
files_to_download = [] # 初始化一个空列表，用于存储要下载的文件
for index, batch_exposure in enumerate(single_group): # 遍历每个单组曝光数据
print('*'*50) # 打印分隔线
print(f"Exposure #{index+1} ({batch_exposure['obs_id']})") # 打印当前曝光的索引和观察ID
# 从列表中提取产品名称以进行过滤
products = Observations.get_product_list(batch_exposure) # 获取当前曝光的产品列表
for filetype, query_dict in file_dict.items(): # 遍历文件类型和查询字典
print('File type:', filetype) # 打印当前文件类型
# 过滤产品以获取符合条件的文件
filtered_products = Observations.filter_products(products,
productType=query_dict['product_type'], # 产品类型
productSubGroupDescription=query_dict['productSubGroupDescription'], # 产品子组描述
calib_level=query_dict['calib_level'], # 校准级别
)
print(filtered_products['productFilename']) # 打印过滤后的产品文件名
files_to_download.extend(filtered_products['productFilename']) # 将文件名添加到下载列表中
print() # 打印空行以便于阅读
print('*'*50) # 打印分隔线
从上面可以看出，在观察列表 (obs_id_table) 中，每个曝光名称都有许多相关文件需要下载。这就是我们需要分批下载的原因。
下载数据#
要实际下载产品，请使用 Observations.download_products() 并提供过滤后的产品列表。
通常情况下，detector1 流水线 不需要进行调整，因此我们可以从 detector1 的输出，即速率文件（rate files）开始，而不是未校准文件（uncal files）。因此，我们只需要下载速率文件和关联文件。如果需要重新运行 detector1 流水线，则需要在 Observations.filter_products 调用中调整 productSubGroupDescription 和 calib_level 以下载未校准文件。
如果数据是专有的，您可能还需要设置您的 API 令牌。绝对不要将您的令牌提交到公共代码库。另一种方法是创建一个单独的配置文件（config_file.yaml），该文件仅对您可读，并包含一个键： ‘mast_token’ : API 令牌
要创建新的 API 令牌，请访问以下链接：
https://auth.mast.stsci.edu/token?suggested_name=Astroquery&suggested_scope=mast:exclusive_access
请注意，下载数据时需要使用版本 >= 0.4.7 的 astroquery，以便调用 flat=True。如果您更喜欢使用早期版本，请从调用中删除该行，下载数据，并将所有下载的子文件夹中的文件移动到由 download_dir 变量定义的单个目录中。
# 检查astroquery库的版本是否高于0.4.7。有关更多信息，请参见上面的说明
import astroquery # 导入astroquery库
astroquery.__version__ # 获取并显示astroquery库的当前版本
import os # 导入os模块，用于处理文件和目录
download_dir = 'data' # 定义下载目录的名称
# 确保下载目录存在；如果不存在，则创建一个新目录
if not os.path.exists(download_dir): # 检查下载目录是否存在
os.mkdir(download_dir) # 创建下载目录
# 现在让我们获取每批观测的产品，并筛选出仅感兴趣的产品。
for index, batch in enumerate(obs_batches): # 遍历每个观测批次
# 进度指示器...
print('\n'+f'Batch #{index+1} / {len(obs_batches)}') # 打印当前批次的索引和总批次数
# 从我们的Astropy观测表中获取`obsid`标识符，以获取产品。
obsids = batch['obsid'] # 获取当前批次的obsid列表
print('Working with the following obsids:') # 打印正在处理的obsid信息
for number, obs_text in zip(obsids, batch['obs_id']): # 遍历每个obsid及其对应的文本
print(f"{number} : {obs_text}") # 打印obsid及其描述
# 获取产品列表
products = Observations.get_product_list(obsids) # 获取与obsid相关的产品列表
# 筛选产品以仅获取感兴趣的产品
filtered_products = Observations.filter_products(products, # 过滤产品列表
productType=["SCIENCE", "INFO"], # 仅选择科学和信息类型的产品
productSubGroupDescription=["RATE", "ASN"], # 不使用"UNCAL"，因为我们可以从速率文件开始
calib_level=[2, 3] # 不使用1，因为不获取UNCAL文件
)
# 下载这些记录的产品。
# `flat=True`选项将所有文件存储在指定的`download_dir`的单个目录中。
manifest = Observations.download_products(filtered_products, # 下载筛选后的产品
download_dir=download_dir, # 指定下载目录
flat=True, # astroquery v0.4.7或更高版本仅适用
)
print('Products downloaded:\n', filtered_products['productFilename']) # 打印下载的产品文件名
如果继续使用WFSS笔记本，我们需要再次确认是否已下载所有剩余笔记本所需的文件。应该已经下载了149个文件。
import glob # 导入glob模块，用于文件路径匹配
import os # 导入os模块，用于操作系统相关的功能
# 获取下载目录中所有以.fits和.json结尾的文件
downloaded_files = glob.glob(os.path.join(download_dir, '*.fits')) + glob.glob(os.path.join(download_dir, '*.json'))
# 打印下载的文件数量和下载目录
print(len(downloaded_files), 'files downloaded to:', download_dir)
检查下载的数据#
此功能的目的是让您更好地了解可用的数据。此外，您将能够使用此数据框（dataframe）选择特定的文件，以便更深入地查看您想要关注的模式（mode）。
ratefile_datadir = 'data/' # 定义存放速率文件的目录
# 首先查找所有已下载的速率文件
rate_files = np.sort(glob.glob(os.path.join(ratefile_datadir, "*rate.fits"))) # 获取所有速率文件并排序
for file_num, ratefile in enumerate(rate_files): # 遍历每个速率文件
rate_hdr = fits.getheader(ratefile) # 获取每个速率文件的主头信息
# 我们想要存储的信息，这些信息可能对我们后续评估数据有用
temp_hdr_dict = {"FILENAME": ratefile, # 文件名
"TARG_RA": [rate_hdr["TARG_RA"]], # 目标的赤经
"TARG_DEC": [rate_hdr["TARG_DEC"]], # 目标的赤纬
"FILTER": [rate_hdr["FILTER"]], # 使用的滤光片; GR150R/GR150C
"PUPIL": [rate_hdr["PUPIL"]], # 使用的光阑; F090W, F115W, F140M, F150W F158M, F200W
"EXPSTART": [rate_hdr['EXPSTART']], # 曝光开始时间 (MJD)
"PATT_NUM": [rate_hdr["PATT_NUM"]], # 在WFSS中抖动模式内的位置编号
"NUMDTHPT": [rate_hdr["NUMDTHPT"]], # 整个抖动模式中的点总数
"XOFFSET": [rate_hdr["XOFFSET"]], # NIRISS的X偏移量（从模式起始位置，单位：角秒）
"YOFFSET": [rate_hdr["YOFFSET"]], # NIRISS的Y偏移量（从模式起始位置，单位：角秒）
}
# 将字典转换为pandas数据框
if file_num == 0: # 如果这是第一个文件，创建初始数据框
rate_df = pd.DataFrame(temp_hdr_dict) # 创建数据框
else: # 否则，为每个文件附加到数据框
new_data_df = pd.DataFrame(temp_hdr_dict) # 创建新的数据框
# 合并两个数据框，创建一个包含所有数据的数据框
rate_df = pd.concat([rate_df, new_data_df], ignore_index=True, axis=0) # 合并数据框
rate_dfsort = rate_df.sort_values('EXPSTART', ignore_index=False) # 按曝光开始时间排序数据框
# 将数据框保存到文件，以便后续读取
outfile = './list_ngdeep_rate.csv' # 输出文件名
rate_dfsort.to_csv(outfile, sep=',', index=False) # 保存数据框为CSV文件
print('Saved:', outfile) # 打印保存信息
# 查看结果数据框
rate_dfsort # 显示排序后的数据框
特别地，让我们来看一下这些速率文件遵循的观测序列。我们已经按曝光时间对上述文件进行了排序，因此它们在数据框中应该已经是按时间顺序排列的。
FILTER = CLEAR 表示直接图像，而 FILTER=GR150R 或 FILTER=GR150C 表示分散图像。PUPIL 是使用的阻挡滤光片。前14个曝光构成了第一组序列：直接图像 -> 光栅 -> 直接图像 -> 光栅。对于分散图像和直接图像，还有多个抖动位置。多个直接图像的抖动将合并在 image3 中，而多个分散图像可以合并在 spec3 中。
# 显示 rate_df 数据框中指定列的前14行
rate_df[['EXPSTART', 'FILTER', 'PUPIL', 'PATT_NUM', 'XOFFSET', 'YOFFSET']].head(14)
下面展示的是前14个速率文件，以便直观了解上述序列。网格线作为视觉指导，帮助识别任何进行的抖动（dithers）。
# plot set up
fig = plt.figure(figsize=(20, 35)) # 创建一个20x35英寸的图形
cols = 3 # 每行的子图数量
rows = int(np.ceil(14 / cols)) # 计算需要的行数，向上取整
# loop over the first 14 rate files and plot them
for plt_num, rf in enumerate(rate_dfsort[0:14]['FILENAME']): # 遍历前14个速率文件
# determine where the subplot should be
xpos = (plt_num % 40) % cols # 计算子图在当前行中的位置
ypos = ((plt_num % 40) // cols) # 计算子图在当前列中的位置，使用整除以获得整数
# make the subplot
ax = plt.subplot2grid((rows, cols), (ypos, xpos)) # 创建子图
# open the data and plot it
with fits.open(rf) as hdu: # 打开FITS文件
data = hdu[1].data # 获取数据
data[np.isnan(data)] = 0 # 将NaN数据填充为0，以帮助matplotlib的颜色刻度
ax.imshow(data, vmin=0, vmax=1.5, origin='lower') # 显示数据图像，设置颜色范围
# adding in grid lines as a visual aid
for gridline in [500, 1000, 1500]: # 添加网格线作为视觉辅助
ax.axhline(gridline, color='black', alpha=0.5) # 添加水平网格线
ax.axvline(gridline, color='black', alpha=0.5) # 添加垂直网格线
ax.set_title(f"#{plt_num+1}: {hdu[0].header['FILTER']} {hdu[0].header['PUPIL']} Dither{hdu[0].header['PATT_NUM']}") # 设置子图标题
