SunPy coordinates

The SunPy coordinates submodule is an implementation of the common solar physics coordinate frames using the Astropy coordinates framework.

Getting Started

The easiest interface to the coordinates module is through the SkyCoord class:

>>> import astropy.units as u
>>> from astropy.coordinates import SkyCoord
>>> from sunpy.coordinates import frames

>>> c = SkyCoord(-100*u.arcsec, 500*u.arcsec, frame=frames.Helioprojective)
>>> c = SkyCoord(x=-72241.0*u.km, y=361206.1*u.km, z=589951.4*u.km, frame=frames.Heliocentric)
>>> c = SkyCoord(70*u.deg, -30*u.deg, frame=frames.HeliographicStonyhurst)
>>> c
<SkyCoord (HeliographicStonyhurst: obstime=None): (lon, lat, rad) in (deg, deg, km)
    (70.0, -30.0, 695508.0)>

It is also possible to use strings to define the frame but in that case make sure to explicitly import sunpy.coordinates as it registers solar coordinate frames with astropy coordinates.:

>>> import astropy.units as u
>>> from astropy.coordinates import SkyCoord

>>> import sunpy.coordinates
>>> c = SkyCoord(-100*u.arcsec, 500*u.arcsec, frame='helioprojective')
>>> c
<SkyCoord (Helioprojective: D0=149597870.7 km, obstime=None, L0=0.0 deg, B0=0.0 deg, rsun=695508.0 km): (Tx, Ty) in arcsec
  (-100.,  500.)>

SunPy implements support for the following solar physics coordinate systems:

for a complete description of these frames see sunpy.coordinates.frames, for a more detailed description of the frames see Thompson (2006)

SkyCoord and all other coordinates objects also support array coordinates. These work the same as single-value coordinates, but they store multiple coordinates in a single object. When you’re going to apply the same operation to many different coordinates, this is a better choice than a list of SkyCoord objects, because it will be much faster than applying the operation to each SkyCoord in a for loop.

>>> c = SkyCoord([-500, 400]*u.arcsec, [100, 200]*u.arcsec, frame=frames.Helioprojective)
>>> c
<SkyCoord (Helioprojective: obstime=None, rsun=695508.0 km, observer=<HeliographicStonyhurst Coordinate (obstime=None): (lon, lat, radius) in (deg, deg, AU)
    ( 0.,  0.,  1.)>): (Tx, Ty) in arcsec
    [(-500.,  100.), ( 400.,  200.)]>
>>> c[0]
<SkyCoord (Helioprojective: obstime=None, rsun=695508.0 km, observer=<HeliographicStonyhurst Coordinate (obstime=None): (lon, lat, radius) in (deg, deg, AU)
    ( 0.,  0.,  1.)>): (Tx, Ty) in arcsec
    (-500.,  100.)>

Accessing Coordinates

Individual coordinates can be accessed via attributes on the SkyCoord object, but the names of the components of the coordinates for each frame differ. For a full description of all the properties of the frames see sunpy.coordinates.frames.

Helioprojective

For the helioprojective frame the coordinates are accessed as Tx and Ty representing theta x and y. These are the same coordinates that are often referred to as ‘solar-x’ and ‘solar-y’.:

>>> c = SkyCoord(-500*u.arcsec, 100*u.arcsec, frame=frames.Helioprojective)
>>> c.Tx
<Longitude180 -500.0 arcsec>
>>> c.Ty
<Latitude 100.0 arcsec>

Heliocentric

Heliocentric normally a Cartesian frame so the coordinates are accessed as x, y, z:

>>> c = SkyCoord(-72241.0*u.km, 361206.1*u.km, 589951.4*u.km, frame=frames.Heliocentric)
>>> c.x
<Quantity -72241.0 km>
>>> c.y
<Quantity 361206.1 km>
>>> c.z
<Quantity 589951.4 km>

HeliographicStonyhurst and HeliographicCarrington

Both the heliographic frames use latitude, longitude and radius which are accessed as follows:

