World coordinate system (WCS) related functions

Data Structures

struct  muse_wcs_object
 WCS object to store data needed while computing the astrometric solution. More...
struct  muse_wcs
 A structure containing a spatial two-axis WCS. More...

Defines

#define MUSE_WCS_KEYS
 regular expression for WCS properties
#define MUSE_WCS_DETIMAGE_EXTNAME   "ASTROMETRY_DETECTION"
 name for special astrometry detection image

Enumerations

enum  muse_wcs_centroid_type { MUSE_WCS_CENTROID_GAUSSIAN = 0, MUSE_WCS_CENTROID_MOFFAT, MUSE_WCS_CENTROID_BOX }
 

Type of centroiding algorithm to use.

More...

Functions

muse_wcs_objectmuse_wcs_object_new (void)
 Allocate memory for a new muse_wcs_object object.
void muse_wcs_object_delete (muse_wcs_object *aWCSObj)
 Deallocate memory associated to a muse_wcs_object.
cpl_table * muse_wcs_centroid_stars (muse_image *aImage, float aSigma, muse_wcs_centroid_type aCentroid)
 Detect and centroid stars on an image of an astrometric exposure.
cpl_error_code muse_wcs_locate_sources (muse_pixtable *aPixtable, float aSigma, muse_wcs_centroid_type aCentroid, muse_wcs_object *aWCSObj)
 Determine the centers of stars on an astrometric exposure.
cpl_error_code muse_wcs_solve (muse_wcs_object *aWCSObj, cpl_table *aReference, float aRadius, float aFAccuracy, int aIter, float aThresh)
 Find the world coordinate solution for a given field using coordinates of detected sources and coordinates of reference values.
cpl_error_code muse_wcs_optimize_solution (muse_wcs_object *aWCSObj, float aDetSigma, muse_wcs_centroid_type aCentroid, cpl_table *aReference, float aRadius, float aFAccuracy, int aIter, float aRejSigma)
 Find an optimal astrometry solution given a detection image.
cpl_propertylist * muse_wcs_create_default (void)
 Create FITS headers containing a default (relative) WCS.
cpl_propertylist * muse_wcs_apply_cd (const cpl_propertylist *aExpHeader, const cpl_propertylist *aCalHeader)
 Apply the CDi_j matrix of an astrometric solution to an observation.
cpl_error_code muse_wcs_project_tan (muse_pixtable *aPixtable, const cpl_propertylist *aWCS)
 Carry out a gnomonic projection of a pixel table into native spherical coordinates.
cpl_error_code muse_wcs_position_celestial (muse_pixtable *aPixtable, double aRA, double aDEC)
 Convert native to celestial spherical coordinates in a pixel table.
cpl_error_code muse_wcs_celestial_from_pixel (cpl_propertylist *aHeader, double aX, double aY, double *aRA, double *aDEC)
 Convert from pixel coordinates to celestial spherical coordinates.
cpl_error_code muse_wcs_pixel_from_celestial (cpl_propertylist *aHeader, double aRA, double aDEC, double *aX, double *aY)
 Convert from celestial spherical coordinates to pixel coordinates.
cpl_error_code muse_wcs_projplane_from_celestial (cpl_propertylist *aHeader, double aRA, double aDEC, double *aX, double *aY)
 Convert from celestial spherical coordinates to projection plane coordinates.
cpl_error_code muse_wcs_projplane_from_pixel (cpl_propertylist *aHeader, double aX, double aY, double *aXOut, double *aYOut)
 Convert from pixel coordinates to projection plane coordinates.
cpl_error_code muse_wcs_pixel_from_projplane (cpl_propertylist *aHeader, double aX, double aY, double *aXOut, double *aYOut)
 Convert from projection plane coordinates to pixel coordinates.
cpl_error_code muse_wcs_get_angles (cpl_propertylist *aHeader, double *aXAngle, double *aYAngle)
 Compute the rotation angles (in degrees) from the FITS header WCS.
cpl_error_code muse_wcs_get_scales (cpl_propertylist *aHeader, double *aXScale, double *aYScale)
 Compute the spatial scales (in degrees) from the FITS header WCS.
muse_wcsmuse_wcs_new (cpl_propertylist *aHeader)
 Create a new WCS structure from a given FITS header.
static void muse_wcs_celestial_from_pixel_fast (muse_wcs *aWCS, double aX, double aY, double *aRA, double *aDEC)
 Convert from pixel coordinates to celestial spherical coordinates.
static void muse_wcs_pixel_from_celestial_fast (muse_wcs *aWCS, double aRA, double aDEC, double *aX, double *aY)
 Convert from celestial spherical coordinates to pixel coordinates.
static void muse_wcs_projplane_from_pixel_fast (muse_wcs *aWCS, double aX, double aY, double *aXOut, double *aYOut)
 Convert from pixel coordinates to projection plane coordinates.
static void muse_wcs_pixel_from_projplane_fast (muse_wcs *aWCS, double aX, double aY, double *aXOut, double *aYOut)
 Convert from projection plane coordinates to pixel coordinates.

