[AstroPy] Inverting SIP distortions

[Leo Singer]

On Nov 1, 2013, at 1:00 PM, astropy-request at scipy.org wrote:

> Hello,
> I'm working on a problem that requires the inversion of SIP
> corrections: that is, I have points in world coordinates that I want
> to project onto an image plane.
> 
> My FITS header contains both the (A_p_q, B_p_q) coefficients for the
> forward transformation, as well as the (AP_p_q, BP_p_q) required for
> the backward transformation, according to the Shupe et al 2005
> reference.
> 
> The astropy WCS object has methods sip_foc2pix() and sip_pix2foc(),
> which appeared to be what I needed, but these don't seem to actually
> invert each other. That is:
> 
>>>> from astropy.wcs import WCS
>>>> wcs = WCS(header)
>>>> Xpix = [[1, 1]]
>>>> Xfoc = wcs.sip_pix2foc(Xpix, 0)
>>>> print wcs.sip_foc2pix(Xfoc, 0)
> [[ 1.2522118   1.74039233]]
> 
> I would have expected this to recover the input pixels, but it's
> clearly not the case.  Is there any functionality in astropy that will
> utilize the available AP/BP coefficients in my header to perform the
> inverse SIP distortion?  Is there anything I'm doing wrong or missing
> here?
> Thanks,
>   Jake

Hi Jake,

One other option to be aware of: the development version of Astropy (what will be version 0.3) contains a new method WCS.all_world2pix that converts world to pixel coordinates. It uses an iterative scheme, so it's not exact but you can set the tolerance to be as accurate as you like.

See the documentation:
http://docs.astropy.org/en/latest/api/astropy.wcs.wcs.WCS.html#astropy.wcs.wcs.WCS.all_world2pix

See also the merged pull request:
https://github.com/astropy/astropy/pull/1281

Best,
Leo Singer
Graduate Student @ LIGO-Caltech