>>> c = SkyCoord(70*u.deg, -30*u.deg, frame=frames.HeliographicStonyhurst)
>>> c.lat
<Latitude -30.0 deg>
>>> c.lon
<Longitude180 70.0 deg>
>>> c.radius
<Distance 695508.0 km>

Transforming Between Coordinate Frames

Both SkyCoord and BaseCoordinateFrame instances have a transform_to method. This can be used to transform the frame to any other frame, either implemented in SunPy or in Astropy the astropy documentation for more details. An example of transforming the center of the solar disk to Carrington coordinates is:

>>> c = SkyCoord(0*u.arcsec, 0*u.arcsec, frame=frames.Helioprojective, obstime="2017-07-26")
>>> c
<SkyCoord (Helioprojective: obstime=2017-07-26 00:00:00, rsun=695508.0 km, observer=<HeliographicStonyhurst Coordinate (obstime=2017-07-26 00:00:00): (lon, lat, radius) in (deg, deg, AU)
    ( 0.,  5.31701821,  1.01567428)>): (Tx, Ty) in arcsec
    ( 0.,  0.)>

>>> c.transform_to(frames.HeliographicCarrington)
<SkyCoord (HeliographicCarrington: obstime=2017-07-26 00:00:00): (lon, lat, radius) in (deg, deg, km)
    (-76.00701638,  5.31701821,  695508.00000058)>

It is also possible to transform to any coordinate system implemented in Astropy. This can be used to find the position of the solar limb in AltAz equatorial coordinates:

>>> from astropy.coordinates import EarthLocation, AltAz

>>> time = '2017-07-11 15:00'
>>> greenbelt = EarthLocation(lat=39.0044*u.deg, lon=-76.8758*u.deg)
>>> greenbelt_frame = AltAz(obstime=time, location=greenbelt)

>>> west_limb = SkyCoord(900*u.arcsec, 0*u.arcsec, frame=frames.Helioprojective, obstime=time)
>>> west_limb.transform_to(greenbelt_frame)
<AltAz Coordinate (obstime=2017-07-11 15:00:00.000, location=(1126916.53031967, -4833386.58391627, 3992696.622115747) m, pressure=0.0 hPa, temperature=0.0 deg_C, relative_humidity=0, obswl=1.0 micron): (az, alt, distance) in (deg, deg, m)
    ( 111.40839171,  57.16645763,   1.51860261e+11)>

Observer Location Information

Both Helioprojective and Heliocentric frames are defined by the location of the observer. For example in Helioprojective the observer is at the origin of the coordinate system. This information is encoded in the Helioprojective and Heliocentric frames as the observer attribute, which is itself an instance of the HeliographicStonyhurst frame. The default observer location is set to the position of the Earth (using get_body_heliographic_stonyhurst) as long as the obstime attribute is specified. If the obstime attribute is not set then you will be unable to transform the frame unless an explicit observer is specified, as the time is required to calculate the location of the Earth. The location of the observer is automatically populated from meta data when coordinate frames are created using map.

It is possible to convert from a Helioprojective frame with one observer location to another Helioprojective frame with a different observer location, by converting through Heliographic, this does involve making an assumption of the radius of the Sun to calculate the position on the solar sphere. The conversion can be performed as follows:

# Input coordinate
>>> hpc1 = SkyCoord(0*u.arcsec, 0*u.arcsec, observer="earth", obstime="2017-07-26", frame=frames.Helioprojective)
# Define a new Helioprojective frame with a different observer.
>>> hpc_out = sunpy.coordinates.Helioprojective(observer="venus", obstime="2017-07-26")
# Perform the transformation from one to the other.
>>> hpc2 = hpc1.transform_to(hpc_out)

An example with two maps, named aia and stereo:

>>> hpc1 = SkyCoord(0*u.arcsec, 0*u.arcsec, frame=aia.coordinate_frame)
>>> hpc2 = hpc1.transform_to(stereo.coordinate_frame)

Design of the Coordinates Module

This module works by defining a collection of Frames (sunpy.coordinates.frames), which exists on a transformation graph, where the transformations between the coordinate frames are then defined and registered with the transformation graph (sunpy.coordinates.transformations). It is also possible to transform SunPy frames to Astropy frames.

