psCoord.h

Go to the documentation of this file.

/** @file  psCoord.h
*
*  @brief Contains basic coordinate transformation definitions and operations
*
*  This file defines the basic types for astronomical coordinate
*  transformation
*
*  @ingroup CoordinateTransform
*
*  @author GLG, MHPCC
*
*  @version $Revision: 1.1.1.1 $ $Name:  $
*  @date $Date: 2005/06/15 21:08:12 $
*
*  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
*/

#ifndef PS_COORD_H
#define PS_COORD_H

#include "psType.h"
#include "psImage.h"
#include "psArray.h"
#include "psList.h"
#include "psFunctions.h"
#include "psTime.h"

/// @addtogroup CoordinateTransform
/// @{

/** Euclidiean Coordinate System.
 *
 *  Both detector and sky positions will be used extensively in the IPP. One
 *  coordinate system to be used is linear coordinates which conform to
 *  Euclidean geometry.
 *
 */
typedef struct
{
    double x;                   ///< x position
    double y;                   ///< y position
    double xErr;                ///< Error in x position
    double yErr;                ///< Error in y position
}
psPlane;

/** Angular Coordinate System
 *
 *  Both detector and sky positions will be used extensively in the IPP. One
 *  coordinate system to be used is angular coordinates for which additional
 *  care must often be taken in comparison to a euclidiean coordinate system.
 *
 */
typedef struct psSphere
{
    double r;                   ///< RA
    double d;                   ///< Dec
    double rErr;                ///< Error in RA
    double dErr;                ///< Error in Dec
}
psSphere;

/** 2D Polynomial Transform
 *
 *  A transform between coordinate systems that consists simply of two 2D
 *  polynomials to transform both components - the output coordinates depend
 *  only on the input coordinates and no other quantities of objects at those
 *  coordinates.
 *
 */
typedef struct
{
    psDPolynomial2D* x;         ///< 2D polynomial transform of X coordinates
    psDPolynomial2D* y;         ///< 2D polynomial transform of Y coordinates
}
psPlaneTransform;

/** 4D Polynomial Transform
 *
 *  A transform between coordinate systems that consists of two 4D polynomials
 *  in which the output coordinates are also specified to be a function of the
 *  magnitude and color of the object with the given coordinates. This type of
 *  coordinate transformation is necessary to represent the (color-dependent)
 *  optical distortions caused by the atmosphere and camera optics, and the
 *  possibly effects of charge transfer inefficiency.
 *
 *  The lowest two terms are the x and y axis of the target system.  The higher
 *  two terms may represent magnitude and color terms.
 */
typedef struct
{
    psDPolynomial4D* x;         ///< 4D polynomial transform of X coordinates
    psDPolynomial4D* y;         ///< 4D polynomial transform of Y coordinates
}
psPlaneDistort;

/** Spherical Transform Definition
 *
 *  We need to be able to convert between ICRS, Galactic and Ecliptic
 *  coordinates, and potentially between arbitrary spherical coordinate
 *  systems. All of these basic spherical transformations represent rotations
 *  of the spherical coordinate reference. We specify a general
 *  transformation function which takes a structure, psSphereTransform,
 *  defining the transformation between two spherical coordinate systems
 *
 */
typedef struct
{
    double alphaP;                    ///< Longitude of the target system pole in the source system
    double cosDeltaP;                 ///< Cosine of target pole latitude in the source system
    double sinDeltaP;                 ///< Sine of target pole latitude in the source system
    double phiP;                      ///< Longitude of the ascending node in the target system
}
psSphereTransform;

/** Projection type for projection/deprojection
 *
 *  @see psProject, psDeproject
 *
 */
typedef enum {
    PS_PROJ_TAN,                ///< Tangent projection
    PS_PROJ_SIN,                ///< Sine projection
    PS_PROJ_AIT,                ///< Aitoff projection
    PS_PROJ_PAR,                ///< Par projection
    //    PS_PROJ_GLS,                ///< GLS projection
    //    PS_PROJ_CAR,                ///< CAR projection
    //    PS_PROJ_MER,                ///< MER projection
    PS_PROJ_NTYPE               ///< Number of types; must be last.
} psProjectionType;

