2 Basic concepts

   
2.1 Spherical projection

As indicated in Fig. 1, the first step in transforming (x,y)coordinates in the plane of projection to celestial coordinates is to convert them to native longitude and latitude, $(\phi,\theta)$. The equations for the transformation depend on the particular projection and this will be specified via the CTYPE ia keyword. Paper I defined "4-3'' form for such purposes; the rightmost three-characters are used as an algorithm code that in this paper will specify the projection. For example, the stereographic projection will be coded as STG. Some projections require additional parameters that will be specified by the FITS keywords PV i_ma for $m = 0,1,2,\ldots$, also introduced in Paper I. These parameters may be associated with the longitude and/or latitude coordinate as specified for each projection. However, definition of the three-letter codes for the projections and the equations, their inverses and the parameters which define them, form a large part of this work and will be discussed in Sect. 5. The leftmost four characters of CTYPE ia are used to identify the celestial coordinate system and will be discussed in Sect. 3.

   
2.2 Reference point of the projection

The last step in the chain of transformations shown in Fig. 1 is the spherical rotation from native coordinates, $(\phi,\theta)$, to celestial coordinates $(\alpha ,\delta )$. Since a spherical rotation is completely specified by three Euler angles it remains only to define them.

In principle, specifying the celestial coordinates of any particular native coordinate pair provides two of the Euler angles (either directly or indirectly). In the AIPS convention, the CRVAL ia keyword values for the celestial axes specify the celestial coordinates of the reference point and this in turn is associated with a particular point on the projection. For zenithal projections that point is the sphere's point of tangency to the plane of projection and this is the pole of the native coordinate system. Thus the AIPS convention links a celestial coordinate pair to a native coordinate pair via the reference point. Note that this association via the reference point is purely conventional; it has benefits which are discussed in Sect. 5 but in principle any other point could have been chosen.

Section 5 presents the projection equations for the transformation of (x,y) to $(\phi,\theta)$. The native coordinates of the reference point would therefore be those obtained for (x,y) = (0,0). However, it may happen that this point lies outside the boundary of the projection, for example as for the ZPN projection of Sect. 5.1.7. Therefore, while this work follows the AIPS approach it must of necessity generalize it.

 

 

Table 1: Summary of important variable names and other symbols used throughout the paper.

Variable(s)	Meaning	Related FITS keywords (if any)
i	Index variable for world coordinates
j	Index variable for pixel coordinates
a	Alternate version code, blank or A to Z
pj	Pixel coordinates
rj	Reference pixel coordinates	CRPIX ja
mij	Linear transformation matrix	CD i_ja or PC i_ja
si	Coordinate scales	CDELT ia
xi	Intermediate world coordinates (in general)
(x,y)	Projection plane coordinates
$(\phi,\theta)$	Native longitude and latitude
$(\alpha ,\delta )$	Celestial longitude and latitude
$(\phi _0,\theta _0)$	Native longitude and latitude of the fiducial point	PV i_1a$^\dag$, PV i_2a$^\dag$
$(\alpha _0,\delta _0)$	Celestial longitude and latitude of the fiducial point	CRVAL ia
$(\phi_{\rm p},\theta_{\rm p})$	Native longitude and latitude of the celestial pole	LONPOLE a ( = PV i_3a$^\dag$), LATPOLE a, ( = PV i_4a$^\dag$)
$(\alpha_{\rm p},\delta_{\rm p})$	Celestial longitude and latitude of the native pole $(\delta_{\rm p} = \theta_{\rm p})$
$\arg()$	Inverse tangent function that returns the correct quadrant

$^\dag$ Associated with longitude axis i.

Accordingly we specify only that a fiducial celestial coordinate pair $(\alpha _0,\delta _0)$ given by the CRVAL ia will be associated with a fiducial native coordinate pair $(\phi _0,\theta _0)$ defined explicitly for each projection. For example, zenithal projections all have $(\phi _0,\theta _0) = (0,90\hbox {$^\circ $ })$, while cylindricals have $(\phi _0,\theta _0) = (0,0)$. The AIPS convention has been honored here as far as practicable by constructing the projection equations so that $(\phi _0,\theta _0)$ transforms to the reference point, (x,y) = (0,0). Thus, apart from the one exception noted, the fiducial celestial and native coordinates are the celestial and native coordinates of the reference point and we will not normally draw a distinction.

It is important to understand why $(\phi _0,\theta _0)$ differs for different projection types. There are two main reasons; the first makes it difficult for it to be the same, while the second makes it desirable that it differs. Of the former, some projections such as Mercator's, diverge at the native pole, therefore they cannot have the reference point there because that would imply infinite values for CRPIX ja. On the other hand, the gnomonic projection diverges at the equator so it can't have the reference point there for the same reason. Possibly $(\phi _0,\theta _0)$ chosen at some mid-latitude could satisfy all projections, but that leads us to the second reason.

