WCS作为使用astropy加载的数据立方体切片的matplotlib投影？答案

【问题标题】：WCS as matplotlib projection for a data-cube slice loaded using astropy?WCS作为使用astropy加载的数据立方体切片的matplotlib投影？
【发布时间】：2016-08-27 08:57:40
【问题描述】：

我有一个名为“my_cube.fits”的带有 WCS 的 FITS 文件。该文件在轴 1 和 2（X 和 Y）上具有空间信息，在轴 3 (Z) 上具有光谱信息。当我使用 astropy.io.fits 加载它时，光谱轴为 0，空间轴为 1 和 2。文件加载如下：

    import astropy.io.fits as pyfits

    filename = 'my_cube.fits'
    my_data = pyfits.getdata(filename)
    my_header = pyfits.getheader(filename)

我一直在使用 matplotlib 显示数据，我想知道如何使用它的 WCS 显示我的数据立方体的单个光谱框架。比方说：

    from astropy.wcs import WCS
    from matplotlib import pyplot as plt

    my_wcs = WCS(my_header)
    fig = plt.figure()
    ax = fig.add_subplot(111, projection=my_wcs)
    ax.imshow(my_data[5, :, :])

    plt.show()

如果我这样做，我有：

    ...
    File "/usr/local/lib/python3.4/dist-packages/matplotlib/figure.py", line 1005, in add_subplot
    a = subplot_class_factory(projection_class)(self, *args, **kwargs)
    File "/usr/local/lib/python3.4/dist-packages/matplotlib/axes/_subplots.py", line 73, in __init__
    self._axes_class.__init__(self, fig, self.figbox, **kwargs)
    File "/usr/local/lib/python3.4/dist-packages/wcsaxes/core.py", line 49, in __init__
    self.patch = self.coords.frame.patch
    File "/usr/local/lib/python3.4/dist-packages/wcsaxes/frame.py", line 129, in patch
    self._update_patch_path()
    File "/usr/local/lib/python3.4/dist-packages/wcsaxes/frame.py", line 115, in _update_patch_path
    self.update_spines()
    File "/usr/local/lib/python3.4/dist-packages/wcsaxes/frame.py", line 192, in update_spines
    self['b'].data = np.array(([xmin, ymin], [xmax, ymin]))
    File "/usr/local/lib/python3.4/dist-packages/wcsaxes/frame.py", line 40, in data
    self._world = self.transform.transform(self._data)
    File "/usr/local/lib/python3.4/dist-packages/wcsaxes/transforms.py", line 166, in transform
    world = self.wcs.wcs_pix2world(pixel_full, 1)
    File "/usr/local/lib/python3.4/dist-packages/astropy/wcs/wcs.py", line 1329, in wcs_pix2world
'output', *args, **kwargs)
    File "/usr/local/lib/python3.4/dist-packages/astropy/wcs/wcs.py", line 1231, in _array_converter
    return _return_single_array(xy, origin)
    File "/usr/local/lib/python3.4/dist-packages/astropy/wcs/wcs.py", line 1212, in _return_single_array
"of shape (N, {0})".format(self.naxis))
    ValueError: When providing two arguments, the array must be of shape (N, 3)

【问题讨论】：

标签： python matplotlib astropy pyfits

【解决方案1】：

您的 WCS 对象包含两个空间轴和光谱轴，因为您使用完整标题对其进行了初始化（我假设您的标题还包含光谱 WCS 解决方案）。因此，仅使用空间创建 2D 图是行不通的，因为您对子图的投影是 3D。

我自己没有这样做，也没有示例文件，但 WCS 文档提到您可以使用 the sub function, or even the celestial functions 获取单个轴或天体（空间）轴；这些函数将返回一个仅包含这些轴的 WCS 对象。

因此，您可能可以使用以下内容：

my_wcs = WCS(my_header).celestial
fig = plt.figure()
ax = fig.add_subplot(111, projection=my_wcs)

【讨论】：

celestial 是属性，不是函数，所以不需要调用
那将需要在文档中进行修复； API 文档是正确的，但我链接到的部分将其称为函数。

【解决方案2】：

这是spectral-cube 的一个很好的用例，它有效地包装了astropy.io.fits 以供多维数据集使用。

假设您安装了wcsaxes，Evert 的解决方案一旦得到纠正，就会起作用。

要为此使用spectral_cube，一个简单的例子是：

cube = SpectralCube.read(filename)
cube[5,:,:].quicklook() # uses aplpy

# using wcsaxes
fig = plt.figure()
ax = fig.add_axes([0.1,0.1,0.8,0.8], projection=cube[5,:,:].wcs)
ax.imshow(cube[5,:,:]) # you may need cube[5,:,:].value depending on mpl version

【讨论】：