- NAME
WCSSPH2XY
- PURPOSE
Convert spherical coordinates to x and y (map) angular coordinates
- EXPLANATION
Convert spherical (longitude and latitude -- sky) coordinates to x
and y intermediate world coordinates (still nominally in degrees) in
the projection plane of the map. This procedure is the inverse of
WCSXY2SPH. See WCS_DEMO for example of use.
This is a lower level procedure -- given a FITS header, the user will
usually use ADXY which will then call WCSSPH2XY with the appropriate
parameters.
- CATEGORY
Mapping and Auxiliary FITS Routine
- CALLING SEQUENCE
wcssph2xy, longitude, latitude, x, y, [ map_type , CTYPE = ,
FACE =, PV1= PV2= , CRVAL = , CRXY = , LONGPOLE = ,
LATPOLE = , PHI0 = , NORTH_OFFSET =, SOUTH_OFFSET =, BADINDEX =]
- INPUT PARAMETERS
longitude - longitude of data, scalar or vector, in degrees
latitude - latitude of data, same number of elements as longitude,
in degrees
map_type - optional positional parameter, numeric scalar (0-26)
corresponding to a particular map projection. This is not a
FITS standard, it is simply put in to allow function similar
to that of less general map projection procedures (eg AITOFF).
The following list gives the map projection types and their
respective numbers.
S Number Name Comments
e code
- ------ ----------------------- -----------------------------------
F 0 Default = Plate Carree
P 1 Zenithal perspective PV2_1 required
N 2 Gnomic AZP w/ mu = 0
N 3 Orthographic PV2_1,PV2_2 optional
G 4 Stereographic AZP w/ mu = 1
C 5 Zenithal Equidistant
N 6 Zenithal polynomial PV2_0, PV2_1....PV2_20 possible
A 7 Zenithal equal area
R 8 Airy PV2_1 required
P 9 Cylindrical perspective PV2_1 and PV2_2 required
R 10 Plate Carree
R 11 Mercator
A 12 Cylindrical equal area PV2_1 required
P 13 Conical perspective PV2_1 and PV2_2 required
D 14 Conical equidistant PV2_1 and PV2_2 required
E 15 Conical equal area PV2_1 and PV2_2 required
O 16 Conical orthomorphic PV2_1 and PV2_2 required
N 17 Bonne's equal area PV2_1 required
O 18 Polyconic
L 19 Sanson-Flamsteed (GLS is allowed as a synonym for SFL)
R 20 Parabolic
T 21 Hammer-Aitoff
L 22 Mollweide
C 23 Cobe Quadrilateralized convergence of inverse is poor
Spherical Cube
C 24 Quadrilateralized
Spherical Cube
C 25 Tangential Spherical Cube
P 26 Slant Zenithal Projection PV2_1,PV2_2, PV2_3 optional
X 27 HealPix
T 28 HealCart (Cartesian approximation of Healpix)
H 29 HEALPix butterfly projection
- OPTIONAL INPUT KEYWORD PARAMETERS
CTYPE - One, two, or three element vector containing 8 character
strings corresponding to the CTYPE1, CTYPE2, and CTYPE3
FITS keywords:
CTYPE[0] - first four characters specify standard system
('RA--','GLON' or 'ELON' for right ascension, Galactic
longitude or ecliptic longitude respectively), second four
letters specify the type of map projection (eg '-AIT' for
Aitoff projection)
CTYPE[1] - first four characters specify standard system
('DEC-','GLAT' or 'ELAT' for declination, galactic latitude
or ecliptic latitude respectively; these must match
the appropriate system of ctype1), second four letters of
ctype2 must match second four letters of ctype1.
CTYPE[2] - if present must be the 8 character string,'CUBEFACE',
only used for spherical cube projections to identify an axis
as containing the face on which each x and y pair of
coordinates lie.
PV2 - Vector of projection parameter associated with latitude axis
PV2 will have up to 21 elements for the ZPN projection, up to 3
for the SIN projection and no more than 2 for any other
projection. The first element corresponds to PV2_1, the
second to PV2_2, etc.
CRXY - 2 element vector giving the x and y coordinates of the
reference point. if this is not set the offset is [0,0].
Used to implement (x0,y0) in Sect 2.5 of Griesen & Calabretta 2002
Do not confuse with CRPIX.
arameters simply passed to WCS_ROTATE:
CRVAL - 2 element vector containing standard system coordinates (the
longitude and latitude) of the reference point
PV1 - Vector of projection parameters associated with longitude
LONGPOLE - native longitude of standard system's North Pole
LATPOLE - "target" native latitude of the standard system's North Pole
arameters intended to enhance invertability:
NORTH_OFFSET - offset (radians) added to input points near north pole.
SOUTH_OFFSET - offset (radians) added to input points near south pole.
- OUTPUT PARAMETERS
x - x coordinate of data, same number of elements as longitude, in
degrees; if CRXY is set, then x will be returned offset by
crxy(0). NOTE: x in all map projections increases to the
left, not the right.
y - y coordinate of data, same number of elements as longitude, in
degrees; if CRXY is set, y will be returned offset by crxy[1]
- OPTIONAL OUTPUT KEYWORD PARAMETERS
FACE - a output variable used for spherical cube projections to
designate the face of the cube on which the x and y
coordinates lie. Will contain the same number of elements as
X and Y. Must contain at least 1 arbitrary element on input
If FACE is NOT defined on input, it is assumed that the
spherical cube projection is laid out over the whole sky
in the "sideways T" configuration.
DINDEX - vector, list of transformed points too close to poles.
- NOTES
The conventions followed here are described in more detail in
"Representations of Celestial Coordinates in FITS" by Calabretta
and Greisen (2002, A&A, 395, 1077; also see
http://fits.gsfc.nasa.gov/fits_wcs.html). The general
scheme outlined in that article is to first use WCS_ROTATE to convert
coordinates in one of three standard systems (celestial, galactic,
or ecliptic) into a "native system" of latitude and longitude. The
latitude and longitude are then converted into x and y coordinates
which depend on the map projection which is performed. The rotation
from standard to native coordinates can be skipped if one so desires.
This procedure necessitates two basic sections. The first converts
"standard" coordinates to "native" coordinates while the second converts
"native" coordinates to x and y coordinates. The first section is
simply a call to WCS_ROTATE, while the second contains the guts of
the code in which all of the map projection is done. This procedure
can be called in a form similar to AITOFF, EQPOLE, or QDCB by calling
wcssph2xy with a fifth parameter specifying the map projection by
number and by not using any of the keywords related to the map
projection type (e.g. CTYPE).
- PROCEDURE
The first task of the procedure is to do general error-checking to
make sure the procedure was called correctly and none of the
parameters or keywords conflict. This is particularly important
because the procedure can be called in two ways (either using
FITS-type keywords or using a number corresponding to a map projection
type). All variables are converted into double precision values and
angular measurements are converted from degrees into radians.
If necessary, longitude values are converted into the range -pi to pi.
Any latitude points close to the of the poles are mapped to a specific
latitude of from the pole so that the map transformations become
completely invertible. The magnitude of this correction is given by
the keywords NORTH_OFFSET and SOUTH_OFFSET and a list of affected
points is optionally returned in the "badindex" output parameter.
The next task of the procedure is to convert the "standard"
coordinates to "native" coordinates by rotating the coordinate system.
This rotation is performed by the procedure WCS_ROTATE and is governed
by the keywords CRVAL and LONGPOLE. The final task of the WCSSPH2XY
is to take "native" latitude and longitude coordinates and convert
them into x and y coordinates. Any map specific error-checking is
done at this time. All of the equations were obtained from
"Representations of Celestial Coordinates in FITS" and cases needing
special attention are handled appropriately (see the comments with
individual map projections for more information on special cases).
Note that a further transformation (using the CD matrix) is required
to convert the (x,y) coordinates to pixel coordinates.
- COMMON BLOCKS
none
- PROCEDURES CALLED
WCS_ROTATE
- AUTHOR
Rick Balsano
- MODIFICATIONS/REVISION LEVEL
1.1 8/31/93
2.3 9/15/93 W. Landsman (HSTX) Update quad cube coords, vectorize
keywords
2.4 12/29/93 I. Freedman (HSTX) Eliminated LU decomposition
2.5 1/5/93 I. Freedman (HSTX) Offset keywords / bad point index
2.6 Dec 94 Compute pole for transformations where the reference
pixel is at the native origin W. Landsman (HSTX)
2.7 May 95 Change internal variable BETA for V4.0 compatibility
2.8 June 95 Change loop indices from integer to long
2.9 3/18/96 Change FACE usage for cube projections to match WCSLIB
C/FORTRAN software library.
2.10 02/18/99 Fixed implementation of ARC algorithm
2.11 June 2003 Update conic projections, add LATPOLE keyword
2.12 Aug 2003, N.Rich - Fix pre-V5.5 bug from previous update
2.13 Sep 2003, W. Landsman CTYPE keywords need not be 8 characters
2.14 Jan 2004, W. Landsman don't modify scalars, fix PARabolic code
2.15 Feb 2004, W. Landsman Fix AZP and AIR algorithms
3.0 May 2004 W. Landsman Support extended SIN (=NCP), slant zenithal
(SZP), and zenithal polynomial (ZPN) projections, use
PV2 keyword vector instead of PROJP1, PROJP2
3.1 Jul 2005 W.Landsman/C. Markwardt Set unprojectable points in
tangent projection to NaN
3.1.1 Jul 2005 Fixed 3.1 mod to work for scalars
3.2 Dec 2005 Fixed Airy projection for latitude centered at 90 deg
3.3 Aug 2007 R. Munoz, W.Landsman Correct treatment of PV1_2 and
PV2_2 parameters
3.4 Oct 2007 Sergey Koposov Support HEALPIX projection
3.4.1 June 2009 Check for range of validity of ZPN polynomial W.L.
3.5 May 2012 Benjamin Alan Weaver, Add nonstandard HEALCART
projection, Allow map_index to be > 25
3.5.1 May 2013 W. Landsman Allow GLS as a synonym for SFL
3.6 Jul 2013 J. P. Leahy added XPH projection, apply polar offsets
only for cylindrical & conic projections.
3.6.1 Dec 2013 W. Landsman Polar offsets done in radians
3.6.2 Jan 2016 W. Landsman Lat and Long can have different size so long
as they have the same number of elements