[AstroPy] Reading in wavelength-calibrated spectra

[Rick White]

I've attached a module that reads at least some versions of these IRAF multispec format files and returns flux and wavelength arrays.  I've had an IDL version of this for a long time now (well-tested), but Kevin's message prompted me to make a Python translation.  This is only lightly tested (finished it this morning!) but it worked on a few different files I had lying around.  I already sent it to Kevin off list and he confirmed that it worked on his file.  So I'm circulating it in case other people find it useful.

I agree with the comment that it will be useful for astropy to be able to read this format.  I don't have any code that writes the format (and I'm not convinced we want any).  This format is a fine supporting example for Perry's argument that we need to replace FITS format with something more general!

Here's a sketch of how to use it:

import readmultispec
rdict = readmultispec.readmultispec('yourfitsfile.fits')
wavelen = rdict['wavelen']
flux = rdict['flux']

The returned dictionary includes a couple of other values (the FITS header, the raw wavelength coefficients from the header).   The header might be useful; the 'wavefields' array is probably not very useful.

The FITS file can be in a variety of formats, including multi-order echelle spectra.  In the general case you'll get back arrays that look like:

wavelen[norders, nwavelengths]
flux[ncomponents, norders, nwavelengths]

If it is a 1-D spectrum (with just one order), these become

wavelen[nwavelengths]
flux[ncomponents, nwavelengths]

The components dimension for IRAF-extracted spectra typically has 4 elements: the extracted spectrum = flux[0,*], an alternative extracted spectrum = flux[1,*], the subtracted sky = flux[2,*], and the error estimates = flux[3,*].  Not all files will have all 4 of these components.  And it's not that easy to figure out from the header exactly what you've got when the number of components is less than 4.  But generally the first component should be the extracted spectrum.

This is not completely general, but it probably could be a reasonable starting point for a complete module.
Cheers,
Rick White

On Aug 14, 2012, at 4:10 PM, Kevin Gullikson wrote:

> Hi,
> 
> I have some spectra that were reduced in iraf, and I would like to read them in to python to do some fitting that I already have code for. However, I am having trouble reading in the wavelength calibration. The spectra are echelle, and the wavelengths are defined with this thing in the header that looks like this:
> 
> "
> WAT2_001= 'wtype=multispec spec1 = "1 32 2 10672.221264362 0.086046403007761 20'
> WAT2_002= '46 0. 6.65 26.40 1. 0. 1 5 1. 2046. 10761.5346117191 87.998952480064'
> WAT2_003= ' -1.32409633589482 -0.0165054046288388 -0.00680394594704411" spec2 ='
> WAT2_004= ' "2 33 2 10348.8582697 0.083445037904806 2046 0. 34.60 58.07 1. 0. 1'
> WAT2_005= ' 5 1. 2046. 10435.4699172677 85.3379030819273 -1.2833101652519 -0.01'
> WAT2_006= '53518242619121 -0.0057861452751027" spec3 = "3 34 2 10043.881872825 '
> ...
> "
> 
> The numbers give the conversion from pixel coordinates to wavelength coordinates.
> I was wondering if there was a python function that could correctly read all that in, parse it, and give me a spectrum in flux vs wavelength? I have tried playing with pywcs but that seems only to work for images?
> 
> Cheers,
> 
> Kevin Gullikson _______________________________________________
> AstroPy mailing list
> AstroPy at scipy.org
> http://mail.scipy.org/mailman/listinfo/astropy
-------------- next part --------------
A non-text attachment was scrubbed...
Name: readmultispec.py
Type: text/x-python-script
Size: 9148 bytes
Desc: readmultispec.py
URL: <http://mail.python.org/pipermail/astropy/attachments/20120815/f5cead48/attachment.bin>