/** Parameter set for projection/deprojection
 *
 *  @see psProject, psDeproject
 *
 */
typedef struct
{
    double R;                   ///< Coordinates of projection center
    double D;                   ///< Coordinates of projection center
    double Xs;                  ///< plate-scale in X direction
    double Ys;                  ///< plate-scale in Y direction
    psProjectionType type;      ///< Projection type
}
psProjection;

/** Mode for Offset calculation between two sky positions
 *
 *  @see  psSphereGetOffset, psSphereSetOffset
 *
 */
typedef enum {
    PS_SPHERICAL,               ///< offset corresponds to an angular offset
    PS_LINEAR                   ///< offset corresponds to a linear offset
} psSphereOffsetMode;

/** The units of the offset
 *
 *  @see  psSphereGetOffset, psSphereSetOffset
 *
 */
typedef enum {
    PS_ARCSEC,                  ///< Arcseconds
    PS_ARCMIN,                  ///< Arcminutes
    PS_DEGREE,                  ///< Degrees
    PS_RADIAN                   ///< Radians
} psSphereOffsetUnit;

/** Allocates a psPlane
 *
 *  @return psPlane*     resulting plane structure.
 */

psPlane* psPlaneAlloc(void);

/** Allocates a psSphere
 *
 *  @return psSphere*     resulting sphere structure.
 */

psSphere* psSphereAlloc(void);


/** Allocates a psPlaneTransform transform.
 *
 *  @return psPlaneTransform*     resulting plane transform
 */

psPlaneTransform* psPlaneTransformAlloc(
    psS32 n1,  ///< The order of the x term in the transform.
    psS32 n2   ///< The order of the y term in the transform.
);

/** Applies the psPlaneTransform transform to a specified coordinate
 *
 *  @return psPlane*     resulting coordinate based on transform
 */
psPlane* psPlaneTransformApply(
    psPlane* out,                      ///< a psPlane to recycle.  If NULL, a new one is generated.
    const psPlaneTransform* transform, ///< the transform to apply
    const psPlane* coords              ///< the coordinate to apply the transform above.
);

/** Allocates a psPlaneDistort transform.
 *
 *  @return psPlaneDistort*     resulting plane distort transform
 */

psPlaneDistort* psPlaneDistortAlloc(
    psS32 n1,  ///< The order of the w term in the transform.
    psS32 n2,  ///< The order of the x term in the transform.
    psS32 n3,  ///< The order of the y term in the transform.
    psS32 n4   ///< The order of the z term in the transform.
);


/** Applies the psPlaneDistort transform to a specified coordinate
 *
 *  @return psPlane*     resulting coordinate based on transform
 */
psPlane* psPlaneDistortApply(
    psPlane* out,                      ///< a psPlane to recycle.  If NULL, a new one is generated.
    const psPlaneDistort* transform,   ///< the transform to apply
    const psPlane* coords,             ///< the coordinate to apply the transform above.
    float term3,                       ///< third term -- maybe magnitude
    float term4                        ///< forth term -- maybe color
);

/** Allocator for psSphereTransform
 *
 *  @return psSphereTransform*         newly allocated struct
 */

psSphereTransform* psSphereTransformAlloc(
    double alphaP,                      ///< north pole latitude
    double deltaP,                      ///< north pole longitude?
    double phiP                         ///< defines the longitude in the input system of the equatorial intersection between the two systems (e.g, the first point of Ares).
);

/** Applies the psSphereTransform transform for a specified coordinate
 *
 *  @return psSphere*      resulting coordinate based on transform
 */
psSphere* psSphereTransformApply(
    psSphere* out,                     ///< a psSphere to recycle.  If NULL, a new one is generated.
    const psSphereTransform* transform,///< the transform to apply
    const psSphere* coord              ///< the coordinate to apply the transform above.x
);