Variables

const muse_cpltable_def muse_wcs_detections_def []
 Definition of the table structure for the astrometric field detections.
const muse_cpltable_def muse_wcs_reference_def []
 Definition of the table structure for the astrometric reference sources.
const muse_cpltable_def muse_wcs_reference_def []
 Definition of the table structure for the astrometric reference sources.

Detailed Description


Define Documentation

#define MUSE_WCS_DETIMAGE_EXTNAME   "ASTROMETRY_DETECTION"

name for special astrometry detection image

Extension name used to save the image that is used to detect the stars in the astrometric field.

Definition at line 56 of file muse_wcs.h.

Referenced by muse_wcs_locate_sources(), and muse_wcs_optimize_solution().

#define MUSE_WCS_KEYS
Value:
"^C(TYPE|UNIT|RPIX|RVAL|DELT|SYER|RDER)|^CD[0-9]+_[0-9]+|" \
                       "^WCSAXES$|^L[OA][NT]POLE$"

regular expression for WCS properties

Regular expression to copy/erase all keywords relevant to world coordinates this list is mostly taken from the FITS Standard v3.0 (Tab. 8.2, p. 77).

Definition at line 48 of file muse_wcs.h.

Referenced by muse_datacube_save(), muse_datacube_save_recimages(), muse_fov_load(), muse_image_save(), muse_pixtable_create(), muse_resampling_cube(), and muse_wcs_project_tan().


Enumeration Type Documentation

Type of centroiding algorithm to use.

Enumerator:
MUSE_WCS_CENTROID_GAUSSIAN 

Gaussian profile fit

MUSE_WCS_CENTROID_MOFFAT 

Moffat profile fit

MUSE_WCS_CENTROID_BOX 

simplest version: marginal centroid in a box

Definition at line 90 of file muse_wcs.h.


Function Documentation

cpl_propertylist* muse_wcs_apply_cd ( const cpl_propertylist *  aExpHeader,
const cpl_propertylist *  aCalHeader 
)

Apply the CDi_j matrix of an astrometric solution to an observation.

Parameters:
aExpHeader the header of the observation
aCalHeader the header giving the astrometric solution
Returns:
the header with CDi_j matrix that merges both rotations or NULL on error.

This function uses matrix multiplication to merge position angle on the sky with the astrometric solution. The main output are the four spatial CDi_j entries in the output propertylist.

Exceptions:
set CPL_ERROR_NULL_INPUT, return NULL aCalHeader and/or aExpHeader are NULL

Definition at line 1234 of file muse_wcs.c.

References muse_astro_posangle(), muse_pfits_get_cd(), muse_wcs_get_angles(), and muse_wcs_get_scales().

Referenced by muse_flux_integrate_cube(), and muse_wcs_project_tan().

cpl_error_code muse_wcs_celestial_from_pixel ( cpl_propertylist *  aHeader,
double  aX,
double  aY,
double *  aRA,
double *  aDEC 
)

Convert from pixel coordinates to celestial spherical coordinates.

Parameters:
aHeader the input header containing the WCS of the exposure
aX the horizontal pixel position
aY the vertical pixel position
aRA the output right-ascension in degrees
aDEC the output declination in degrees
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader is used to do the transformation, only the gnomonic (TAN) is supported.

Uses Eqns (14), (15), (55), and (2) from Calabretta & Greisen 2002 A&A 395, 1077 (Paper II). We use that phi_p = 180 deg for zenithal projections (like TAN). The equations are actually implemented in muse_wcs_celestial_from_pixel_fast() that is called by this function.

Exceptions:
return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL
return CPL_ERROR_UNSUPPORTED_MODE aHeader does not contain gnomonic (TAN) WCS

Definition at line 1556 of file muse_wcs.c.

References muse_pfits_get_ctype(), muse_wcs_celestial_from_pixel_fast(), and muse_wcs_new().

Referenced by muse_flux_integrate_cube().

static void muse_wcs_celestial_from_pixel_fast ( muse_wcs aWCS,
double  aX,
double  aY,
double *  aRA,
double *  aDEC 
) [static]

Convert from pixel coordinates to celestial spherical coordinates.

Parameters:
aWCS the input WCS structure
aX the horizontal pixel position
aY the vertical pixel position
aRA the output right-ascension in degrees
aDEC the output declination in degrees
Warning:
This function does not do any safety checks. If you are not sure, if the WCS does not contain a proper gnomonic (TAN) setup, use muse_wcs_celestial_from_pixel() instead!
Note:
This function is defined as static inline, i.e. should get integrated into the caller by the compiler, to maximise execution speed.

See muse_wcs_celestial_from_pixel() for the origin of the equations.

Definition at line 138 of file muse_wcs.h.

References cd22, crpix2, and crval2.

Referenced by muse_resampling_collapse_pixgrid(), muse_resampling_euro3d(), and muse_wcs_celestial_from_pixel().

