适合文件处理 (astropy.io.fits
)¶
介绍¶
这个 astropy.io.fits
包提供对FITS文件的访问。FITS(flexibleimagetransportsystem)是一种便携式文件标准,在天文学界被广泛用于存储图像和表格。
入门¶
本节简要介绍使用 astropy.io.fits
. 我们的目标是演示软件包的基本功能,而不需要太多细节。如果您是第一次使用或从未使用过 astropy
或者PyFITS,这是你应该开始的地方。另请参见 FAQ 对于常见问题的答案。
备注
如果要以FITS格式读取或写入单个表,建议的方法是通过高级 统一文件读写接口 . 具体请参见 Unified I/O FITS 部分。
读取和更新现有的FITS文件¶
打开FITS文件¶
备注
这个 astropy.io.fits.util.get_testdata_filepath()
函数用于访问随附的数据 astropy
. 若要使用您自己的数据,请使用 astropy.io.fits.open()
,它采用相对路径或绝对路径。
一旦 astropy.io.fits
使用标准约定加载包 1, 我们可以打开现有文件:
>>> from astropy.io import fits
>>> fits_image_filename = fits.util.get_testdata_filepath('test0.fits')
>>> hdul = fits.open(fits_image_filename)
这个 open()
函数有几个可选参数,这些参数将在后面的章节中讨论。在上面的例子中,默认模式是“readonly”。open函数返回一个名为 HDUList
哪个是 list
-就像HDU对象的集合。HDU(头数据单元)是FITS文件结构的最高级别组件,由头和(通常)数据数组或表组成。
在上述公开电话之后, hdul[0]
是主HDU, hdul[1]
是第一个扩展HDU等(如果有任何扩展),依此类推。应该指出的是 astropy
在引用hdu和头卡时使用基于零的索引,尽管FITS标准(设计时考虑到Fortran)使用基于一的索引。
这个 HDUList
有一个有用的方法 HDUList.info()
,它总结了打开的FITS文件的内容:
>>> hdul.info()
Filename: ...test0.fits
No. Name Ver Type Cards Dimensions Format
0 PRIMARY 1 PrimaryHDU 138 ()
1 SCI 1 ImageHDU 61 (40, 40) int16
2 SCI 2 ImageHDU 61 (40, 40) int16
3 SCI 3 ImageHDU 61 (40, 40) int16
4 SCI 4 ImageHDU 61 (40, 40) int16
完成打开的文件后,使用 HDUList.close()
方法:
>>> hdul.close()
使用可以避免手动关闭文件 open()
作为上下文管理器:
>>> with fits.open(fits_image_filename) as hdul:
... hdul.info()
Filename: ...test0.fits
No. Name Ver Type Cards Dimensions Format
0 PRIMARY 1 PrimaryHDU 138 ()
1 SCI 1 ImageHDU 61 (40, 40) int16
2 SCI 2 ImageHDU 61 (40, 40) int16
3 SCI 3 ImageHDU 61 (40, 40) int16
4 SCI 4 ImageHDU 61 (40, 40) int16
退出后 with
文件范围将自动关闭。这是(通常)在Python中打开文件的首选方法,因为即使发生异常,它也会关闭文件。
如果文件是用 lazy_load_hdus=False
,关闭HDUList后,仍可以访问所有标头。头和数据可能被访问,也可能不被访问,这取决于数据是否被触摸以及它们是否被内存映射;有关详细信息,请参阅后面的章节。
使用大文件¶
这个 open()
函数支持 memmap=True
参数,该参数允许用mmap访问每个HDU的数组数据,而不是一次全部读入内存。这对于处理无法完全放入物理内存的非常大的阵列特别有用。在这里 memmap=True
默认情况下,此值从配置项获取 astropy.io.fits.Conf.use_memmap
.
这对较小的文件也有最小的影响,尽管一些操作(如按顺序读取数组数据)可能会产生一些额外的开销。在32位系统上,大于2到3 GB的阵列不能被mmap'd(这很好,因为到那时你很可能会用尽物理内存),但64位系统在这方面的限制要小得多。
警告
打开文件时使用 memmap=True
,因为mmap是如何工作的,这意味着当访问HDU数据时(即。, hdul[0].data
)FITS文件的另一个句柄由mmap打开。这意味着即使打过电话 hdul.close()
mmap仍然持有一个打开的数据句柄,这样就可以由一些粗心的程序访问它,这些程序是在假定.data属性的内存中包含所有数据的情况下构建的。
要强制关闭mmap,请等待包含 HDUList
对象以超出范围,或手动调用 del hdul[0].data
. (只要没有其他对数据数组的引用,这种方法就有效。)
无符号整数¶
由于FITS格式源自Fortran,FITS本身不支持图像或表格中的无符号整数数据。但是,有一个共同的约定是将无符号整数存储为有符号整数,以及 转移 说明(a BZERO
带值关键字 2 ** (BITPIX - 1)
)将所有有符号整数上移到无符号整数。例如,在写入值时 0
作为无符号32位整数,它存储在FITS文件中 -32768
,以及header关键字 BZERO = 32768
.
astropy
默认情况下识别并应用此约定,以便所有看起来应该被解释为无符号整数的数据都会自动转换(这适用于图像和表)。在 astropy
v1.1.0之前的版本 not 自动应用,需要传递参数 uint=True
到 open()
. 在v1.1.0或更高版本中,这是默认设置。
即使 uint=False
, the BZERO
shift仍然被应用,但是返回的数组是“float64”类型。要完全禁用缩放/移动,请使用 do_not_scale_image_data=True
(见 为什么包含整数数据的图像会意外地转换为浮点数? 更多详细信息,请参阅常见问题解答)。
使用压缩文件¶
备注
在文件中讨论了在文件中使用的压缩文件 Compressed Image Data .
这个 open()
函数将无缝地打开用gzip、bzip2或pkzip压缩的FITS文件。注意,在本文中,我们讨论的是一个FITS文件,它是用这些实用程序(例如,a。适合gz文件)。
使用压缩文件时有一些限制。例如,对于包含多个压缩文件的Zip文件,只能访问第一个文件。另外,bzip2不支持append或update访问模式。
编写文件时(例如,使用 writeto()
函数),将根据给定的文件扩展名或正在写入的现有文件中使用的压缩来确定压缩。
使用非标准文件¶
什么时候? astropy.io.fits
读取不符合FITS标准的FITS文件,它将尝试对不符合标准的字段进行有根据的解释。这可能并不总是成功的,并且可能在访问头时触发警告,或者在写入文件时触发异常。可以使用控制写入输出文件的字段的验证 output_verify
参数 open()
. 打开供读取的文件可以用方法进行验证和修复 HDUList.verify
. 此方法通常在打开文件后但在访问任何头或数据之前调用:
>>> with fits.open(fits_image_filename) as hdul:
... hdul.verify('fix')
... data = hdul[1].data
在上面的示例中,调用 hdul.verify("fix")
要求 astropy.io.fits
修复不兼容的字段并打印信息性消息。除了 "fix"
在FITS下描述 验证
参见
FITS 验证 .
使用配合标题¶
如前所述 HDUList
是一个带有 .header
和 .data
属性,可用于访问HDU的标头和数据部分。
对于那些不熟悉FITS头的人,它们由一个80字节的“卡片”组成,其中卡片包含关键字、值和注释。关键字和注释都必须是字符串,而值可以是字符串或整数、浮点数、复数或 True
/False
. 关键字通常在标题中是唯一的,除了少数特殊情况。
header属性是一个header实例,另一个 astropy
对象。要获得与header关键字相关联的值,请执行以下操作(a la Python dicts)::
>>> hdul = fits.open(fits_image_filename)
>>> hdul[0].header['DATE']
'01/04/99'
获取关键字“DATE”的值,该关键字是字符串“01/04/99”。
尽管在FITS文件中关键字名称总是大写,但使用 astropy
为方便用户,不区分大小写。如果指定的关键字名称不存在,它将引发 KeyError
例外。
我们还可以通过索引(a la Python list)来获得关键字值:
>>> hdul[0].header[7]
32768.0
这个例子返回第八个关键字(就像Python列表一样,它是0索引的)的值-a float-32768.0。
类似地,可以在中更新关键字的值 astropy
,通过关键字名称或索引::
>>> hdr = hdul[0].header
>>> hdr['targname'] = 'NGC121-a'
>>> hdr[27] = 99
但是请注意,几乎所有的应用程序代码都应该通过关键字名称而不是位置索引来更新头值。这是因为大多数FITS关键字可能出现在标题的任何位置。
也可以通过将关键字指定为元组来更新与关键字相关联的值和注释:
>>> hdr = hdul[0].header
>>> hdr['targname'] = ('NGC121-a', 'the observation target')
>>> hdr['targname']
'NGC121-a'
>>> hdr.comments['targname']
'the observation target'
像dict一样,您也可以使用上面的语法添加新的关键字/值对(也可以选择添加注释)。在这种情况下,新的卡片被追加到标题的末尾(除非它是注释或历史之类的注释关键字,在这种情况下,它被附加到带有该关键字的最后一张卡片之后)。
更新现有卡或附加新卡的另一种方法是使用 Header.set()
方法:
>>> hdr.set('observer', 'Edwin Hubble')
注释或历史记录的添加方式与普通卡片相同,但在这种情况下,始终会创建新卡片,而不是更新现有的历史记录或注释卡片:
>>> hdr['history'] = 'I updated this file 2/26/09'
>>> hdr['comment'] = 'Edwin Hubble really knew his stuff'
>>> hdr['comment'] = 'I like using HST observations'
>>> hdr['history']
I updated this file 2/26/09
>>> hdr['comment']
Edwin Hubble really knew his stuff
I like using HST observations
注意:注意不要混淆评论卡和普通卡片的评论值。
要更新现有的评论或历史记录卡,请按索引引用它们:
>>> hdr['history'][0] = 'I updated this file on 2/27/09'
>>> hdr['history']
I updated this file on 2/27/09
>>> hdr['comment'][1] = 'I like using JWST observations'
>>> hdr['comment']
Edwin Hubble really knew his stuff
I like using JWST observations
若要查看FITS文件中显示的整个标题(结束卡和填充被去除),请单独输入header对象,或 print(repr(hdr))
::
>>> hdr
SIMPLE = T / file does conform to FITS standard
BITPIX = 16 / number of bits per data pixel
NAXIS = 0 / number of data axes
...
>>> print(repr(hdr))
SIMPLE = T / file does conform to FITS standard
BITPIX = 16 / number of bits per data pixel
NAXIS = 0 / number of data axes
...
只进入 print(hdr)
也可以工作,但在大多数显示器上可能不是很清晰,因为它显示的标题是写入FITS文件本身,这意味着卡片之间没有换行符。对于新用户来说,这是一个常见的混淆源。
也可以查看标题的一部分:
>>> hdr[:2]
SIMPLE = T / file does conform to FITS standard
BITPIX = 16 / number of bits per data pixel
上面只显示了前两张牌。
要获取所有关键字的列表,请使用 Header.keys()
方法就像你用dict::
>>> list(hdr.keys())
['SIMPLE', 'BITPIX', 'NAXIS', ...]
实例:
也见 编辑配合标题 .
结构关键词¶
FITS关键字混合元数据和解析文件所需的文件结构的关键信息。这些 结构的 关键字由内部管理 astropy.io.fits
而且,一般来说,用户不应触碰。相反,应该使用 astropy.io.fits
类(见下面的示例)。
FITS标准使用的特定结构关键字集因HDU类型而异。下表列出了与每种HDU类型关联的关键字:
HDU类型 |
结构关键词 |
---|---|
所有 |
|
|
|
|
|
|
|
|
还有许多其他保留关键字,例如数据缩放或表的列属性,如中所述 FITS Standard . 其中大多数可以通过 Column
或者HDU对象 hdu.name
设置 EXTNAME
或 hdu.ver
对于 EXTVER
. 作为常见操作的结果,检查和/或更新结构关键字。例如,当:
设置数据。这个
NAXIS*
关键字是从数据形状设置的 (.data.shape
)BITPIX
从数据类型 (.data.dtype
)设置收割台。它的关键字根据数据属性更新(如上所述)。
正在写入文件。所有必要的关键字都会被删除、更新或添加到标题中。
调用HDU的验证方法(例如。,
PrimaryHDU.verify()
). 有些关键字可以自动修复。
在这些情况下,用户可能分配给这些关键字的任何手写值都将被覆盖。
使用图像数据¶
如果HDU的数据是图像,HDU对象的data属性将返回一个 numpy
ndarray
对象。参见 numpy
有关操作这些数值数组的详细信息,请参阅文档:
>>> data = hdul[1].data
在这里, data
指向第二个HDU中的数据对象(第一个HDU, hdul[0]
,它对应于“SCI”扩展。或者,可以通过扩展名(在EXTNAME关键字中指定)访问扩展名:
>>> data = hdul['SCI'].data
如果有多个扩展名具有相同的EXTNAME,则需要将EXTVER值与EXTNAME一起指定为tuple;例如:
>>> data = hdul['sci',2].data
注意EXTNAME也不区分大小写。
归还的人 numpy
对象有许多属性和方法,供用户获取有关数组的信息,例如:
>>> data.shape
(40, 40)
>>> data.dtype.name
'int16'
因为图像数据是 numpy
对象,我们可以对其进行切片、查看和执行数学运算。要查看x=5,y=2处的像素值:
>>> print(data[1, 4])
348
请注意,与C(与Fortran不同)一样,Python是0索引的,索引的第一个轴是最慢的,变化最快的轴是最后一个;也就是说,对于二维图像,对应于FITS NAXIS1关键字的快速轴(X轴)是第二个索引。类似地,在Python中,x=11到20(含)和y=31到40(含)的1索引子段将给出如下:
>>> data[30:40, 10:20]
array([[350, 349, 349, 348, 349, 348, 349, 347, 350, 348],
[348, 348, 348, 349, 348, 349, 347, 348, 348, 349],
[348, 348, 347, 349, 348, 348, 349, 349, 349, 349],
[349, 348, 349, 349, 350, 349, 349, 347, 348, 348],
[348, 348, 348, 348, 349, 348, 350, 349, 348, 349],
[348, 347, 349, 349, 350, 348, 349, 348, 349, 347],
[347, 348, 347, 348, 349, 349, 350, 349, 348, 348],
[349, 349, 350, 348, 350, 347, 349, 349, 349, 348],
[349, 348, 348, 348, 348, 348, 349, 347, 349, 348],
[349, 349, 349, 348, 350, 349, 349, 350, 348, 350]], dtype=int16)
要更新像素或子区域的值,请执行以下操作:
>>> data[30:40, 10:20] = data[1, 4] = 999
此示例更改两个像素的值[1, 4] 还有那部分[30:40,10:20] 新值999。见 Numpy documentation 有关Python样式数组索引和切片的详细信息。
下一个数组操作示例是将图像数据从计数转换为通量:
>>> photflam = hdul[1].header['photflam']
>>> exptime = hdr['exptime']
>>> data = data * photflam / exptime
>>> hdul.close()
请注意,对整个图像执行这样的操作需要将整个图像保存在内存中。这个例子在适当的地方执行乘法,这样就不会生成副本,但是原始图像必须首先能够放入主内存中。对于大多数观察来说,这在现代个人电脑上不应该是一个问题。
如果此时要保留所做的所有更改并将其写入新文件,则可以使用 HDUList.writeto()
方法(见下文)。
实例:
也见 从FITS文件读取并打印图像 .
使用表数据¶
本节介绍如何使用 fits
直接打包。但在某些情况下,高层 统一文件读写接口 通常就足够了,而且使用起来更方便一些。见 Unified I/O FITS 详细信息。
与图像一样,FITS表扩展的数据部分位于 .data
属性:
>>> fits_table_filename = fits.util.get_testdata_filepath('tb.fits')
>>> hdul = fits.open(fits_table_filename)
>>> data = hdul[1].data # assuming the first extension is a table
如果你熟悉 numpy
recarray
对象时,您会发现表数据基本上是一个具有一些额外属性的记录数组。但熟悉记录数组并不是本指南的先决条件。
要查看表的第一行:
>>> print(data[0])
(1, 'abc', 3.7000000715255736, False)
表中的每一行都是 FITS_record
对象,它看起来像一个包含异构数据类型元素的(Python)元组。在本例中:一个整数、一个字符串、一个浮点数和一个布尔值。所以表数据只是这样的记录的一个数组。更常见的是,用户可能以列方式访问数据。这可以通过使用 field()
方法。要获取表的第一列(或NumPy术语中的“field”,在这里与“column”互换使用),请使用:
>>> data.field(0)
array([1, 2]...)
A numpy
返回具有指定字段数据类型的对象。
与标题关键字一样,列可以按索引引用,如上所述,也可以按名称引用:
>>> data.field('c1')
array([1, 2]...)
当按名称访问列时,也可以进行类似dict的访问(甚至更可取):
>>> data['c1']
array([1, 2]...)
在大多数情况下,最好按列的名称访问列,因为列名完全独立于它在表中的物理顺序。与头关键字一样,列名不区分大小写。
但是我们怎么知道表中有哪些列呢?首先,我们将介绍表HDU的另一个属性: columns
属性:
>>> cols = hdul[1].columns
此属性是 ColDefs
(列定义)对象。如果我们使用 ColDefs.info()
方法从交互提示::
>>> cols.info()
name:
['c1', 'c2', 'c3', 'c4']
format:
['1J', '3A', '1E', '1L']
unit:
['', '', '', '']
null:
[-2147483647, '', '', '']
bscale:
['', '', 3, '']
bzero:
['', '', 0.4, '']
disp:
['I11', 'A3', 'G15.7', 'L6']
start:
['', '', '', '']
dim:
['', '', '', '']
coord_type:
['', '', '', '']
coord_unit:
['', '', '', '']
coord_ref_point:
['', '', '', '']
coord_ref_value:
['', '', '', '']
coord_inc:
['', '', '', '']
time_ref_pos:
['', '', '', '']
它将显示表中所有列的属性,例如它们的名称、格式、bscales、bzero等。一个类似的输出将显示列名及其格式,可以从脚本中使用以下命令打印:
>>> hdul[1].columns
ColDefs(
name = 'c1'; format = '1J'; null = -2147483647; disp = 'I11'
name = 'c2'; format = '3A'; disp = 'A3'
name = 'c3'; format = '1E'; bscale = 3; bzero = 0.4; disp = 'G15.7'
name = 'c4'; format = '1L'; disp = 'L6'
)
我们也可以单独获取这些属性;例如:
>>> cols.names
['c1', 'c2', 'c3', 'c4']
返回字段名的(Python)列表。
因为每个字段都是 numpy
反对,我们将拥有 numpy
使用的工具。我们可以重新分配(更新)值:
>>> data['c4'][:] = 0
取一列的平均值:
>>> data['c3'].mean()
5.19999989271164
等等。
实例:
也见 访问存储在MEF扩展名中的数据 .
保存文件更改¶
如前所述,当用户打开一个文件,对头或数据进行一些更改后,用户可以使用 HDUList.writeto()
保存更改。这将获取内存中头和数据的版本,并将它们写入磁盘上的新FITS文件。可以对内存中的数据执行后续操作,并将其写入另一个不同的文件,而无需将原始数据重新复制到(更多)内存中:
hdul.writeto('newtable.fits')
将写入的当前内容 hdulist
到一个新的磁盘文件新文件.fits. 如果文件是以更新模式打开的,则 HDUList.flush()
方法也可用于写入自 open()
,返回原始文件。这个 close()
方法将对以更新模式打开的FITS文件执行相同的操作:
with fits.open('original.fits', mode='update') as hdul:
# Change something in hdul.
hdul.flush() # changes are written back to original.fits
# closing the file will also flush any changes and prevent further writing
创建新的配合文件¶
创建新图像文件¶
到目前为止,我们已经演示了如何读取和更新现有的FITS文件。但是从头开始创建一个新的FITS文件怎么样?这样的任务在 astropy
对于图像HDU。我们将首先演示如何创建一个FITS文件,该文件仅由主HDU和图像数据组成。
首先,我们创建一个 numpy
数据部分的对象:
>>> import numpy as np
>>> n = np.arange(100.0) # a simple sequence of floats from 0.0 to 99.9
接下来,我们创建一个 PrimaryHDU
对象来封装数据:
>>> hdu = fits.PrimaryHDU(n)
然后我们创建一个 HDUList
要包含新创建的主HDU,并写入新文件:
>>> hdul = fits.HDUList([hdu])
>>> hdul.writeto('new1.fits')
就这样!事实上, astropy
甚至还提供了最后两行完成相同行为的快捷方式:
>>> hdu.writeto('new2.fits')
这将向FITS文件写入单个HDU,而不必手动将其封装到 HDUList
对象优先。
创建新表文件¶
备注
如果要创建 二元的 适合没有其他HDU的桌子,您可以使用 Table
而是写信给FITS。这比“低层”配合界面复杂:
>>> from astropy.table import Table
>>> t = Table([[1, 2], [4, 5], [7, 8]], names=('a', 'b', 'c'))
>>> t.write('table1.fits', format='fits')
等效代码使用 astropy.io.fits
看起来像这样:
>>> from astropy.io import fits
>>> import numpy as np
>>> c1 = fits.Column(name='a', array=np.array([1, 2]), format='K')
>>> c2 = fits.Column(name='b', array=np.array([4, 5]), format='K')
>>> c3 = fits.Column(name='c', array=np.array([7, 8]), format='K')
>>> t = fits.BinTableHDU.from_columns([c1, c2, c3])
>>> t.writeto('table2.fits')
HDU需要更多的信息,因为HDU需要更多的信息。首先,表只能是扩展HDU,而不是主扩展。有两种FITS表扩展名:ASCII和binary。我们将在这里使用二进制表示例。
要从头开始创建表,我们需要首先定义列,通过构造 Column
对象及其数据。假设我们有两列,第一列包含字符串,第二列包含浮点数:
>>> import numpy as np
>>> a1 = np.array(['NGC1001', 'NGC1002', 'NGC1003'])
>>> a2 = np.array([11.1, 12.3, 15.2])
>>> col1 = fits.Column(name='target', format='20A', array=a1)
>>> col2 = fits.Column(name='V_mag', format='E', array=a2)
备注
没有必要创建 Column
如果数据存储在 structured array .
下一步,创建一个 ColDefs
(列定义)所有列的对象::
>>> cols = fits.ColDefs([col1, col2])
现在,使用 BinTableHDU.from_columns()
功能:
>>> hdu = fits.BinTableHDU.from_columns(cols)
此函数返回(在本例中)a BinTableHDU
.
用于表示FITS表的数据结构称为 FITS_rec
并且是从 numpy.recarray
接口。当创建一个新的表HDU时,各个列数组将被组合成一个 FITS_rec
数组。
您可以创建 BinTableHDU
更简洁,不为单个列创建中间变量,也不需要手动创建 ColDefs
对象:
>>> hdu = fits.BinTableHDU.from_columns(
... [fits.Column(name='target', format='20A', array=a1),
... fits.Column(name='V_mag', format='E', array=a2)])
现在您可以将这个新的表HDU直接写入FITS文件,如下所示:
>>> hdu.writeto('table3.fits')
此快捷方式将自动创建一个没有数据的最小主HDU,并将其前置到表HDU以创建有效的FITS文件。如果您在主HDU中需要额外的数据或头关键字,您仍然可以创建一个 PrimaryHDU
对象并使用 HDUList
,如下一节所述。
创建具有多个扩展名的文件¶
在前面的例子中,我们创建了一个有意义的扩展名(a PrimaryHDU
或 BinTableHDU
). 要创建具有多个扩展名的文件,我们需要创建扩展名hdu并将它们附加到 HDUList
.
首先,我们为图像扩展创建一些数据:
>>> import numpy as np
>>> n = np.ones((3, 3))
>>> n2 = np.ones((100, 100))
>>> n3 = np.ones((10, 10, 10))
请注意,不同扩展的数据形状不需要相同。接下来,将数据单独放入 PrimaryHDU
和 ImageHDU
物体::
>>> primary_hdu = fits.PrimaryHDU(n)
>>> image_hdu = fits.ImageHDU(n2)
>>> image_hdu2 = fits.ImageHDU(n3)
一个多扩展名FITS文件不被限制为只有成像或表格数据,我们可以混合它们。为了说明这一点,我们将使用上一节中的示例来生成一个 BinTableHDU
::
>>> c1 = fits.Column(name='a', array=np.array([1, 2]), format='K')
>>> c2 = fits.Column(name='b', array=np.array([4, 5]), format='K')
>>> c3 = fits.Column(name='c', array=np.array([7, 8]), format='K')
>>> table_hdu = fits.BinTableHDU.from_columns([c1, c2, c3])
现在当我们创建 HDUList
我们列出所有要包括的扩展名:
>>> hdul = fits.HDUList([primary_hdu, image_hdu, table_hdu])
因为 HDUList
表现得像 list
例如,我们还可以附加 ImageHDU
已经存在的 HDUList
::
>>> hdul.append(image_hdu2)
多重延伸 HDUList
就像那些只有 PrimaryHDU
,因此要保存文件,请使用 HDUList.writeto()
如上图所示。
备注
FITS标准强制所有文件只有一个 PrimaryHDU
这是文件中的第一个HDU。本标准在 HDUList.writeto()
如果不满足,将引发错误。见 output_verify
选择权 HDUList.writeto()
找到修复或忽略这些警告的方法。
在前面的示例中 PrimaryHDU
包含实际数据。在某些情况下,最好将 PrimaryHDU
只有基本的标题信息。为此,首先创建一个新的 Header
对象封装要包含在主HDU中的任何关键字,然后像以前一样创建 PrimaryHDU
::
>>> hdr = fits.Header()
>>> hdr['OBSERVER'] = 'Edwin Hubble'
>>> hdr['COMMENT'] = "Here's some commentary about this FITS file."
>>> empty_primary = fits.PrimaryHDU(header=hdr)
当我们创建一个新的主HDU时,如上面的例子所示,它将自动包含任何额外的头关键字 必修的 按FITS格式(关键字如 SIMPLE
和 NAXIS
例如)。一般来说,用户不必手动管理这些关键字,而应该只创建和修改观察特定的信息关键字。
然后我们创建一个HDUList,它包含主HDU和任何其他HDU想要的HDU::
>>> hdul = fits.HDUList([empty_primary, image_hdu2, table_hdu])
实例:
也见 从头开始创建多扩展配合(MEF)文件 .
便利功能¶
astropy.io.fits
还提供了一些高级(“便利”)功能。这样一个方便的功能是一个“屏蔽”操作来实现一个任务。通过使用这些“便利”功能,用户不必担心打开或关闭文件;所有的内务处理都是隐式的。
警告
这些函数对于交互式Python会话和不太复杂的分析脚本很有用,但不应用于应用程序代码,因为它们效率很低。例如,每次调用 getval()
需要重新分析整个FITS文件。重复使用这些函数的代码应改为使用 open()
直接访问数据结构。
第一个功能是 getheader()
,以获取HDU的标头。下面是几个获取标题的示例。此函数只需要文件名。其余参数是可选的,可以灵活地指定用户想要访问哪个HDU::
>>> from astropy.io.fits import getheader
>>> hdr = getheader(fits_image_filename) # get default HDU (=0), i.e. primary HDU's header
>>> hdr = getheader(fits_image_filename, 0) # get primary HDU's header
>>> hdr = getheader(fits_image_filename, 2) # the second extension
>>> hdr = getheader(fits_image_filename, 'sci') # the first HDU with EXTNAME='SCI'
>>> hdr = getheader(fits_image_filename, 'sci', 2) # HDU with EXTNAME='SCI' and EXTVER=2
>>> hdr = getheader(fits_image_filename, ('sci', 2)) # use a tuple to do the same
>>> hdr = getheader(fits_image_filename, ext=2) # the second extension
>>> hdr = getheader(fits_image_filename, extname='sci') # first HDU with EXTNAME='SCI'
>>> hdr = getheader(fits_image_filename, extname='sci', extver=2)
不明确的规范将引发异常:
>>> getheader(fits_image_filename, ext=('sci', 1), extname='err', extver=2)
Traceback (most recent call last):
...
TypeError: Redundant/conflicting extension arguments(s): ...
获取头之后,可以访问其中的信息,例如获取和修改关键字值::
>>> fits_image_2_filename = fits.util.get_testdata_filepath('o4sp040b0_raw.fits')
>>> hdr = getheader(fits_image_2_filename, 0) # get primary hdu's header
>>> filter = hdr['filter'] # get the value of the keyword "filter'
>>> val = hdr[10] # get the 11th keyword's value
>>> hdr['filter'] = 'FW555' # change the keyword value
对于header关键字,header就像一个字典,以及一个列表。用户可以通过名称或数字索引访问关键字,如本章前面所述。
如果用户只需要读取一个关键字,则 getval()
函数可以进一步简化为一个调用,而不是上面示例中所示的两个调用:
>>> from astropy.io.fits import getval
>>> # get 0th extension's keyword FILTER's value
>>> flt = getval(fits_image_2_filename, 'filter', 0)
>>> flt
'Clear'
>>> # get the 2nd sci extension's 11th keyword's value
>>> val = getval(fits_image_2_filename, 10, 'sci', 2)
>>> val
False
函数 getdata()
获取HDU的数据。类似 getheader()
,它只需要输入FITS文件名,而扩展名是通过可选参数指定的。它有一个额外的可选参数头。如果header设置为True,此函数将同时返回data和header,否则只返回data::
>>> from astropy.io.fits import getdata
>>> # get 3rd sci extension's data:
>>> data = getdata(fits_image_filename, 'sci', 3)
>>> # get 1st extension's data AND header:
>>> data, hdr = getdata(fits_image_filename, 1, header=True)
上面介绍的函数是供阅读的。下面几个函数演示了方便编写的函数:
>>> fits.writeto('out.fits', data, hdr)
这个 writeto()
函数使用提供的数据和可选的头来写入输出FITS文件。
>>> fits.append('out.fits', data, hdr)
这个 append()
函数将使用提供的数据和可选的头附加到现有的FITS文件。如果指定的输出文件不存在,它将创建一个。
from astropy.io.fits import update
update(filename, dat, hdr, 'sci') # update the 'sci' extension
update(filename, dat, 3) # update the 3rd extension
update(filename, dat, hdr, 3) # update the 3rd extension
update(filename, dat, 'sci', 2) # update the 2nd SCI extension
update(filename, dat, 3, header=hdr) # update the 3rd extension
update(filename, dat, header=hdr, ext=5) # update the 5th extension
这个 update()
函数将使用输入数据/标头更新指定的扩展名。第三个参数可以是与数据关联的头。如果第三个参数不是头,则假定它(和其他位置参数)是扩展规范。头和扩展规范也可以是关键字参数。
这个 printdiff()
函数将打印两个FITS文件的差异报告,包括标题和数据。前两个参数必须是two-FITS文件名或FITS-file对象与匹配的数据类型(即,如果使用字符串指定文件名,则两个输入都必须是字符串)。第三个参数是可选的扩展规范,具有相同的调用格式 getheader()
和 getdata()
. 此外,您可以添加 FITSDiff
班级。
from astropy.io.fits import printdiff
# get a difference report of ext 2 of inA and inB
printdiff('inA.fits', 'inB.fits', ext=2)
# ignore HISTORY and COMMMENT keywords
printdiff('inA.fits', 'inB.fits', ignore_keywords=('HISTORY','COMMENT')
最后, info()
函数将打印指定的FITS文件的信息:
>>> fits.info(fits_image_filename)
Filename: ...test0.fits
No. Name Ver Type Cards Dimensions Format
0 PRIMARY 1 PrimaryHDU 138 ()
1 SCI 1 ImageHDU 61 (40, 40) int16
2 SCI 2 ImageHDU 61 (40, 40) int16
3 SCI 3 ImageHDU 61 (40, 40) int16
4 SCI 4 ImageHDU 61 (40, 40) int16
这是最有用的方便函数之一,它可以在不查看任何细节的情况下获得给定文件包含的内容的概述。
命令行实用程序¶
为了方便起见 astropy
的子包在您的系统上安装实用程序,这些程序允许在不必打开Python解释器的情况下执行常见任务。这些实用程序包括:
性能提示¶
可以为设置数据数组 PrimaryHDU
和 ImageHDU
到A dask 数组。如果将其写入磁盘,将在写入时计算dask阵列,这将避免使用过多内存:
>>> import dask.array as da
>>> array = da.random.random((1000, 1000))
>>> from astropy.io import fits
>>> hdu = fits.PrimaryHDU(data=array)
>>> hdu.writeto('test_dask.fits')