pyfits backend fails to read in a SWIFT RMF

[DougBurke]

From #354 (comment)

If you use the pyfits/astropy IO backend, you can't load in the SWIFT RMF swxpc0to12s6_20130101v014.rmf - one location for this file is http://darts.jaxa.jp/pub/legacy.gsfc.nasa.gov/caldb/data/swift/xrt/cpf/rmf/swxpc0s6_20130101v014.rmf

With pycrates (in CIAO 4.9) I can load it in:

sherpa-1> rmf = unpack_rmf('swxpc0to12s6_20130101v014.rmf')
sherpa-2> print(rmf)
name     = swxpc0to12s6_20130101v014.rmf
detchans = 1024
energ_lo = Float64[2400]
energ_hi = Float64[2400]
n_grp    = UInt64[2400]
f_chan   = UInt32[2400]
n_chan   = UInt32[2400]
matrix   = Float64[2457600]
offset   = 0
e_min    = Float64[1024]
e_max    = Float64[1024]

but with astropy 1.3.2 installed

In [2]: rmf = sherpa.astro.ui.unpack_rmf('swxpc0to12s6_20130101v014.rmf')
---------------------------------------------------------------------------
IOErr                                     Traceback (most recent call last)
<ipython-input-2-8f4410910367> in <module>()
----> 1 rmf = sherpa.astro.ui.unpack_rmf('swxpc0to12s6_20130101v014.rmf')

<string> in unpack_rmf(arg)

/home/djburke/sherpa/sherpa-rgs/sherpa/astro/ui/utils.py in unpack_rmf(self, arg)
   5363 
   5364         """
-> 5365         return sherpa.astro.instrument.RMF1D(sherpa.astro.io.read_rmf(arg))
   5366 
   5367     # DOC-TODO: add an example of a grating/multiple response

/home/djburke/sherpa/sherpa-rgs/sherpa/astro/io/__init__.py in read_rmf(arg)
    181     read_rmf( RMFCrate )
    182     """
--> 183     data, filename = backend.get_rmf_data(arg)
    184     return DataRMF(filename, **data)
    185 

/home/djburke/sherpa/sherpa-rgs/sherpa/astro/io/pyfits_backend.py in get_rmf_data(arg, make_copy)
    747         data['matrix'] = None
    748         raise IOErr('badfile', filename,
--> 749                     " MATRIX column not a variable length field")
    750 
    751     good = (data['n_grp'] > 0)

IOErr: 'swxpc0to12s6_20130101v014.rmf' is not a filename or  MATRIX column not a variable length field

I haven't looked to see if this is due to a change in the backend code since this was written or, more likely, if it's just a bug in our code.

[olaurino]

It would be good if at some point we could align the backends so that they only differ in the way they invoke the different libraries to get data and metadata from the files. At least as much as possible. In other terms, the actual backend should be the implementation of the abstraction that extracts information (e.g. get the value for such column from such HDU) rather than the full logic. Admittedly that's much broader than the scope of the bug itself.

[DougBurke]

So, the pyfits_backend code explicitly wants a FITS variable-length array for the MATRIX column: line 746 is

if not isinstance(data['matrix'], _VLF):

where _VLF is astropy.io.fits.column._VLF.

However, the MATRIX column does not have to be a variable-length array, as discussed in the OGIP standard:

https://heasarc.gsfc.nasa.gov/docs/heasarc/caldb/docs/memos/cal_gen_92_002/cal_gen_92_002.html#tth_sEc3.1.2

Unsurprisingly, the SWIFT RMF just has a simple 2D matrix rather than taking advantage of the space-saving horror that is the variable-length-array approach:

% dmlist swxpc0to12s6_20130101v014.rmf cols
 
--------------------------------------------------------------------------------
Columns for Table Block MATRIX
--------------------------------------------------------------------------------
 
ColNo  Name                 Unit        Type             Range
   1   ENERG_LO             keV          Real4          -Inf:+Inf            
   2   ENERG_HI             keV          Real4          -Inf:+Inf            
   3   N_GRP                             Int2           -                    
   4   F_CHAN                            Int2           0:32767              
   5   N_CHAN                            Int2           -                    
   6   MATRIX[1024]                      Real4(1024)    -Inf:+Inf            

The crates_backend delegates the RMF reading to code in pycrates that handles this case, whereas we have to do it ourselves with astropy.io.fits, and we have just missed this possibility.

[DougBurke]

@olaurino - while I agree with you that a lot of consolidation of logic would be good - indeed, just looking at the RMF reading code they create subtly-different error messages when files can't be interpreted as an RMF file - it's also shows how hard it can be. The pycrates code is simpler here since it can rely on the logic and code in pycrates.rmfcratedataset.RMFCrateDataset to do all the hard work, whereas we don't have a similar abstraction in astropy.io.fits.

[dhuppenkothen]

I just found this issue. I'm having the same problem with simulations for both Athena X-IFU and eXTP SFA. They're a long way off, but it'd be useful to be able to work with the simulations in sherpa!

[DougBurke]

@dhuppenkothen do you have a response that you can share that shows this problem, or see if PR #358 works for you?

[dhuppenkothen]

PR #358 indeed works for me! Thanks for that! It'd be really useful to get that into master with the next release.

[DougBurke]

Thanks for checking. We're a bit resource constrained at the moment, so development of Sherpa is going slower than we'd like.