Merge pull request #68 from Deuchnord/occultation

Add support for occultations

kosmorrolib/data.py

@@ -20,6 +20,8 @@ from abc import ABC, abstractmethod
 from typing import Union
 from datetime import datetime
 
+from numpy import pi, arcsin
+
 from skyfield.api import Topos, Time
 from skyfield.vectorlib import VectorSum as SkfPlanet
 
@@ -40,6 +42,7 @@ MOON_PHASES = {
 EVENTS = {
     'OPPOSITION': {'message': _('%s is in opposition')},
     'CONJUNCTION': {'message': _('%s and %s are in conjunction')},
+    'OCCULTATION': {'message': _('%s occults %s')},
     'MAXIMAL_ELONGATION': {'message': _("%s's largest elongation")}
 }
 
@@ -106,16 +109,19 @@ class Object(ABC):
     def __init__(self,
                  name: str,
                  skyfield_name: str,
-                 ephemerides: AsterEphemerides or None = None):
+                 ephemerides: AsterEphemerides or None = None,
+                 radius: float = None):
         """
         Initialize an astronomical object
 
         :param str name: the official name of the object (may be internationalized)
         :param str skyfield_name: the internal name of the object in Skyfield library
+        :param float radius: the radius (in km) of the object
         :param AsterEphemerides ephemerides: the ephemerides associated to the object
         """
         self.name = name
         self.skyfield_name = skyfield_name
+        self.radius = radius
         self.ephemerides = ephemerides
 
     def get_skyfield_object(self) -> SkfPlanet:
@@ -125,6 +131,18 @@ class Object(ABC):
     def get_type(self) -> str:
         pass
 
+    def get_apparent_radius(self, time: Time, from_place) -> float:
+        """
+        Calculate the apparent radius, in degrees, of the object from the given place at a given time.
+        :param time:
+        :param from_place:
+        :return:
+        """
+        if self.radius is None:
+            raise ValueError('Missing radius for %s object' % self.name)
+
+        return 360 / pi * arcsin(self.radius / from_place.at(time).observe(self.get_skyfield_object()).distance().km)
+
 
 class Star(Object):
     def get_type(self) -> str:
