NIRCam Wisp 移除#
作者: Ben Sunnquist (bsunnquist@stsci.edu)
最新更新: 2024年7月24日
使用案例: NIRCam 成像探测器 A3, A4, B3, 和 B4。
数据: PID 1063 观测 6 成像平场
测试管道版本: 1.15.1
引言#
本笔记本演示了如何从 NIRCam 成像数据中去除 光晕。光晕是影响 A3、A4、B3 和 B4 探测器的散射光特征。对于给定的滤光片，光晕在同一探测器位置出现，仅其亮度在不同曝光之间变化；因此，可以通过缩放和减去光晕模板（即所有光晕出现的中位数组合）来从科学数据中去除它们。
本笔记本中使用的光晕模板将在 数据 部分下载，但它们也可以在 NIRCam 光晕模板 Box 文件夹 的 “version3” 文件夹中找到。未来对这些模板的更新将添加到同一 Box 文件夹中，鼓励用户定期检查此文件夹，以确保他们拥有最新版本。
本笔记本使用 subtract_wisp.py 代码来缩放和减去光晕。该代码可以单独在 Python 中使用，如果需要并行校准大量文件，建议使用该代码，但本笔记本将用于演示可用的各种参数，以优化光晕去除。对于每个笔记本单元，我们还将展示在 Python 中运行等效命令。
目录#
导入#
如何创建一个环境来运行这个笔记本：
conda create -n jwst1140 python=3.11 notebook
source activate jwst1140
pip install -r requirements.txt
在安装之前，可以在 requirements.txt 文件中更新此笔记本使用的管道版本。
# 从subtract_wisp.py代码中加载函数
from subtract_wisp import make_segmap, process_file, subtract_wisp
# 加载重要的包
from astropy.io import fits # 用于处理FITS文件
from astroquery.mast import Mast, Observations # 用于查询和获取JWST数据
import glob # 用于文件路径操作
import os # 用于操作系统相关功能
import urllib # 用于处理URL
import matplotlib.pyplot as plt # 用于绘图
%matplotlib inline # 在Jupyter Notebook中内联显示图形
数据#
我们将使用来自程序1063第6次观测的F200W成像数据，这些数据来自探测器A3、A4、B3和B4。这是6幅稀疏场的图像，受到了明亮光晕的影响。我们还将使用相应的长波图像进行源检测，因为长波不受光晕的影响。
请使用astroquery从MAST下载相关数据。
警告：如果在此步骤中该笔记本被中断，下载的文件可能不完整，可能会导致后续崩溃！
# 从MAST选择相关数据
params = {"columns": "dataURI, filename, exp_type", # 请求的列，包括数据URI、文件名和曝光类型
"filters": [{"paramName": "program", "values": ['1063']}, # 过滤条件：程序编号为1063
{"paramName": "observtn", "values": ['6']}, # 过滤条件：观测编号为6
{"paramName": "exposure", "values": ['00005']}, # 过滤条件：曝光编号为00005
{"paramName": "visit", "values": ['004']}, # 过滤条件：访问编号为004
{"paramName": "detector", "values": ['NRCA3', 'NRCA4', 'NRCB3', 'NRCB4', 'NRCALONG', 'NRCBLONG']}, # 过滤条件：探测器类型
{"paramName": "productLevel", "values": ['2b']}] # 过滤条件：产品级别为2b
}
t = Mast().service_request('Mast.Jwst.Filtered.Nircam', params) # 向MAST发送请求以获取数据
for row in t: # 遍历返回的数据行
if '_cal' in row['filename']: # 只选择包含'_cal'的文件名，即校准文件
result = Observations().download_file(row['dataURI'], cache=False) # 下载文件，不使用缓存
接下来，我们将下载wisp模板。这些模板相当于NIRCam wisp模板Box文件夹中的“version3”文件夹。建议用户定期检查此文件夹以获取任何更新。在本笔记本中，仅需要F200W模板：
WISP_NRCA3_F200W_CLEAR.fits
WISP_NRCA4_F200W_CLEAR.fits
WISP_NRCB3_F200W_CLEAR.fits
WISP_NRCB4_F200W_CLEAR.fits
# 下载WISP模板
boxlink = 'https://data.science.stsci.edu/redirect/JWST/jwst-data_analysis_tools/nircam_wisp_templates/' # WISP模板的链接
for detector in ['NRCA3', 'NRCA4', 'NRCB3', 'NRCB4']: # 遍历每个探测器
boxfile = os.path.join(boxlink, 'WISP_{}_F200W_CLEAR.fits'.format(detector)) # 构建每个探测器的文件链接
urllib.request.urlretrieve(boxfile, os.path.basename(boxfile)) # 下载文件并保存为基本文件名
# 检查模板文件是否已下载
template_files = ['WISP_NRCA3_F200W_CLEAR.fits', 'WISP_NRCA4_F200W_CLEAR.fits',
'WISP_NRCB3_F200W_CLEAR.fits', 'WISP_NRCB4_F200W_CLEAR.fits']
# 遍历每个模板文件
for tem_file in template_files:
# 如果文件不存在
if not os.path.isfile(tem_file):
# 打印提示信息，告知用户文件缺失
print(f'file {tem_file} does not exist, please download it from Box')
运行wisp减法代码#
在本节中，我们将展示如何使用默认参数在所有短波图像上运行wisp减法代码。唯一的例外是我们将show_plot=True设置为允许在笔记本中显示诊断图。
主要处理函数process_file()将segmap创建和wisp缩放/减法步骤结合在一起。代码使用的基本过程如下：
使用相应的长波图像生成源掩模，因为长波不受wisp的影响。接下来的步骤中忽略这些像素。此步骤在
make_segmap()函数中进行。加载相关的wisp模板。
对模板应用一系列缩放因子。对于每个缩放因子，将其乘以wisp模板，并从输入科学数据中减去结果。记录在每个测试的缩放因子下wisp区域内的残余噪声。在此步骤中应用对残余1/f（水平）噪声的修正，以帮助拟合。对残余进行多项式拟合，并根据此拟合选择噪声最低的缩放因子。此步骤在
subtract_wisp()函数中进行。将选择的wisp因子乘以wisp模板，以生成最终的wisp模型。此步骤在
subtract_wisp()函数中进行。从输入数据中减去最终的wisp模型，并输出wisp减法后的数据、wisp模型以及总结结果的诊断图。此步骤在
subtract_wisp()函数中进行。
# 收集所有科学数据文件
files = glob.glob('*_cal.fits') # 使用glob模块获取所有以'_cal.fits'结尾的文件
# 过滤出短波文件，排除包含'long'的文件
files = [cal_file for cal_file in files if 'long' not in cal_file] # 只保留不包含'long'的文件
运行每个数据集中的 process_file() 函数。代码生成以下内容：
输入文件的去雾（wisp）校正版本，文件名与输入文件相同，只是在末尾添加了 _wisp 后缀（注意：这里使用的后缀可以通过
suffix参数更改，也可以设置为空字符串以覆盖输入文件）。从输入图像中减去的去雾模型，文件名与输入文件相同，只是在末尾添加了 _wisp_model 后缀。
一张诊断图，显示原始图像、去雾校正图像、减去的去雾模型，以及一张显示所有测试过的去雾缩放因子及其对应残差噪声的图。噪声最低的因子是应用于去雾模板以生成去雾模型的因子。
注意：调用 process_files() 而不是 process_file() 将会并行处理文件。这将节省时间，但图表不会在笔记本中内联显示。
# 从科学数据中移除光晕
for file in files: # 遍历文件列表中的每个文件
results = process_file(file, show_plot=True) # 处理当前文件并显示图表
Python 等效命令：
python subtract_wisp.py –files *_cal.fits
使用自定义设置运行wisp减法代码#
在本节中，我们将使用自定义设置运行wisp去除代码，这些设置可能有助于优化wisp的去除。这里显示的不同设置可以与其他示例中的设置结合使用，也可以单独使用。
示例 1：新的细丝缩放方法#
默认情况下，wisp因子是基于产生最低残余噪声的结果来决定的。然而，缩放也可以基于最小化wisp影响区域内外的整体信号水平差异，通过设置scale_method='median'来实现。对于这种缩放方法，建议将残差的多项式拟合关闭（poly_degree=0），以简单选择产生最低差异的因子。
在缩放过程中，默认情况下会纠正残余的1/f噪声（correct_rows=True）；然而，对于具有强烈的例如放大器偏移或奇偶列残差的数据集，这种相同的纠正也可以在垂直方向上应用（correct_cols=True）。请注意，这些纠正不会在输出文件中被减去，它们仅用于帮助wisp缩放过程。
下面，我们将使用这些设置进行新的运行。
# 使用缩放方法去除光晕
# 处理JWST数据文件，指定文件名和参数
results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', # 输入文件名
scale_method='median', # 使用中位数作为缩放方法
poly_degree=0, # 多项式的度数为0（即不使用多项式拟合）
correct_cols=True, # 是否校正列
show_plot=True) # 是否显示处理结果的图形
Python 等效命令：python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --scale_method median --poly_degree 0 --correct_cols
有时，例如在拥挤的场域中，可能没有足够的空白区域来正确缩放 wisp。在这种情况下，简单地减去未缩放的 wisp 模板可能会提供最佳结果（scale_wisp=False，或者在 Python 中类似地使用 --no-scale_wisp）。
示例 2：数据质量标记#
另一种选择是标记图像数据质量数组中任何具有高于某个阈值的 wisp 信号的像素，但不从科学数据中减去 wisp 本身。如果您希望在将多个图像合并为马赛克时忽略这些像素，或在进行光度测量时标记它们，这可能会很有用。这可以通过将 flag_wisp_thresh 设置为所选择的 wisp 信号阈值（单位假定与输入文件相同）并将 sub_wisp=False 来实现。这里使用的数据质量值由 dq_val 参数设置。
dq_val 的默认值为 1，相当于 JWST 流水线的 DO_NOT_USE 标志，并且在 image3 的 drizzling 过程中默认被忽略。如果您希望这些像素被标记，但在该步骤中不被忽略，另一种选择是使用 dq_val=1073741824，这相当于 JWST 流水线的 OTHER_BAD_PIXEL 标志。有关各种数据质量标志的更多详细信息，请参见 这里。
# 不移除光晕，只标记它们
results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', # 处理JWST的FITS文件
flag_wisp_thresh=0.03, # 设置光晕标记阈值为0.03
dq_val=1073741824, # 设置数据质量标志值
sub_wisp=False, # 不从数据中减去光晕
show_plot=True) # 显示处理结果的图形
Python 等效命令：python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --flag_wisp_thresh 0.03 --dq_val 1073741824 --no-sub_wisp
检查此文件中的数据质量数组是否已被适当更新。
# 检查数据质量 (DQ) 数组
dq = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_wisp.fits', 'DQ') # 从 FITS 文件中获取 DQ 数据
# 显示被标记为 OTHER_BAD_PIXEL 的像素
dq = (dq & 1073741824 != 0).astype(int) # 仅选择 DQ 值为 1073741824 的像素，并转换为整数类型
plt.imshow(dq, cmap='gray', origin='lower', vmin=0, vmax=0.1) # 使用灰度图显示 DQ 数据，设置原点在下方，vmin 和 vmax 控制显示范围
示例 3：其他杂项设置#
这个示例将展示一些其他设置，这些设置可能单独使用或与上述其他设置组合使用都很有用。完整的可选设置列表可以在 make_segmap() 和 subtract_wisp() 的文档字符串中查看。
在这个示例中，我们将在输入图像本身上进行源查找（而不是对应的长波图像），通过设置 seg_from_lw=False 来实现。我们将增加源检测的 sigma 值，以帮助避免标记出光晕本身，并通过设置 save_segmap=True 来保存结果分割掩模以供检查。分割掩模的文件名与输入文件相同，只是在末尾添加了 _seg 后缀。
我们还将仅减去信号强度超过 0.01 MJy/sr 的光晕（min_wisp=0.01）。为了节省一些时间，我们将仅测试在 0.5 到 1.5 之间以 0.05 为步长的光晕尺度因子。
# 获取make_segmap()的所有可选设置列表，使用help(make_segmap)
#
# 使用特殊设置
results = process_file('jw01063006004_02101_00005_nrcb4_cal.fits', # 处理指定的JWST数据文件
seg_from_lw=False, # 不从长波段生成分段图
sigma=1.5, # 设置信号噪声比的阈值
save_segmap=True, # 保存生成的分段图
min_wisp=0.01, # 设置最小的光斑强度
factor_min=0.5, # 设置最小因子
factor_max=1.5, # 设置最大因子
factor_step=0.05, # 设置因子的步长
show_plot=True) # 显示生成的图表
Python 等效命令：python subtract_wisp.py --files jw01063006004_02101_00005_nrcb4_cal.fits --no-seg_from_lw --sigma 1.5 --save_segmap --min_wisp 0.01 --factor_min 0.5 --factor_max 1.5 --factor_step 0.05
如上所示，这些自定义设置导致了过度减法。这可能是由于在wisp区域（例如，wisp本身被标记为一个源）中源查找不佳，因此建议使用长波长图像进行源查找。我们可以通过检查生成的分割图来确认这一点。
# 导入必要的库
from astropy.io import fits # 用于处理FITS文件
import matplotlib.pyplot as plt # 用于绘图
# 读取JWST数据文件
data = fits.getdata('jw01063006004_02101_00005_nrcb4_cal.fits') # 获取JWST图像数据
# 读取分割图文件
segmap = fits.getdata('jw01063006004_02101_00005_nrcb4_cal_seg.fits') # 获取分割图数据
# 显示图像
fig, axes = plt.subplots(1, 2, figsize=(20, 10)) # 创建一个1行2列的子图
# 显示原始图像数据
axes[0].imshow(data, origin='lower', cmap='gray', vmin=0.18, vmax=0.3) # 显示图像，设置灰度色图和显示范围
# 显示分割图
axes[1].imshow(segmap, origin='lower', cmap='gray', vmin=0, vmax=0.1) # 显示分割图，设置灰度色图和显示范围
使用自定义 WISP 模板#
本节展示了如何在使用主代码中包含的wisp减法程序的同时，结合自定义的wisp模板。如果您有自己的wisp模板，或者想对现有模板应用手动缩放因子，这可能会很有用。
在这个例子中，我们将使用现有模板，但将其乘以我们自己的缩放因子，并关闭缩放程序（scale_wisp=False），这样代码就会直接应用给定的wisp模型。请注意，在使用自定义模板时，关闭缩放程序并不是必需的；用户可以提供自定义模板，同时保持缩放程序开启。由于在这个例子中缩放程序已关闭，源查找就不再相关，因此我们可以跳过分割图创建函数make_segmap()，直接运行subtract_wisp()。
# 获取WISP模板并应用手动缩放因子
wisp_data = fits.getdata('WISP_NRCB4_F200W_CLEAR.fits') * 1.05 # 从FITS文件中读取WISP数据并乘以1.05进行缩放
# 使用自定义WISP模板处理文件
results = subtract_wisp('jw01063006004_02101_00005_nrcb4_cal.fits', wisp_data, scale_wisp=False, show_plot=True) # 调用subtract_wisp函数处理指定的FITS文件，传入WISP数据，禁用缩放并显示图像
分段运行代码#
在本节中，我们将展示如何分步骤运行完整的WISP去除过程。这可能有助于提供更多的灵活性，或者在寻找适合您数据的最佳设置时使用。我们将使用多种自定义设置——有关更多详细信息，请参见函数文档字符串。
# 要处理的文件
cal_file = 'jw01063006004_02101_00005_nrcb4_cal.fits'
# 要使用的WISP模板。我们将再次使用“自定义”模板。
wisp_data = fits.getdata('WISP_NRCB4_F200W_CLEAR.fits') * 1.05 # 从文件中获取WISP数据并进行缩放
# 生成分割图
segmap_data = make_segmap(cal_file, sigma=1.0, npixels=8, dilate_segmap=10) # 创建分割图，使用指定的参数
# 缩放并减去WISP模板
results = subtract_wisp(cal_file, wisp_data=wisp_data, segmap_data=segmap_data, show_plot=True) # 从校准文件中减去WISP模板，并显示结果图