Coordinates_module

More...

Modules
	Coordinates_module_internal_classes
	Internal Coordinates_module classes and functions.

Classes
class	casacore::Coordinate
	More...
class	casacore::CoordinateSystem
	More...
class	casacore::CoordinateUtil
	More...
class	casacore::DirectionCoordinate
	More...
class	casacore::FITSCoordinateUtil
	More...
class	casacore::FrequencyAligner< T >
	More...
class	casacore::GaussianConvert
	More...
class	casacore::LinearCoordinate
	More...
class	casacore::ObsInfo
	More...
class	casacore::Projection
	More...
class	casacore::QualityCoordinate
	More...
class	casacore::SpectralCoordinate
	More...
class	casacore::StokesCoordinate
	More...
class	casacore::TabularCoordinate
	More...

Functions
ostream &	casacore::operator<< (ostream &os, const ObsInfo &info)
	Global functions More...

Detailed Description

Classes to interconvert pixel and world (physical) coordinates

See below for an overview of the classes in this module.

Prerequisite

• Knowledge of astronomical coordinate conversions in general. Probably the best documents are the papers by Mark Calabretta and Eric Greisen. The initial draft from 1996 can be found at http://www.atnf.csiro.au/~mcalabre. It is this draft that the Coordinate classes are based upon. Since then, this paper has evolved into three which can be found at the above address, and will be published in the Astronomy and Astrophysics Supplement Series (probably in 2000). The design has changed since the initial draft. When these papers are finalized, and the IAU has ratified the new standards, WCSLIB (Mark Calabretta's implementation of these conventions) will be revised for the new designs. At that time, the Coordinate classes may also be revised.
• Generic Casacore classes; especially those in the Arrays module.
• The Measures module.

Review Status

Reviewed By:

Peter Barnes

Date Reviewed:

1999/12/24

Synopsis

The primary notion is that a Coordinate can interconvert between a length "n" Vector<Double> (the pixel coordinate) and a length "m" Vector<Double> (the "world" coordinate). Note that "m" and "n" do not in principle have to be the same (so that one can get both the RA and DEC from an image slice, for example), however in practice they currently always are. Each Coordinate has the full mapping from pixel to world coordinates, i.e. it has a reference value, reference pixel, increments, and an arbitrary transformation matrix. To go from a pixel to a world coordinate the following steps are applied:

1. The reference pixel is subtracted from the pixel position.
2. The result is multiplied by the transformation matrix.
3. The result is multiplied by an increment per output axis to convert it into physical coordinates.
4. For some coordinate types (e.g., Direction), a non-linear function is applied to this result.

The classes are arranged as follows. The base class is Coordinate which defines the interface. Classes derived from it are

1. DirectionCoordinate
   
2. LinearCoordinate
   
3. SpectralCoordinate
   
4. TabularCoordinate
   
5. StokesCoordinate
   
6. CoordinateSystem
   

Other classes are Projection
which is used to specify an astronomical projection for DirectionCoordinates, and LinearXform
a helper class which the application programmer will not interact with.

CoordinateSystem is the class that application programmers will usually interact with. A CoordinateSystem consists of a collection of the other classes derived from Coordinate. Normally one group will be for RA/DEC, another for a Stokes axis, and another group for the spectral axis. The axes may be transposed arbitrarily, for example RA could be the first axis, and DEC the third.

Normally the CoordinateSystem being manipulated will be embedded in a PagedImage or other object. Note that the axes of the PagedImage do not necessarily map directly to the axes in the CoordinateSystem. Functionality is provided to determine this mapping.

One or more axes from the CoordinateSystem may be removed. Pixel axes and/or world axes may be removed. You are encouraged to leave all the world axes when you remove pixel axes.
If a world axis is removed, the corresponding pixel axis is also removed. This means that one can be sure that a pixel axis always has a corresponding world axis (it makes no sense otherwise). The opposite is not necessarily true: a world axis can exist without a pixel axis.

The linear transformation and sky projection computations are carried out in an underlying library – WCSLIB – written by Mark Calabretta of the ATNF.

Caution: All pixels coordinates are zero relative;

Example

First, let's make a DirectionCoordinate — used to represent a direction, usually an RA/DEC, but it could also be, e.g., an AZ/EL pair.