cpl_table* muse_wcs_centroid_stars ( muse_image aImage,
float  aSigma,
muse_wcs_centroid_type  aCentroid 
)

Detect and centroid stars on an image of an astrometric exposure.

Parameters:
aImage reconstructed image of the astrometric exposure
aSigma the sigma level used for object detection
aCentroid type of centroiding algorithm to use
Returns:
the table with object properties on success or NULL on failure
Remarks:
aCentroid be either MUSE_WCS_CENTROID_GAUSSIAN, MUSE_WCS_CENTROID_MOFFAT, or MUSE_WCS_CENTROID_BOX.

Find objects in this image and determine the centroid of each object by fitting a two-dimensional function of the given type to it. Record the resulting central pixel coordinates, an estimate of their positional error, and their fluxes and FWHMs in the output table. Exclude objects where FWHM estimates cannot be determined from the image data, since they likely lie on the edge of the image.

Exceptions:
set CPL_ERROR_NULL_INPUT, return NULL aImage or its data or stat components are NULL
set CPL_ERROR_ILLEGAL_INPUT, return NULL aSigma is not positive
set CPL_ERROR_ILLEGAL_INPUT, return NULL an unknown centroiding type was requested in aCentroid
set CPL_ERROR_DATA_NOT_FOUND, return NULL no sources were detected at the given sigma level

Definition at line 178 of file muse_wcs.c.

References muse_image::data, muse_image::dq, muse_image::header, muse_cpltable_new(), muse_image_delete(), muse_image_new(), muse_pfits_get_fwhm_end(), muse_pfits_get_fwhm_start(), muse_pfits_get_mode(), muse_quality_image_reject_using_dq(), MUSE_UTILS_CENTROID_MEDIAN, muse_utils_fit_moffat_2d(), muse_utils_image_get_centroid_window(), MUSE_WCS_CENTROID_BOX, MUSE_WCS_CENTROID_GAUSSIAN, MUSE_WCS_CENTROID_MOFFAT, and muse_image::stat.

Referenced by muse_wcs_locate_sources(), and muse_wcs_optimize_solution().

cpl_propertylist* muse_wcs_create_default ( void   ) 

Create FITS headers containing a default (relative) WCS.

Returns:
A cpl_propertylist * on success, failures are not expected to occur.

Units in the WCS are in deg (as required by the FITS standard, ยง8.2), approximately relative to a MUSE field center. Rotation is not known in this function, so assumed to be zero. This needs to be adapted when applying this WCS to data.

Definition at line 1171 of file muse_wcs.c.

Referenced by muse_flux_integrate_cube(), and muse_wcs_solve().

cpl_error_code muse_wcs_get_angles ( cpl_propertylist *  aHeader,
double *  aXAngle,
double *  aYAngle 
)

Compute the rotation angles (in degrees) from the FITS header WCS.

Parameters:
aHeader the input header containing the WCS of the exposure
aXAngle the output angle in x-direction
aYAngle the output angle in y-direction
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader, i.e. the CDi_j matrix, is used to compute the angles.

References:

Exceptions:
return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL
propagate the error of cpl_propertylist_get_double() some elements of the CDi_j matrix don't exist in aHeader

Definition at line 1775 of file muse_wcs.c.

References muse_pfits_get_cd().

Referenced by muse_wcs_apply_cd(), and muse_wcs_solve().

cpl_error_code muse_wcs_get_scales ( cpl_propertylist *  aHeader,
double *  aXScale,
double *  aYScale 
)

Compute the spatial scales (in degrees) from the FITS header WCS.

Parameters:
aHeader the input header containing the WCS of the exposure
aXScale the output scale in x-direction
aYScale the output scale in y-direction
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader, i.e. the CDi_j matrix, is used to compute the scales.

References:

Exceptions:
return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL
propagate the error of cpl_propertylist_get_double() some elements of the CDi_j matrix don't exist in aHeader

Definition at line 1822 of file muse_wcs.c.

References muse_pfits_get_cd().

Referenced by muse_dar_correct(), muse_resampling_collapse_pixgrid(), muse_resampling_image(), muse_wcs_apply_cd(), and muse_wcs_solve().

cpl_error_code muse_wcs_locate_sources ( muse_pixtable aPixtable,
float  aSigma,
muse_wcs_centroid_type  aCentroid,
muse_wcs_object aWCSObj 
)

Determine the centers of stars on an astrometric exposure.

Parameters:
aPixtable the input pixel table containing the astrometric exposure
aSigma the sigma level used for object detection
aCentroid type of centroiding algorithm to use
aWCSObj the MUSE WCS object
Returns:
CPL_ERROR_NONE on success, another CPL error code on failure
Remarks:
The table returned in aWCSObj->detected contains the x and y positions of all detected sources. This function also stores the center (components xcenter, ycenter and ra, dec) of the dataset in aWCSObj as well as the datacube used for the detection (component cube).
This function assumes that the differential atmospheric refraction was already corrected in the data beforehand, this will be checked by searching for the respective FITS headers (MUSE_HDR_PT_DAR_NAME or MUSE_HDR_PT_DAR_CORR). If this is not the case, a warning is printed and an error is set, but the data is fully processed and the return code is not changed.
aCentroid be either MUSE_WCS_CENTROID_GAUSSIAN, MUSE_WCS_CENTROID_MOFFAT, or MUSE_WCS_CENTROID_BOX.

