如何使用Python裁剪fits.gz文件并查看标注中心点的裁剪后图像

依赖库 #

• astropy：用于FITS文件读写、Cutout2D裁剪实现
• matplotlib：用于裁剪结果可视化与中心点标记
• numpy：数值处理基础依赖

完整实现代码 #

from astropy.io import fits
from astropy.nddata import Cutout2D
import matplotlib.pyplot as plt
from matplotlib.patches import Circle
import numpy as np

# 配置参数，替换为你的实际值
filename = "your_file.fits.gz"
x_center = 1024  # 替换为原始图像下的中心点X坐标
y_center = 2048  # 替换为原始图像下的中心点Y坐标
cutout_size = (2048, 2048)  # 裁剪尺寸（高，宽）
save_path = "cutout_result.fits"

# 1. 读取原始FITS文件
with fits.open(filename) as hdu_list:
    raw_data = hdu_list[0].data
    raw_header = hdu_list[0].header

# 2. 执行裁剪
# 若中心点靠近边缘需要允许不完整裁剪，可添加参数 mode='partial', fill_value=np.nan
cutout = Cutout2D(raw_data, (x_center, y_center), cutout_size)

# 3. 保存裁剪结果为FITS格式
# 绑定裁剪后的头信息，保证坐标体系正确性
new_hdu = fits.PrimaryHDU(data=cutout.data, header=cutout.to_header(raw_header))
new_hdu.writeto(save_path, overwrite=True)

# 4. 可视化核验，标注中心点
fig, ax = plt.subplots(figsize=(8,8))
# 对数拉伸适配天文图像动态范围，普通图像可直接用cutout.data替代np.log10(cutout.data)
ax.imshow(np.log10(cutout.data), cmap='gray', origin='lower')
# 绘制半径10像素的红色空心圆作为中心点标记，可按需调整尺寸、颜色
center_cutout = cutout.center_cutout
circle = Circle((center_cutout.x, center_cutout.y), radius=10, edgecolor='red', facecolor='none', linewidth=2)
ax.add_patch(circle)
plt.axis('off')
plt.show()

批量处理提示 #

如需处理多个文件，只需将上述逻辑放入遍历所有fits.gz文件的循环中，逐一传入对应文件路径和中心点坐标即可。

内容的提问来源于stack exchange，提问作者J. O'Connor