使用astropy从旧数据编写一个新的FITS文件

Question

使用astropy从旧数据编写一个新的FITS文件

3

我有一个非常简单的操作需要在一个FITS文件上执行（数据为numpy数组格式），但我无法让astropy将其保存为新文件或覆盖现有文件。

我正在重写一些旧代码，该代码使用numpy pyfits模块处理天文FITS文件 - 我想更新它以使用astropy.io fits模块。具体来说，我处理的某些数据是3D的，而另一些则是4D的。 4D的东西只是一种约定 - 第四个轴不包含任何有用的信息（可以在此处找到数据示例：http://www.mpia.de/THINGS/Data_files/NGC_628_NA_CUBE_THINGS.FITS）。所以我更喜欢删除额外的轴，然后我的余下的代码就不需要任何特殊的要求了。

这是我用的旧pyfits-based代码，它完全正常：

import numpy
import pyfits

filename = 'NGC628.fits'
outfile  = 'NGC628_reshaped.fits'

# Get the shape of the file
fitsfile=pyfits.open(filename)
image = fitsfile[0].data
header =fitsfile[0].header
z = image.shape[1]  # No. channels                  
y = image.shape[2]  # No. x pixels
x = image.shape[3]  # No. y pixels

newimage = numpy.reshape(image,[z,y,x])

pyfits.core.writeto(outfile,newimage,header, clobber=True)

这里没有什么复杂的东西，它只是重新整理一个数组并将其保存到新文件中。非常好。现在我想用astropy fits模块来替换它。如果我执行以下操作：

import numpy
from astropy.io import fits as pyfits

fitsfile = pyfits.open('NGC628.fits', mode='update')
image  = fitsfile[0].data
header = fitsfile[0].header
z = image.shape[1]                  
y = image.shape[2]
x = image.shape[3]
image = numpy.reshape(image,[z,y,x])

如果目前为止一切顺利，那么"image"数组已经正确地被重塑，如image.shape所证实。但我无法想象如何将其保存到新的（或旧的）FITS文件中。使用旧语法：

fitsfile.writeto('NGC628_2.fits',image,header)

...会出现错误信息："AttributeError: 'numpy.ndarray' object has no attribute 'lower'。如果按照astropy文档，只是省略图像和头部并尝试保存到原始文件中：

fitsfile.writeto('NGC628.fits')

然后我收到一个错误，说文件已经存在。如果我提供“overwrite=True”关键字，则会抱怨“WinError 32：另一个进程正在使用该文件，因此无法访问文件：NGCC628.fits”。该文件绝对没有在任何其他程序中打开。

如果我指定新的文件名NGC628_2.fits，那么Python就会崩溃（将我发送回命令提示符），并且不会出现任何错误。只写入了包含标头数据但不包含图像数据的非常小的文件。编辑：如果我使用正确的命令使用图像和标头数据编写新的FITS文件，例如pyfits.writeto（'NGC628_2.fits'，image，header），则完全相同的事情会发生。

为了使事情更加混乱，如果我执行稍微简单的操作，例如将所有图像数据设置为恒定值，然后关闭文件：

import numpy
from astropy.io import fits as pyfits
fitsfile = pyfits.open('NGC628.fits', mode='update')
image  = fitsfile[0].data
header = fitsfile[0].header
image[:,:,:,:] = 5.0
fitsfile.close()

然后这个方法是可行的 - 原始文件现在是一个数组，每个值都等于5。从astropy文档中我了解到，简单地以更新模式打开文件并关闭它应该就足够了，在这种情况下确实如此。但是当重塑图像时，这个技巧却不起作用 - FITS文件保持不变。

那么我究竟做错了什么？无论是更新原始文件还是保存到新文件都可以（最好是后者），但我都无法正确地执行任何一项操作。

编辑：我使用的是Python版本3.5.3，numpy版本1.17.3，astropy版本3.2.3，并且运行在Windows 10上。

- Rhysy

1个回答

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Iguananaut · Accepted Answer

好的，我认为您一开始犯了一个小错误，而您一直简单地忽略了它，这导致您陷入了寻找不同解决方案的困境。您的其他尝试也很可能会成功，只是因为一个小细节而失败。我将先讲述一开始出了什么问题，然后再解释一些后续情况出了什么问题，所以有点长...

首先，在您最初的代码中，您写了：

pyfits.core.writeto(outfile,newimage,header, clobber=True)

这可以更简单地写成：

pyfits.writeto(outfile, newimage, header, clobber=True)

pyfits.core模块是一个实现细节，几乎没有直接引用它的理由。您可以直接调用此函数，如pyfits.writeto()。

如果您这样做，那么只需更改导入即可使您现有的代码完全像以前一样工作。

from astropy.io import fits as pyfits

唯一的注意事项是参数被重命名为，但我认为仍然可以使用，并带有弃用警告。

第一个错误

你改变了它的名称，似乎忽略了它应该改成：

fitsfile.writeto('NGC628_2.fits',image,header)

这并不是同一件事。在第一个案例中，它是高级writeto()“便捷函数”。但现在你正在调用fitsfile.writeto。这里的fitsfile不是astropy.io.fits模块，而是文件对象本身。

>>> type(fitsfile)
<class 'astropy.io.fits.hdu.hdulist.HDUList'>

HDUList类还有一个方法HDUList.writeto，它执行类似的功能，但接受不同的参数。您可以通过多种方式进行检查，例如：

