问题描述
嗨，我想在fits文件的python上绘制光谱。 ESO提供了有关如何显示一维光谱的指南，以及可以正常工作的代码，如下所示：
import sys
from astropy.io import fits as pyfits
import numpy as np
#
hdulist = pyfits.open( "your_1d_spectrum_here.fits" )
# print column information
hdulist[1].columns
# get to the data part (in extension 1)
scidata = hdulist[1].data
wave = scidata[0][0]
arr1 = scidata[0][1]
arr2 = scidata[0][2]
# etc.
# where arr1 will contain the data corresponding to the column named: hdulist[1].columns[1]
# where arr2 will contain the data corresponding to the column named: hdulist[1].columns[2]
# etc.
# To plot using maptplotlib:
import matplotlib.pyplot as plt
plt.plot(wave,arr1)
plt.show()
但是当将这段代码放入并使用s1d文件时，我得到IndexError：列表索引超出范围
HDUlist.info（）包含一个索引为0的主HDU
No. Name Ver Type Cards Dimensions Format
0 PRIMARY 1 PrimaryHDU 3051 (313088,) float32
并打印此数据表明它是一维数组： [45.182148 71.06376 3.1499128 ... 631.7653 628.2799 632.91364]
所以我对如何绘制光谱图一无所知，因为它只有一维。
Eso还有一个更复杂的代码，可绘制https://archive.eso.org/cms/eso-data/help/1dspectra.html#Python网站上的光谱图，文件为1dspectrum.py，可以下载。在此输入文件后，也会出现索引错误，这种错误可以在gyazo https://gyazo.com/e30a1f39e5f33821b5a5c54f3939a028
中看到尝试过Iguananaut的解决方案 如下所示：
import numpy as np
from astropy.io import fits
from astropy.units import u
import matplotlib.pyplot as plt
import specutils
f = fits.open('C:/Users/Rp199/Desktop/J0608_59_harps_2018/HARPS.2018-05-23T23_44_45.005_s1d_A.fits')
f.info()
No. Name Ver Type Cards Dimensions Format
0 PRIMARY 1 PrimaryHDU 3055 (313093,) float32
检查标题，它们相同，所以我删除了BUNIT并将CTYPE1更改为CUNIT1
del header['BUNIT']
header['CTYPE1'] = header['CUNIT1']
del header['CTYPE1']
然后我尝试使用Spectrum1D.read
spec = specutils.Spectrum1D.read(f)
spec
这给我以下错误
IORegistryError Traceback (most recent call last)
<ipython-input-28-62aceed26d07> in <module>
----> 1 spec = specutils.Spectrum1D.read(f)
2 spec
~\anaconda3\lib\site-packages\astropy\nddata\mixins\ndio.py in __call__(self,*args,**kwargs)
54
55 def __call__(self,**kwargs):
---> 56 return registry.read(self._cls,**kwargs)
57
58
~\anaconda3\lib\site-packages\astropy\io\registry.py in read(cls,format,**kwargs)
517 fileobj = args[0]
518
--> 519 format = _get_valid_format(
520 'read',cls,path,fileobj,args,kwargs)
521
~\anaconda3\lib\site-packages\astropy\io\registry.py in _get_valid_format(mode,kwargs)
597 if len(valid_formats) == 0:
598 format_table_str = _get_format_table_str(cls,mode.capitalize())
--> 599 raise IORegistryError("Format could not be identified based on the"
600 " file name or contents,please provide a"
601 " 'format' argument.\n"
IORegistryError: Format could not be identified based on the file name or contents,please provide a 'format' argument.
The available formats are:
Format Read Write Auto-identify
----------------- ---- ----- -------------
APOGEE apStar Yes No Yes
APOGEE apVisit Yes No Yes
APOGEE aspcapStar Yes No Yes
ASCII Yes No Yes
ECSV Yes No Yes
HST/COS Yes No Yes
HST/STIS Yes No Yes
IPAC Yes No Yes
JWST s2d Yes No Yes
JWST s3d Yes No Yes
JWST x1d Yes No Yes
SDSS-I/II spSpec Yes No Yes
SDSS-III/IV spec Yes No Yes
Subaru-pfsObject Yes No Yes
iraf Yes No Yes
muscles-sed Yes No Yes
tabular-fits Yes Yes Yes
wcs1d-fits Yes No Yes
我假设由于我的文件是s1d.fits文件，因此格式将为wcs1d-fits？我尝试浏览Fits查看器中的标题以识别格式，但找不到任何内容。将wcs1d-fits设置为以下格式时：
spec = specutils.Spectrum1D.read(f,format="wcs1d-fits")
spec
我已经退回了，
TypeError Traceback (most recent call last)
<ipython-input-29-c038104effad> in <module>
----> 1 spec = specutils.Spectrum1D.read(f,format="wcs1d-fits")
2 spec
~\anaconda3\lib\site-packages\astropy\nddata\mixins\ndio.py in __call__(self,**kwargs)
521
522 reader = get_reader(format,cls)
--> 523 data = reader(*args,**kwargs)
524
525 if not isinstance(data,cls):
~\anaconda3\lib\site-packages\specutils\io\default_loaders\wcs_fits.py in wcs1d_fits_loader(file_name,spectral_axis_unit,flux_unit,hdu_idx,**kwargs)
63 logging.info("Spectrum file looks like wcs1d-fits")
64
---> 65 with fits.open(file_name,**kwargs) as hdulist:
66 header = hdulist[hdu_idx].header
67 wcs = WCS(header)
~\anaconda3\lib\site-packages\astropy\io\fits\hdu\hdulist.py in fitsopen(name,mode,memmap,save_backup,cache,lazy_load_hdus,**kwargs)
162 raise ValueError(f'Empty filename: {name!r}')
163
--> 164 return HDUList.fromfile(name,165 lazy_load_hdus,**kwargs)
166
~\anaconda3\lib\site-packages\astropy\io\fits\hdu\hdulist.py in fromfile(cls,**kwargs)
401 """
402
--> 403 return cls._readfrom(fileobj=fileobj,mode=mode,memmap=memmap,404 save_backup=save_backup,cache=cache,405 lazy_load_hdus=lazy_load_hdus,**kwargs)
~\anaconda3\lib\site-packages\astropy\io\fits\hdu\hdulist.py in _readfrom(cls,data,**kwargs)
1052 if not isinstance(fileobj,_File):
1053 # instantiate a FITS file object (ffo)
-> 1054 fileobj = _File(fileobj,cache=cache)
1055 # The Astropy mode is determined by the _File initializer if the
1056 # supplied mode was None
~\anaconda3\lib\site-packages\astropy\utils\decorators.py in wrapper(*args,**kwargs)
533 warnings.warn(message,warning_type,stacklevel=2)
534
--> 535 return function(*args,**kwargs)
536
537 return wrapper
~\anaconda3\lib\site-packages\astropy\io\fits\file.py in __init__(self,overwrite,cache)
193 self._open_filename(fileobj,overwrite)
194 else:
--> 195 self._open_filelike(fileobj,overwrite)
196
197 self.fileobj_mode = fileobj_mode(self._file)
~\anaconda3\lib\site-packages\astropy\io\fits\file.py in _open_filelike(self,overwrite)
544
545 if mode == 'ostream':
--> 546 self._overwrite_existing(overwrite,False)
547
548 # Any "writeable" mode requires a write() method on the file object
~\anaconda3\lib\site-packages\astropy\io\fits\file.py in _overwrite_existing(self,closed)
442 # The file will be overwritten...
443 if ((self.file_like and hasattr(fileobj,'len') and fileobj.len > 0) or
--> 444 (os.path.exists(self.name) and os.path.getsize(self.name) != 0)):
445 if overwrite:
446 if self.file_like and hasattr(fileobj,'truncate'):
~\anaconda3\lib\genericpath.py in exists(path)
17 """Test whether a path exists. Returns False for broken symbolic links"""
18 try:
---> 19 os.stat(path)
20 except (OSError,ValueError):
21 return False
TypeError: stat: path should be string,bytes,os.PathLike or integer,not method
解决方法
我看了从this page下载的随机频谱。实际上，它仅包含一个HDU，因此，如果尝试从不存在的HDU中访问数据，则会得到IndexError。
您链接到的代码可能很旧（毕竟它位于“存档”页面上），并且仅是示例，不一定适合于从任何FITS文件中绘制数据。
specutils软件包具有许多实用程序，可用于分析和绘制一维光谱。它也可以读取许多常见FITS格式的光谱，而无需手动进行太多操作。
这就是我所做的。首先，我打开了文件，尤其是这个文件：
>>> f = fits.open('https://www.eso.org/sci/facilities/lasilla/instruments/harps/inst/monitoring/sundata/ceres_2006-05-22_s1d.fits')
>>> f.info()
f = fits.open('https://www.eso.org/sci/facilities/lasilla/instruments/harps/inst/monitoring/sundata/ceres_2006-05-22_s1d.fits')
看看标题，尤其是头十几个标题关键字，可以为我们提供有关如何解释文件中数据的许多线索：
>>> header = f[0].header
>>> header[:12]
SIMPLE = T / file does conform to FITS standard
BITPIX = -32 / number of bits per data pixel
NAXIS = 1 / number of data axes
NAXIS1 = 313237 / length of data axis 1
EXTEND = T / FITS dataset may contain extensions
COMMENT FITS (Flexible Image Transport System) format is defined in 'Astronomy
COMMENT and Astrophysics',volume 376,page 359; bibcode: 2001A&A...376..359H
CRPIX1 = 1. / Reference pixel
CRVAL1 = 3781.19 / Coordinate at reference pixel
CDELT1 = 0.01 / Coordinate increment par pixel
CTYPE1 = 'Angstrom' / Units of coordinate
BUNIT = 'Relative Flux' / Units of data values
我们可以从BUNIT = 'Relative Flux'推断出，文件中的单个数组包含一些未指定单位的通量值（该文件也以几种方式格式错误，但我会在后面介绍）。
它不包含光谱阵列，而是光谱从第3781.19Å开始，在第一个像素处（我们从CRVAL1和CRPIX1获得；请记住CRPIXn像FITS文件中的其他索引一样，使用以1开头的Fortran样式索引，而不是以0开头的C / Python，因此这表示数组的第0个元素是在3781.19Å处的光谱通量）。
它也说CDELT1 = 0.01 / Coordinate increment par pixel，所以每个像素都是+0.01Å。
所以我们可以这样做：
>>> import numpy as np
>>> from astropy.units import u
>>> flux = f[0].data
>>> ref = header['CRVAL1']
>>> step = header['CDELT1']
>>> wave = np.arange(ref,ref + (len(flux) * step),step) * u.AA
>>> plt.plot(wave,flux)
Specutils原则上可以简化此过程。要从FITS文件中读取Spectrum1D（如果格式受支持），您可以执行以下操作：
>>> from specutils import Spectrum1D
>>> spec = Spectrum1D.read(f)
不幸的是，在这种情况下，Specutils无法读取文件，我得到了：
ValueError: 'Relative Flux' did not parse as unit: At col 0,Relative is not a valid unit. If this is meant to be a custom unit,define it with 'u.def_unit'. To have it recognized inside a file reader or other code,enable it with 'u.add_enabled_units'. For details,see https://docs.astropy.org/en/latest/units/combining_and_defining.html
这是我早些时候写的文件格式错误的意思的一部分。对于将出现在FITS文件的BUNIT标头中的任何单位，“相对通量”不是有效的已知名称。由于我不知道这种情况下的通量单位应该是什么，所以我删除了它：
>>> del header['BUNIT']
此文件的另一个问题是它没有正确使用CTYPE1关键字。这应该是CUNIT1。我不知道他们为什么这样做。我再次更新了标题：
>>> header['CUNIT1'] = header['CTYPE1]
>>> del header['CTYPE1']
现在可以使用了：
>>> spec = Spectrum1D.read(f)
>>> spec
<Spectrum1D(flux=<Quantity [0.0043928,0.00500706,0.00511265,...,0.00435138,0.00433926,0.00433523]>,spectral_axis=<SpectralAxis [3781.19,3781.2,3781.21,6913.53,6913.54,6913.55] Angstrom>)>
并可以按照如下方式绘制（仅遵循specutils文档）：
>>> from astropy.visualization import quantity_support
>>> quantity_support()
>>> fig,ax = plt.subplots()
>>> ax.step(spec.spectral_axis,spec.flux)
>>> fig.show()
设置时间 控制面板
错误1：Request method ‘DELETE‘ not supported 错误还原：...