如何叠加两个具有不同 WCS/分辨率的 .fits 图像？

Question

我正在尝试从一个物体的一张图像（比如说在波段 A 中）绘制轮廓，在同一物体的较低分辨率图像之上（比如说在波段 Z 中）。

两个图像都是大尺寸文件，因此我必须为每个图像创建一个 2D 剪切图。然而，带 Z 中的图像的像素远少于带 A 的像素。因此，当我尝试在同一个图中绘制它们时，带 Z 中的图像根据像素尺寸绘制，因此左下角非常小该图（参见下面链接中的图）。我需要根据天文角度尺寸而不是像素来绘制这两个图像。

import numpy as np
import matplotlib
%matplotlib inline
import matplotlib.pyplot as plt
from astropy.io import fits
from astropy import units as u
from astropy.coordinates import SkyCoord
from astropy.wcs import WCS
from astropy.nddata import Cutout2D
from reproject import reproject_interp

# Read in .fits images and WCS information
image_a = fits.open("image_a.fits")
image_z = fits.open("image_z.fits")
wcs_a = WCS(image_a.header).celestial
wcs_z = WCS(image_z.header).celestial

# Define object RA and DEC
obj_coords = SkyCoord(ra=135.19081*u.degree, dec=0.68991393*u.degree, frame='fk5')

# Create 2D cutouts of the object in each band in a 6 by 6 arcsec box
size = u.Quantity((6, 6), u.arcsec)
stamp_a = Cutout2D(image_a.data, obj_coords, size, wcs=wcs_a)
stamp_z = Cutout2D(image_z.data, obj_coords, size, wcs=wcs_z)

# Plot the the stamp in band Z with contours from band A
ax = plt.subplot(projection = wcs_z)
plt.imshow(stamp_z.data, cmap='gray', origin='lower', alpha=1)
ax.contour(stamp_a.data, levels = [2*sigma,3*sigma,4*sigma,5*sigma], 
colors='green')

结果图（点我）

我还看到了一些关于使用重新投影的内容。我尝试过，但不认为我完全理解它，因为它只会产生一个很小的像素。这是我尝试过的：

array, footprint = reproject_interp((stamp_a.data, wcs_a), image_z.header)
plt.imshow(array, origin='lower')

如果我缺少任何信息来回答这个问题，请告诉我。在使用 IDL 几年后，我上周刚刚切换到 python，所以如果它有点粗糙，我深表歉意。

Answer 1

如何使用reproject 的问题是您通过了(stamp_a.data, wcs_a)，但wcs_aWCS 是来自原始图像，而不是来自图章。

您可以从Cutout2D图像中获取与您的图章相匹配的 WCS 对象。我认为更改为(stamp_a.data, stamp_a.wcs)会给你一个正确的结果。

你看过astropy.visualization.wcsaxes了吗？你可以在这里看到一个例子查看如何绘制图像以及如何使用不同的 WCS 覆盖第二张图像的轮廓的如果这对您有用，它将比使用Cutout2D和reproject手动更简单。

我不确定性能；如果您的图像很大，您可能仍然想像现在一样进行剪切。如果内存/执行时间不是问题，您可以只绘制整个图像，然后设置您想要的范围wcsaxes（请参阅他们的教程文档）。

是自己进行重新投影还是让幕后matplotlib进行重新投影也是一个品味问题。astropy.visualisation.wcsaxes直接使用reproject代码会更多一些，也更复杂，但可以让您更好地控制所使用的精确重投影方法（例如插值或较慢的通量守恒reproject_exact），并且可能更容易完全理解正在发生的事情。