Resample the input pixel table to a datacube with very coarse wavelength sampling (50 Angstrom/pixel) and use the central three wavelength planes to reconstruct a median-combined image. This special image is also appended to the list of reconstructed image in the created datacube, under the (extension) name MUSE_WCS_DETIMAGE_EXTNAME, as 2nd image after the white-light image. The special detection image is the passed on to muse_wcs_centroid_stars() to compute the exact centroids using the given method aCentroid in a table.

Exceptions:
return CPL_ERROR_NULL_INPUT inputs aPixtable, its header component, or aWCSObj are NULL
return CPL_ERROR_ILLEGAL_INPUT aSigma was non-positive
return CPL_ERROR_ILLEGAL_INPUT an unknown centroiding type was requested
output warning, set CPL_ERROR_UNSUPPORTED_MODE, but continue processing input pixel table was not corrected for DAR
propagate error of muse_resampling_cube() cube could not be created
propagate error image could not be created from cube
propagate error from muse_wcs_centroid_stars() sources detection/centroiding failed

Definition at line 493 of file muse_wcs.c.

References muse_resampling_params::crsigma, muse_resampling_params::crtype, muse_image::data, muse_datacube::data, muse_datacube::dq, muse_image::dq, muse_datacube::header, muse_image::header, muse_pixtable::header, muse_combine_median_create(), muse_datacube_collapse(), muse_datacube_save(), MUSE_HDR_PT_DAR_CORR, MUSE_HDR_PT_DAR_NAME, muse_image_new(), muse_image_save(), muse_imagelist_delete(), muse_imagelist_new(), muse_imagelist_set(), muse_pfits_get_dec(), muse_pfits_get_ra(), muse_pixtable_reset_dq(), muse_pixtable_save(), MUSE_RESAMPLE_WEIGHTED_DRIZZLE, MUSE_RESAMPLING_CRSTATS_MEDIAN, muse_resampling_cube(), muse_resampling_params_delete(), muse_resampling_params_new(), muse_table_load_filter(), MUSE_WCS_CENTROID_BOX, MUSE_WCS_CENTROID_GAUSSIAN, MUSE_WCS_CENTROID_MOFFAT, muse_wcs_centroid_stars(), MUSE_WCS_DETIMAGE_EXTNAME, muse_resampling_params::pfx, muse_datacube::recimages, muse_datacube::recnames, muse_datacube::stat, and muse_image::stat.

Referenced by muse_postproc_process_exposure().

muse_wcs* muse_wcs_new ( cpl_propertylist *  aHeader  ) 

Create a new WCS structure from a given FITS header.

Parameters:
aHeader the input header containing the input WCS
Returns:
A newly allocated muse_wcs *.

The world coordinate system from aHeader, i.e. the CDi_j matrix, the CRPIXi, CRVALi, are used to fill the structure. The iscelsph element of the structure is set to FALSE, it can be easily set using other means.

The returned pointer has to be deallocated using cpl_free().

Exceptions:
return the structure initialized with zeros, except cd11 = cd22 = cddet = 1. aHeader is NULL or does not contain the needed values
set CPL_ERROR_SINGULAR_MATRIX, but continue normally the determinant of the CDi_j matrix (the cddet element) is zero

Definition at line 1867 of file muse_wcs.c.

References cd22, cddet, crpix2, crval2, muse_pfits_get_cd(), muse_pfits_get_crpix(), and muse_pfits_get_crval().

Referenced by muse_euro3dcube_collapse(), muse_pixgrid_create(), muse_resampling_collapse_pixgrid(), muse_resampling_euro3d(), muse_wcs_celestial_from_pixel(), muse_wcs_pixel_from_celestial(), muse_wcs_pixel_from_projplane(), and muse_wcs_projplane_from_pixel().

void muse_wcs_object_delete ( muse_wcs_object aWCSObj  ) 

Deallocate memory associated to a muse_wcs_object.

Parameters:
aWCSObj input MUSE wcs object

Just calls the required *_delete() functions for each component before freeing the memory for the pointer itself. As a safeguard, it checks if a valid pointer was passed, so that crashes cannot occur.

Definition at line 89 of file muse_wcs.c.

References muse_datacube_delete().

muse_wcs_object* muse_wcs_object_new ( void   ) 

Allocate memory for a new muse_wcs_object object.

Returns:
a new muse_wcs_object * or NULL on error
Remarks:
The returned object has to be deallocated using muse_wcs_object_delete().
This function does not allocate the contents of the elements, these have to be allocated with the respective *_new() functions.

Simply allocate memory to store the pointers of the muse_wcs_object structure.

