1"""Vector functions and their composition."""
2
3from jplephem.names import target_names as _jpl_code_name_dict
4from numpy import max
5from .constants import C_AUDAY
6from .descriptorlib import reify
7from .errors import DeprecationError
8from .functions import length_of
9from .positionlib import build_position
10from .timelib import Time
11
12class VectorFunction(object):
13    """Given a time, computes a corresponding position."""
14
15    ephemeris = None
16
17    @reify
18    def vector_name(self):
19        return type(self).__name__
20
21    @reify
22    def center_name(self):
23        return _jpl_name(self.center)
24
25    @reify
26    def target_name(self):
27        return _jpl_name(self.target)
28
29    def __repr__(self):
30        return '<{0} {1}>'.format(type(self).__name__, str(self))
31
32    def __str__(self):
33        if self.target is self:
34            return self.target_name
35        return self.arrow_str()
36
37    def arrow_str(self):
38        return '{0} {1} -> {2}'.format(
39            self.vector_name, self.center_name, self.target_name,
40        )
41
42    def __add__(self, other):
43        if self.target != other.center:
44            if other.target == self.center:
45                self, other = other, self
46            else:
47                raise ValueError(
48                    "you can only add two vectors"
49                    " if the target where one of the vectors ends"
50                    " is the center where the other vector starts"
51                )
52
53        self_vfs = getattr(self, 'vector_functions', None) or (self,)
54        other_vfs = getattr(other, 'vector_functions', None) or (other,)
55
56        return VectorSum(self.center, other.target, self_vfs + other_vfs)
57
58    def __neg__(self):
59        return ReversedVector(self)
60
61    def __sub__(self, other):
62        if self.center != other.center:
63            raise ValueError(
64                "you can only subtract two vectors"
65                " if they both start at the same center"
66            )
67
68        self_vfs = getattr(self, 'vector_functions', None) or (self,)
69        other_vfs = getattr(other, 'vector_functions', None) or (other,)
70        other_vfs = tuple(reversed([-vf for vf in other_vfs]))
71
72        return VectorSum(other.target, self.target, other_vfs + self_vfs)
73
74    def at(self, t):
75        """At time ``t``, compute the target's position relative to the center.
76
77        If ``t`` is an array of times, then the returned position object
78        will specify as many positions as there were times.  The kind of
79        position returned depends on the value of the ``center``
80        attribute:
81
82        * Solar System Barycenter: :class:`~skyfield.positionlib.Barycentric`
83        * Center of the Earth: :class:`~skyfield.positionlib.Geocentric`
84        * Difference: :class:`~skyfield.positionlib.Geometric`
85        * Anything else: :class:`~skyfield.positionlib.ICRF`
86
87        """
88        if not isinstance(t, Time):
89            raise ValueError('please provide the at() method with a Time'
90                             ' instance as its argument, instead of the'
91                             ' value {0!r}'.format(t))
92        p, v, gcrs_position, message = self._at(t)
93        center = self.center
94        position = build_position(p, v, t, center, self.target)
95        position._ephemeris = self.ephemeris
96        position._observer_gcrs_au = gcrs_position
97        position.message = message
98        return position
99
100    def _observe_from_bcrs(self, observer):
101        if self.center != 0:
102            raise ValueError('you can only observe() a body whose vector'
103                             ' center is the Solar System Barycenter,'
104                             ' but this vector has the center {0}'
105                             .format(self.center_name))
106        return _correct_for_light_travel_time(observer, self)
107
108    def geometry_of(self, other):
109        raise DeprecationError(
110"""the geometry_of() method has, alas, been deprecated
111
112This old method has been replaced by an improved interface.  If you just
113need your software working again, install Skyfield 0.9.1 for a quick fix:
114
115    pip install skyfield==0.9.1
116
117Or, to update your old code, replace each operation that looks like:
118
119    position = boston.geometry_of(satellite).at(t)
120
121with the vector math that was previously hiding inside the old method:
122
123    position = (satellite - boston).at(t)""")
124
125    def topos(self, latitude=None, longitude=None, latitude_degrees=None,
126              longitude_degrees=None, elevation_m=0.0, x=0.0, y=0.0):
127        raise DeprecationError(
128"""the topos() method has, alas, been deprecated
129
130This old method has been replaced by an improved interface.  If you just
131need your software working again, install Skyfield 0.9.1 for a quick fix:
132
133    pip install skyfield==0.9.1
134
135Or, to update your old code, replace each operation that looks like:
136
137    boston = earth.topos(...)
138
139with the vector math that was previously hiding inside the old method:
140
141    from skyfield.api import Topos
142    boston = earth + Topos(...)""")
143
144    def satellite(self, text):
145        raise DeprecationError(
146"""the satellite() method has, alas, been deprecated
147
148This old method has been replaced by an improved interface.  If you just
149need your software working again, install Skyfield 0.9.1 for a quick fix:
150
151    pip install skyfield==0.9.1
152
153Or, to update your old code, replace each operation that looks like:
154
155    sat = earth.satellite(tle_text)
156
157with the vector math (and the little bit of text manipulation) that was
158previously hiding inside the old method:
159
160    from skyfield.api import EarthSatellite
161    line1, line2 = tle_text.splitlines()[-2:]
162    sat = earth + EarthSatellite(line1, line2)""")
163
164class ReversedVector(VectorFunction):
165    def __init__(self, vector_function):
166        self.center = vector_function.target
167        self.target = vector_function.center
168        self.vector_function = vector_function
169
170    @reify
171    def vector_name(self):
172        return 'Reversed ' + self.vector_function.vector_name
173
174    @reify
175    def center_name(self):
176        return self.vector_function.target_name
177
178    @reify
179    def target_name(self):
180        return self.vector_function.center_name
181
182    def __neg__(self):
183        return self.vector_function
184
185    def _at(self, t):
186        p, v, _, message = self.vector_function._at(t)
187        return -p, -v, None, message
188
189class VectorSum(VectorFunction):
190    def __init__(self, center, target, vector_functions):
191        self.center = center
192        self.target = target
193        self.vector_functions = vector_functions
194
195        # For now, just grab the first ephemeris we can find.
196        ephemerides = (segment.ephemeris for segment in vector_functions
197                       if segment.ephemeris)
198        self.ephemeris = next(ephemerides, None)
199
200    def __str__(self):
201        vector_functions = self.vector_functions
202        lines = [' ' + segment.arrow_str() for segment in vector_functions]
203        return 'Sum of {0} vectors:\n{1}'.format(
204            len(vector_functions),
205            '\n'.join(lines),
206        )
207
208    def __repr__(self):
209        return '<Vector{0}>'.format(self)
210
211    def _at(self, t):
212        p, v = 0.0, 0.0
213        gcrs_position = None
214        vfs = self.vector_functions
215        for vf in vfs:
216            p2, v2, _, message = vf._at(t)
217            if vf.center == 399:
218                gcrs_position = -p
219            p += p2
220            v += v2
221        if vfs[0].center == 0 and vf.center == 399:
222            gcrs_position = p2
223        return p, v, gcrs_position, message
224
225def _correct_for_light_travel_time(observer, target):
226    """Return a light-time corrected astrometric position and velocity.
227
228    Given an `observer` that is a `Barycentric` position somewhere in
229    the solar system, compute where in the sky they will see the body
230    `target`, by computing the light-time between them and figuring out
231    where `target` was back when the light was leaving it that is now
232    reaching the eyes or instruments of the `observer`.
233
234    """
235    t = observer.t
236    ts = t.ts
237    whole = t.whole
238    tdb_fraction = t.tdb_fraction
239
240    cposition = observer.position.au
241    cvelocity = observer.velocity.au_per_d
242
243    tposition, tvelocity, gcrs_position, message = target._at(t)
244
245    distance = length_of(tposition - cposition)
246    light_time0 = 0.0
247    for i in range(10):
248        light_time = distance / C_AUDAY
249        delta = light_time - light_time0
250        if max(abs(delta)) < 1e-12:
251            break
252
253        # We assume a light travel time of at most a couple of days.  A
254        # longer light travel time would best be split into a whole and
255        # fraction, for adding to the whole and fraction of TDB.
256        t2 = ts.tdb_jd(whole, tdb_fraction - light_time)
257
258        tposition, tvelocity, gcrs_position, message = target._at(t2)
259        distance = length_of(tposition - cposition)
260        light_time0 = light_time
261    else:
262        raise ValueError('light-travel time failed to converge')
263    return tposition - cposition, tvelocity - cvelocity, t, light_time
264
265def _jpl_name(target):
266    if not isinstance(target, int):
267        return type(target).__name__
268    name = _jpl_code_name_dict.get(target)
269    if name is None:
270        return str(target)
271    return '{0} {1}'.format(target, name)
272