Positions within these Frames are stored as a Representation of a coordinate, a representation being a description of a point in a Cartesian, spherical or cylindrical system (sunpy.coordinates.representation). A frame that contains a representation of one or many points is said to have been ‘realized’.

For a more in depth look at the design and concepts of the Astropy coordinates system see Overview of astropy.coordinates concepts

Frames and SkyCoord

The SkyCoord class is a high level wrapper around the astropy.coordinates package. It provides an easier way to create and transform coordinates, by using string representations for frames rather than the classes themselves and some other usability improvements, for more information see the SkyCoord documentation.

The main advantage provided by SkyCoord is the support it provides for caching Frame attributes. Frame attributes are extra data specified with a frame, some examples in sunpy.coordinates are obstime or observer for observer location. Only the frames where this data is meaningful have these attributes, i.e. only the Helioprojective frames have observer. However, when you transform into another frame and then back to a projective frame using SkyCoord it will remember the attributes previously provided, and repopulate the final frame with them. If you were to do transformations using the Frames alone this would not happen.

The most important implication for this in sunpy.coordinates is the rsun parameter in the projective frames. If you create a projective frame with a rsun attribute, if you convert back to a projective frame it will be set correctly. It should also be noted that, if you create a Heliographic frame and then transform to a projective frame with an rsun attribute, it will not match the radius coordinate in the Heliographic frame. This is because you may mean to be describing a point above the defined ‘surface’ of the Sun.

Coordinates and WCS

The sunpy.coordinates package provides a mapping between FITS-WCS CTYPE convention and the coordinate frames as defined in sunpy.coordinates. This is used via the astropy.wcs.utils.wcs_to_celestial_frame function, with which the SunPy frames are registered upon being imported. This list is used by packages such as wcsaxes to convert from astropy.wcs.WCS objects to coordinate frames.

The sunpy.map.GenricMap class creates astropy.wcs.WCS objects as amap.wcs, however, it adds some extra attributes to the WCS object to be able to fully specify the coordinate frame. It adds heliographic_observer and rsun.

If you want to obtain a un-realized coordinate frame corresponding to a GenericMap object you can do the following:

>>> import sunpy.map
>>> from sunpy.data.sample import AIA_171_IMAGE

>>> amap = sunpy.map.Map(AIA_171_IMAGE)
>>> amap.observer_coordinate
<Helioprojective Frame (obstime=2011-06-07 06:33:02.770000, rsun=696000000.0 m, observer=<HeliographicStonyhurst Coordinate (obstime=None): (lon, lat, radius) in (deg, deg, m)
    ( 0.,  0.048591,   1.51846026e+11)>)>

which is equivalent to:

>>> from astropy.wcs.utils import wcs_to_celestial_frame
>>> wcs_to_celestial_frame(amap.wcs)
<Helioprojective Frame (obstime=2011-06-07 06:33:02.770000, rsun=696000000.0 m, observer=<HeliographicStonyhurst Coordinate (obstime=None): (lon, lat, radius) in (deg, deg, m)
    ( 0.,  0.048591,   1.51846026e+11)>)>

sunpy.coordinates Package

Functions

get_body_heliographic_stonyhurst(body[, time]) Return a HeliographicStonyhurst frame for the location of a solar-system body at a specified time.
get_earth([time]) Return a SkyCoord for the location of the Earth at a specified time in the HeliographicStonyhurst frame.
get_sun_B0([time]) Return the B0 angle for the Sun at a specified time, which is the heliographic latitude of the Sun-disk center as seen from Earth.
get_sun_L0([time]) Return the L0 angle for the Sun at a specified time, which is the Carrington longitude of the Sun-disk center as seen from Earth.
get_sun_P([time]) Return the position (P) angle for the Sun at a specified time, which is the angle between geocentric north and solar north as seen from Earth, measured eastward from geocentric north.
get_sun_orientation(location[, time]) Return the orientation angle for the Sun from a specified Earth location and time.
get_sunearth_distance([time]) Return the distance between the Sun and the Earth at a specified time.
solar_wcs_frame_mapping(wcs) This function registers the coordinates frames to their FITS-WCS coordinate type values in the astropy.wcs.utils.wcs_to_celestial_frame registry.

