1# -*- coding: utf-8 -*-
2# Licensed under a 3-clause BSD style license - see LICENSE.rst
3
4from astropy import units as u
5from astropy.utils.decorators import format_doc
6from astropy.coordinates import representation as r
7from astropy.coordinates.baseframe import BaseCoordinateFrame, base_doc
8from astropy.coordinates.attributes import TimeAttribute, QuantityAttribute
9from .utils import EQUINOX_J2000, DEFAULT_OBSTIME
10
11__all__ = ['GeocentricMeanEcliptic', 'BarycentricMeanEcliptic',
12           'HeliocentricMeanEcliptic', 'BaseEclipticFrame',
13           'GeocentricTrueEcliptic', 'BarycentricTrueEcliptic',
14           'HeliocentricTrueEcliptic',
15           'HeliocentricEclipticIAU76', 'CustomBarycentricEcliptic']
16
17
18doc_components_ecl = """
19    lon : `~astropy.coordinates.Angle`, optional, keyword-only
20        The ecliptic longitude for this object (``lat`` must also be given and
21        ``representation`` must be None).
22    lat : `~astropy.coordinates.Angle`, optional, keyword-only
23        The ecliptic latitude for this object (``lon`` must also be given and
24        ``representation`` must be None).
25    distance : `~astropy.units.Quantity` ['length'], optional, keyword-only
26        The distance for this object from the {0}.
27        (``representation`` must be None).
28
29    pm_lon_coslat : `~astropy.units.Quantity` ['angualar speed'], optional, keyword-only
30        The proper motion in the ecliptic longitude (including the ``cos(lat)``
31        factor) for this object (``pm_lat`` must also be given).
32    pm_lat : `~astropy.units.Quantity` ['angular speed'], optional, keyword-only
33        The proper motion in the ecliptic latitude for this object
34        (``pm_lon_coslat`` must also be given).
35    radial_velocity : `~astropy.units.Quantity` ['speed'], optional, keyword-only
36        The radial velocity of this object.
37"""
38
39
40@format_doc(base_doc,
41            components=doc_components_ecl.format('specified location'),
42            footer="")
43class BaseEclipticFrame(BaseCoordinateFrame):
44    """
45    A base class for frames that have names and conventions like that of
46    ecliptic frames.
47
48    .. warning::
49            In the current version of astropy, the ecliptic frames do not yet have
50            stringent accuracy tests.  We recommend you test to "known-good" cases
51            to ensure this frames are what you are looking for. (and then ideally
52            you would contribute these tests to Astropy!)
53    """
54
55    default_representation = r.SphericalRepresentation
56    default_differential = r.SphericalCosLatDifferential
57
58
59doc_footer_geo = """
60    Other parameters
61    ----------------
62    equinox : `~astropy.time.Time`, optional
63        The date to assume for this frame.  Determines the location of the
64        x-axis and the location of the Earth (necessary for transformation to
65        non-geocentric systems). Defaults to the 'J2000' equinox.
66    obstime : `~astropy.time.Time`, optional
67        The time at which the observation is taken.  Used for determining the
68        position of the Earth. Defaults to J2000.
69"""
70
71
72@format_doc(base_doc, components=doc_components_ecl.format('geocenter'),
73            footer=doc_footer_geo)
74class GeocentricMeanEcliptic(BaseEclipticFrame):
75    """
76    Geocentric mean ecliptic coordinates.  These origin of the coordinates are the
77    geocenter (Earth), with the x axis pointing to the *mean* (not true) equinox
78    at the time specified by the ``equinox`` attribute, and the xy-plane in the
79    plane of the ecliptic for that date.
80
81    Be aware that the definition of "geocentric" here means that this frame
82    *includes* light deflection from the sun, aberration, etc when transforming
83    to/from e.g. ICRS.
84
85    The frame attributes are listed under **Other Parameters**.
86    """
87
88    equinox = TimeAttribute(default=EQUINOX_J2000)
89    obstime = TimeAttribute(default=DEFAULT_OBSTIME)
90
91
92@format_doc(base_doc, components=doc_components_ecl.format('geocenter'),
93            footer=doc_footer_geo)
94class GeocentricTrueEcliptic(BaseEclipticFrame):
95    """
96    Geocentric true ecliptic coordinates.  These origin of the coordinates are the
97    geocenter (Earth), with the x axis pointing to the *true* (not mean) equinox
98    at the time specified by the ``equinox`` attribute, and the xy-plane in the
99    plane of the ecliptic for that date.
100
101    Be aware that the definition of "geocentric" here means that this frame
102    *includes* light deflection from the sun, aberration, etc when transforming
103    to/from e.g. ICRS.
104
105    The frame attributes are listed under **Other Parameters**.
106    """
107
108    equinox = TimeAttribute(default=EQUINOX_J2000)
109    obstime = TimeAttribute(default=DEFAULT_OBSTIME)
110
111
112doc_footer_bary = """
113    Other parameters
114    ----------------
115    equinox : `~astropy.time.Time`, optional
116        The date to assume for this frame.  Determines the location of the
117        x-axis and the location of the Earth and Sun.
118        Defaults to the 'J2000' equinox.
119"""
120
121
122@format_doc(base_doc, components=doc_components_ecl.format("barycenter"),
123            footer=doc_footer_bary)
124class BarycentricMeanEcliptic(BaseEclipticFrame):
125    """
126    Barycentric mean ecliptic coordinates.  These origin of the coordinates are the
127    barycenter of the solar system, with the x axis pointing in the direction of
128    the *mean* (not true) equinox as at the time specified by the ``equinox``
129    attribute (as seen from Earth), and the xy-plane in the plane of the
130    ecliptic for that date.
131
132    The frame attributes are listed under **Other Parameters**.
133    """
134
135    equinox = TimeAttribute(default=EQUINOX_J2000)
136
137
138@format_doc(base_doc, components=doc_components_ecl.format("barycenter"),
139            footer=doc_footer_bary)
140class BarycentricTrueEcliptic(BaseEclipticFrame):
141    """
142    Barycentric true ecliptic coordinates.  These origin of the coordinates are the
143    barycenter of the solar system, with the x axis pointing in the direction of
144    the *true* (not mean) equinox as at the time specified by the ``equinox``
145    attribute (as seen from Earth), and the xy-plane in the plane of the
146    ecliptic for that date.
147
148    The frame attributes are listed under **Other Parameters**.
149    """
150
151    equinox = TimeAttribute(default=EQUINOX_J2000)
152
153
154doc_footer_helio = """
155    Other parameters
156    ----------------
157    equinox : `~astropy.time.Time`, optional
158        The date to assume for this frame.  Determines the location of the
159        x-axis and the location of the Earth and Sun.
160        Defaults to the 'J2000' equinox.
161    obstime : `~astropy.time.Time`, optional
162        The time at which the observation is taken.  Used for determining the
163        position of the Sun. Defaults to J2000.
164"""
165
166
167@format_doc(base_doc, components=doc_components_ecl.format("sun's center"),
168            footer=doc_footer_helio)
169class HeliocentricMeanEcliptic(BaseEclipticFrame):
170    """
171    Heliocentric mean ecliptic coordinates.  These origin of the coordinates are the
172    center of the sun, with the x axis pointing in the direction of
173    the *mean* (not true) equinox as at the time specified by the ``equinox``
174    attribute (as seen from Earth), and the xy-plane in the plane of the
175    ecliptic for that date.
176
177    The frame attributes are listed under **Other Parameters**.
178
179    {params}
180
181
182    """
183
184    equinox = TimeAttribute(default=EQUINOX_J2000)
185    obstime = TimeAttribute(default=DEFAULT_OBSTIME)
186
187
188@format_doc(base_doc, components=doc_components_ecl.format("sun's center"),
189            footer=doc_footer_helio)
190class HeliocentricTrueEcliptic(BaseEclipticFrame):
191    """
192    Heliocentric true ecliptic coordinates.  These origin of the coordinates are the
193    center of the sun, with the x axis pointing in the direction of
194    the *true* (not mean) equinox as at the time specified by the ``equinox``
195    attribute (as seen from Earth), and the xy-plane in the plane of the
196    ecliptic for that date.
197
198    The frame attributes are listed under **Other Parameters**.
199
200    {params}
201
202
203    """
204
205    equinox = TimeAttribute(default=EQUINOX_J2000)
206    obstime = TimeAttribute(default=DEFAULT_OBSTIME)
207
208
209@format_doc(base_doc, components=doc_components_ecl.format("sun's center"),
210            footer="")
211class HeliocentricEclipticIAU76(BaseEclipticFrame):
212    """
213    Heliocentric mean (IAU 1976) ecliptic coordinates.  These origin of the coordinates are the
214    center of the sun, with the x axis pointing in the direction of
215    the *mean* (not true) equinox of J2000, and the xy-plane in the plane of the
216    ecliptic of J2000 (according to the IAU 1976/1980 obliquity model).
217    It has, therefore, a fixed equinox and an older obliquity value
218    than the rest of the frames.
219
220    The frame attributes are listed under **Other Parameters**.
221
222    {params}
223
224
225    """
226
227    obstime = TimeAttribute(default=DEFAULT_OBSTIME)
228
229
230@format_doc(base_doc, components=doc_components_ecl.format("barycenter"),
231            footer="")
232class CustomBarycentricEcliptic(BaseEclipticFrame):
233    """
234    Barycentric ecliptic coordinates with custom obliquity.
235    These origin of the coordinates are the
236    barycenter of the solar system, with the x axis pointing in the direction of
237    the *mean* (not true) equinox of J2000, and the xy-plane in the plane of the
238    ecliptic tilted a custom obliquity angle.
239
240    The frame attributes are listed under **Other Parameters**.
241    """
242
243    obliquity = QuantityAttribute(default=84381.448 * u.arcsec, unit=u.arcsec)
244