Working With Regions¶
Starlink AST is primarily a library for working with world coordinate systems (WCS) for astronomical data. It has extensive functionality related to coordinate system translations, platting, and more. One of the primary tools of the library (and the initial focus of the Cornish Python interface) is the handling of regions. Regions can be Cartesian, e.g. a pixel grid on a CCD, or else in a celestial sphere frame, e.g. a region on an ASTSkyFrame
. For the latter, all lines connecting vertices are great circles.
All region objects are all subclassed from ASTRegion
. See the C library reference documentation for the full details of the objects. Region classes include:
When creating a region, note that the frame the region to be defined in must be specified. The class ASTICRSFrame
is provided as a convenience, defined as sky frame in the ICRS system with epoch 2000.0 and equinox 2000.0. Note that regions default to the ASTICRSFrame
if one is not provided.
Region objects are defined in the top level of the cornish
namespace.
All Regions¶
The ASTRegion
class is the superclass for all Cornish region objects. It contains a great deal of functionality. The API reference is a worth looking at to see what is available. A short listing:
points
: Return a list of points that describe the region. The shape is dependent on the type of region.boundingCircle()
: Returns a newASTCircle
object that bounds the given region.overlaps()
: Check whether two regions overlap.isIdenticalTo()
: Check whether two regions are identical.isFullyWithin()
: Check whether one region is fully within another.fullyEncloses()
: Check whether one region is fully encloses another.boundaryPointMesh()
: Returns an array of evenly distributed points that cover the boundary of the region.interiorPointMesh()
: Returns an array of evenly distributed points that cover the surface of the region.containsPoint()
: Determine whether a given point lies within a given region.
See the API reference for more methods and properties.
Circles¶
Creating Circles¶
Circles can be defined as either a center point and a radius or else a center point and another on the circumference. Coordinates can be specified an astropy.coodinates.SkyCoord
object or pairs of values in degrees. The examples below demonstrate various ways to create circle regions.
from cornish import ASTCircle, ASTICRSFrame, ASTSkyFrame
from cornish.constants import SYSTEM_GALACTIC, EQUINOX_J2010
from astropy.coordinates import SkyCoord
import astropy.units as u
# note that the default frame is ICRS, epoch=2000.0, equinox=2000.0
# defined as center + radius
# --------------------------
# using Astropy objects
center = SkyCoord(ra="12d42m22s", dec="-32d18m58s")
circle = ASTCircle(center=center, radius=2.0*u.deg)
# using float values, defaults to degrees
circle = ASTCircle(center=[12.7061, -31.6839], radius=2.0) # assumes degrees
circle = ASTCircle(center=[12.7061*u.deg, -31.6839*u.deg], radius=2.0*u.deg) # Quanitites also accepted
# defined as center + circumference point
# ---------------------------------------
circle = ASTCircle(center=center, edge_point=[12.7061, -32.6839]) # edge_point also takes SkyCoord
# define the circle in another frame
# ----------------------------------
gal_frame = ASTSkyFrame(system=SYSTEM_GALACTIC)
gal_frame.equinox = EQUINOX_J2010
ASTCircle(frame=gal_frame, center=center, radius=2.0*u.deg)
Circle Properties¶
Circles have radius
and centre
properties as one might expect, and both can be directly modified:
circle.radius
>>> <Quantity 2. deg>
circle.centre # or "center" if you prefer...
>>> array([ 12.70611111, -32.31611111]) # output in degrees
New circles can be created by a scale factor or increased by addition from an existing circle.
scaled_circle = circle * 2.0
scaled_circle.radius
>>> <Quantity 4. deg>
larger_circle = circle + 6*u.deg
larger_circle.radius
>>> <Quantity 8. deg>
Converting to Polygons¶
For code that requires a polygon region as an input the method toPolygon()
will convert a circle to an ASTPolygon
. The default is to sample 200 points for the polygon but this can be customized by using the npoints parameter (often even 20 are sufficient). Note that all of the polygon points fall on the circle’s circumference, so the resulting region is fully inscribed by the original circle.
polygon = circle.toPolygon()
finer_polygon = circle(toPolygon(npoints=200))
All regions have a boundingCircle()
property that returns an ASTCircle
that bounds the region. In the case of ASTCircle
objects, this method returns the original circle.
Polygons¶
A polygon is a collection of vertices that lie in a specific frame. The default frame ASTICRSFrame
is used if none is specified.
from cornish import ASTPolygon, ASTICRSFrame
import numpy as np
points = np.array([[ 12.70611111, -30.31611111],
[ 13.42262189, -30.41196836],
[ 14.07300863, -30.69069244],
[ 14.59623325, -31.12642801],
[ 14.94134955, -31.67835614],
[ 15.07227821, -32.29403528],
[ 14.97204342, -32.91392471],
[ 14.6459242 , -33.47688136],
[ 14.12273328, -33.92626054],
[ 13.4533703 , -34.21603194],
[ 12.70611111, -34.31611111],
[ 11.95885193, -34.21603194],
[ 11.28948894, -33.92626054],
[ 10.76629802, -33.47688136],
[ 10.4401788 , -32.91392471],
[ 10.33994401, -32.29403528],
[ 10.47087267, -31.67835614],
[ 10.81598897, -31.12642801],
[ 11.3392136 , -30.69069244],
[ 11.98960033, -30.41196836]])
polygon = ASTPolygon(frame=ASTICRSFrame(), points=points)
Points can be specified as an array of coordinate points (as above) or as parallel arrays of each dimension (which is just points.T
from above):
points = np.array([[ 12.70611111, 13.42262189, 14.07300863, 14.59623325, ...],
[-30.31611111, -30.41196836, -30.69069244, -31.12642801, ...]])
Todo
Provide example of how to convert a region from one frame to another.
Boxes¶
Todo
Box section coming soon! (But it’s pretty straightforward from the ASTBox
API.)
Compound Regions¶
Todo
Compound regions section coming soon! (But it’s pretty straightforward from the ASTBox
API.)
From FITS Files¶
Cornish is able to create regions based on image FITS headers alone. The example below shows how to create a region object based on the area covered by a FITS image from the header. The example file below can be downloaded here.
from cornish import ASTPolygon
from astropy.io import fits
filename = "frame-g-006174-2-0094.fits.bz2"
with fits.open(filename) as hdu_list:
hdu1 = hdu_list[0]
polygon = ASTPolygon.fromFITSHeader(hdu1.header)