Definition at line 71 of file muse_wcs.c.

Referenced by muse_postproc_process_exposure().

cpl_error_code muse_wcs_optimize_solution ( muse_wcs_object aWCSObj,
float  aDetSigma,
muse_wcs_centroid_type  aCentroid,
cpl_table *  aReference,
float  aRadius,
float  aFAccuracy,
int  aIter,
float  aRejSigma 
)

Find an optimal astrometry solution given a detection image.

Parameters:
aWCSObj the MUSE WCS object
aDetSigma the (negative!) sigma level used for object detection
aCentroid type of centroiding algorithm to use
aReference reference coordinates for the same field
aRadius initial radius for pattern matching in pixels
aFAccuracy factor to set initial accuracy relative to mean position accuracy in aWCSObj->detected
aIter number of iterations
aRejSigma the rejection threshold
Returns:
CPL_ERROR_NONE on success, another CPL error code on failure
Note:
aWCSObj has to contain the cube with the associated detection image, as produced by muse_wcs_locate_sources().
Remarks:
The propertylist returned in aWCSObj->wcs contains the world coordinate solution.

Repeatedly call muse_wcs_centroid_stars() and muse_wcs_solve() to find the best detection sigma level for the most accurate astrometric solution. The repeated measurements start with a detection sigma of |aDetSigma| and loop down to the 1.0-sigma level.

The best solution is identified, as the one that gives the highest value of

nid / 50. * min(medrms_x) / medrms_x * min(medrms_y) / medrms_y

where nid is the identified number of stars, 50 is something of a lower limit for the useful number of stars, and medrms_[xy] are the median scatter left in the solution in x and y.

The final detection sigma value is recorded in the header keyword MUSE_HDR_WCS_DETSIGMA.

Note:
This function changes the CPL messaging level and is hence not supposed to be used from parallel code.
Exceptions:
return CPL_ERROR_NULL_INPUT input WCS object or its cube component are NULL
return CPL_ERROR_DATA_NOT_FOUND the astrometric detection image is missing from the cube
return CPL_ERROR_ILLEGAL_INPUT aDetSigma is not negative
return CPL_ERROR_ILLEGAL_INPUT an unknown centroiding type was requested
return CPL_ERROR_ILLEGAL_INPUT table of reference coordinates is empty or NULL
return CPL_ERROR_BAD_FILE_FORMAT table of detected or reference coordinates don't contain the expected columns
return CPL_ERROR_ILLEGAL_INPUT aRadius and/or aFAccuracy are non-positive
return CPL_ERROR_ILLEGAL_OUTPUT no objects could be identified at any sigma level
propagate return code of muse_wcs_solve() the solution with the finally selected detection sigma failed

Definition at line 985 of file muse_wcs.c.

References muse_datacube::header, muse_cpltable_check(), muse_imagelist_get(), muse_pfits_get_mode(), MUSE_WCS_CENTROID_BOX, MUSE_WCS_CENTROID_GAUSSIAN, MUSE_WCS_CENTROID_MOFFAT, muse_wcs_centroid_stars(), MUSE_WCS_DETIMAGE_EXTNAME, muse_wcs_solve(), muse_datacube::recimages, and muse_datacube::recnames.

Referenced by muse_postproc_process_exposure().

cpl_error_code muse_wcs_pixel_from_celestial ( cpl_propertylist *  aHeader,
double  aRA,
double  aDEC,
double *  aX,
double *  aY 
)

Convert from celestial spherical coordinates to pixel coordinates.

Parameters:
aHeader the input header containing the WCS of the exposure
aRA the right-ascension in degrees
aDEC the declination in degrees
aX the output horizontal pixel position
aY the output vertical pixel position
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader is used to do the transformation, only the gnomonic (TAN) is supported.

Uses Eqns (5), (12), (13), and (54) from Calabretta & Greisen 2002 A&A 395, 1077 (Paper II). We use that phi_p = 180 deg for zenithal projections (like TAN). The equations are actually implemented in muse_wcs_pixel_from_celestial_fast() that is called by this function.

Exceptions:
set and return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL
set and return CPL_ERROR_UNSUPPORTED_MODE aHeader does not contain gnomonic (TAN) WCS
set and return CPL_ERROR_SINGULAR_MATRIX, set contents of aX and aY to NAN determinant of the CDi_j matrix in aHeader is zero

Definition at line 1600 of file muse_wcs.c.

References cddet, crval2, muse_pfits_get_ctype(), muse_wcs_new(), and muse_wcs_pixel_from_celestial_fast().

Referenced by muse_wcs_solve().

static void muse_wcs_pixel_from_celestial_fast ( muse_wcs aWCS,
double  aRA,
double  aDEC,
double *  aX,
double *  aY 
) [static]

Convert from celestial spherical coordinates to pixel coordinates.

