JuliaAstro / FITSIO.jl

Flexible Image Transport System (FITS) file support for Julia
http://juliaastro.org/FITSIO.jl/
MIT License
55 stars 29 forks source link

Add PrimaryHDU #111

Open giordano opened 5 years ago

giordano commented 5 years ago

Currently, we have the following HDUs:

julia> subtypes(HDU)
3-element Array{Any,1}:
 ASCIITableHDU
 ImageHDU     
 TableHDU

However the FITS expects at least one other type, "Primary HDU", see https://fits.gsfc.nasa.gov/fits_primer.html.

One example of FITS file with a Primary HDU has been reported at https://github.com/JuliaAstro/AstroImages.jl/pull/9: EUVEngc4151imgx.fits. FITSIO.jl shows the following table:

julia> using FITSIO

julia> file = download("https://fits.gsfc.nasa.gov/samples/EUVEngc4151imgx.fits")
"/tmp/julia6zY8AX"

julia> FITS(file)
File: /tmp/julia6zY8AX
Mode: "r" (read-only)
HDUs: Num  Name             Type   
      1                     Image  
      2    ds               Image  
      3    sw_night         Image  
      4    mw               Image  
      5    lw               Image  
      6    ds_limits        Table  
      7    sw_night_limits  Table  
      8    mw_limits        Table  
      9    lw_limits        Table

instead Astropy reports:

>>> from astropy.io import fits
>>> fits.open("/tmp/julia6zY8AX").info()
Filename: /tmp/julia6zY8AX
No.    Name      Ver    Type      Cards   Dimensions   Format
  0  PRIMARY       1 PrimaryHDU      62   ()      
  1  ds            1 ImageHDU        72   (512, 512)   int16   
  2  sw_night      1 ImageHDU        78   (2048, 300)   int16   
  3  mw            1 ImageHDU        78   (2048, 300)   int16   
  4  lw            1 ImageHDU        78   (2048, 300)   int16   
  5  ds_limits     1 BinTableHDU     71   3R x 3C   [8A, 1E, 1E]   
  6  sw_night_limits    1 BinTableHDU     71   2R x 3C   [12A, 1E, 1E]   
  7  mw_limits     1 BinTableHDU     71   2R x 3C   [12A, 1E, 1E]   
  8  lw_limits     1 BinTableHDU     71   2R x 3C   [12A, 1E, 1E]
kbarbary commented 5 years ago

Multi-extension FITS files often use this style of "header-only" primary HDU, for metadata that applies to all the extensions in the file. It is also used for FITS files containing only table HDUs, as the first HDU cannot be a table HDU.

It's not clear that a separate subtype is necessary for the "primary HDU." As prior art, the Python package fitsio (on which the FITSIO.jl API is based) does not have a separate type for the primary HDU. The primary HDU is simply the first HDU (with the additional restriction that it can't be a table HDU). So it can only be an ImageHDU (with or without a data segment). Personally, I find this helpful: when you see that it is an ImageHDU, you immediately know that there's nothing special about this first HDU (other than the fact that it perhaps contains no data; but this can actually be the case for any ImageHDU!).

I do think it would be helpful to have a quick method for testing whether an ImageHDU contains data, and also helpful to print the dimension information in show(::FITS), as astropy.io.fits does.

Without looking at the issue, I suspect that the right approach in JuliaAstro/AstroImages.jl#9 should be to look for the first ImageHDU with a non-empty data array. (And perhaps specifically a 2-d data array at that!).

giordano commented 5 years ago

Without looking at the issue, I suspect that the right approach in JuliaAstro/AstroImages.jl#9 should be to look for the first ImageHDU with a non-empty data array. (And perhaps specifically a 2-d data array at that!).

I didn't look into how Astropy detects a PrimaryHDU, but I was thinking about something like this.

giordano commented 5 years ago

For the record, also fv indicates the first extension as "Image", however it is named "Primary", instead FITSIO.jl leaves it unnamed, Astropy calls it "PRIMARY".

image

ziotom78 commented 5 years ago

There is another usage of the Primary HDU that I have found sometimes in astronomical and cosmological codes: using it to encode some complex metadata in a high-level format (XML, JSON, YAML) that are associated with the file as a whole.

I heard this technique during a discussion at ADASS2016 («The future of FITS», if I remember correctly). I have used it several times, the idea is to use something like this:

import JSON
import FITSIO

# Encode some complex structure
s = JSON.json(Dict(:a => 10, :b => ["foo", "bar"]))
f = FITSIO.FITS("newfile.fits", "w");
# Convert the string into an array of bytes
write(f, convert.(UInt8, collect(s)))
close(f)

Don't know how much this trick is widespread, however. Maybe AstroImages.jl could just check that the image in the PrimaryHDU has more than one row.

giordano commented 5 years ago

A side note about "the future of FITS": there is already a julia package to work with ASDF files: https://github.com/eschnett/ASDF.jl