/** Creates the appropriate transform for converting from ICRS to Ecliptic
 *  coordinate systems.
 *
 *  @return psSphereTransform*     transform for ICRS->Ecliptic coordinate systems
 */
psSphereTransform* psSphereTransformICRSToEcliptic(
    psTime *time                        ///< the time for which the resulting transform will be valid
);

/** Creates the appropriate transform for converting from Ecliptic to ICRS
 *  coordinate systems.
 *
 *  @return psSphereTransform*     transform for Ecliptic->ICRS coordinate systems
 */
psSphereTransform* psSphereTransformEclipticToICRS(
    psTime *time                        ///< the time for which the resulting transform will be valid
);

/** Creates the appropriate transform for converting from ICRS to Galactic
 *  coordinate systems.
 *
 */
psSphereTransform* psSphereTransformICRSToGalactic(void);

/** Creates the appropriate transform for converting from Galactic to ICRS
 *  coordinate systems.
 *
 */
psSphereTransform* psSphereTransformGalacticToICRS(void);

/** Allocates memory for a psProjection structure
 *
 *  @return psProjection*    psProjection structure
 */
psProjection* psProjectionAlloc(
    psF64 R,                   ///< Right-ascension of projection center.
    psF64 D,                   ///< Declination of projection center.
    psF64 Xs,                  ///< Scale in x-dimension
    psF64 Ys,                  ///< Scale in y-dimension
    psProjectionType type
);

/** Projects a spherical coordinate to a linear coordinate system
 *
 *  @return psPlane*    projected coordinate
 */
psPlane* psProject(
    const psSphere* coord,             ///< coordinate to project
    const psProjection* projection     ///< parameters of the projection
);

/** Reverse projection of a linear coordinate to a spherical coordinate system
 *
 *  @return psPlane*    projected coordinate
 */
psSphere* psDeproject(
    const psPlane* coord,              ///< coordinate to project
    const psProjection* projection     ///< parameters of the projection
);

/** Determines the offset (RA,Dec) on the sky between two positions.
 *
 *  Both an offset mode and an offset unit may be defined. The mode may be
 *  either PS_SPHERICAL, in which case the specified offset corresponds to an
 *  offset in angles, or it may be PS_LINEAR, in which case the offset
 *  corresponds to a linear offset in a local projection. The offset unit may
 *  be in one of PS_ARCSEC, PS_ARCMIN, PS_DEGREE, and PS_RADIAN, which
 *  specifies the units of the offset only.
 *
 *  @return psSphere*        the offset between position1 and position2
 */
psSphere* psSphereGetOffset(
    const psSphere* position1,         ///< first position for calculating offset
    const psSphere* position2,         ///< second position for calculating offset
    psSphereOffsetMode mode,           ///< type of offset can be PS_SPHERICAL or PS_LINEAR
    psSphereOffsetUnit unit            ///< specifies the units of offset only
);

/** Applies the given offset to a coordinate.
 *
 *  Both an offset mode and an offset unit may be defined. The mode may be
 *  either PS_SPHERICAL, in which case the specified offset corresponds to an
 *  offset in angles, or it may be PS_LINEAR, in which case the offset
 *  corresponds to a linear offset in a local projection. The offset unit may
 *  be in one of PS_ARCSEC, PS_ARCMIN, PS_DEGREE, and PS_RADIAN, which
 *  specifies the units of the offset only.
 *
 *  @return psSphere*              the original position with the given offset applied.
 */
psSphere* psSphereSetOffset(
    const psSphere* position,          ///< coordinate of origin
    const psSphere* offset,            ///< coordinate of offset to apply
    psSphereOffsetMode mode,           ///< corresponds to an offset in angles or local projection
    psSphereOffsetUnit unit            ///< specifies the units of offset only
);

/** Generates the complete spherical rotation to account for precession
 *  between two times.  The equinoxes shall be Julian equinoxes.
 *
 *  @return psSphere* the resulting spherical rotation
 */