Matrix<Double> xform(2,2);                                    // 1
xform = 0.0; xform.diagonal() = 1.0;                          // 2
Quantum<Double> refLon(135.0, "deg");
Quantum<Double> refLat(60.0, "deg");
Quantum<Double> incLon(-1.0, "deg");
Quantum<Double> incLat(1.0, "deg");
DirectionCoordinate radec(MDirection::J2000,                  // 3
Projection(Projection::SIN),          // 4
refLon, refLat,                       // 5
incLon, incLat,                       // 6
xform,                                // 7
128, 128);                            // 8

• 1-2:Here we set up a diagonal transformation matrix. Normally this matrix should be diagonal, however if you wanted to introduce a rotation or skew, you would do it through this matrix.
• 3:This defines the astronomical type of the world coordinate. Most of the time it will probably be J2000 or B1950, but there are many other possibilities as listed in the MDirection class header.
• 4:The Projection class defines the "geometry" that is used to map xy<–>world. SIN is the most common projection for radio interferometers. Note that SIN can optionally take parameters as defined in Calabretta and Greisen. If not provided, they default to 0.0, which is the "old" SIN convention.
• 5:Set the reference position to RA=135, DEC=60 degrees. Note that the native units of a DirectionCoordinate is radians.
• 6: Set the increments to -1 degree in RA, and +1 degree in DEC.
• 7: Set the previously defined transformation matrix.
• 8: Set the zero-relative reference pixel. Note that it does not have to be incremental. At the reference pixel, the world coordinate has the reference value.

Although we happeend to create our DirectionCoordinate with Quanta in degrees, these have been converted to radians by the constructor. We can set the native units to degrees if we wish as follows:

Vector<String> units(2); units = "deg";                       //  9
radec.setWorldAxisUnits(units);                               // 10

The increment and reference value are updated appropriately.

Set up a couple of vectors to use the world and pixel coordinate values.

Vector<Double> world(2), pixel(2);                            // 11
pixel = 138.0;                                                // 12

We use 138 as an abitrary pixel position which is near the reference pixel so we can tell if the answers look foolish or not.

We can actually perform a transformation like this as follows. If it succeeds we print the value of the world coordinate.

Bool ok = radec.toWorld(world, pixel);                        // 13
if (!ok) {                                                    // 14
cout << "Error: " << radec.errorMessage() << endl;          // 15
return 1;                                                   // 16
}                                                             // 17
cout << world << " <--- " << pixel << endl;         // 18

There is an overloaded "toWorld" function that produces an MDirection in case you want to, e.g., find out what the position in B1950 coordinates would be.

The reverse transformation takes place similarly:

ok = radec.toPixel(pixel, world);                             // 19

Suppose we have an image with a Stokes axis. It can be set up as follows:

Vector<Int> iquv(4);
iquv(0) = Stokes::I; iquv(1) = Stokes::Q;                    // 21
iquv(2) = Stokes::U; iquv(3) = Stokes::V;                    // 22
StokesCoordinate stokes(iquv);                               // 23

We create an integer array the same length as the Stokes axis, and place the corresponding Stokes enum into each element of the array. The values must be unique, e.g. there can only be one "I" plane. Besides the generic Vector<Double> toWorld/toPixel interface, you can also directly interconvert between Stokes enum and and (zero-relative) plane number:

Int plane;                                                   // 24
ok = stokes.toPixel(plane, Stokes::Q);                       // 25

Here it will return True and set plane to 1. On the other hand, it would return False for:

ok = stokes.toPixel(plane, Stokes::XX);                      // 26

since "XX" is not one of the Stokes enumerations we used to create this coordinate.

A Spectral ("frequency") coordinate may be created as follows:

SpectralCoordinate spectral(MFrequency::TOPO,               // 27
1.4E+9,                         // 28
2.0E+4,                         // 29
0,                              // 30
1420.40575E+6);                 // 31

The default frequency units of a spectral coordinate are Hz, although they may be changed to whatever is convenient. The first line (27) defines the type of frequency we have – topocentric here. The second (28) line defines the frequency at the reference pixel, 0 (28) here. The channel increment is defined on line 29. A rest frequency of the spectral may be provided. It is useful in calculating doppler velocities. These calculations are carried out by the MFrequency and MDoppler classes of the Measures system.

Motivation

The primary motivation is to provide support for converting pixel locations in an image to physical ("world") positions.

To Do

• Add measures interfaces that handle reference frame conversions
• offset coordinates

Function Documentation

operator<<()

ostream& casacore::operator<<	(	ostream &	os,
		const ObsInfo &	info
	)

Global functions

Output declaration - useful for debugging.