Subaru/FOCAS support

[monodera]

There is an IFU mode for FOCAS at Subaru Telescope (https://www2.nao.ac.jp/~shinobuozaki/focasifu/) and I wanted to use ZAP for data taken with it. I've just copied and pasted the code blocks for KCWI so that ZAP runs for my FOCAS IFU data. The current FOCAS IFU pipeline () produces a FITS file with the hdu[0] for the cube data.

The modification I made (https://github.com/musevlt/zap/commit/9d7c596cb15648c278ed1b88f971615b770589db) apparently works nicely with my science cubes (though I haven't done extensive test).

However, it failed to produce the explained variance curves with the following error message. Do you have any ideas how to fix it?

I also think that it would be nice to have options to specify which extension (and axis configuration) ZAP should look at rather than coding instrument-by-instrument basis.

Python 3.6.9 |Anaconda, Inc.| (default, Jul 30 2019, 13:42:17)
Type 'copyright', 'credits' or 'license' for more information
IPython 7.8.0 -- An enhanced Interactive Python. Type '?' for help.

In [1]: import zap
zobj = zap.
In [2]: zobj = zap.process('cube01.fits', outcubefits='output01.fits')
[INFO] Running ZAP 2.1 !
WARNING: FITSFixedWarning: RADECSYS= 'FK5 ' / The equatorial coordinate system
the RADECSYS keyword is deprecated, use RADESYSa. [astropy.wcs.wcs]
[WARNING] FITSFixedWarning: RADECSYS= 'FK5 ' / The equatorial coordinate system
the RADECSYS keyword is deprecated, use RADESYSa.
WARNING: FITSFixedWarning: 'datfix' made the change 'Set MJD-OBS to 58664.000000 from DATE-OBS.
Set DATE-END to '2019-06-30T06:35:48.104' from MJD-END'. [astropy.wcs.wcs]
[WARNING] FITSFixedWarning: 'datfix' made the change 'Set MJD-OBS to 58664.000000 from DATE-OBS.
Set DATE-END to '2019-06-30T06:35:48.104' from MJD-END'.
[INFO] Cleaning NaN values in the cube
[INFO] Rejected 0 spaxels with more than 25.0% NaN pixels
[INFO] Fixing 0 remaining NaN pixels
[INFO] _nanclean - Time: 0.36 sec.
[INFO] Extract to 2D, 1536 valid spaxels (100%)
[INFO] _extract - Time: 0.08 sec.
[INFO] Median zlevel subtraction
[INFO] _zlevel - Time: 0.12 sec.
[INFO] Applying Continuum Filter, cftype=median
[INFO] Using cfwidth=300
[INFO] _continuumfilter - Time: 2.55 sec.
[INFO] _prepare - Time: 3.14 sec.
[INFO] Calculating SVD on 2 segments ([[   0 3561]
 [3561 3781]])
[INFO] _msvd - Time: 0.60 sec.
[INFO] Compute number of components
[INFO] Choosing [33  6] eigenspectra for segments
[INFO] Reconstructing Sky Residuals
[INFO] reconstruct - Time: 0.05 sec.
[INFO] Applying correction and reshaping data product
[INFO] Cube file saved to output01.fits
[INFO] Zapped! (took 4.60 sec.)

In [3]: zobj = zap.process('cube01.fits', outcubefits='output01-2.fits', varcurvefits='output01_varcurve.fits')
[INFO] Running ZAP 2.1 !
WARNING: FITSFixedWarning: RADECSYS= 'FK5 ' / The equatorial coordinate system
the RADECSYS keyword is deprecated, use RADESYSa. [astropy.wcs.wcs]
[WARNING] FITSFixedWarning: RADECSYS= 'FK5 ' / The equatorial coordinate system
the RADECSYS keyword is deprecated, use RADESYSa.
WARNING: FITSFixedWarning: 'datfix' made the change 'Set MJD-OBS to 58664.000000 from DATE-OBS.
Set DATE-END to '2019-06-30T06:35:48.104' from MJD-END'. [astropy.wcs.wcs]
[WARNING] FITSFixedWarning: 'datfix' made the change 'Set MJD-OBS to 58664.000000 from DATE-OBS.
Set DATE-END to '2019-06-30T06:35:48.104' from MJD-END'.
[INFO] Cleaning NaN values in the cube
[INFO] Rejected 0 spaxels with more than 25.0% NaN pixels
[INFO] Fixing 0 remaining NaN pixels
[INFO] _nanclean - Time: 0.02 sec.
[INFO] Extract to 2D, 1536 valid spaxels (100%)
[INFO] _extract - Time: 0.08 sec.
[INFO] Median zlevel subtraction
[INFO] _zlevel - Time: 0.06 sec.
[INFO] Applying Continuum Filter, cftype=median
[INFO] Using cfwidth=300
[INFO] _continuumfilter - Time: 2.55 sec.
[INFO] _prepare - Time: 2.75 sec.
[INFO] Calculating SVD on 2 segments ([[   0 3561]
 [3561 3781]])
[INFO] _msvd - Time: 0.48 sec.
[INFO] Compute number of components
[INFO] Choosing [33  6] eigenspectra for segments
[INFO] Reconstructing Sky Residuals
[INFO] reconstruct - Time: 0.05 sec.
[INFO] Applying correction and reshaping data product
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-3-a0bd10bfa0c3> in <module>
----> 1 zobj = zap.process('cube01.fits', outcubefits='output01-2.fits', varcurvefits='output01_varcurve.fits')

~/anaconda3/envs/astroconda/lib/python3.6/site-packages/zap/zap.py in process(cubefits, outcubefits, clean, zlevel, cftype, cfwidthSVD, cfwidthSP, nevals, extSVD, skycubefits, mask, interactive, ncpu, pca_class, n_components, overwrite, varcurvefits)
    186
    187     if varcurvefits is not None:
--> 188         zobj.writevarcurve(varcurvefits=varcurvefits, overwrite=overwrite)
    189
    190     zobj.mergefits(outcubefits, overwrite=overwrite)

~/anaconda3/envs/astroconda/lib/python3.6/site-packages/zap/zap.py in writevarcurve(self, varcurvefits, overwrite)
    804         """Write the explained variance curves to an individual fits file."""
    805         from astropy.table import Table
--> 806         table = Table([m.explained_variance_ for m in self.models])
    807         hdu = fits.table_to_hdu(table)
    808         _newheader(self, hdu.header)

~/anaconda3/envs/astroconda/lib/python3.6/site-packages/astropy/table/table.py in __init__(self, data, masked, names, dtype, meta, copy, rows, copy_indices, **kwargs)
    536
    537         # Finally do the real initialization
--> 538         init_func(data, names, dtype, n_cols, copy)
    539
    540         # Set table meta.  If copy=True then deepcopy meta otherwise use the

~/anaconda3/envs/astroconda/lib/python3.6/site-packages/astropy/table/table.py in _init_from_list(self, data, names, dtype, n_cols, copy)
    836             cols.append(col)
    837
--> 838         self._init_from_cols(cols)
    839
    840     def _init_from_ndarray(self, data, names, dtype, n_cols, copy):

~/anaconda3/envs/astroconda/lib/python3.6/site-packages/astropy/table/table.py in _init_from_cols(self, cols)
    886         if len(lengths) > 1:
    887             raise ValueError('Inconsistent data column lengths: {0}'
--> 888                              .format(lengths))
    889
    890         # Set the table masking

ValueError: Inconsistent data column lengths: {1536, 220}

In [4]:

[saimn]

Hi,

Great to know that ZAP works for another instrument!

I also think that it would be nice to have options to specify which extension (and axis configuration) ZAP should look at rather than coding instrument-by-instrument basis.

Yes I agree that having an option to specify the extension would it be useful. For KCWI I went the easy way and just added the special case on INSTRUMEN, which is also interesting because then you don't have to specify the extension number. So, I would be happy to integrate your commit adding FOCAS with the INSTRUMEN check (please make a pull request with this commit), and you or someone else wants to have a more flexible way to select the extension I will be happy to review. I'm no more working on MUSE so I cannot spend too much time on this.

About your bug with writevarcurve, ZAP used to use "segments", dividing the wavelength axis in smaller pieces to optimize the processing. I changed this to use only one segment by default (though it is still possible to use more if needed), but with hard coded values of (0, 10'000) Angstroms, which is not very clever... What is the wavelength coverage of your data ? The log above suggest that something happens with this (Calculating SVD on 2 segments ([[ 0 3561] [3561 3781]])). And then, maybe I did not test well enough writevarcurve with more than one segment.

[monodera]

Thanks a lot for the reply. I made a pull request with the commit.

The wavelength coverage of my data cube is 5600A to 10980A, extending outside of 10000A limit (is this the reason?).

As shown in the web site (https://www.naoj.org/Observing/Instruments/FOCAS/ifu/), the wavelength coverage depends on the choice of grism and filter configuration.

[saimn]

Ok, so yes this is the reason. You can modify the SKYSEG to check, and we should set it to something higher by default. This is not ideal but as I explain it is because historically zap was using these segments (this is in the paper), and when I switched to one segment I wanted to preserve this possibility, just in case someone wants to play with it.

[saimn]

I removed the hard-coded limits, without any additional change since the code can actually deal with it directly and will use the wavelength limits from the data. Please check that it works for you, thanks.

[monodera]

Thanks for the modification. ZAP runs successfully to produce the explained variance curve. Because of the different segmentation, there are some differences about an order of <~10%. Is this expected?

The varcurve looks as attached (x is the row index and the y is numbers stored in the table). Does this look reasonable?

[saimn]

Yes I would expect some differences but mostly for the second segment ([3561 3781]) since it was small and in my experience from MUSE it was giving much better results with only one segment. The varcurve looks reasonable, I guess zap chose to use something like 4 or 5 eigenvalues ? That's not a lot but it depends on your type of data, it's probably enough to explain the sky components in your case.