Different projection types are best suited to different purposes. For example, zenithal projections are best for mapping the region in the vicinity of a point, often a pole; cylindrical projections are appropriate for the neighborhood of a great circle, usually an equator; and the conics are suitable for small circles such as parallels of latitude. Thus, it would be awkward if a cylindrical used to map, say, a few degrees on either side of the galactic plane, had its reference point, and thus CRPIX ja and CRVAL ia, at the native pole, way outside the map boundary. In formulating the projection equations themselves the native coordinate system is chosen to simplify the geometry as much as possible. For the zenithals the natural formulation has (x,y) = (0,0) at the native pole, whereas for the cylindricals the equations are simplest if (x,y) = (0,0) at a point on the equator.

As discussed above, a third Euler angle must be specified and this will be given by the native longitude of the celestial pole, $\phi_{\rm p}$, specified by the new FITS keyword

LONPOLE a (floating-valued).

The default value of LONPOLE a will be 0 for $\delta_0\ge\theta_0$ or $180\hbox{$^\circ$ }$ for $\delta_0<\theta_0$. This is the condition for the celestial latitude to increase in the same direction as the native latitude at the reference point. Thus, for example, in zenithal projections the default is always $180\hbox{$^\circ$ }$ (unless $\delta_0 = 90\hbox{$^\circ$ }$) since $\theta_0 = 90\hbox{$^\circ$ }$. In cylindrical projections, where $\theta_0 = 0$, the default value for LONPOLE a is 0 for $\delta_0 \ge 0$, but it is $180\hbox{$^\circ$ }$ for $\delta_0 < 0$.

2.3 Spherical coordinate rotation

Since $(\phi _0,\theta _0)$ differs for different projections it is apparent that the relationship between $(\alpha _0,\delta _0)$ and the required Euler angles also differs.

For zenithal projections, $(\phi _0,\theta _0) = (0,90\hbox {$^\circ $ })$ so the CRVAL ia specify the celestial coordinates of the native pole, i.e. $(\alpha_0,\delta_0) = (\alpha_{\rm p},\delta_{\rm p})$. There is a simple relationship between the Euler angles for consecutive rotations about the Z-, X-, and Z-axes and $\alpha_{\rm p}$, $\delta_{\rm p}$ and $\phi_{\rm p}$; the ZXZ Euler angles are $(\alpha_{\rm p}+90\hbox{$^\circ$ }, 90\hbox{$^\circ$ }-\delta_{\rm p}, \phi_{\rm p}-90\hbox{$^\circ$ })$. Given this close correspondence it is convenient to write the Euler angle transformation formulæ directly in terms of $\alpha_{\rm p}$, $\delta_{\rm p}$ and $\phi_{\rm p}$:

$$\begin{array}{rcl} \alpha & = & \alpha_{\rm p} + \arg~(\sin\... ...eta \cos\delta_{\rm p}\cos(\phi-\phi_{\rm p}) ) , \end{array}$$ (2)

where $\arg~()$ is an inverse tangent function that returns the correct quadrant, i.e. if $(x,y) = (r\cos\beta,r\sin\beta)$ with r > 0 then $\arg~(x,y) = \beta$. Note that, if $\delta_{\rm p} = 90\hbox{$^\circ$ }$, Eqs. (2) become

$$\begin{array}{rcl} \alpha & = & \alpha_{\rm p} + \phi - \phi... ... - 180\hbox{$^\circ$ }, \\ \delta & = & \theta , \end{array}$$ (3)

which may be used to define a simple change in the origin of longitude. Likewise for $\delta_{\rm p} = -90\hbox{$^\circ$ }$

$$\begin{array}{rcl} \alpha & = & \alpha_{\rm p} - \phi + \phi_{\rm p} , \\ \delta & = & -\theta . \end{array}$$ (4)

The inverse equations are

$$\begin{array}{rcl} \phi & = & \phi_{\rm p} + \arg~(\sin\delt... ... \cos\delta_{\rm p}\cos(\alpha-\alpha_{\rm p})) . \end{array}$$ (5)

Useful relations derived from Eqs. (2) and (5) are

 

$$\cos\delta \cos(\alpha - \alpha_{\rm p}) \sin\theta\cos\delta_{\rm p} - \cos\theta \sin\delta_{\rm p} \cos(\phi-\phi_{\rm p}) ,$$ =

$$\cos\delta \sin(\alpha - \alpha_{\rm p}) -\cos\theta \sin(\phi - \phi_{\rm p}) ,$$ = (6)

 

$$\cos\theta \cos(\phi - \phi_{\rm p}) \sin\delta \cos\delta_{\rm p} -\cos\delta \sin\delta _{\rm p} \cos(\alpha - \alpha_{\rm p}) ,$$ =