@@ -212,13 +230,13 @@ def skyfield_to_moon_phase(times: [Time], vals: [int], now: Time) -> Union[MoonP
 
 MONTHS = ['JAN', 'FEB', 'MAR', 'APR', 'MAY', 'JUN', 'JUL', 'AUG', 'SEP', 'OCT', 'NOV', 'DEC']
 
-ASTERS = [Star(_('Sun'), 'SUN'),
-          Satellite(_('Moon'), 'MOON'),
-          Planet(_('Mercury'), 'MERCURY'),
-          Planet(_('Venus'), 'VENUS'),
-          Planet(_('Mars'), 'MARS'),
-          Planet(_('Jupiter'), 'JUPITER BARYCENTER'),
-          Planet(_('Saturn'), 'SATURN BARYCENTER'),
-          Planet(_('Uranus'), 'URANUS BARYCENTER'),
-          Planet(_('Neptune'), 'NEPTUNE BARYCENTER'),
-          Planet(_('Pluto'), 'PLUTO BARYCENTER')]
+ASTERS = [Star(_('Sun'), 'SUN', radius=696342),
+          Satellite(_('Moon'), 'MOON', radius=1737.4),
+          Planet(_('Mercury'), 'MERCURY', radius=2439.7),
+          Planet(_('Venus'), 'VENUS', radius=6051.8),
+          Planet(_('Mars'), 'MARS', radius=3396.2),
+          Planet(_('Jupiter'), 'JUPITER BARYCENTER', radius=71492),
+          Planet(_('Saturn'), 'SATURN BARYCENTER', radius=60268),
+          Planet(_('Uranus'), 'URANUS BARYCENTER', radius=25559),
+          Planet(_('Neptune'), 'NEPTUNE BARYCENTER', radius=24764),
+          Planet(_('Pluto'), 'PLUTO BARYCENTER', radius=1185)]

kosmorrolib/dumper.py

@@ -116,6 +116,7 @@ class JsonDumper(Dumper):
         if isinstance(obj, Object):
             obj = obj.__dict__
             obj.pop('skyfield_name')
+            obj.pop('radius')
             obj['object'] = obj.pop('name')
             obj['details'] = obj.pop('ephemerides')
             return obj

kosmorrolib/events.py

@@ -56,7 +56,18 @@ def _search_conjunction(start_time: Time, end_time: Time) -> [Event]:
 
             for i, time in enumerate(times):
                 if is_conjs[i]:
-                    conjunctions.append(Event('CONJUNCTION', [aster1, aster2], time.utc_datetime()))
+                    aster1_pos = (aster1.get_skyfield_object() - earth).at(time)
+                    aster2_pos = (aster2.get_skyfield_object() - earth).at(time)
+                    distance = aster1_pos.separation_from(aster2_pos).degrees
+
+                    if distance - aster2.get_apparent_radius(time, earth) < aster1.get_apparent_radius(time, earth):
+                        occulting_aster = [aster1,
+                                           aster2] if aster1_pos.distance().km < aster2_pos.distance().km else [aster2,
+                                                                                                                aster1]
+
+                        conjunctions.append(Event('OCCULTATION', occulting_aster, time.utc_datetime()))
+                    else:
+                        conjunctions.append(Event('CONJUNCTION', [aster1, aster2], time.utc_datetime()))
 
         computed.append(aster1)
 

kosmorrolib/locales/messages.pot

@@ -8,7 +8,7 @@ msgid ""
 msgstr ""
 "Project-Id-Version: kosmorro 0.6.2\n"
 "Report-Msgid-Bugs-To: EMAIL@ADDRESS\n"
-"POT-Creation-Date: 2020-03-24 13:44+0100\n"
+"POT-Creation-Date: 2020-03-29 14:06+0200\n"
 "PO-Revision-Date: YEAR-MO-DA HO:MI+ZONE\n"
 "Last-Translator: FULL NAME <EMAIL@ADDRESS>\n"
 "Language-Team: LANGUAGE <LL@li.org>\n"
@@ -17,90 +17,95 @@ msgstr ""
 "Content-Transfer-Encoding: 8bit\n"
 "Generated-By: Babel 2.8.0\n"
 
-#: kosmorrolib/data.py:30
+#: kosmorrolib/data.py:32
 msgid "New Moon"
 msgstr ""
 
-#: kosmorrolib/data.py:31
+#: kosmorrolib/data.py:33
 msgid "Waxing crescent"
 msgstr ""
 
-#: kosmorrolib/data.py:32
+#: kosmorrolib/data.py:34
 msgid "First Quarter"
 msgstr ""
 
-#: kosmorrolib/data.py:33
+#: kosmorrolib/data.py:35
 msgid "Waxing gibbous"
 msgstr ""
 
-#: kosmorrolib/data.py:34
+#: kosmorrolib/data.py:36
 msgid "Full Moon"
 msgstr ""
 
-#: kosmorrolib/data.py:35
+#: kosmorrolib/data.py:37
 msgid "Waning gibbous"
 msgstr ""
 
-#: kosmorrolib/data.py:36
+#: kosmorrolib/data.py:38
 msgid "Last Quarter"
 msgstr ""
 
-#: kosmorrolib/data.py:37
+#: kosmorrolib/data.py:39
 msgid "Waning crescent"
 msgstr ""
 
-#: kosmorrolib/data.py:41
+#: kosmorrolib/data.py:43
 #, python-format
 msgid "%s is in opposition"
 msgstr ""
 
-#: kosmorrolib/data.py:42
+#: kosmorrolib/data.py:44
 #, python-format
 msgid "%s and %s are in conjunction"
 msgstr ""
 
-#: kosmorrolib/data.py:43
+#: kosmorrolib/data.py:45
+#, python-format
+msgid "%s occults %s"
+msgstr ""
+
+#: kosmorrolib/data.py:46
 #, python-format
 msgid "%s's largest elongation"
 msgstr ""
 
-#: kosmorrolib/data.py:215
+#: kosmorrolib/data.py:233
 msgid "Sun"
 msgstr ""
 
-#: kosmorrolib/data.py:216
+#: kosmorrolib/data.py:234
 msgid "Moon"
 msgstr ""
 
-#: kosmorrolib/data.py:217
+#: kosmorrolib/data.py:235
 msgid "Mercury"
 msgstr ""
 
-#: kosmorrolib/data.py:218
+#: kosmorrolib/data.py:236
 msgid "Venus"
 msgstr ""
 
-#: kosmorrolib/data.py:219
+#: kosmorrolib/data.py:237
 msgid "Mars"
 msgstr ""
 
-#: kosmorrolib/data.py:220
+#: kosmorrolib/data.py:238
 msgid "Jupiter"
 msgstr ""
 
-#: kosmorrolib/data.py:221
+#: kosmorrolib/data.py:239
 msgid "Saturn"
 msgstr ""
 
-#: kosmorrolib/data.py:222
+#: kosmorrolib/data.py:240
 msgid "Uranus"
 msgstr ""
 
-#: kosmorrolib/data.py:223
+#: kosmorrolib/data.py:241
 msgid "Neptune"
 msgstr ""
 
-#: kosmorrolib/data.py:224
+#: kosmorrolib/data.py:242
 msgid "Pluto"
 msgstr ""
 
@@ -116,68 +121,68 @@ msgstr ""
 msgid "{hours}:{minutes}"
 msgstr ""
 
-#: kosmorrolib/dumper.py:147
+#: kosmorrolib/dumper.py:148
 msgid "Expected events:"
 msgstr ""
 
-#: kosmorrolib/dumper.py:151
+#: kosmorrolib/dumper.py:152
 msgid "Note: All the hours are given in UTC."
 msgstr ""
 
-#: kosmorrolib/dumper.py:156
+#: kosmorrolib/dumper.py:157
 msgid "Note: All the hours are given in the UTC{offset} timezone."
 msgstr ""
 
-#: kosmorrolib/dumper.py:202 kosmorrolib/dumper.py:271
+#: kosmorrolib/dumper.py:203 kosmorrolib/dumper.py:272
 msgid "Object"
 msgstr ""
 
-#: kosmorrolib/dumper.py:203 kosmorrolib/dumper.py:272
+#: kosmorrolib/dumper.py:204 kosmorrolib/dumper.py:273
 msgid "Rise time"
 msgstr ""
 
-#: kosmorrolib/dumper.py:204 kosmorrolib/dumper.py:273
+#: kosmorrolib/dumper.py:205 kosmorrolib/dumper.py:274
 msgid "Culmination time"
 msgstr ""
 
-#: kosmorrolib/dumper.py:205 kosmorrolib/dumper.py:274
+#: kosmorrolib/dumper.py:206 kosmorrolib/dumper.py:275
 msgid "Set time"
 msgstr ""
 
-#: kosmorrolib/dumper.py:219 kosmorrolib/dumper.py:277
+#: kosmorrolib/dumper.py:220 kosmorrolib/dumper.py:278
 msgid "Moon phase:"
 msgstr ""
 
-#: kosmorrolib/dumper.py:220
+#: kosmorrolib/dumper.py:221
 msgid "{next_moon_phase} on {next_moon_phase_date} at {next_moon_phase_time}"
 msgstr ""
 
-#: kosmorrolib/dumper.py:258
+#: kosmorrolib/dumper.py:259
 msgid "A Summary of your Sky"
 msgstr ""
 
-#: kosmorrolib/dumper.py:262
+#: kosmorrolib/dumper.py:263
 msgid ""
 "This document summarizes the ephemerides and the events of {date}. It "
 "aims to help you to prepare your observation session. All the hours are "
 "given in {timezone}."
 msgstr ""
 
-#: kosmorrolib/dumper.py:268
+#: kosmorrolib/dumper.py:269
 msgid ""
 "Don't forget to check the weather forecast before you go out with your "
 "material."
 msgstr ""
 
-#: kosmorrolib/dumper.py:270
+#: kosmorrolib/dumper.py:271
 msgid "Ephemerides of the day"
 msgstr ""
 
-#: kosmorrolib/dumper.py:279
+#: kosmorrolib/dumper.py:280
 msgid "Expected events"
 msgstr ""
 
-#: kosmorrolib/dumper.py:354
+#: kosmorrolib/dumper.py:355
 msgid ""
 "Building PDFs was not possible, because some dependencies are not "
 "installed.\n"

test/__init__.py

@@ -1,5 +1,5 @@
+from .core import *
+from .data import *
 from .dumper import *
 from .ephemerides import *
 from .events import *
-from .core import *
-

test/data.py

@@ -0,0 +1,17 @@
+import unittest
+
+from kosmorrolib import data, core
+
+
+class DataTestCase(unittest.TestCase):
+    def test_object_radius_must_be_set_to_get_apparent_radius(self):
+        o = data.Planet('Saturn', 'SATURN')
+
+        with self.assertRaises(ValueError) as context:
+            o.get_apparent_radius(core.get_timescale().now(), core.get_skf_objects()['earth'])
+
+        self.assertEqual(('Missing radius for Saturn object',), context.exception.args)
+
+
+if __name__ == '__main__':
+    unittest.main()

test/events.py

@@ -38,6 +38,8 @@ class EventTestCase(unittest.TestCase):
 
         (date(2020, 3, 24), [Event('MAXIMAL_ELONGATION', [ASTERS[2]], datetime(2020, 3, 24, 1, 56), details='27.8°'),
                              Event('MAXIMAL_ELONGATION', [ASTERS[3]], datetime(2020, 3, 24, 21, 58), details='46.1°')]),
+
+        (date(2005, 6, 16), [Event('OCCULTATION', [ASTERS[1], ASTERS[5]], datetime(2005, 6, 16, 6, 31))])
     )
 
     @data_provider(expected_events_provider)