Upgrade to wcslib 5.3 - #3751

So I take it this still includes the version number improvements?

@embray: Yes -- if you mean that wcslib now has a function to get the runtime version of the library, yes. It also supports a couple of new projection types.

@mdboom - this seems fine to me as a temporary measure. Do you have a sense of how long it will be for compatibility with the SIP stuff? (I notice it also include TPV support, which is a useful addition as well!)

@eteq: I haven't really done much investigation into the effort involved. The tricky bit will be maintaining the current API which allows for running with or without SIP, only SIP etc. (since wcslib, fairly reasonably, always does everything).

@stargaser - I noticed in #3669 that you said you might be able to help with testing - do you plan to do that on this version?

Also, @olebole, does the changelog look good to you here? (you asked about it in #3669)

Also, @mdboom, I just tested this on Dark Energy Survey data from the NOAO science archive (which is in the TPV projection), and it appears to load correctly. So that's a very important use case that now works!

@eteq changelog is fine with me.
@mdboom For Debian, I would upgrade to wcslib-5.5 once astropy works with it; so for me there is in thefirst order no need to support both. But: I am planning to do an astropy backport to Debian Jessie (this currently has still astropy version 0.4.2). Jessie has only wcslib-4.24, and backporting wcslib-5.5 is complicated since its ABI changed (and so this would break existing programs). So, the perfect world for me would be to support both, at least for one version or so.

@olebole: I have no plans to drop support for wcslib 4.24.

@eteq Yes I am planning to help with testing on this version. I would like to test on Spitzer, WISE and PTF images -- the latter include both TPV and SIP.

@mdboom I've done some testing on this branch on an image from the Palomar Transient Factory. The header contains distortion keywords for both SIP and TPV conventions, with CTYPEn ending in 'TAN-SIP'.

When the header is set up for TPV, I get very different answers from all_pix2world from SIP. DS9 gives me the same results for the TPV form and SIP form.

Further, for the TPV WCS, I get the same results for all_pix2world and wcs_pix2world. This makes sense when I think about it, since wcslib is handling the whole transformation in both cases.

@stargaser: I'm not surprised at your results. With this branch, both astropy and wcslib are performing the SIP transformation. There is ongoing work to make that happen, but there turn out to be a great deal of changes required, and wcslib is not yet doing the SIP transformation correctly in all cases.

Did you get a warning emitted for the SIP header? If not, I'd consider that a bug -- that's all this PR tries to do at this point, since it knowingly doesn't do the right thing.

@mdboom I didn't get a warning for the SIP header.

I made a version of the FITS file with all the SIP keywords removed, and with CTYPEn set to end with 'TPV', and I got the same results as the 'TPV' case with SIP keywords still present.

@mdboom Let me make sure I've got my branch properly synced. It's wcs/upgrade-5.3 (like at the top of this pull request), right?

Yes -- that is the correct branch.

@mdboom I've completely cleaned out and rebuilt this branch. I'm getting the same results.

For simplicity, I'll just show the case of SIP-keywords only (no PVi_j), and TPV-keywords only.

In [14]: w_siponly.all_pix2world(np.array([[1935.0,160.0]]),1)
Out[14]: array([[ 276.62615374,  -24.78997499]])

In [13]: w_tpvonly.all_pix2world(np.array([[1935.0,160.0]]),1)
Out[13]: array([[ 274.80694571,  -25.9746477 ]])

In [16]: w_siponly.wcs_pix2world(np.array([[1935.0,160.0]]),1)
Out[16]: array([[ 276.65616811,  -24.81478424]])

In [17]: w_tpvonly.wcs_pix2world(np.array([[1935.0,160.0]]),1)
Out[17]: array([[ 274.80694571,  -25.9746477 ]])

The TPV cases are over a degree away from the SIP results. In DS9, I get the same coordinates for both TPV and SIP.

Would it work to pickle the headers and send them to you, somehow? (please specify the how!) The data themselves are not public, but the pixel data is not needed for this test.

I see now what you're doing -- and the results are surprising. You can just load your FITS files in astropy.io.fits, set all the data to zeros and write it back out, and post as a gist or something somewhere.

@mdboom Here are links to FITS files:

SIP-only: https://caltech.box.com/shared/static/7czst4o0whjhc39dxrwqh93gbpw4kik0.fits

TPV-only: https://caltech.box.com/shared/static/5zfm3633qs17r6xjpwoxuwjxg13yajls.fits

(sorry, Box won't let me rename the direct links...)

Thanks. I've looked at this a bit -- I'll let you know if I find anything.

@mdboom - I just tested this now with a DECam image (basically anything you download at http://www.portal-nvo.noao.edu/search/query called "CTIO 4m + DECam"), which is TPV, and I'm getting something weird: when I do pix2world on any of the individual CCD HDUs, I always get the same coordinates (the center of the field). And world2pix always gives infs and nans. Any idea what might be going on there?

@eteq: Can you provide a self-contained example? There's lots of ways to shoot one's self in the foot along the way with what you're describing, and I just want to make sure it's not any of the obvious ones.

Superceded by #3845.

@stargaser: The new branch against wcslib 5.5 in #3845 is better (and at least is consistent in the wcs_pix2world case), but there is still a difference between SIP and TPV. I'm going to forward your question along to Mark Calabretta to see if he has any insights.

@mdboom Just to tie up the loose ends, I've tested the files I linked above with the current astropy master e15cbb6 and now it works as I would expect.

In [22]: w_siponly.all_pix2world(np.array([[1935.0,160.0]]),1)
Out[22]: array([[ 276.62615374,  -24.78997499]])

In [23]: w_tpvonly.all_pix2world(np.array([[1935.0,160.0]]),1)
Out[23]: array([[ 276.62615374,  -24.78997499]])

In [24]: w_siponly.wcs_pix2world(np.array([[1935.0,160.0]]),1)
Out[24]: array([[ 276.65616811,  -24.81478424]])

In [25]: w_tpvonly.wcs_pix2world(np.array([[1935.0,160.0]]),1)
Out[25]: array([[ 276.62615374,  -24.78997499]])

Thanks for the report! Glad to hear it.