Implement ephemeris calculations for leap years

.gitignore

@@ -1,5 +1,3 @@
-# See https://help.github.com/articles/ignoring-files/ for more about ignoring files.
-
 # dependencies
 /node_modules
 /.pnp
@@ -21,3 +19,136 @@
 npm-debug.log*
 yarn-debug.log*
 yarn-error.log*
+
+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+# C extensions
+*.so
+
+# Distribution / packaging
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+share/python-wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+MANIFEST
+
+# PyInstaller
+#  Usually these files are written by a python script from a template
+#  before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.nox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*.cover
+*.py,cover
+.hypothesis/
+.pytest_cache/
+cover/
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+db.sqlite3
+db.sqlite3-journal
+
+# Flask stuff:
+instance/
+.webassets-cache
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+.pybuilder/
+target/
+
+# Jupyter Notebook
+.ipynb_checkpoints
+
+# IPython
+profile_default/
+ipython_config.py
+
+# PEP 582; used by e.g. github.com/David-OConnor/pyflow
+__pypackages__/
+
+# Celery stuff
+celerybeat-schedule
+celerybeat.pid
+
+# SageMath parsed files
+*.sage.py
+
+# Environments
+.env
+.venv
+env/
+venv/
+ENV/
+env.bak/
+venv.bak/
+
+# Spyder project settings
+.spyderproject
+.spyproject
+
+# Rope project settings
+.ropeproject
+
+# mkdocs documentation
+/site
+
+# mypy
+.mypy_cache/
+.dmypy.json
+dmypy.json
+
+# Pyre type checker
+.pyre/
+
+# pytype static type analyzer
+.pytype/
+
+# Cython debug symbols
+cython_debug/
+
+# JetBrains
+.idea/
+
+# Ephemeris downloads
+*.bsp

gen_calendar.py

@@ -0,0 +1,91 @@
+import json
+import os
+import traceback
+from functools import cached_property
+from multiprocessing import Pool
+
+from skyfield.almanac import seasons
+from skyfield.api import load
+from skyfield.jpllib import SpiceKernel
+from skyfield.searchlib import find_discrete
+from skyfield.timelib import Time, julian_day
+
+DE441_PART_1 = 'https://naif.jpl.nasa.gov/pub/naif/generic_kernels/spk/planets/de441_part-1.bsp'
+DE441_PART_2 = 'https://naif.jpl.nasa.gov/pub/naif/generic_kernels/spk/planets/de441_part-2.bsp'
+ts = load.timescale()
+
+YEAR_RANGE = range(-13200, 17191)
+
+
+class SeasonEphemeris:
+    @cached_property
+    def _de440(self) -> SpiceKernel:
+        return load('de440.bsp')
+
+    @cached_property
+    def _de441a(self) -> SpiceKernel:
+        return load(DE441_PART_1)
+
+    @cached_property
+    def _de441b(self) -> SpiceKernel:
+        return load(DE441_PART_2)
+
+    # noinspection PyStatementEffect
+    def initialize(self):
+        self._de440
+        self._de441a
+        self._de441b
+
+    def __getstate__(self):
+        state = self.__dict__.copy()
+        state.pop('_de440', None)
+        state.pop('_de441a', None)
+        state.pop('_de441b', None)
+        return state
+
+    def __setstate__(self, state):
+        self.__dict__.update(state)
+
+    def _eph_for(self, year: int) -> SpiceKernel:
+        if 1550 < year < 2650:
+            return self._de440
+        elif year < 1969:
+            return self._de441a
+        else:
+            return self._de441b
+
+    def autumnal_equinox(self, year: int) -> Time:
+        try:
+            t, _ = find_discrete(ts.utc(year, 9, 1), ts.utc(year, 9, 31), seasons(self._eph_for(year)))
+            return t[0]
+        except Exception:
+            traceback.print_exc()
+            raise
+
+
+def paris_jdn(t: Time) -> int:
+    y, mo, d, h, m, s = t.ut1_calendar()
+    ds = h * 3600 + m * 60 + int(round(s)) + 561
+    return int(julian_day(y, mo, d) + ds // 86400)
+
+
+def fr_year(gregorian_year: int) -> int:
+    return gregorian_year - 1791
+
+
+def main():
+    ephemeris = SeasonEphemeris()
+    pool = Pool(os.cpu_count())
+    equinoxes = pool.map(ephemeris.autumnal_equinox, YEAR_RANGE)
+
+    result = {
+        'start_jd': paris_jdn(equinoxes[0]),
+        'start_year': fr_year(int(equinoxes[0].ut1_calendar()[0])),
+        'leap': [int(paris_jdn(b) - paris_jdn(a) == 366) for a, b in zip(equinoxes, equinoxes[1:])]
+    }
+    with open('cal.json', 'w') as f:
+        json.dump(result, f, separators=(',', ':'))
+
+
+if __name__ == '__main__':
+    main()