Parameters:
aWCS the input WCS structure
aRA the right-ascension in radians!
aDEC the declination in radians!
aX the output horizontal pixel position
aY the output vertical pixel position
Warning:
This function does not do any safety checks. If you are not sure, if the WCS does not contain a proper gnomonic (TAN) setup, use muse_wcs_pixel_from_celestial() instead!
Note:
This function is defined as static inline, i.e. should get integrated into the caller by the compiler, to maximise execution speed.

See muse_wcs_pixel_from_celestial() for the origin of the equations.

Definition at line 176 of file muse_wcs.h.

References cd22, cddet, crpix2, and crval2.

Referenced by muse_euro3dcube_collapse(), muse_pixgrid_create(), and muse_wcs_pixel_from_celestial().

cpl_error_code muse_wcs_pixel_from_projplane ( cpl_propertylist *  aHeader,
double  aX,
double  aY,
double *  aXOut,
double *  aYOut 
)

Convert from projection plane coordinates to pixel coordinates.

Parameters:
aHeader the input header containing the WCS of the exposure
aX the x position
aY the y position
aXOut the output horizontal pixel position
aYOut the output vertical pixel position
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader is used to do the transformation.

Note:
It is assumed that the header only contains a linear transformation, but no check is carried out for that.
Exceptions:
set and return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL
set and return CPL_ERROR_SINGULAR_MATRIX, set contents of aX and aY to NAN determinant of the CDi_j matrix in aHeader is zero

Definition at line 1736 of file muse_wcs.c.

References cddet, muse_wcs_new(), and muse_wcs_pixel_from_projplane_fast().

static void muse_wcs_pixel_from_projplane_fast ( muse_wcs aWCS,
double  aX,
double  aY,
double *  aXOut,
double *  aYOut 
) [static]

Convert from projection plane coordinates to pixel coordinates.

Parameters:
aWCS the input WCS structure
aX the x position
aY the y position
aXOut the output horizontal pixel position
aYOut the output vertical pixel position
Warning:
This function does not do any safety checks. If you are not sure, if the WCS does not contain a linear setup, use muse_wcs_pixel_from_projplane() instead!
Note:
This function is defined as static inline, i.e. should get integrated into the caller by the compiler, to maximise execution speed.

See muse_wcs_pixel_from_projplane() for the origin of the equations.

Definition at line 245 of file muse_wcs.h.

References cd22, cddet, crpix2, and crval2.

Referenced by muse_euro3dcube_collapse(), muse_pixgrid_create(), and muse_wcs_pixel_from_projplane().

cpl_error_code muse_wcs_position_celestial ( muse_pixtable aPixtable,
double  aRA,
double  aDEC 
)

Convert native to celestial spherical coordinates in a pixel table.

Parameters:
aPixtable the input pixel table containing the exposure to be converted
aRA the reference right-ascension in degrees
aDEC the reference declination in degrees
Returns:
CPL_ERROR_NONE for success, any other value for failure
Remarks:
The resulting WCS transform is directly applied and returned in the input pixel table, changing units of the spatial columns from "rad" to "deg". MUSE_PIXTABLE_XPOS is transformed to RA and MUSE_PIXTABLE_YPOS to DEC.
This function updates a FITS header (MUSE_HDR_PT_WCS) with a string value (MUSE_HDR_PT_WCS_POSI) to the pixel table, for information.

Loop through the input pixel table and evaluate the input WCS at every pixel to find the correct celestial coordinate for a given input position. This uses standard coordinate transforms as defined in Calabretta & Greisen 2002 A&A 395, 1077 (Paper II) for the exposure in the pixel table. The output contains a full spatial WCS with "celestial spherical" coordinates.

Exceptions:
return CPL_ERROR_NULL_INPUT the input pixel table or its header are NULL
return CPL_ERROR_INCOMPATIBLE_INPUT the input pixel table WCS type is not native spherical
return CPL_ERROR_UNSUPPORTED_MODE the input WCS is not for a gnomonic (TAN) projection

Definition at line 1463 of file muse_wcs.c.

References muse_pixtable::header, MUSE_HDR_PT_WCS, MUSE_HDR_PT_WCS_POSI, muse_pfits_get_ctype(), muse_pixtable_compute_limits(), muse_pixtable_get_nrow(), muse_pixtable_wcs_check(), MUSE_PIXTABLE_WCS_NATSPH, and muse_pixtable::table.

Referenced by muse_xcombine_tables().

cpl_error_code muse_wcs_project_tan ( muse_pixtable aPixtable,
const cpl_propertylist *  aWCS 
)

Carry out a gnomonic projection of a pixel table into native spherical coordinates.

Parameters:
aPixtable the input pixel table containing the exposure to be corrected
aWCS FITS headers containing the WCS solution
Returns:
CPL_ERROR_NONE for success, any other value for failure
Remarks:
The resulting WCS transform is directly applied and returned in the input pixel table, changing units of the spatial columns from "pix" to "rad". The table column MUSE_PIXTABLE_XPOS is transformed to phi and MUSE_PIXTABLE_YPOS to theta - pi/2.
This function adds a FITS header (MUSE_HDR_PT_WCS) with a string value (MUSE_HDR_PT_WCS_PROJ) to the pixel table, for information.

