From: fred at manono.span (Fred Patt (301)286-4569) To: fitsbits-request at fits.CV.NRAO.EDU Subject: Comments on Draft WCS Standard Date: Thu, 18 Feb 93 11:12:32 EST
I was offered a chance by Richard White to review your draft World Coordinate System standard for FITS, due to my interest in map projections in general and the quad sphere in particular. My recent involvement has been with Earth science data sets but I was also involved with the COBE mission for several years and am familiar with their skymap implementation. I have advocated the adoption of the quad sphere by local Earth science projects and was responsible for the use of this scheme in the Advanced Very High Resolution Radiometer (AVHRR) Pathfinder Land prototype processing system here at Goddard. I provided details of the projection equations to the AVHRR group in the form of pseudocode, which has been implemented in C by the developers. The data system manager, Mary James of Goddard Code 902, has given me encouragement to publish my work on this project.
My point here is that an exact form of the projection was derived in the second of the two references (O'Neill and Laubscher, 1976) in terms of polar coordinates, and this is the form which was implemented by AVHRR. The reference only derived the projection in the forward direction, from sky to cube coordinates, but I found that the derivation of the reverse projection was very straightforward. I also found that the equations could be significantly simplified by direct use of cartesian coordinates, thereby eliminating the polar angles as intermediaries and, as a side benefit, reducing the number of trig function evaluations required.
The information regarding the COBE implementation by Immanuel Freedman is completely accurate. COBE has accepted the non-exact form of the transformation, but I suggest that the description of the exact form is at least as appropriate for your work. The combination of the overall complexity of the projection equations and the errors in the projection described for COBE would seem to discourage most potential users.
I have attached a set of summary notes which I provided for the AVHRR group and the derivation of the reverse projection and cartesian form equations for the exact transformation (in computer-language syntax; it would be reasonable for me to put these into a more appropriate notation if desired). I hope that this information may be of use to you.
The quadrilaterlized spherical cube, or quad sphere for short, is an equal-area mapping and binning scheme for data collected on a spherical surface (either Earth data or the celestial sphere). It was first proposed by Chan and O'Neill in 1975 for the Naval Environmental Prediction Research Facility (Reference 1).
There are two key elements to the quad sphere:
o The mapping consists of projecting the sphere onto the faces of an inscribed cube using a curvilinear projection which preserves area. The sphere is divided into six equal areas which correspond to the faces of the cube. The vertices of the cube correspond to the cartesian coordinates defined by |x|=|y|=|z| on the unit sphere. For an Earth projection, the cube is normally oriented with one face normal to the North Pole and one face centered on the Greenwich meridian (although any definition of pole and meridian could be used). The faces of the cube are divided into square bins, where the number of bins along each edge is a power of 2, selected to produce the desired bin size. Thus the number of bins on each face is 2**(2*N), where N is the binning level, and the total number of bins is 6*2**(2*N). For example, a level of 10 gives 1024x1024 bins on each face and 6291456 (6*2**20) total bins, which are 23.605 square arcminutes (1.99737E-6 steradians) in size.
o The bins are numbered serially, rather than being rastered as for an image. The bin numbers are determined as follows. The total number of bits required for the bin numbers at level N is 2*N+3, where the 3 MSBs are used for the face numbers and the remaining bits are used to number the bins within each face. The faces are numbered 0-5 with 0 being the North face, 1 through 4 being equatorial with 1 corresponding to Greenwich, and 5 being South. Thus at level 10, face 0 has bin numbers 0-1048577, face 1 has numbers 1048576-2097151, etc. Within each face the bins are numbered serially from one corner (the convention is to start at the "lower left") to the opposite corner, with the ordering such that each pair of bits corresponds to a level of bin resolution. This ordering in effect is a two-dimensional binary tree, which is referred to as the quad-tree The conversion between bin numbers and coordinates is straightforward. If 4-byte integers are used for the bin numbers the maximum practical bin level is 14, which uses 31 of the 32 bits and results in a bin size of 0.0922 square arcminutes (7.80223E-9 steradians). The advantages of this numbering approach are stated below.
In principal the mapping and numbering schemes are separable; the projection onto the cube could be used with another bin numbering scheme, and the numbering scheme itself could be used with any arrangement of bins which can be partitioned as a set of square arrays. Used together, they comprise a flexible and efficient system for storing map data.
The quad sphere projection does not produce singularities at the poles or elsewhere, as do some other equal-area mapping schemes. The distortion is moderate over the entire sphere, so that at no point are shapes distorted beyond recognition.
Individual faces can be used to display parts of the sphere without remapping by converting from pixel order to raster order; coordinate grids can be overlaid on the faces to assist in registration.
The bin numbering allows the resolution to be changed (by factors of two) simply by adding or deleting LSBs. This is particularly useful if it is desired to increase the bin size, for example to compare maps with different resolutions. This is performed simply by dividing the bin numbers by four for each factor of two increase in bin size.
The bin numbers tend be closer together in sequence for neighboring bins.
The quadtree bin numbering allows maps to be stored as FITS BINTABLES with no wasted storage needed for blank areas (unfilled bins). Only bins containing valid data need be included in a file, as the bin numbers themselves identify the locations of the bins.
The serial numbering frees up one array dimension for another use (e.g., time, a third coordinate, multichannel data, etc.).
The equations for the polar form of the transformation from the sphere to the cube are given in reference 2 as equations (3-21) and (3-38). I have taken this approach further in two areas: 1) inversion of these equations to derive the cube-to-sphere transformation; and 2) simplification of these equations to perform the transformation directly in terms of cartesian coordinates. Both of these are straightforward as shown below.
1. DERIVATION OF THE CUBE-TO-SPHERE TRANSFORMATION
The exact form of the sphere-to-cube tranformation in reference 2 is derived in terms of the co-elevation and azimuth angles (theta, phi) on the sphere and an equivalent set of angles (mu,nu) on the cube. The equal-area tranformation is performed using the following equations:
TAN(mu) = (12/pi)*(theta+ACOS(SIN(theta)*SIN(pi/4))-pi/2) (3-21)
TAN(nu) = SQRT((SEC(mu))**2*(1-COS(theta))/(1-COS(ATAN(SEC(theta))))) (3-38)
The inversion of these equations is as follows. Using the trig identity ASIN(x)+ACOS(x)=pi/2, and the knowledge that SIN(pi/4) = 1/SQRT(2), equation (3-21) can be written as
TAN(mu) = (12/pi)*(theta-ASIN(SIN(theta)/SQRT(2))) (3-21a)
and then rearranged as
theta - (pi/12)*TAN(mu) = ASIN(SIN(theta)/SQRT(2))
Taking the sine of both sides gives
SIN(theta)*COS((pi/12)*TAN(mu)) - COS(theta)*SIN((pi/12)*TAN(mu)) = SIN(theta)/SQRT(2)
Dividing through by COS(theta) and solving for TAN(theta) yields the desired inverse relationship:
TAN(theta) = SIN((pi/12)*TAN(mu))/(COS((pi/12)*TAN(mu))-1/SQRT(2)) (1)
The inversion of equation (3-38) is performed by squaring both sides and solving for COS(phi):
COS(phi) = 1 - (TAN(nu))**2*(1-COS(ATAN(SEC(theta))))/(SEC(mu)**2 (2)
As in the original derivation, this equation used the prior determination of theta as a function of mu.
2. EXPRESSION OF THE TRANSFORMATIONS IN CARTESIAN COORDINATES
Both sets of transformation equations can be expressed in cartesian coordinates (three-dimensional coordinates of the unit vector on the sphere and two-dimensional coordinates on the cube face). Assume the unit vector has coordinates (q,r,s) and the coordinates on the cube face are (x,y), where q is normal to the plane of the cube face and (r,s) cooresponds to (x,y). Then the functions of angles can be expressed as follows:
TAN(theta) = s/r SIN(theta) = s/SQRT(r*r+s*s) COS(phi) = q TAN(mu) = y/x TAN(nu) = SQRT(x*x+y*y)
Substituting these into the transformation equations results in the following revised forms for the equations.
y/x = (12/pi)*(ATAN(s/r)- ASIN(s/SQRT(2*r*r+2*s*s))) (3)
Equation (3-38) can first be simplified by the following substitution:
COS(ATAN(SEC(theta))) = 1/SQRT(1+(SEC(theta))**2)
Substituting this and the appropriate cartesian coordinate expressions gives:
x*x+y*y = (1+(y/x)**2)*(1-q)/(1-1/SQRT(2+(s/r)**2))
which can be re-written to eliminate y:
x = sqrt((1-q)/(1-1/SQRT(2+(s/r)**2))) (4)
Note that there is no ambiquity in the sign of the square root since the derivation of (3-38) assumes that x is positive. This and equation (3) for y/x above give the solution for both y and x, noting also that (3-21) assumes that x is not less than y.
Equation (1) above can be rewritten also by direct substitution:
s/r = SIN((pi/12)*(y/x))/(COS((pi/12)*(y/x))-1/SQRT(2)) (5)
Equation (2) above can be rewritten by solving (4) for q:
q = 1 - x*x*(1-1/SQRT(2+(s/r)**2))) (6)
Thus equations (3) through (6) can be used to convert directly from coordinates on the sphere to the cube and back without the need to convert to angles as intermediate variables, and the number of trigonometric function evaluations required is greatly reduced.
1. Chan, F.K. and O'Neill (1975), Feasibility Study of a Quadrilateralized Spherical Cube Earth Data Base, Computer Sciences Corporation, EPRF Technical Report 2-75 (CSC). Prepared for the Environmental Prediction Research Facility, Monterey, California.
2. O'Neill, I.M. and Laubscher, R.E. (1976), Extended Studies of a Quadrilateralized Spherical Cube Earth Data Base, Computer Sciences Corporation, NEPRF Technical Report 3-76 (CSC). Prepared for the Naval Environmental Prediction Research Facility, Monterey, California.
End of forwarded message ——-