>>> help(fitsfile.writeto)
Help on method writeto in module astropy.io.fits.hdu.hdulist:

writeto(fileobj, output_verify='exception', overwrite=False, checksum=False) method of astropy.io.fits.hdu.hdulist.HDUList instance
    Write the `HDUList` to a new file.

它唯一需要的参数是要写入文件的文件名。与其他的writeto不同，它不接受数据或标头参数，因为HDUList已经是一个包含相关数据和标头的HDUs集合。实际上，如果您只想对现有的FITS文件进行一些更改并将这些更改写入新文件，则不应像最初一样使用高级的writeto() - 它主要用于从头开始创建数组并快速将其写入新的FITS文件。所以您最初的代码也可以这样写：

import numpy as np
import astropy.io.fits as pyfits

filename = 'NGC628.fits'
outfile  = 'NGC628_reshaped.fits'

# Get the shape of the file
fitsfile = pyfits.open(filename)
image = fitsfile[0].data
z = image.shape[1]  # No. channels                  
y = image.shape[2]  # No. x pixels
x = image.shape[3]  # No. y pixels

# Replace the primary HDU's data in-place; this just manipulates the
# in-memory HDU data structure; it does not change anything on disk
fitsfile[0].data = np.reshape(image, [z,y,x])

# Write the new HDU structure to outfile
fitsfile.writeto(outfile, overwrite=True)

更新模式

接下来，您尝试使用 mode='update' 打开文件进行修改，但这并不是 update 模式的正确使用方式。 mode='update' 更类似于以写入模式打开普通文本文件，例如 open('foo.txt', 'w')：它旨在原地对磁盘上的现有文件进行修改。当文件关闭时，您所做的任何修改都会被刷新到磁盘上，无需使用 writeto() 或类似方法。

您写道：

使用旧语法：

fitsfile.writeto('NGC628_2.fits',image,header)

但正如我之前所写的那样，在这里没有“旧语法”，您只是试图使用 HDUList.writeto() 方法而不是 writeto() 函数。

如果我提供关键字 "overwrite=True"，则会出现错误 "WinError 32: 由于正在被另一个进程使用，因此无法访问文件：NGCC628.fits"。该文件肯定没有在任何其他程序中打开。

我看到您正在使用Windows操作系统，这是Windows的一个特定限制：如果文件已经有任何打开的句柄，Windows不允许对其进行写入、移动或删除操作。我认为这里的错误信息有误导性（此消息直接来自Windows本身）：“它正在被另一个进程使用”。它应该是类似于“它正在被此进程或其他进程使用”。当您执行以下操作时：

fitsfile = pyfits.open('NGC628.fits', mode='update')

你现在已经打开了该文件的句柄，因此在关闭它之前，Windows不允许你覆盖该文件（例如fitsfile.close()）。

如果我指定新的文件名NGC628_2.fits，那么Python会崩溃（将我送回命令提示符），没有任何错误。只写入一个非常小的文件，其中包含头数据但没有图像数据。

现在听起来像是一个真正的错误。它听起来像是你的Python解释器发生了段错误。我无法解释这个问题。但是，我尝试按照你描述的相同步骤（在Windows上）重现它，但我没能成功。最终的fitsfile.writeto('NGC628_2.fits')没有出现问题。我可以想象的一件事是，当你用mode='update'打开一个文件时，内部状态管理可能变得非常复杂，因为它必须跟踪你对FITS数据结构所做的所有更改，以便能够正确地在现有文件中移动东西。在尝试许多轻微病态操作的过程中（比如在更新文件时尝试覆盖文件），可能会导致某些未定义状态的出现。不过，我无法找到相同的方法。

“正确”使用mode='update'来原地修改文件的用法可能如下：

with pyfits.open(filename, mode='update') as fitsfile:
    image = fitsfile[0].data
    z, y, x = image.shape[1:]
    # Replace the original HDU data in-place
    fitsfile[0].data = image.reshape((z, y, x))

就是这样！我会习惯使用with语句：当你在with语句的主体内执行需要打开文件的任何操作时，当退出with语句时，它将确保将你所做的所有更改（如果有）刷新到磁盘并关闭文件。即使只是以只读模式打开文件，我也会使用它。尤其是在Windows上，因为像你发现的那样，Windows对于你是否正确关闭了文件尤为挑剔。

类似地，你的最后一个示例可以重写为：

with pyfits.open('NGC628.fits', mode='update') as fitsfile:
    image  = fitsfile[0].data
    header = fitsfile[0].header
    image[:,:,:,:] = 5.0

最后，你写道：

我从astropy文档中了解到，仅在更新模式下打开文件并关闭即可，而在这种情况下确实如此。但是，当重新调整图像大小时，这个技巧不起作用 - FITS文件不会改变。

但我不确定你尝试了什么；您没有展示那个尝试。如果我要猜测，您可能编写了类似以下内容的内容：

    image = np.reshape(image, (z, y, x))

但是这只是用一个新的数组替换了image变量所指向的值。它没有更新原始的fitsfile [0] .data。之前的image [:] = 0.5示例之所以有效，是因为在这里image指向与fitsfile [0] .data相同的数组，并且您正在就地修改其内容。但是np.reshape（...）不是就地操作；它返回一个新数组。请参见我的先前示例，了解正确的方法。