Loop through the input pixel table and evaluate the input WCS at every pixel to find the correct position for a given input position. This uses standard WCS transforms to produce "native spherical" coordinates in gnomonic (TAN) projection as defined in Calabretta & Greisen 2002 A&A 395, 1077 (Paper II) for the exposure in the pixel table.

Exceptions:
return CPL_ERROR_NULL_INPUT the input pixel table, its header, or the input WCS headers object are NULL
return CPL_ERROR_INCOMPATIBLE_INPUT the input pixel table WCS type is not pixel-based
return CPL_ERROR_UNSUPPORTED_MODE the input WCS is not for a gnomonic (TAN) projection

Definition at line 1323 of file muse_wcs.c.

References muse_pixtable::header, MUSE_HDR_PT_PREDAR_XHI, MUSE_HDR_PT_PREDAR_XLO, MUSE_HDR_PT_PREDAR_YHI, MUSE_HDR_PT_PREDAR_YLO, MUSE_HDR_PT_WCS, MUSE_HDR_PT_WCS_PROJ, MUSE_HDR_PT_XHI, MUSE_HDR_PT_XLO, MUSE_HDR_PT_YHI, MUSE_HDR_PT_YLO, muse_pfits_get_cd(), muse_pfits_get_crpix(), muse_pfits_get_ctype(), muse_pixtable_compute_limits(), muse_pixtable_get_nrow(), muse_pixtable_wcs_check(), MUSE_PIXTABLE_WCS_PIXEL, muse_wcs_apply_cd(), MUSE_WCS_KEYS, and muse_pixtable::table.

Referenced by muse_postproc_process_exposure().

cpl_error_code muse_wcs_projplane_from_celestial ( cpl_propertylist *  aHeader,
double  aRA,
double  aDEC,
double *  aX,
double *  aY 
)

Convert from celestial spherical coordinates to projection plane coordinates.

Parameters:
aHeader the input header containing the WCS of the exposure
aRA the right-ascension in degrees
aDEC the declination in degrees
aX the output horizontal proj. plane position
aY the output vertical proj. plane position
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader is used to do the transformation, only the gnomonic (TAN) is supported.

Uses Eqns (5), (12), (13), and (54) from Calabretta & Greisen 2002 A&A 395, 1077 (Paper II). We use that phi_p = 180 deg for zenithal projections (like TAN).

XXX this duplicates most of the code of muse_wcs_pixel_from_celestial()

Exceptions:
return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL
return CPL_ERROR_UNSUPPORTED_MODE aHeader does not contain gnomonic (TAN) WCS
return CPL_ERROR_ILLEGAL_INPUT determinant of the CDi_j matrix in aHeader is zero

Definition at line 1655 of file muse_wcs.c.

References muse_pfits_get_crval(), and muse_pfits_get_ctype().

cpl_error_code muse_wcs_projplane_from_pixel ( cpl_propertylist *  aHeader,
double  aX,
double  aY,
double *  aXOut,
double *  aYOut 
)

Convert from pixel coordinates to projection plane coordinates.

Parameters:
aHeader the input header containing the WCS of the exposure
aX the horizontal pixel position
aY the vertical pixel position
aXOut the output x position
aYOut the output y position
Returns:
CPL_ERROR_NONE for success, any other value for failure

The world coordinate system from aHeader is used to do the transformation.

Note:
It is assumed that the header only contains a linear transformation, but no check is carried out for that.
Exceptions:
return CPL_ERROR_NULL_INPUT aHeader, aX, and/or aY are NULL

Definition at line 1703 of file muse_wcs.c.

References muse_wcs_new(), and muse_wcs_projplane_from_pixel_fast().

static void muse_wcs_projplane_from_pixel_fast ( muse_wcs aWCS,
double  aX,
double  aY,
double *  aXOut,
double *  aYOut 
) [static]

Convert from pixel coordinates to projection plane coordinates.

Parameters:
aWCS the input WCS structure
aX the horizontal pixel position
aY the vertical pixel position
aXOut the output x position
aYOut the output y position
Warning:
This function does not do any safety checks. If you are not sure, if the WCS does not contain a linear setup, use muse_wcs_projplane_from_pixel() instead!
Note:
This function is defined as static inline, i.e. should get integrated into the caller by the compiler, to maximise execution speed.

See muse_wcs_projplane_from_pixel() for the origin of the equations.

Definition at line 216 of file muse_wcs.h.

References cd22, crpix2, and crval2.

Referenced by muse_resampling_collapse_pixgrid(), muse_resampling_euro3d(), and muse_wcs_projplane_from_pixel().

cpl_error_code muse_wcs_solve ( muse_wcs_object aWCSObj,
cpl_table *  aReference,
float  aRadius,
float  aFAccuracy,
int  aIter,
float  aThresh 
)

Find the world coordinate solution for a given field using coordinates of detected sources and coordinates of reference values.

