File SkyWcs.h¶
-
namespace
lsst
Class for a simple mapping implementing a generic AstrometryTransform.
Remove all non-astronomical counts from the Chunk Exposure’s pixels.
Forward declarations for lsst::utils::Cache
For details on the Cache class, see the Cache.h file.
It uses a template rather than a pointer so that the derived classes can use the specifics of the transform. The class simplePolyMapping overloads a few routines.
A base class for image defects
Numeric constants used by the Integrate.h integrator routines.
Compute Image Statistics
- Note
Gauss-Kronrod-Patterson quadrature coefficients for use in quadpack routine qng. These coefficients were calculated with 101 decimal digit arithmetic by L. W. Fullerton, Bell Labs, Nov 1981.
- Note
The Statistics class itself can only handle lsst::afw::image::MaskedImage() types. The philosophy has been to handle other types by making them look like lsst::afw::image::MaskedImage() and reusing that code. Users should have no need to instantiate a Statistics object directly, but should use the overloaded makeStatistics() factory functions.
-
namespace
afw
-
namespace
geom
Functions
-
Eigen::Matrix2d
makeCdMatrix
(lsst::geom::Angle const &scale, lsst::geom::Angle const &orientation = 0 * lsst::geom::degrees, bool flipX = false)¶ Make a WCS CD matrix
- Return
the CD matrix, where element (i-1, j-1) corresponds to FITS keyword CDi_j and i, j have range [1, 2]
- Parameters
[in] scale
: Pixel scale as an angle on sky/pixels[in] orientation
: Position angle of pixel +Y, measured from N through E. At 0 degrees, +Y is along N and +X is along W/E if flipX false/true At 90 degrees, +Y is along E and +X is along N/S if flipX false/true[in] flipX
: Flip x axis? See orientation for details.
-
std::shared_ptr<SkyWcs>
makeFlippedWcs
(SkyWcs const &wcs, bool flipLR, bool flipTB, lsst::geom::Point2D const ¢er)¶ Return a copy of a FITS-WCS with pixel positions flipped around a specified center
- Parameters
[in] wcs
: The initial WCS[in] flipLR
: Flip pixel positions left/right aboutcenter
[in] flipTB
: Flip pixel positions top/bottom aboutcenter
[in] center
: Center pixel position of flip
-
std::shared_ptr<SkyWcs>
makeModifiedWcs
(TransformPoint2ToPoint2 const &pixelTransform, SkyWcs const &wcs, bool modifyActualPixels)¶ Create a new SkyWcs whose pixels are transformed by pixelTransform, as described below.
If modifyActualPixels is true and the cameraGeom::ACTUAL_PIXELS frame exists then pixelTransform is inserted just after the cameraGeom::ACTUAL_PIXELS frame:
newActualPixelsToPixels = pixelTransform -> oldActualPixelsToPixels
This is appropriate for shifting a WCS, e.g. when writing FITS metadata for a subimage.
If modifyActualPixels is false or the cameraGeom::ACTUAL_PIXELS frame does not exist then pixelTransform is inserted just after the cameraGeom::PIXELS frame:
newPixelsToIwc = pixelTransform -> oldPixelsToIwc
This is appropriate for inserting a model for optical distortion.
Other than the change described above, the new SkyWcs will be just like the old one.
- Return
the new WCS
- Parameters
[in] pixelTransform
: Transform to insert[in] wcs
: Input WCS[in] modifyActualPixels
: Location at which to insert the transform; if true and the cameraGeom::ACTUAL_PIXELS frame is present then insert just after the cameraGeom::ACTUAL_PIXELS frame, else insert just after the cameraGeom::PIXELS frame.
-
std::shared_ptr<SkyWcs>
makeSkyWcs
(daf::base::PropertySet &metadata, bool strip = false)¶ Construct a SkyWcs from FITS keywords
This function is preferred over calling the SkyWcs metadata constructor directly because it allows us to change SkyWcs to an abstract base class in the future, without affecting code that constructs a WCS from FITS metadata.
- Parameters
[in] metadata
: FITS header metadata[in] strip
: If true: strip items frommetadata
used to create the WCS, such as RADESYS, EQUINOX, CTYPE12, CRPIX12, CRVAL12, etc. Always keep keywords that might be wanted for other purpposes, including NAXIS12 and date-related keywords such as “DATE-OBS” and “TIMESYS” (but not “EQUINOX”).
- Exceptions
lsst::pex::exceptions::TypeError
: if the metadata does not describe a celestial WCS.
-
std::shared_ptr<SkyWcs>
makeSkyWcs
(lsst::geom::Point2D const &crpix, lsst::geom::SpherePoint const &crval, Eigen::Matrix2d const &cdMatrix, std::string const &projection = "TAN")¶ Construct a simple FITS SkyWcs with no distortion
- Parameters
[in] crpix
: Center of projection on the CCD using the LSST convention: 0, 0 is the lower left pixel of the parent image[in] crval
: Center of projection on the sky[in] cdMatrix
: CD matrix, where element (i-1, j-1) corresponds to FITS keyword CDi_j and i, j have range [1, 2]. May be computed by calling makeCdMatrix.[in] projection
: The name of the FITS WCS projection, e.g. “TAN” or “STG”
-
std::shared_ptr<SkyWcs>
makeSkyWcs
(TransformPoint2ToPoint2 const &pixelsToFieldAngle, lsst::geom::Angle const &orientation, bool flipX, lsst::geom::SpherePoint const &boresight, std::string const &projection = "TAN")¶ Construct a FITS SkyWcs from camera geometry
- Return
a SkyWcs whose sky origin is the boresight and pixel origin is focal plane (0, 0).
- Note
Unlike makeCdMatrix,
orientation
is with respect to the focal plane axes, not the CCD axes. This is because field angle is defined with respect to focal plane axes.- Parameters
[in] pixelsToFieldAngle
: Transformation from pixels to field angle (in radians).[in] orientation
: Position angle of focal plane +Y, measured from N through E at crval. At 0 degrees, +Y is along N and +X is along W/E if flipX false/true. At 90 degrees, +Y is along E and +X is along N/S if flipX false/true.[in] flipX
: Flip x axis? See orientation for details.[in] boresight
: ICRS sky position at the boresight (field angle (0, 0)).[in] projection
: The name of the FITS WCS projection, e.g. “TAN” or “STG”.
-
std::shared_ptr<SkyWcs>
makeTanSipWcs
(lsst::geom::Point2D const &crpix, lsst::geom::SpherePoint const &crval, Eigen::Matrix2d const &cdMatrix, Eigen::MatrixXd const &sipA, Eigen::MatrixXd const &sipB)¶ Construct a TAN-SIP SkyWcs with forward SIP distortion terms and an iterative inverse.
- Parameters
[in] crpix
: Center of projection on the CCD using the LSST convention: 0, 0 is the lower left pixel of the parent image[in] crval
: Center of projection on the sky[in] cdMatrix
: CD matrix, where element (i-1, j-1) corresponds to FITS keyword CDi_j and i, j have range [1, 2]. May be computed by calling makeCdMatrix.[in] sipA
: Forward distortion matrix for axis 1[in] sipB
: Forward distortion matrix for axis 2
-
std::shared_ptr<SkyWcs>
makeTanSipWcs
(lsst::geom::Point2D const &crpix, lsst::geom::SpherePoint const &crval, Eigen::Matrix2d const &cdMatrix, Eigen::MatrixXd const &sipA, Eigen::MatrixXd const &sipB, Eigen::MatrixXd const &sipAp, Eigen::MatrixXd const &sipBp)¶ Construct a TAN WCS with forward and inverse SIP distortion terms.
- Parameters
[in] crpix
: Center of projection on the CCD using the LSST convention: 0, 0 is the lower left pixel of the parent image[in] crval
: Center of projection on the sky[in] cdMatrix
: CD matrix, where element (i-1, j-1) corresponds to FITS keyword CDi_j and i, j have range [1, 2]. May be computed by calling makeCdMatrix.[in] sipA
: Forward distortion matrix for axis 1[in] sipB
: Forward distortion matrix for axis 2[in] sipAp
: Reverse distortion matrix for axis 1[in] sipBp
: Reverse distortion matrix for axis 2
-
std::shared_ptr<TransformPoint2ToPoint2>
makeWcsPairTransform
(SkyWcs const &src, SkyWcs const &dst)¶ A Transform obtained by putting two SkyWcs objects “back to back”.
Provides basic exception safety.
- Return
a Transform whose forward transformation converts from
src
pixels todst
pixels, and whose inverse transformation converts in the opposite direction.- Parameters
src
: the WCS for the source pixelsdst
: the WCS for the destination pixels
-
std::shared_ptr<TransformPoint2ToSpherePoint>
getIntermediateWorldCoordsToSky
(SkyWcs const &wcs, bool simplify = true)¶ Return a transform from intermediate world coordinates to sky
-
std::shared_ptr<TransformPoint2ToPoint2>
getPixelToIntermediateWorldCoords
(SkyWcs const &wcs, bool simplify = true)¶ Return a transform from pixel coordinates to intermediate world coordinates
The pixel frame is is the base frame: cameraGeom::ACTUAL_PIXELS, if present, else cameraGeom::PIXELS.
-
class
SkyWcs
: public lsst::afw::table::io::PersistableFacade<SkyWcs>, public Storable - #include <SkyWcs.h>
A 2-dimensional celestial WCS that transform pixels to ICRS RA/Dec, using the LSST standard for pixels.
SkyWcs is an immutable object that can not only represent any standard FITS WCS, but can also contain arbitrary Transforms, e.g. to model optical distortion or pixel imperfections.
In order to make a SkyWcs that models optical distortion, say, it is usually simplest to start with a standard FITS WCS (such as a TAN WCS) as an approximation, then insert a transform that models optical distortion by calling makeModifiedWcs. However, it is also possible to build a SkyWcs entirely from transforms, if you prefer, by building an ast::FrameDict and constructing the SkyWcs from that.
Frames of reference
SkyWcs internally keeps track of the following frames of reference:
cameraGeom::ACTUAL_PIXELS (optional): “actual” pixel position using the LSST standard. This has the same meaning as the cameraGeom::ACTUAL_PIXELS frame: actual pixels include effects such as “tree ring” distortions and electrical effects at the edge of CCDs. This frame should only be provided if there is a reasonable model for these imperfections.
cameraGeom::PIXELS: nominal pixel position, using the LSST standard. This has the same meaning as the cameraGeom::PIXELS frame: nominal pixels may be rectangular, but are uniform in size and spacing.
IWC: intermediate world coordinates (the FITS WCS concept).
SKY: position on the sky as ICRS, with standard RA, Dec axis order.
Pixel position standards
The LSST standard for pixel position is: 0,0 is the center of the lower left image pixel. The FITS standard for pixel position is: 1,1 is the center of the lower left image pixel.
LSST and FITS also use a different origin for subimages:
LSST pixel position is in the frame of reference of the parent image
FITS pixel position is in the frame of reference of the subimage However, SkyWcs does not keep track of this difference. Code that persists and unpersists SkyWcs using FITS-WCS header cards must handle the offset, e.g. by calling
copyAtShiftedPixelOrigin
Internal details: the contained ast::FrameDict
SkyWcs contains an ast::FrameDict which transforms from pixels to sky (in radians) in the forward direction.
This FrameDict contains the named frames described in frames of reference, e.g. “SKY”, “IWC”, cameraGeom::PIXELS and possibly cameraGeom::ACTUAL_PIXELS. “SKY” is the current frame. If cameraGeom::ACTUAL_PIXELS is present then it is the base frame, otherwise cameraGeom::PIXELS is the base frame.
The “SKY” frame is of type ast::SkyFrame and has the following attributes:
SkyRef
is set to the sky origin of the WCS (ICRS RA, Dec) in radians.SkyRefIs
is set to “Ignored” so that SkyRef is not used in transformations.
The other frames are of type ast::Frame and have 2 axes.
Unnamed Group
-
lsst::geom::SpherePoint
pixelToSky
(lsst::geom::Point2D const &pixel) const Compute sky position(s) from pixel position(s)
-
lsst::geom::SpherePoint
pixelToSky
(double x, double y) const
Unnamed Group
Public Functions
-
SkyWcs
(SkyWcs const&)
-
SkyWcs
(SkyWcs&&)
-
~SkyWcs
()
-
bool
operator==
(SkyWcs const &other) const Equality is based on the string representations being equal
Two SkyWcs constructed the same way will be equal, and a SkyWcs that has been saved and restored will be equal to the original. However, it is possible to construct two SkyWcs that behave identically as far as transforming points go, but will compare as unequal due to subtle internal differences, such as a contained ast::Mapping that has a different ID in one SkyWcs than another.
Thus equality is primarily useful for testing persistence.
-
bool
operator!=
(SkyWcs const &other) const
-
SkyWcs
(daf::base::PropertySet &metadata, bool strip = false) Construct a SkyWcs from FITS keywords
- Parameters
[in] metadata
: FITS header metadata[in] strip
: If true: strip items frommetadata
used to create the WCS, such as RADESYS, EQUINOX, CTYPE12, CRPIX12, CRVAL12, etc. Always keep keywords that might be wanted for other purpposes, including NAXIS12 and date-related keywords such as “DATE-OBS” and “TIMESYS” (but not “EQUINOX”).
- Exceptions
lsst::pex::exceptions::TypeError
: if the metadata does not describe a celestial WCS.
-
SkyWcs
(ast::FrameDict const &frameDict) Construct a SkyWcs from an ast::FrameDict
This is the most general constructor; it can be used to define any celestial WCS. Note that in many cases the result will not be exactly representable as a FITS WCS.
- Parameters
[in] frameDict
: An ast::FrameDict that describes the transformation from pixels to sky. It must meet the requirements outlined in the contained ast::FrameDict.
- Exceptions
lsst::pex::exceptions::TypeError
: ifframeDict
is missing any of the required frames of reference.
-
std::shared_ptr<SkyWcs>
copyAtShiftedPixelOrigin
(lsst::geom::Extent2D const &shift) const Return a copy of this SkyWcs with the pixel origin shifted by the specified amount.
new pixel origin = the old pixel origin + shift
- Parameters
[in] shift
: The amount by which to shift the pixel origin (pixels)
-
std::shared_ptr<daf::base::PropertyList>
getFitsMetadata
(bool precise = false) const Return the WCS as FITS WCS metadata
FITS representations of WCS are described in “Representations of World Coordinates in FITS” by Greisen and Calabretta and several related papers.
- Parameters
[in] precise
: Fail if the WCS cannot be accurately represented as FITS metadata? If False then return an approximation. For now that approximation is pure TAN but as of DM-13170 it will be a fit TAN-SIP. The required precision is set by constant TIGHT_FITS_TOL in SkyWcs.cc
The required precision is hard-coded as constant TIGHT_FITS_TOL in SkyWcs.cc
- Exceptions
lsst::pex::exceptions::RuntimeError
: if precise is true and AST cannot represent this WCS as a FITS WCS to sufficient precision.
-
lsst::geom::Angle
getPixelScale
(lsst::geom::Point2D const &pixel) const Get the pixel scale at the specified pixel position
The scale is the square root of the area of the specified pixel on the sky.
- Warning
Unlike getPixelScale() the value is not cached, even if pixel = pixel origin.
-
lsst::geom::Angle
getPixelScale
() const Get the pixel scale at the pixel origin
The scale is the square root of the area of the specified pixel on the sky.
The value is cached, so this is a cheap call.
-
lsst::geom::Point2D
getPixelOrigin
() const Get the pixel origin, in pixels, using the LSST convention
This is CRPIX1 - 1, CRPIX2 -1 in FITS terminology
-
lsst::geom::SpherePoint
getSkyOrigin
() const Get the sky origin, the celestial fiducial point
This is CRVAL1, CRVAL2 in FITS terminology
-
Eigen::Matrix2d
getCdMatrix
(lsst::geom::Point2D const &pixel) const Get the 2x2 CD matrix at the specified pixel position
The elements are in degrees
-
Eigen::Matrix2d
getCdMatrix
() const Get the 2x2 CD matrix at the pixel origin
The elements are in degrees
-
std::shared_ptr<SkyWcs>
getTanWcs
(lsst::geom::Point2D const &pixel) const Get a local TAN WCS approximation to this WCS at the specified pixel position
-
std::shared_ptr<const ast::FrameDict>
getFrameDict
() const Get the contained FrameDict
The base frame will be cameraGeom::PIXELS or cameraGeom::ACTUAL_PIXELS and the current frame will be SKY, so the forward transform goes from pixels (using the LSST zero convention) to sky ICRS RA, Dec (in radians). For more details see the contained ast::FrameDict
-
std::shared_ptr<const TransformPoint2ToSpherePoint>
getTransform
() const Get a TransformPoint2ToSpherePoint that transforms pixels to sky in the forward direction and sky to pixels in the inverse direction.
-
bool
isFlipped
() const Does the WCS follow the convention of North=Up, East=Left?
- Return
False/true if E is along -X/+X when the N/E axes are rotated so that N is along image +Y.
- Exceptions
lsst::pex::exceptions::RuntimeError
: if the parity cannot be determined because the CD matrix is singular.
-
bool
isPersistable
() const
-
lsst::geom::AffineTransform
linearizePixelToSky
(lsst::geom::SpherePoint const &coord, lsst::geom::AngleUnit const &skyUnit) const Return the local linear approximation to pixelToSky at a point given in sky coordinates.
The local linear approximation is defined such the following is true (ignoring floating-point errors):
wcs.linearizePixelToSky(sky, skyUnit)(wcs.skyToPixel(sky)) == sky.getPosition(skyUnit);
(lsst::geom::AffineTransform::operator() applies the transform in the forward direction)
- Parameters
[in] coord
: Position in sky coordinates where transform is desired.[in] skyUnit
: Units to use for sky coordinates; units of matrix elements will be skyUnits/pixel.
-
lsst::geom::AffineTransform
linearizePixelToSky
(lsst::geom::Point2D const &pixel, lsst::geom::AngleUnit const &skyUnit) const Return the local linear approximation to pixelToSky at a point given in pixel coordinates.
The local linear approximation is defined such the following is true (ignoring floating-point errors):
wcs.linearizePixelToSky(pixel, skyUnit)(pixel) == wcs.pixelToSky(pixel).getPosition(skyUnit)
(lsst::geom::AffineTransform::operator() applies the transform in the forward direction)
- Parameters
[in] pixel
: Position in pixel coordinates where transform is desired.[in] skyUnit
: Units to use for sky coordinates; units of matrix elements will be skyUnits/pixel.
-
lsst::geom::AffineTransform
linearizeSkyToPixel
(lsst::geom::SpherePoint const &coord, lsst::geom::AngleUnit const &skyUnit) const Return the local linear approximation to skyToPixel at a point given in sky coordinates.
The local linear approximation is defined such the following is true (ignoring floating-point errors):
wcs.linearizeSkyToPixel(sky, skyUnit)(sky.getPosition(skyUnit)) == wcs.skyToPixel(sky)
(lsst::geom::AffineTransform::operator() applies the transform in the forward direction)
- Parameters
[in] coord
: Position in sky coordinates where transform is desired.[in] skyUnit
: Units to use for sky coordinates; units of matrix elements will be pixels/skyUnit.
-
lsst::geom::AffineTransform
linearizeSkyToPixel
(lsst::geom::Point2D const &pixel, lsst::geom::AngleUnit const &skyUnit) const Return the local linear approximation to skyToPixel at a point given in pixel coordinates.
The local linear approximation is defined such the following is true (ignoring floating-point errors):
wcs.linearizeSkyToPixel(pixel, skyUnit)(wcs.pixelToSky(pixel).getPosition(skyUnit)) == pixel
(lsst::geom::AffineTransform::operator() applies the transform in the forward direction)
- Parameters
[in] pixel
: Position in pixel coordinates where transform is desired.[in] skyUnit
: Units to use for sky coordinates; units of matrix elements will be pixels/skyUnit.
-
bool
hasFitsApproximation
() const Does this SkyWcs have an approximate SkyWcs that can be represented as standard FITS WCS?
This feature is not yet implemented, so hasFitsApproximation is always false.
-
bool
isFits
() const Return true getFitsMetadata(true) will succeed, false if not.
In other words, true indicates that the WCS can be accurately represented using FITS WCS metadata.
-
void
writeStream
(std::ostream &os) const Serialize this SkyWcs to an output stream
Version 1 format is as follows:
The version number (an integer)
A space
getShortClassName()
A space
A space
The contained ast::FrameDict written using FrameDict.show(os, false)
If and when fits approximation is supported, the approximate WCS will be written as a second FrameDict immediately following the first.
- Parameters
[out] os
: output stream to which to serialize this SkyWcs
-
std::string
writeString
() const Serialize this SkyWcs to a string, using the same format as writeStream.
-
std::shared_ptr<typehandling::Storable>
cloneStorable
() const Create a new SkyWcs that is a copy of this one.
-
std::string
toString
() const Create a string representation of this object.
An example of the output might look like this (note the newlines):
FITS standard SkyWcs: Sky Origin: (0.000000, +45.000000) Pixel Origin: (100, 100) Pixel Scale: 1 arcsec/pixel
-
bool
equals
(typehandling::Storable const &other) const Compare this object to another Storable.
- Return
*this == other
ifother
is a SkyWcs; otherwisefalse
.
Public Static Functions
-
static std::string
getShortClassName
()
Protected Functions
-
std::string
getPersistenceName
() const
-
std::string
getPythonModule
() const
-
void
write
(OutputArchiveHandle &handle) const
Private Functions
-
lsst::geom::AffineTransform
_linearizePixelToSky
(lsst::geom::Point2D const &pix, lsst::geom::SpherePoint const &coord, lsst::geom::AngleUnit const &skyUnit) const¶
-
lsst::geom::AffineTransform
_linearizeSkyToPixel
(lsst::geom::Point2D const &pix, lsst::geom::SpherePoint const &coord, lsst::geom::AngleUnit const &skyUnit) const¶
-
void
_computeCache
()¶ Compute _pixelOrigin and _pixelScaleAtOrigin.
-
Eigen::Matrix2d
-
namespace