astropy.io.适合常见问题解答#

一般性问题#

什么是PyFITS?它与什么相关 astropy 是吗?#

PyFITS 是一个编写的库,并与 Python 用于读、写和操作的编程语言 FITS 格式化的文件。它包括一个高级接口,使标题符合对标题进行高级别和低级别操作的能力,并支持将图像和表数据读取为 NumPy 数组。它还支持在某些FITS文件中发现的更模糊和非标准的格式。

这个 astropy.io.fits 模块与PyFITS相同,但名称有所更改。当发展到 astropy 一开始,很明显核心需求之一就是FITS阅读器。PyFITS不是从头开始,而是被移植到Python中最灵活的FITS阅读器 astropy . 有计划逐步淘汰PyFITS作为一个独立的模块,并将其弃用,取而代之的是 astropy.io.fits . 请看下一个问题。

尽管PyFITS主要是用Python编写的,但它包含一个用C编写的可选模块,该模块是读/写压缩图像数据所必需的。但是,其余的PyFITS函数没有这个扩展模块。

PyFITS的发展状况如何?#

该软件由PyFITS负责维护和编写 Space Telescope Science Institute ,并由授权 AURA 在一个 3-clause BSD license .

It is now exclusively developed as a component of astropy (astropy.io.fits) rather than as a stand-alone module. There are a few reasons for this: The first is simply to reduce development effort; the overhead of maintaining both PyFITS and astropy.io.fits in separate code bases is nontrivial. The second is that there are many features of astropy (units, tables, etc.) from which the astropy.io.fits module can benefit greatly. Since PyFITS is already integrated into astropy, it makes more sense to continue development there rather than make astropy a dependency of PyFITS.

PyFITS过去的主要开发人员和活跃的维护者是Erik Bray。有一个 GitHub project 对于PyFITS,但是PyFITS不再积极开发,所以补丁和问题报告应该发布在Astropy问题跟踪程序上。

当前(也是最后一个)稳定版本是3.4.0。

使用问题#

有些东西没有像我预期的那样起作用。我做错什么了吗?#

可能吧。但是,如果您遵循了文档,而事情仍然没有按预期工作,那么完全可能是文档中有错误,代码中有错误,或者两者都有。所以请随意报告为bug。在FITS文件中也有很多很多角落的案例,几乎每周都有新的案例被发现。 astropy.io.fits 总是在改进,但并不完全支持所有案例。FITS格式(例如,缩放数据)的某些特性难以正确支持,有时可能会导致意外行为。

但是,对于最常见的情况,例如读取和更新FITS标题、图像和表, astropy.io.fits 非常稳定且经过良好测试。在每一个 astropy release它确保它的所有测试都能在各种平台上通过,这些测试覆盖了大多数用例(直到发现新的角落案例)。

astropy 崩溃并输出一长串代码。我该怎么办?#

这个代码列表就是所谓的 stack trace (或者用Python的说法是“回溯”)。当代码中发生未经处理的异常导致程序结束时,这是一种显示异常发生的位置和导致异常的代码的路径的方法。

AS astropy 是为了在其他软件项目中用作一个部分,由 astropy 都是故意的。例如,最常见的异常之一是 KeyError 当试图读取标头中不存在的关键字的值时:

>>> from astropy.io import fits
>>> h = fits.Header()
>>> h['NAXIS']
Traceback (most recent call last):
    ...
KeyError: "Keyword 'NAXIS' not found."

这表明有东西在寻找一个不存在的名为“NAXIS”的关键字。如果在使用 astropy ,它可能表示该软件中存在错误,因为它期望找到一个文件中不存在的关键字。

大多数“预期的”异常都会在回溯结束时输出一条消息,给出异常发生的原因和处理方法。异常中的错误消息越模糊和神秘,则越有可能是由中的错误引起的 astropy . 所以,如果你得到了一个异常,而你真的不知道为什么或者怎么做,你可以把它当作一个bug来报告。