Classes

Heliocentric(*args, **kwargs) A coordinate or frame in the Heliocentric system.
HeliographicCarrington(*args, **kwargs) A coordinate or frame in the Carrington Heliographic system.
HeliographicStonyhurst(*args, **kwargs) A coordinate or frame in the Stonyhurst Heliographic system.
Helioprojective(*args, **kwargs) A coordinate or frame in the Helioprojective (Cartesian) system.
NorthOffsetFrame A frame which is relative to some position and another frame.

Class Inheritance Diagram

Inheritance diagram of sunpy.coordinates.frames.Heliocentric, sunpy.coordinates.frames.HeliographicCarrington, sunpy.coordinates.frames.HeliographicStonyhurst, sunpy.coordinates.frames.Helioprojective, sunpy.coordinates.offset_frame.NorthOffsetFrame

sunpy.coordinates.frames Module

Common solar physics coordinate systems.

This submodule implements various solar physics coordinate frames for use with the astropy.coordinates module.

Classes

HeliographicStonyhurst(*args, **kwargs) A coordinate or frame in the Stonyhurst Heliographic system.
HeliographicCarrington(*args, **kwargs) A coordinate or frame in the Carrington Heliographic system.
Heliocentric(*args, **kwargs) A coordinate or frame in the Heliocentric system.
Helioprojective(*args, **kwargs) A coordinate or frame in the Helioprojective (Cartesian) system.

Class Inheritance Diagram

Inheritance diagram of sunpy.coordinates.frames.HeliographicStonyhurst, sunpy.coordinates.frames.HeliographicCarrington, sunpy.coordinates.frames.Heliocentric, sunpy.coordinates.frames.Helioprojective

sunpy.coordinates.transformations Module

Coordinate Transformation Functions

This module contains the functions for converting one sunpy.coordinates.frames object to another.

Warning

The functions in this submodule should never be called directly, transforming between coordinate frames should be done using the .transform_to methods on BaseCoordinateFrame or SkyCoord instances.

The diagram below shows all of the coordinate systems built into the coordinates package, their aliases (useful for converting other coordinates to them using attribute-style access) and the pre-defined transformations between them. The user is free to override any of these transformations by defining new transformations between these systems, but the pre-defined transformations should be sufficient for typical usage.

The color of an edge in the graph (i.e. the transformations between two frames) is set by the type of transformation; the legend box defines the mapping from transform class name to color.