psSphere* psSpherePrecess(
    psSphere *coords,                  ///< coordinates (modified in-place)
    const psTime *fromTime,            ///< equinox of coords input
    const psTime *toTime               ///< equinox of coords output
);

/** Takes a given transform and inverts it linearly if possible.
 *
 *  @return psPlaneTransform
 *  the linearly inverted transform
*/
psPlaneTransform *p_psPlaneTransformLinearInvert(
    psPlaneTransform *transform        ///<    transform to invert
);


/** Takes a transform and tests whether or not it is a linear projection.
 *
 *  @return psS32
 *  the order of the projection
*/
psS32 p_psIsProjectionLinear(
    psPlaneTransform *transform        ///<     transform to test for linearity
);

/** inverts a given transformation.
 *
 *  It may assume that the input transformation is one-to-one, and that the
 *  inverse transformation may be specified through using polynomials of the
 *  same type and order as the forward transformation. In the event that the
 *  input transformation is linear, an exact solution may be calculated;
 *  otherwise nSamples samples in each axis, covering the region specified by
 *  region shall be used as a grid to fit the best inverse transformation. The
 *  function shall return NULL if it was unable to generate the inverse
 *  transformation; otherwise it shall return the inverse transformation. In
 *  the event that out is NULL, a new psPlaneTransform shall be allocated and
 *  returned.
 *
 *  @return psPlaneTransform*  the resulting inverted transform
 */
psPlaneTransform* psPlaneTransformInvert(
    psPlaneTransform *out,             ///< a transform to recycle, or NULL if one is to be created.
    const psPlaneTransform *in,        ///< transform to invert
    psRegion region,                   ///< region to fit for non-linear transform inversion
    int nSamples                       ///< number of samples in each axis for fit
);

/** Creates a single transformation that has the effect of performing trans1
 *  followed by trans2.
 *
 *  psPlaneTransformCombine takes two transformations (trans1 and trans2) and
 *  returns a single transformation that has the effect of performing trans1
 *  followed by trans2. In the event that the input transformation is linear,
 *  an exact solution may be calculated; otherwise nSamples samples in each
 *  axis, covering the region specified by region shall be used as a grid to
 *  fit the best inverse transformation. The function shall return NULL if it
 *  was unable to generate the transformation; otherwise it shall return the
 *  transformation.
 *
 *  @return psPlaneTransform*    resulting transformation
 */
psPlaneTransform* psPlaneTransformCombine(
    psPlaneTransform *out,             ///< a transform to recycle, or NULL if one is to be created.
    const psPlaneTransform *trans1,    ///< first transform to combine
    const psPlaneTransform *trans2,    ///< first transform to combine
    psRegion region,                   ///< region to cover (for non-linear transforms)
    int nSamples                       ///< number of samples on each axis (for non-linear transforms)
);


/** takes two arrays containing matched coordinates and returns the
 *  best-fitting transformation.
 *
 *  psPlaneTransformFit takes two arrays containing matched coordinates (i.e.,
 *  coordinates in the source array correspond to the coordinates in the dest
 *  array) and returns the best-fitting transformation. The source and dest
 *  will contain psCoords. In the event that the number of coordinates in each
 *  is not identical, the function shall generate a warning, and extra
 *  coordinates in the longer of the two shall be ignored. The trans transform
 *  may not be NULL, since it specifies the desired order, polynomial type and
 *  any polynomial terms to mask. nRejIter rejection iterations shall be
 *  performed, wherein coordinates lying more than sigmaClip standard
 *  deviations from the fit shall be rejected.
 *
 *  @return bool        TRUE if successful, otherwise FALSE.
 */
bool psPlaneTransformFit(
    psPlaneTransform *trans,
    const psArray *source,
    const psArray *dest,
    int nRejIter,
    float sigmaClip
);
//XXX: need to add doxygen comments on the parameters above. -rdd

/// @}

#endif // #ifndef PS_COORD_H

---