0

I am trying to update a FITS file with a new column of data. My file has a Primary HDU, and two other HDUs, each one including a table.

Since adding a new column to the table of an already existing FITS file is a pain (unsolvable, see here and here), I changed my mind and try to focus on creating a new file with a modified table.

This means I have to copy all the rest from the original file (Primary HDU, other HDUs, etc.). Is there a standard way to do this? Or, what is the best (fastest?) way, possibly avoiding to copy each element one by one "by hand"?

Community
  • 1
  • 1
Py-ser
  • 1,475
  • 6
  • 24
  • 46

1 Answers1

0

On the topic of adding a new columns, have you seen this documentation? This is the most straightforward way to create a new table with the new column added. This necessarily involves creating a new binary table HDU since it describes different data.

Or have you looked into the Astropy table interface? It supports reading and writing FITS tables. See here. It basically works the same way but goes to some more efforts to the hide the details. This is the interface the PyFITS/astropy.io.fits interface is gradually being replaced with since it actually provides a good table interface.

Adding a new HDU or replacing an existing HDU in an existing FITS file is simply a matter of opening that file and updating the HDUList data structure (which works like a normal Python list) and writing the updated HDUList to a new file.

A full example might look something like:

try:
    from astropy.io import fits
except ImportError:
    import pyfits as fits

with fits.open('path/to/file.fits') as hdul:
    table_hdu = hdul[1]  # If the table is the first extension HDU
    new_column = fits.Column(name='NEWCOL', format='D', array=np.zeros(len(table_hdu.data)))
    new_columns = fits.ColDefs([new_column])
    new_table_hdu = fits.new_table(table_hdu.columns + new_columns)

    # Replace the original table HDU with the new one
    hdul[1] = new_table_hdu
    hdul.writeto('path/to/new_file.fits')

Something roughly like that should work. This will be easier in Astropy once the new Table interface is fully integrated but for now that's what it involves. There is no reason to do anything "by hand" so to speak.

Iguananaut
  • 15,675
  • 4
  • 43
  • 50