[AstroPy] Running wcs_pix2world for all image pixels

[Michael Droettboom]

On 10/18/2013 04:54 AM, David Berry wrote:
>
>
>
> On 17 October 2013 16:28, Maik Riechert <maik.riechert at arcor.de 
> <mailto:maik.riechert at arcor.de>> wrote:
>
>     I didn't use tranGrid here, just tran! So no tolerance involved.
>
>
> Doh! Sorry - didn't spot that.
>
>
>     NAXIS1  =                 1468 / image width
>     NAXIS2  =                  916 / image height
>     CTYPE1  = 'RA---TAN-SIP' / TAN (gnomic) projection + SIP distortions
>     CTYPE2  = 'DEC--TAN-SIP' / TAN (gnomic) projection + SIP distortions
>     EQUINOX =               2000.0 / Equatorial coordinates definition
>     (yr)
>     LONPOLE =                180.0 / no comment
>     LATPOLE =                  0.0 / no comment
>     CRVAL1  =         302.88771538 / RA  of reference point
>     CRVAL2  =        41.8941029213 / DEC of reference point
>     CRPIX1  =                  735 / X reference pixel
>     CRPIX2  =                  459 / Y reference pixel
>     CUNIT1  = 'deg     ' / X pixel scale units
>     CUNIT2  = 'deg     ' / Y pixel scale units
>     CD1_1   =     -0.0157184611387 / Transformation matrix
>     CD1_2   =     0.00022032143241 / no comment
>     CD2_1   =   -0.000284856284919 / no comment
>     CD2_2   =     -0.0159943489114 / no comment
>     IMAGEW  =                 1468 / Image width,  in pixels.
>     IMAGEH  =                  916 / Image height, in pixels.
>     A_ORDER =                    2 / Polynomial order, axis 1
>     A_0_2   =   -5.10167355829E-05 / no comment
>     A_1_1   =   -5.54267287337E-05 / no comment
>     A_2_0   =    3.68095680847E-05 / no comment
>     B_ORDER =                    2 / Polynomial order, axis 2
>     B_0_2   =    2.49853669312E-05 / no comment
>     B_1_1   =    6.07756179822E-05 / no comment
>     B_2_0   =   -1.70350503136E-05 / no comment
>     AP_ORDER=                    2 / Inv polynomial order, axis 1
>     AP_0_1  =     -0.0014471978634 / no comment
>     AP_0_2  =     5.0935282476E-05 / no comment
>     AP_1_0  =    0.000641040653934 / no comment
>     AP_1_1  =    5.54267224252E-05 / no comment
>     AP_2_0  =     -3.674588665E-05 / no comment
>     BP_ORDER=                    2 / Inv polynomial order, axis 2
>     BP_0_1  =      0.0010860421788 / no comment
>     BP_0_2  =   -2.49999638422E-05 / no comment
>     BP_1_0  =   -0.000569392740582 / no comment
>     BP_1_1  =   -6.07473585161E-05 / no comment
>     BP_2_0  =     1.7000827442E-05 / no comment
>
>
> OK - so I suspect the problems may actually be in astropy rather than 
> AST.  Maybe the maintainers of astropy.wcs could comment.
>
> The above header includes a polynomial distortion correction using the 
> Spitzer SIP conventions. The "A_" and "B_" headers give the 
> coefficients of the forward ( (x,y) to (ra,dec) ) polynomial, and the 
> "AP_" and "BP_" headers give the coefficients of the inverse ( 
> (ra,dec) to (x,y) ) polynomial.
>
> Using the header as supplied to transform the following pixel positions:
>
> x = -196, 735, 4059
> y = -504, 459, 2328
>
> I get, using AST (converted to degrees):
>
> ra = 328.090902, 302.887715, 265.839543
> dec = 53.570723, 41.894103, 8.257963
>
> whereas astropy gives:
>
> ra = 326.992376, 302.866900, 263.545072
> dec = 54.749009, 41.877821, 10.589047
>
> For an independent check I also used SAO WCSTools which gave:
>
> ra =  328.09090, 302.88772, 265.83954
> dec = 53.57072, 41.89410, 8.25796
>
> So AST and WCSTools agree, but astropy is different. This is not 
> conclusive in itself, but I then modified the header by removing all 
> the SIP keywords and removing "-SIP" from the CTYPE1/2 headers. With 
> these changes the header represents an undistorted TAN projection. I 
> then transformed the same pixel positions again using the changed header:
>
> AST gave:
>
> ra = 327.022569, 302.887715, 263.548641
> dec = 54.758313, 41.894103, 10.599558
>
> astropy gave:
>
> ra = 326.992376, 302.866900, 263.545072
> dec = 54.749009, 41.877822, 10.589047
>
> and WCSTools gave:
>
> ra = 327.02257, 302.88772, 263.54864
> dec = 54.75831, 41.89410, 10.59956
>
> So again AST and WCSTools agree, with astropy as the odd one out. But 
> more importantly, the astropy results are exactly the same as when 
> using the header that includes the SIP distortion. This seems to 
> indicate that in fact astropy is ignoring the SIP distortion. You can 
> see that the astropy values are much closer to the AST/WCSTools values 
> when no distortion is included.

You need to use wcs.all_pix2world to include the SIP distortion. SIP is 
not a part of the WCS standard, so we only perform it when explicitly 
asked for.  See:

https://astropy.readthedocs.org/en/v0.2.4/_generated/astropy.wcs.wcs.WCS.html#astropy.wcs.wcs.WCS.wcs_pix2world

Mike

-- 
                    _
|\/|o _|_  _. _ | | \.__  __|__|_|_  _  _ ._ _
|  ||(_| |(_|(/_| |_/|(_)(/_|_ |_|_)(_)(_)| | |

http://www.droettboom.com

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/astropy/attachments/20131018/42ccd517/attachment.html>