为什么在CFITSIO、ds9等中打开文件有效,而在 astropy 是吗?#

正如本常见问题解答中提到的,在处理FITS文件时有许多不寻常的拐角情况。文件可能可以正常工作,但由于错误而不受支持。有时甚至可以让文件在旧版本中工作 astropy ,但由于尚未测试的回归,因此不是更新版本。

FITS格式的另一个问题是,尽管它是旧的,但在某些源文件中出现的许多约定不符合FITS标准。然而,它们是如此的普遍,因此有必要在任何适合的读者中支持它们。连续卡就是这样一个例子。支持的非标准约定 astropy CFITSIO不支持,可能反之亦然。你可能碰到了其中一个案子。

如果 astropy 打开文件时遇到问题,是排除问题是否与 astropy 将文件运行到 fitsverify 程序。对于较小的文件,您甚至可以使用 online FITS verifier . 这些文件在幕后使用了CFITSIO,应该能很好地指示文件是否有错误。如果文件格式不正确,fitsverify将输出错误和警告。

如果fitsverify确认文件没有问题,并且 astropy 仍然无法打开它(特别是如果它生成回溯),则可能中存在错误 astropy .

如何关闭警告消息 astropy 输出到我的控制台?#

astropy 使用Python的内置 warnings 子系统,用于通知代码中可恢复但用户可能希望被通知的异常情况。中最常见的警告之一 astropy.io.fits 更新标头值时发生,必须截断注释以保留空间:

Card is too long, comment is truncated.

由生成的任何控制台输出 astropy 可以假设来自警告子系统。参见Astropy的文档 Python 警告系统 有关如何控制和安静警告的更多信息。

什么惯例 astropy 用于索引,例如图像坐标?#

所有数组和序列 astropy 使用从零开始的索引方案。例如,头中的第一个关键字是 header[0] 不是 header[1] . 这与Python本身以及Python所基于的C语言是一致的。

这可能会让来自IRAF的资深FITS用户感到意外,IRAF通常使用基于1的索引,因为它起源于Fortran。

同样,nxn数组中左上角的像素是 data[0,0] . 二维数组的索引是行主顺序,第一个索引是行号,第二个索引是列号。或者用轴来表示,第一个轴是y轴,第二个轴是x轴。这与Fortran使用的列major order相反,因此FITS。例如,第二个索引引用FITS标题中NAXIS1指定的轴。

一般来说,对于N维数组,行主序意味着最右边的轴在数组数据上线性移动时变化最快。例如,三维数组:

[[[1, 2],
  [3, 4]],
 [[5, 6],
  [7, 8]]]

按行主要顺序线性表示为:

[1, 2, 3, 4, 5, 6, 7, 8]

因为2紧跟在1之后,您可以看到最右边(或最内层)的轴变化最快。

轴心顺序的差异可能需要一些时间来适应,但这是一种必然的罪恶。由于大多数其他Python和C软件都采用行主顺序,所以尝试在 astropy 可能会造成比实际价值更多的困难。

如何打开无法放入内存的非常大的图像?#

astropy.io.fits.open 有一个选项,可以使用内存映射访问HDU的数据部分 mmap . 在 astropy 默认情况下使用此选项。

这意味着访问上面例子中的数据只会根据需要将部分数据读入内存。例如,如果我们只请求图像的一部分,例如 hdul[0].data[100:200] ,则只有第100-200行将被读入内存。这是透明的,就像整个图像已经在内存中一样。这对于表也是一样的。在大多数情况下,这是处理大文件的最佳选择。

若要确保使用内存映射,请添加 memmap=True 参数为 fits.open 。同样,使用 memmap=False 将强制将数据完全读取到内存中。

也可以通过名为 USE_MEMMAP . 将此设置为 0 默认情况下将禁用mmap。

不幸的是,内存映射目前不能很好地处理缩放图像数据,需要对数据应用BSCALE和BZERO因子以生成物理值。目前,这需要足够的内存来容纳整个数组,尽管这是一个将来会得到改进的领域。