Parameters:
aWCSObj the MUSE WCS object
aReference reference coordinates for the same field
aRadius initial radius for pattern matching in pixels
aFAccuracy factor to set initial accuracy relative to mean position accuracy in aWCSObj->detected
aIter number of iterations
aThresh the rejection threshold
Returns:
CPL_ERROR_NONE on success, another CPL error code on failure
Remarks:
The propertylist returned in aWCSObj->wcs contains the world coordinate solution.
Both input object tables (aWCSObj->detected and aReference) are modified in that they are sorted by increasing source flux.

Identify reference objects for the detected objects using pattern matching. Start using a search radius of aRadius pixels, and iteratively decrease it, until no duplicate detections are identified any more. Similarly, iterate the data accuracy (decrease it downwards from the mean positioning error aFAccuracy) until matches are found. Remove the remaining unidentified objects.

The final values used for relative and absolute accuracy and the radius are recorded in the header keywords MUSE_HDR_WCS_RADIUS, MUSE_HDR_WCS_ACCURACY, and MUSE_HDR_WCS_FACCURACY.

Derive the world coordinate plate solution using cpl_wcs_platesol() (taking into account shift, 2 scales, and rotation) and add it as cpl_propertylist element wcs to aWCSObj.

Exceptions:
return CPL_ERROR_NULL_INPUT input WCS object is NULL
return CPL_ERROR_ILLEGAL_INPUT table of detected or reference coordinates are empty or NULL
return CPL_ERROR_BAD_FILE_FORMAT table of detected or reference coordinates don't contain the expected columns
return CPL_ERROR_ILLEGAL_INPUT aRadius and/or aFAccuracy are non-positive
return CPL_ERROR_DATA_NOT_FOUND no objects could be identified
propagate return code cpl_wcs_platesol failed

Definition at line 668 of file muse_wcs.c.

References muse_datacube::header, muse_cplarray_has_duplicate(), muse_cpltable_check(), muse_pfits_get_cd(), muse_wcs_create_default(), muse_wcs_get_angles(), muse_wcs_get_scales(), and muse_wcs_pixel_from_celestial().

Referenced by muse_postproc_process_exposure(), and muse_wcs_optimize_solution().


Variable Documentation

Initial value:
 {
  { "XPOS", CPL_TYPE_DOUBLE, "pix", "%f", "horizontal position", CPL_TRUE },
  { "YPOS", CPL_TYPE_DOUBLE, "pix", "%f", "vertical position", CPL_TRUE },
  { "XERR", CPL_TYPE_DOUBLE, "pix", "%f",
    "error estimate of the horizontal position", CPL_TRUE },
  { "YERR", CPL_TYPE_DOUBLE, "pix", "%f",
    "error estimate of the vertical position", CPL_TRUE },
  { "FLUX", CPL_TYPE_DOUBLE, "count", "%e", "(relative) flux of the source", CPL_TRUE },
  { "XFWHM", CPL_TYPE_DOUBLE, "pix", "%f", "horizontal FWHM", CPL_TRUE },
  { "YFWHM", CPL_TYPE_DOUBLE, "pix", "%f", "vertical FWHM", CPL_TRUE },
  { NULL, 0, NULL, NULL, NULL, CPL_FALSE }
}

Definition of the table structure for the astrometric field detections.

  • XPOS: horizontal position
  • YPOS: vertical position
  • XERR: error estimate of the horizontal position
  • YERR: error estimate of the vertical position
  • FLUX: (relative) flux of the source
  • XFWHM: horizontal FWHM
  • YFWHM: vertical FWHM

Definition at line 116 of file muse_wcs.c.

Definition of the table structure for the astrometric reference sources.

  • SourceID: the source identification
  • RA: right ascension (in decimal degrees)
  • DEC: declination (in decimal degrees)
  • filter: the filter used for the magnitude
  • mag: the object (Vega) magnitude

Definition at line 141 of file muse_wcs.c.

Initial value:
 {
  { "SourceID", CPL_TYPE_STRING, "", "%s", "the source identification", CPL_TRUE },
  { "RA", CPL_TYPE_DOUBLE, "deg", "%f", "right ascension (decimal degrees)", CPL_TRUE },
  { "DEC", CPL_TYPE_DOUBLE, "deg", "%f", "declination (decimal degrees)", CPL_TRUE },
  { "filter", CPL_TYPE_STRING, "", "%s", "the filter used for the magnitude", CPL_TRUE },
  { "mag", CPL_TYPE_DOUBLE, "mag", "%f", "the object (Vega) magnitude", CPL_TRUE },
  { NULL, 0, NULL, NULL, NULL, CPL_FALSE }
}

Definition of the table structure for the astrometric reference sources.

  • SourceID: the source identification
  • RA: right ascension (in decimal degrees)
  • DEC: declination (in decimal degrees)
  • filter: the filter used for the magnitude
  • mag: the object (Vega) magnitude

Definition at line 141 of file muse_wcs.c.


Generated on 26 Jan 2017 for MUSE Pipeline Reference Manual by  doxygen 1.6.1