digraph AstropyCoordinateTransformGraph { GCRS [shape=oval label="GCRS\n`gcrs`"]; GeocentricTrueEcliptic [shape=oval label="GeocentricTrueEcliptic\n`geocentrictrueecliptic`"]; ICRS [shape=oval label="ICRS\n`icrs`"]; HCRS [shape=oval label="HCRS\n`hcrs`"]; CIRS [shape=oval label="CIRS\n`cirs`"]; PrecessedGeocentric [shape=oval label="PrecessedGeocentric\n`precessedgeocentric`"]; HeliocentricTrueEcliptic [shape=oval label="HeliocentricTrueEcliptic\n`heliocentrictrueecliptic`"]; FK5 [shape=oval label="FK5\n`fk5`"]; Galactocentric [shape=oval label="Galactocentric\n`galactocentric`"]; BarycentricTrueEcliptic [shape=oval label="BarycentricTrueEcliptic\n`barycentrictrueecliptic`"]; LSR [shape=oval label="LSR\n`lsr`"]; ITRS [shape=oval label="ITRS\n`itrs`"]; HeliographicStonyhurst [shape=oval label="HeliographicStonyhurst\n`heliographic_stonyhurst`"]; Heliocentric [shape=oval label="Heliocentric\n`heliocentric`"]; HeliographicCarrington [shape=oval label="HeliographicCarrington\n`heliographic_carrington`"]; FK4NoETerms [shape=oval label="FK4NoETerms\n`fk4noeterms`"]; Galactic [shape=oval label="Galactic\n`galactic`"]; FK4 [shape=oval label="FK4\n`fk4`"]; Helioprojective [shape=oval label="Helioprojective\n`helioprojective`"]; GalacticLSR [shape=oval label="GalacticLSR\n`galacticlsr`"]; Supergalactic [shape=oval label="Supergalactic\n`supergalactic`"]; AltAz [shape=oval label="AltAz\n`altaz`"]; SkyOffsetFrame[ shape=oval ]; BaseEclipticFrame[ shape=oval ]; BaseRADecFrame[ shape=oval ]; GCRS -> GCRS[ color = "#d95f02" ]; GCRS -> GeocentricTrueEcliptic[ color = "#d95f02" ]; GCRS -> ICRS[ color = "#d95f02" ]; GCRS -> HCRS[ color = "#d95f02" ]; GCRS -> CIRS[ color = "#d95f02" ]; GCRS -> PrecessedGeocentric[ color = "#d95f02" ]; HeliocentricTrueEcliptic -> ICRS[ color = "#d95f02" ]; ICRS -> GCRS[ color = "#d95f02" ]; ICRS -> FK5[ color = "#1b9e77" ]; ICRS -> HeliocentricTrueEcliptic[ color = "#d95f02" ]; ICRS -> Galactocentric[ color = "#555555" ]; ICRS -> CIRS[ color = "#d95f02" ]; ICRS -> BarycentricTrueEcliptic[ color = "#1b9e77" ]; ICRS -> LSR[ color = "#555555" ]; ICRS -> HCRS[ color = "#555555" ]; Galactocentric -> ICRS[ color = "#555555" ]; BarycentricTrueEcliptic -> ICRS[ color = "#1b9e77" ]; ITRS -> ITRS[ color = "#d95f02" ]; ITRS -> CIRS[ color = "#d95f02" ]; HeliographicStonyhurst -> Heliocentric[ color = "#783001" ]; HeliographicStonyhurst -> HeliographicCarrington[ color = "#783001" ]; HeliographicStonyhurst -> HCRS[ color = "#1b9e77" ]; HeliographicCarrington -> HeliographicStonyhurst[ color = "#783001" ]; FK5 -> FK5[ color = "#1b9e77" ]; FK5 -> FK4NoETerms[ color = "#1b9e77" ]; FK5 -> Galactic[ color = "#1b9e77" ]; FK5 -> ICRS[ color = "#1b9e77" ]; FK4 -> FK4[ color = "#d95f02" ]; FK4 -> FK4NoETerms[ color = "#d95f02" ]; Helioprojective -> Helioprojective[ color = "#783001" ]; Helioprojective -> Heliocentric[ color = "#783001" ]; HCRS -> HCRS[ color = "#d95f02" ]; HCRS -> HeliographicStonyhurst[ color = "#1b9e77" ]; HCRS -> ICRS[ color = "#555555" ]; Heliocentric -> Helioprojective[ color = "#783001" ]; Heliocentric -> HeliographicStonyhurst[ color = "#783001" ]; GalacticLSR -> Galactic[ color = "#555555" ]; FK4NoETerms -> FK5[ color = "#1b9e77" ]; FK4NoETerms -> FK4[ color = "#d95f02" ]; FK4NoETerms -> FK4NoETerms[ color = "#1b9e77" ]; FK4NoETerms -> Galactic[ color = "#1b9e77" ]; Supergalactic -> Galactic[ color = "#7570b3" ]; CIRS -> GCRS[ color = "#d95f02" ]; CIRS -> AltAz[ color = "#d95f02" ]; CIRS -> ITRS[ color = "#d95f02" ]; CIRS -> CIRS[ color = "#d95f02" ]; CIRS -> ICRS[ color = "#d95f02" ]; LSR -> ICRS[ color = "#555555" ]; GeocentricTrueEcliptic -> GCRS[ color = "#d95f02" ]; Galactic -> FK5[ color = "#1b9e77" ]; Galactic -> Supergalactic[ color = "#7570b3" ]; Galactic -> FK4NoETerms[ color = "#1b9e77" ]; Galactic -> GalacticLSR[ color = "#555555" ]; AltAz -> AltAz[ color = "#d95f02" ]; AltAz -> CIRS[ color = "#d95f02" ]; PrecessedGeocentric -> GCRS[ color = "#d95f02" ]; overlap=false }

  • AffineTransform:

  • FunctionTransform:

  • FunctionTransformWithFiniteDifference:

  • StaticMatrixTransform:

  • DynamicMatrixTransform:

Functions

hgs_to_hgc(hgscoord, hgcframe) Transform from Heliographic Stonyhurst to Heliograpic Carrington.
hgc_to_hgs(hgccoord, hgsframe) Convert from Heliograpic Carrington to Heliographic Stonyhurst.
hcc_to_hpc(helioccoord, heliopframe) Convert from Heliocentic Cartesian to Helioprojective Cartesian.
hpc_to_hcc(heliopcoord, heliocframe) Convert from Helioprojective Cartesian to Heliocentric Cartesian.
hcc_to_hgs(helioccoord, heliogframe) Convert from Heliocentric Cartesian to Heliographic Stonyhurst.
hgs_to_hcc(heliogcoord, heliocframe) Convert from Heliographic Stonyhurst to Heliograpic Carrington.
hpc_to_hpc(heliopcoord, heliopframe) This converts from HPC to HPC, with different observer location parameters.
hcrs_to_hgs(hcrscoord, hgsframe) Convert from HCRS to Heliographic Stonyhurst (HGS).
hgs_to_hcrs(hgscoord, hcrsframe) Convert from Heliographic Stonyhurst to HCRS.

sunpy.coordinates.ephemeris Module

Ephemeris calculations using SunPy coordinate frames

Functions

get_body_heliographic_stonyhurst(body[, time]) Return a HeliographicStonyhurst frame for the location of a solar-system body at a specified time.
get_earth([time]) Return a SkyCoord for the location of the Earth at a specified time in the HeliographicStonyhurst frame.
get_sun_B0([time]) Return the B0 angle for the Sun at a specified time, which is the heliographic latitude of the Sun-disk center as seen from Earth.
get_sun_L0([time]) Return the L0 angle for the Sun at a specified time, which is the Carrington longitude of the Sun-disk center as seen from Earth.
get_sun_P([time]) Return the position (P) angle for the Sun at a specified time, which is the angle between geocentric north and solar north as seen from Earth, measured eastward from geocentric north.
get_sunearth_distance([time]) Return the distance between the Sun and the Earth at a specified time.
get_sun_orientation(location[, time]) Return the orientation angle for the Sun from a specified Earth location and time.

sunpy.coordinates.representation Module

SunPy specific representations.

This submodule extends astropy.coordinates.representations with two Spherical representation classes, primarily for use with the Helioprojective coordinate system, due to the convention of Longitude going from -180 to 180 degrees.

Classes

Longitude180 Quantity class that represents Longitude.
SphericalWrap180Representation(lon, lat, …) Representation of points in 3D Spherical coordinates.
UnitSphericalWrap180Representation(lon, lat) Representation of points in 3D Spherical coordinates.

Class Inheritance Diagram

Inheritance diagram of sunpy.coordinates.representation.Longitude180, sunpy.coordinates.representation.SphericalWrap180Representation, sunpy.coordinates.representation.UnitSphericalWrap180Representation

sunpy.coordinates.offset_frame Module

Classes

NorthOffsetFrame A frame which is relative to some position and another frame.

Class Inheritance Diagram

Inheritance diagram of sunpy.coordinates.offset_frame.NorthOffsetFrame

sunpy.coordinates.wcs_utils Module

Functions

solar_wcs_frame_mapping(wcs) This function registers the coordinates frames to their FITS-WCS coordinate type values in the astropy.wcs.utils.wcs_to_celestial_frame registry.

Attribution

Some of this documentation was adapted from Astropy under the terms of the BSD License.

This package was initially developed by Pritish Chakraborty as part of GSOC 2014 and Stuart Mumford.