$$\cos\theta \sin(\phi - \phi_{\rm p}) -\cos\delta \sin(\alpha - \alpha_{\rm p}) .$$ = (7)

A matrix method of handling the spherical coordinate rotation is described in Appendix B.

   
2.4 Non-polar $(\phi _0,\theta _0)$

Projections such as the cylindricals and conics for which $(\phi_0,\theta_0) \neq (0,90\hbox{$^\circ$ })$ are handled by providing formulae to compute $(\alpha_{\rm p},\delta_{\rm p})$ from $(\alpha _0,\delta _0)$ whence the above equations may be used.

Given that $(\alpha _0,\delta _0)$ are the celestial coordinates of the point with native coordinates $(\phi _0,\theta _0)$, Eqs. (6) and (7) may be inverted to obtain

 

$$\delta_{\rm p} \arg~(\cos\theta_0 \cos(\phi_{\rm p}-\phi_0) , \sin\theta_0) \;$$

$$\pm \cos^{-1} \left( \frac{\sin\delta_0} {\sqrt{1 - \cos^2\theta_0\sin^2(\phi_{\rm p}-\phi_0)}} \right) ,$$ (8)

  

$$\sin(\alpha_0-\alpha_{\rm p}) \sin(\phi_{\rm p}-\phi_0)\cos\theta_0 / \cos\delta_0 ,$$ (9)

$$\cos(\alpha_0-\alpha_{\rm p}) \frac{\sin\theta_0 - \sin\delta_{\rm p}\sin\delta_0} {\cos\delta_{\rm p}\cos\delta_0} ,$$ (10)

whence Eqs. (2) may be used to determine the celestial coordinates. Note that Eq. (8) contains an ambiguity in the sign of the inverse cosine and that all three indicate that some combinations of $\phi _0$, $\theta _0$, $\alpha_0$, $\delta_0$, and $\phi_{\rm p}$ are not allowed. For these projections, we must therefore adopt additional conventions:

1.

Equations (9) and (10) indicate that $\alpha_{\rm p}$ is undefined when $\delta_0 = \pm90\hbox{$^\circ$ }$. This simply represents the longitude singularity at the pole and forces us to define $\alpha_{\rm p} = \alpha_0$ in this case.

2.

If $\delta_{\rm p} = \pm90\hbox{$^\circ$ }$ then the longitude at the native pole is $\alpha_{\rm p} = \alpha_0 + \phi_{\rm p} - \phi_0 - 180\hbox{$^\circ$ }$ for $\delta_{\rm p} = 90\hbox{$^\circ$ }$ and $\alpha_{\rm p} = \alpha_0 - \phi_{\rm p} + \phi_0$ for $\delta_{\rm p} = -90\hbox{$^\circ$ }$.

3.

Some combinations of $\phi _0$, $\theta _0$, $\delta_0$, and $\phi_{\rm p}$ produce an invalid argument for the $\cos^{-1}()$ in Eq. (8). This is indicative of an inconsistency for which there is no solution for $\delta_{\rm p}$. Otherwise Eq. (8) produces two solutions for $\delta_{\rm p}$. Valid solutions are ones that lie in the range $-90\hbox{$^\circ$ }$ to $+90\hbox{$^\circ$ }$, and it is possible in some cases that neither solution is valid.

Note, however, that if $\phi_0 = 0$, as is usual, then when LONPOLE a ( $\equiv \phi_{\rm p}$) takes its default value of 0 or $180\hbox{$^\circ$ }$ (depending on $\theta _0$) then any combination of $\delta_0$ and $\theta _0$ produces a valid argument to the $\cos^{-1}()$ in Eq. (8), and at least one of the solutions is valid.

4.

Where Eq. (8) has two valid solutions the one closest to the value of the new FITS keyword

LATPOLE a (floating-valued)

is chosen. It is acceptable to set LATPOLE a to a number greater than $+90\hbox{$^\circ$ }$ to choose the northerly solution (the default if LATPOLE a is omitted), or a number less than $-90\hbox{$^\circ$ }$ to select the southern solution.

5.

Equation (8) often only has one valid solution (because the other lies outside the range $-90\hbox{$^\circ$ }$ to $+90\hbox{$^\circ$ }$). In this case LATPOLE a is ignored.

6.

For the special case where $\theta_0 = 0$, $\delta_0 = 0$, and $\phi_{\rm p} - \phi_0 =\pm90\hbox{$^\circ$ }$ then $\delta_{\rm p}$ is not determined and LATPOLE a specifies it completely. LATPOLE a has no default value in this case.

These rules governing the application of Eqs. (8-10) are certainly the most complex of this formalism. FITS writers are well advised to check the values of $\phi _0$, $\theta _0$, $\alpha_0$, $\delta_0$, and $\phi_{\rm p}$ against them to ensure their validity.

   
User-specified $(\phi _0,\theta _0)$

