Re: Astropy Fits issue

[Erik Bray]

Hi Jeremy (and CCing Astropy-dev since this should be opened as an
issue but I'm too lazy),

I see the problem. This is a bug in astropy.io.fits (née PyFITS),
which to my horror it seems has been around just about forever. You
have a FITS file with data in the ASCII table format which used to be
more common, but isn't used a lot today which is maybe why it hasn't
turned up.

As you wrote, pyfits seems to think the format code 'E16.6' means the
string data should be read into single precision floats, and there are
many values in that table that cannot be represented in
single-precision floating point (not just in the 'Flux' column
either). I can see how this bug came about--in the *binary* table
FITS format there is a column format code 'E' which definitely is for
single-precision floats. But for ASCII tables, to quote the FITS
standard:

"""
There is no difference between the F, D, and E formats; the content of
the string determines its interpretation. Numbers requiring more
precision and/or range than the local computer can support may be
represented. It is good form to specify a D format in TFORMn for a
column of an ASCII table when that column will contain numbers that
cannot be accurately represented in 32-bit IEEE binary format (see
Appendix E).
"""

If your FITS file were in "good form" it would use "D" instead of "E"
since it contains values that cannot be represented in 32-bit floating
point. But apparently that is just "good form" and not a requirement
(I "love" FITS). In any case it's possible that was misinterpreted at
some point to mean that floats in an ASCII table should *only* be read
into double-precision machine floats if the format is "D", and that is
not the case. In fact for modern use there's no reason not to just
default to doubles.

There is an easy workaround you can employ right now. Just tweak the
table that pyfits uses to map ASCII table format codes to Numpy type
codes:

>>> from astropy.io import fits
>>> fits.column.ASCII2NUMPY['E'] = fits.column.ASCII2NUMPY['F8'] = 'f8'

Then proceed as normal. The code should be fixed so that that's the default.

Thanks,
Erik


On Thu, Sep 1, 2016 at 8:39 PM, CHASTENET Jérémy (OBS)
<jeremy.c...@astro.unistra.fr> wrote:
> Good afternoon,
>
> I apologize for using your personal email address; I couldn't manage to find your professional one.
> Please let me know if you'd rather have me use a different one.
>
> I'm a PhD student at the Space Telescope Institute. I'm encountering an issue while using the Fits function of Astropy. I've been advised to reach out to you since you've been contributing to that part.
>
> As I am reading a fits file, some of the columns output 'Inf' numbers. I've identified this as a precision definition issue. One of the keyword for that column is format: 'E16.6'.
> I was wondering if there is a way to go around the default precision and force the reading in double precision.
> Please find attached the said fits file. (The column of interested in 'Flux', column 4)
> (As a check, I know the number are not corrupted since IDL can read them out.)
>
> I hope you'll be able to help.
> Best regards,
> Jeremy Chastenet
>

[Erik Bray]

On Fri, Sep 2, 2016 at 10:08 AM, Erik Bray <erik....@gmail.com> wrote:
>>>> from astropy.io import fits
>>>> fits.column.ASCII2NUMPY['E'] = fits.column.ASCII2NUMPY['F8'] = 'f8'

To be clear, there's a typo here. This should be:

>>> fits.column.ASCII2NUMPY['E'] = fits.column.ASCII2NUMPY['F'] = 'f8'