另一种方法是sections接口,它目前只适用于图像数据(即非表)。它在很大程度上被对mmap的更好的支持所取代,但是在虚拟内存空间更有限的系统上,例如在32位系统上,它可能仍然有用。对缩放图像数据的支持也不完整,尽管这将被修复。请参阅上的文档 image sections 有关使用此接口的详细信息。

如何从头开始创建非常大的FITS文件?#

此示例演示如何使用以下命令从头开始创建大文件(大于内存容量) astropy.io.fits

通常情况下,要创建单个图像FITS文件,需要执行如下操作:

import os
import numpy as np
from astropy.io import fits

data = np.zeros((40000, 40000), dtype=np.float64)
hdu = fits.PrimaryHDU(data=data)

然后使用 astropy.io.fits.writeto() 方法将新文件写出到磁盘:

hdu.writeto("large.fits")

然而,一个40000 x 40000的双精度数组几乎是12 GB!大多数系统不能仅仅为了写到磁盘而在内存中创建它。为了高效地创建如此大的文件,需要一些额外的工作和一些假设。

首先,预测标题的大小(例如,有多少关键字)是很有帮助的。FITS标头必须以2880字节块的形式写入,每个块的大小足以容纳36个关键字(包括最终块中的END关键字)。典型的报头有1到4个块,但有时会更多。

因为我们写入FITS文件的第一件事是头,所以我们想要写足够的头块,以便有足够的填充来添加新的关键字,而不必调整整个文件的大小。假设您希望标头默认使用4个数据块。然后,不包括Astropy将自动添加的末尾卡片,创建标题并将其填充到36*4张卡片中。

创建一个存根数组来初始化HDU;它的确切大小无关紧要,只要它具有所需的维数:

data = np.zeros((100, 100), dtype=np.float64)
hdu = fits.PrimaryHDU(data=data)
header = hdu.header
while len(header) < (36 * 4 - 1):
    header.append()  # Adds a blank card to the end

现在将NAXISn关键字调整为所需的数组大小,并仅将头文件写出到文件中。使用 hdu.writeto() 方法将导致AstPy“有帮助地”重置NAXISn关键字以匹配伪数组的大小。这是因为它努力确保只写入有效的FITS文件。相反,我们可以只将头写入使用 astropy.io.fits.Header.tofile 方法:

header["NAXIS1"] = 40000
header["NAXIS2"] = 40000
header.tofile("large.fits")

最后,将文件末尾加长以匹配数据的长度(加上标题的长度)。在大多数系统上,通过查找文件末尾并写入单个字节,可以非常高效地完成此操作,如下所示:

with open("large.fits", "rb+") as fobj:
    # Seek past the length of the header, plus the length of the
    # Data we want to write.
    # 8 is the number of bytes per value, i.e. abs(header['BITPIX'])/8
    # (this example is assuming a 64-bit float)
    # The -1 is to account for the final byte that we are about to
    # write:
    fobj.seek(len(header.tostring()) + (40000 * 40000 * 8) - 1)
    fobj.write(b"\0")

更广泛地说,可以这样写:

shape = tuple(header[f"NAXIS{ii}"] for ii in range(1, header["NAXIS"] + 1))
with open("large.fits", "rb+") as fobj:
    fobj.seek(
        len(header.tostring()) + (np.prod(shape) * np.abs(header["BITPIX"] // 8)) - 1
    )
    fobj.write(b"\0")

在现代操作系统上,这将导致文件(在标题之后)被零填充,直到保存40000 x 40000映像所需的~12 GB。在支持稀疏文件创建的文件系统(大多数Linux文件系统,但不是大多数Mac使用的HFS+文件系统)上,这是一个非常快速、高效的操作。在其他系统上,您的里程数可能会有所不同。

这不是构建大文件的唯一方法,但可能是最安全的方法之一。这种方法还可以用于创建大型多扩展FIT文件,只需稍加小心。

对于创建非常大的表,也可以使用这种方法,尽管提前确定表需要多少行是很困难的。一般来说,使用 astropy.io.fits module is currently discouraged for the creation and manipulation of large tables. The FITS format itself is not designed for efficient on-disk or in-memory manipulation of table structures. For large, heavy-duty table data it might be better too look into using HDF5 通过 PyTables 类库。这个 Astropy Table 接口还可以在不同的磁盘表格式之间提供一个抽象层(例如,用于在FITS和HDF5之间转换表)。

PyTables在幕后利用了NumPy,可以用FITS要求的相同格式将二进制表数据写入磁盘。然后就可以将表序列化为FITS格式进行分发。在某种程度上,这个常见问题解答可能会提供一个如何做到这一点的例子。

为什么包含整数数据的图像会意外地转换为浮点数?#

如果图像的标题包含可选的BSCALE和/或BZERO关键字(即BSCALE!=1和/或b0!=0),则必须根据以下公式将文件中的原始数据重新缩放为其物理值:

physical_value = BZERO + BSCALE * array_value

由于BZERO和BSCALE是浮点值,因此结果值也必须是浮点值。如果原始值是16位整数,则结果值是单精度(32位)浮点。如果原始值是32位整数,则结果值为双精度(64位浮点)。

这种自动缩放很容易让您措手不及,如果您没有预料到它,因为它不会发生,直到HDU的数据部分被访问(以允许像更新头部一样的事情,而不重新缩放数据)。例如::

>>> fits_scaledimage_filename = fits.util.get_testdata_filepath('scale.fits')

>>> hdul = fits.open(fits_scaledimage_filename)
>>> image = hdul[0]
>>> image.header['BITPIX']
16
>>> image.header['BSCALE']
0.045777764213996
>>> data = image.data  # Read the data into memory
>>> data.dtype.name    # Got float32 despite BITPIX = 16 (16-bit int)
'float32'
>>> image.header['BITPIX']  # The BITPIX will automatically update too
-32
>>> 'BSCALE' in image.header  # And the BSCALE keyword removed
False

这样做的原因是,一旦用户访问数据,他们也可能对其进行操作和计算。如果数据被强制保留为整数,则会损失大量精度。因此,最好是在不丢失数据的情况下出错,代价是在开始时造成一些混乱。

如果要保存数据,必须先返回整数 scale 方法:

>>> image.scale('int32')
>>> image.header['BITPIX']
32
>>> hdul.close()

或者,如果使用 mode='update' 随着 scale_back=True 参数,则在保存之前,将自动对数据重新应用原始的BSCALE和BZERO缩放。通常这是不可取的,尤其是当从浮点值转换回无符号整数值时。但在需要根据物理值的变化修改原始数据的情况下,这可能很有用。

为了防止重新缩放发生(这有利于更新标头-即使您不想让代码访问数据,在这里还是要小心出错),请使用 do_not_scale_image_data 打开文件时的参数::

>>> hdul = fits.open(fits_scaledimage_filename, do_not_scale_image_data=True)
>>> image = hdul[0]
>>> image.data.dtype.name
'int16'
>>> hdul.close()

为什么在头中分配浮点值时会丢失精度?#

FITS标准允许使用两种格式在标题值中存储浮点数。“固定”格式要求数字的ASCII表示形式为头卡的字节11到30,并且要右对齐。这将为任何注释字符串留下标准的字符数。

固定格式不够宽,无法表示可完全精确存储在64位浮点中的值的全部范围。所以FITS还支持一种“自由”格式,在这种格式中,ASCII表示可以存储在任何地方,使用卡的全部70个字节(在关键字之后)。

目前 astropy 只支持写入固定格式(它可以读取两种格式),因此所有分配给标头的浮点值都以固定格式存储。有计划增加对更灵活格式的支持。

同时,可以通过从字符串手动格式化卡片图像来添加或更新卡片,因为它应该出现在FITS文件中:

>>> c = fits.Card.fromstring('FOO     = 1234567890.123456789')
>>> h = fits.Header()
>>> h.append(c)
>>> h
FOO     = 1234567890.123456789

只要你不把新值赋给'FOO'通过 h['FOO'] = 123 ,将保持标题值与您格式化时完全相同(只要它根据FITS标准有效)。

为什么读取FITS表中的行这么慢?#

返回的每个表数据数组的基础 astropy.io.fits is a numpy recarray which is a numpy array type specifically for representing structured array data (i.e., a table). As with normal image arrays, astropy accesses the underlying binary data from the FITS file via mmap (see the question "What performance differences are there between astropy.io.fits and fitsio? “为了更深入地解释这一点)。然后底层mmap作为 recarray 一般来说,这是一种非常有效的读取数据的方法。

然而,对于许多(如果不是大多数)FITS表来说,这并不是那么简单。对于许多列,必须在“磁盘上”(在FITS文件中)的实际数据和返回给用户的数据值之间进行转换。例如,FITS二进制表表示布尔值的方式与 numpy 期望它们被表示,“逻辑”列需要动态转换为一种格式 numpy (因此用户)可以理解。此问题也适用于通过 TSCALnTZEROn 标题关键字。

支持所有这些“FITS-ism”会带来很多开销,这些开销可能不是所有表都需要的,但是仍然很常见。这并不是说即使支持FITS的特性,它也不能更快——例如,CFITSIO支持所有相同的特性,但速度要快几个数量级。 astropy 在这里也可以做得更好,还有很多已知的问题导致经济放缓。有很多加速的机会,补丁是受欢迎的。同时,对于带有FITS表的高性能应用程序,一些用户可能会发现 fitsio 他们更喜欢类库。

我在一个循环中打开了许多FITS文件,结果发现OSError:打开的文件太多了#

假设你有一些代码,比如:

from astropy.io import fits

for filename in filenames:
    with fits.open(filename) as hdul:
        for hdu in hdul:
            hdu_data = hdul.data
            # Do some stuff with the data

细节可能有所不同,但定性的观点是,许多hdu和/或FITS文件的数据是在一个循环中访问的。这可能会导致以下异常:

Traceback (most recent call last):
  File "<stdin>", line 2, in <module>
OSError: [Errno 24] Too many open files: 'my_data.fits'

中所解释的 note on working with large files ,因为 astropy 默认情况下使用mmap读取FITS文件中的数据,即使您使用 HDUList.close 保持对该文件的句柄打开,以便仍然可以继续透明地读取内存映射的数据数组。

numpy 支持mmap是这样的文件映射直到覆盖时才关闭 ndarray 对象没有对它的引用,因此释放了内存。然而,当在大量文件(甚至hdu)上快速循环时,这可能不会立即发生。或者在某些情况下,如果HDU对象持续存在,则附加到它的数据数组也可能会持续存在。建议的解决方法是 手动 删除 .data 属性,以便 ndarray 引用被释放,并且可以关闭mmap:

from astropy.io import fits

for filename in filenames:
    with fits.open(filename) as hdul:
        for hdu in hdul:
            hdu_data = hdul.data
            # Do some stuff with the data
            # ...
            # Don't need the data anymore; delete all references to it
            # so that it can be garbage collected
            del hdu_data
            del hdu.data

在某些极端情况下,文件的打开和关闭速度足够快,以至于Python的垃圾回收器不能经常释放它们(因此释放文件句柄)。为了缓解这种情况,您的代码可以通过调用 gc.collect() 在循环的末尾。

在将来的版本中,如果需要,在关闭FITS文件时自动执行这种清理将更加方便。

使用页眉 ['纳西2'] +=1不会将另一行添加到我的表中#

NAXIS 类似的关键词是适合的 结构的 关键字,用户不应修改。它们由自动更新 astropy.io.fits 当检查数据和标头的有效性时。看到了吗 结构关键词 更多信息。

要向表中添加行,可以修改实际数据。

与其他适合读者的比较#

有什么区别astropy.io.适合菲西奥呢?#

这个 astropy.io.fits module (originally PyFITS) is a "pure Python" FITS reader in that all of the code for parsing the FITS file format is in Python, though numpy is used to provide access to the FITS data via the ndarray interface. astropy.io.fits currently also accesses the CFITSIO 支持FITS Tile压缩约定,但此功能是可选的。除了读取压缩图像外,它不使用CFITSIO。

fitsio _另一方面,是CFITSIO库的Python包装器。读取FITS格式的所有繁重工作都由CFITSIO处理,而 fitsio 提供了使用面向对象API的更好方法,包括提供 numpy 接口连接到从CFITSIO读取的FITS文件。其中大部分是用C编写的(以提供Python和CFITSIO之间的接口),其余的是用Python编写的。Python端主要提供文档和用户级API。

因为 fitsio 它继承了CFITSIO的大部分优点和缺点,尽管它还有一个优点,即提供比直接使用CFITSIO更方便的API。

为什么Astropy采用PyFITS作为FITS阅读器而不是fitsio?#

当Astropy项目刚开始时,很明显它的核心组件之一应该是读写FITS文件的子模块,因为很多其他组件都可能依赖于这个功能。当时 fitsio 一揽子计划刚刚起步(大约可以追溯到2011年),而PyFITS已经建立(可以追溯到2000年之前)。它已经是一个成熟的软件包,支持在野外发现的绝大多数FITS文件,包括过时的格式,如“随机组”FITS文件仍然在射电天文学界广泛使用。

尽管PyFITS接口的许多方面都经过了多年的发展,但其中的大部分仍然保持不变,并且对于在Python中处理FITS文件的天文学家来说已经很熟悉了。大多数现有的培训材料都是基于PyFITS的。PyFITS是STScI开发的,它也为开发Astropy提供了大量的资源,着眼于将Astropy集成到STScI自己的软件栈中。由于STScI的大多数Python软件都使用PyFITS,因此它是实现这种转换的唯一实际选择。

最后,尽管CFITSIO fitsio )可以读取任何符合FITS标准的FITS文件,它不支持在野外添加到FITS文件中的所有非标准约定。虽然它确实支持其中的一些约定(例如CONTINUE cards,以及在一定程度上,HIERARCH cards),但是要在大型复杂的C代码库中添加对其他约定的支持并不容易。

PyFITS面向对象的设计使得支持非标准约定在大多数情况下变得更加容易,因此PyFITS可以在读取和返回的FITS文件类型上更加灵活 有用的 数据来自。这包括更好地支持不能满足FITS标准的文件,但仍然包含有用的数据,这些数据应该足够可读,以纠正任何违反FITS标准的行为。例如,在非英语区域中,一个常见错误是将非ASCII字符插入到FITS头中。这不是一个有效的FITS文件,但在某种意义上仍然是可读的。在CFITSIO中,支持这样的结构误差更为困难,因为CFITSIO假定结构更为刚性。

两者之间有什么性能差异astropy.io.适合菲西奥呢?#

有两个主要的性能领域需要关注:读取/解析FITS头和读取FITS数据(像图像一样的数组以及表)。

在头部区域, fitsio 在大多数情况下都要快得多。这在很大程度上是由于(几乎)纯C实现(由于CFITSIO的使用),但也因为它更严格,并且不支持许多本地约定和其他特殊情况 astropy.io.fits 尝试在其纯Python实现中支持。

也就是说,这种差异很小,只有在打开包含数千个hdu的文件时,或者在连续读取数千个FITS文件的头文件时才会成为瓶颈(无论哪种情况下,差异都不是一个数量级)。

在数据方面,情况稍微复杂一些,需要了解如何 astropy.io.fits 与CFITSIO和 fitsio . 首先,了解它们在内存管理方面的差异是很重要的。

astropy.io.fits 默认情况下,使用mmap提供对FITS文件中原始二进制数据的访问。Mmap是一个系统调用(或者现在在大多数情况下是libc中用于较低级别系统调用的包装器),它允许用户空间应用程序执行与操作系统使用页面文件(交换空间)作为虚拟内存时所做的相同的事情:它允许将磁盘上文件中的数据分页到物理内存中一页(或者在实践中通常是这样)几页)一次根据需要。文件的这些缓存页也可以从系统上的所有进程访问,因此多个进程可以从同一个文件中读取,而不会增加额外的开销。在读取文件中的所有数据的情况下,使用mmap与一次将整个数据读入物理内存之间的性能差异可能会因系统、硬件以及当前系统上发生的其他情况而有很大差异,但mmap几乎总是会更好。

原则上,它需要更多的开销,因为访问每个页面将导致页面错误,并且系统需要更多的磁盘请求。但实际上,操作系统会非常积极地优化这一点,尤其是对于最常见的顺序访问的情况,同样在现实中,将整个内容读入内存仍然会导致大量的页面错误。对于随机访问,将所有数据存储在物理内存中总是最好的,不过使用mmap通常也会很好。(大多数用户通常不会以完全随机的顺序访问一个文件中的所有数据—通常会最频繁地访问其中的一些部分,因此操作系统将尽可能地将这些页保存在物理内存中。)对于读取FITS文件(或磁盘上大多数大数据)的最常见情况,这是最佳选择,尤其是临时用户,因此默认启用。

另一方面,CFITSIO/`fitsio``并不假定存在诸如mmap和页面缓存之类的技术。因此,它实现了自己的I/O缓冲区的LRU缓存,以FITS著名的2880字节块大小将从磁盘读取的FITS文件部分存储在内存中。I/O缓冲区被大量使用,特别是在内存中保存头。但是对于大数据读取(例如,从文件读取整个图像),它 does 绕过缓存,而是直接从磁盘读取到用户提供的内存缓冲区。

然而,即使CFITSIO直接从文件中读取数据,这在很大程度上仍然不如使用mmap。通常,当您的操作系统从磁盘读取文件时,它会尽可能多地将读取的内容缓存在物理内存中(在其页缓存中),这样以后访问这些相同的页面就不需要再进行昂贵的磁盘读取了。使用mmap时也会发生这种情况,因为数据必须在某个时刻从磁盘复制到RAM中。不同之处在于,当使用mmap访问数据时,程序能够读取该数据 直接地 从操作系统的页面缓存中取出(只要它只被读取)。另一方面,当从文件中读取数据到本地缓冲区(如with fread())时,首先将数据读入页缓存(如果尚未存在),然后从页缓存复制到本地缓冲区。因此,每次读取每个页面至少执行一个额外的内存拷贝(需要两倍的物理内存,如果文件很大,需要从缓存中删除页面,则可能需要大量分页)。

CFITSIO的用户API通常是通过让用户分配足够大的内存缓冲区来保存他们想要读取的图像/表(或者至少是他们感兴趣的部分)。有一些helper函数用于确定要分配的适当空间量。然后传入一个指向缓冲区的指针,CFITSIO处理所有读取操作(通常使用上述过程),并将结果复制到用户缓冲区中。对于大的读取,它直接从文件读取到缓冲区,但是如果需要缩放数据,它会先在CFITSIO自己的缓冲区中停止,然后将重新缩放的值写入用户缓冲区(如果已请求重新缩放)。不管怎样,这意味着如果你的程序希望一次在内存中保存一个完整的图像,它将使用与数据大小相同的RAM。对于大多数应用程序来说,处理较小的数据部分更好(也足够),但这需要额外的复杂性。另一方面,使用mmap可以更有效地管理这种复杂性。

一个非正式的测试证明了这一点。该测试在四个尺寸为256x256、1024x1024、4096x4096和256x1024x1024的简单拟合图像上进行。每个图像都是在测试前生成的,并填充随机的64位浮点值。使用这两种方法进行了类似的试验 astropy.io.fitsfitsio . 使用每个库的基本语义打开FITS文件的句柄,然后将文件的整个数据数组复制到内存中的临时数组中(例如,如果我们要将图像blitting到视频缓冲区)。为 astropy 测试内容如下:

def read_test_astropy(filename):
    with fits.open(filename, memmap=True) as hdul:
        data = hdul[0].data
        c = data.copy()

测试是在IPython中进行的,该系统的内核版本为2.6.32,采用6核Intel Xeon X5650 CPU,每核2.67 GHz,内存11.6 GB,使用:

for filename in filenames:
    print(filename)
    %timeit read_test_astropy(filename)

在哪里? filenames 只是上述生成的示例文件的列表。结果是:

256x256.fits
1000 loops, best of 3: 1.28 ms per loop
1024x1024.fits
100 loops, best of 3: 4.24 ms per loop
4096x4096.fits
10 loops, best of 3: 60.6 ms per loop
256x1024x1024.fits
1 loops, best of 3: 1.15 s per loop

为了 fitsio 测试是:

def read_test_fitsio(filename):
    with fitsio.FITS(filename) as f:
        data = f[0].read()
        c = data.copy()

这也在所有示例文件上循环运行,产生的结果是:

256x256.fits
1000 loops, best of 3: 476 µs per loop
1024x1024.fits
100 loops, best of 3: 12.2 ms per loop
4096x4096.fits
10 loops, best of 3: 136 ms per loop
256x1024x1024.fits
1 loops, best of 3: 3.65 s per loop

应该明确的是,样本文件是用新的随机数据重写的 astropy 测试和fitsio测试,所以他们没有从操作系统的页面缓存中读取相同的数据。Fitsio在小的(256x256)图像上要快得多,因为在这种情况下,时间是由解析头决定的。如前所述,这在CFITSIO中要快得多。但是,随着数据量的增加,头解析不再占主导地位, astropy.io.fits 使用mmap的速度大约是前者的两倍。这种差异几乎完全是由于它需要大约一半的内存副本来读取数据,如前所述。也就是说,更广泛的基准测试可能会非常有趣。

这也不是说 astropy.io.fits 在所有情况下都做得更好。在某些情况下,它目前被fitsio吹走了。见后面的问题。

为什么fitsio比 astropy 在阅读桌上?#

在许多情况下,情况并非如此:要么没有区别,要么在 astropy 取决于您要对表执行的操作以及表的列类型或列数。但是,有些情况下 fitsio 可以大大加快速度,主要是因为上面的“为什么从FITS表中读取行这么慢?”?`_"

原则上,一个表和一个像素数组没有区别。但是数组的每个元素不是像素而是某种记录结构(例如,两个浮点、一个布尔值和一个20个字符的字符串字段)。正如64位浮点是数组中的8字节记录一样,这样一个表中的一行可以被认为是一维行数组中的37字节(在前面的例子中)记录。所以,原则上,在回答问题时所解释的一切“在性能上有什么区别astropy.io.适合菲西奥呢?`_“与其他数组一样,它同样适用于表。

然而,FITS表有许多额外的复杂性,有时会妨碍数据直接从磁盘流式传输,而是需要从磁盘上的FITS格式转换为对用户更直接有用的格式。一个常见的例子是FITS如何表示二进制表中的布尔值。另一个更复杂的例子是可变长度数组。

正如“为什么从FITS表中读取行这么慢?”?_", `astropy.io.fits 当前不能尽可能有效地处理其中一些情况,特别是在用户只希望从表中读取几行的情况下。另一方面,Fitsio有一个更好的接口,可以一次从表中复制一行并对该行执行必要的转换 only ,而不是在从中提取行的一个或多个整列上。因此,在许多情况下 fitsio 获得更好的性能,对于许多性能关键的表操作,应该是首选。

Fitsio还公开了一种微语言(在CFITSIO中实现),用于对表进行高效的类似SQL的查询(只有单个表-没有连接或类似的东西)。此格式,在 CFITSIO documentation 在某些情况下,可以执行比使用 numpy 它需要创建一个中间掩码数组来执行行选择。