In Sect. 2.2 we formally decoupled $(\alpha _0,\delta _0)$ from the reference point and associated it with $(\phi _0,\theta _0)$. One implication of this is that it should be possible to allow $(\phi _0,\theta _0)$to be user-specifiable. This may be useful in some circumstances, mainly to allow CRVAL ia to match a point of interest rather than some predefined point which may lie well outside the image and be of no particular interest. We therefore reserve keywords PV i_1a and PV i_2a attached to longitude coordinate i to specify $\phi _0$ and $\theta _0$ respectively.

By itself, this prescription discards the AIPS convention and lacks utility because it breaks the connection between CRVAL ia and any point whose pixel coordinates are given in the FITS header. New keywords could be invented to define these pixel coordinates but this would introduce additional complexity and still not satisfy the AIPS convention. The solution adopted here is to provide an option to force (x,y) = (0,0) at $(\phi _0,\theta _0)$ by introducing an implied offset (x₀,y₀) which is computed for $(\phi _0,\theta _0)$ from the relevant projection equations given in Sect. 5. This is to be applied to the (x,y) coordinates when converting to or from pixel coordinates. The operation is controlled by the value of PV i_0a attached to longitude coordinate i; the offset is to be applied only when this differs from its default value of zero.

This construct should be considered advanced usage, of which Figs. 33 and 34 provide an example. Normally we expect that PV i_1a and PV i_2a will either be omitted or set to the projection-specific defaults given in Sect. 5.

   
2.6 Encapsulation

So that all required transformation parameters can be contained completely within the recognized world coordinate system (WCS) header cards, the values of LONPOLE a and LATPOLE a may be recorded as PV i_3a and PV i_4a, attached to longitude coordinate i, and these take precedence where a conflict arises.

We recommend that FITS writers include the PV i_1a, PV i_2a, PV i_3a, and PV i_4a cards in the header, even if only to denote the correct use of the default values.

Note carefully that these are associated with the longitude coordinate, whereas the projection parameters defined later are all associated with the latitude coordinate.

Figure 2: A linear equatorial coordinate system (top) defined via the methods of Wells et al. (1981), and the corresponding oblique system constructed using the methods of this paper. The reference coordinate $(\alpha _0,\delta _0)$ for each is at right ascension $8^{\rm hr}$, declination $+60\hbox {$^\circ $ }$ (marked). The two sets of FITS header cards differ only in their CTYPE ia keyword values. The non-oblique graticule could be obtained in the current system by setting $(\alpha _0,\delta _0,\phi _{\rm p})=(120\hbox {$^\circ $ },0,0)$.

2.7 Change of coordinate system

A change of coordinate system may be effected in a straightforward way if the transformation from the original system, $(\alpha ,\delta )$, to the new system, $(\alpha',\delta')$, and its inverse are known. The new coordinates of the $(\phi _0,\theta _0)$, namely $(\alpha'_0,\delta'_0)$, are obtained simply by transforming $(\alpha _0,\delta _0)$. To obtain $\phi'_{\rm p}$, first transform the coordinates of the pole of the new system to the original celestial system and then transform the result to native coordinates via Eq. (5) to obtain $(\phi'_{\rm p},\theta'_{\rm p})$. As a check, compute $\delta'_{\rm p}$ via Eq. (8) and verify that $\theta'_{\rm p} = \delta'_{\rm p}$.

   
2.8 Comparison with linear coordinate systems

It must be stressed that the coordinate transformation described here differs from the linear transformation defined by Wells et al. (1981) even for some simple projections where at first glance they may appear to be the same. Consider the plate carrée projection defined in Sect. 5.2.3 with $\phi = x$, $\theta = y$ and illustrated in Fig. 18. Figure 1 shows that while the transformation from (x,y) to $(\phi,\theta)$ may be linear (in fact identical), there still remains the non-linear transformation from $(\phi,\theta)$ to $(\alpha ,\delta )$. Hence the linear coordinate description defined by the unqualified CTYPE ia pair of RA, DEC which uses the Wells et al. prescription will generally differ from that of RA--CAR, DEC--CAR with the same CRVAL ia, etc. If LONPOLE a assumes its default value then they will agree to first order at points near the reference point but gradually diverge at points away from it.

Figure 2 illustrates this point for a plate carrée projection with reference coordinates of $8^{\rm hr}$ right ascension and $+60\hbox {$^\circ $ }$ declination and with $\phi_{\rm p} = 0$. It is evident that since the plate carrée has $(\phi _0,\theta _0) = (0,0)$, a non-oblique graticule may only be obtained by setting $\delta_0 = 0$ with $\phi_{\rm p} = 0$. It should also be noted that where a larger map is to be composed of tiled submaps the coordinate description of a submap should only differ in the value of its reference pixel coordinate.