{ "metadata": { "name": "", "signature": "sha256:f571c46c9bbc4f79b7735c9948e57326aa5a7c0f4de17c74c97cc1ee04c1290f" }, "nbformat": 3, "nbformat_minor": 0, "worksheets": [ { "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "In the last two notebooks about proper motion$^{1,2}$ I saw that the big dipper was placed in really strange positions in the night sky in the far-away future, and in the ancient past. We resolved that the proper way to compute the positions of stars while mimicking an observer on earth was to vary *both* year and epoch together to get the proper celestial orientation.\n", "\n", "I filed an issue in the PyEphem issue tracker$^3$ just to get confirmation, that varying time *and* epoch together was the best route. It was, but there is an intrinsic limitation of the IAU precession model that [PyEphem](http://rhodesmill.org/pyephem/) inherits from libastro. [Skyfield](http://rhodesmill.org/skyfield/) has a similar limitation. As Brandon Rhodes says$^4$\n", "\n", "> Far-past and far-future epochs are broken\n", "> \n", "> If we try computing the location of the Big Dipper in the skies of the far past or future, we get nonsense: 50,000 years ago PyEphem thinks the Dipper was in the southern sky (!), and 50,000 years in the future it places it up at the north celestial pole\n", "\n", "This particular notebook wound up being more of personal research on the subject and less code, because it is my first time reading about precession matrices. I originally intended to utilize data from HORIZONS to correct for precession, because it uses a long-term precession model that can be used for many hundreds of thousands of years. I found out this was a misconception of mine, because the precession model is \"baked into\" this JPL ephemeris. I include this anyway, because I think the notes are helpful to anyone who may be in the same place as me.\n", "\n", "The follow-up [vondrak.ipynb](https://github.com/digitalvapor/asterisms/blob/master/notebooks/vondrak.ipynb) ([nbviewer](http://nbviewer.ipython.org/github/digitalvapor/asterisms/blob/master/notebooks/vondrak.ipynb)) does implement a long-term precession matrix from beginning to end.\n", "\n", "1. [proper-motion.ipynb](https://github.com/digitalvapor/asterisms/blob/master/notebooks/proper-motion.ipynb) ([nbviewer](http://nbviewer.ipython.org/github/digitalvapor/asterisms/blob/master/notebooks/proper-motion.ipynb))\n", "2. [big-dipper.ipynb](https://github.com/digitalvapor/asterisms/blob/master/notebooks/big-dipper.ipynb) ([nbviewer](http://nbviewer.ipython.org/github/digitalvapor/asterisms/blob/master/notebooks/big-dipper.ipynb))\n", "3. [PyEphem issue #61](https://github.com/brandon-rhodes/pyephem/issues/61)\n", "4. [github-issue-61.ipynb](https://github.com/brandon-rhodes/pyephem/blob/master/issues/github-issue-61.ipynb) ([nbviewer](http://nbviewer.ipython.org/github/brandon-rhodes/pyephem/blob/master/issues/github-issue-61.ipynb))\n", "\n", "##What to do?\n", "Well, I'm going to have to find a different model or formula to intercede when past the bounds of J2000 \u00b11000 years.\n", "\n", "Summary of conclusions:\n", "\n", "* The answer to my original question *should I vary epoch and time together for an ancient observer on earth or an observer in the far-flung future?* is yes.\n", "* The IAU precession formulas the code is currently based on gives an accuracy of 0.05 milliarcseconds per cycle from AD 1000 to AD 3000. Past that point we aren't sure (yet). See *[Improvement of the IAU 2000 precession model](http://www.aanda.org/articles/aa/full/2005/10/aa1908/aa1908.html)* (2004).\n", "* A quick eyeball shows that anything past J2000 \u00b1 20000 years is completely erratic.\n", "* Take stellar depictions involving precession with a grain of salt. Apparently a lot of popular programs don't account for our earth's changing orientation.\n", "* Basemap is *much easier* to use than plotting using an equirectangular projection and provides less distortion.\n", "\n", "I looked into Barry's tip. Thanks! We're getting into stuff that is way beyond me calculation-wise, but I found a paper by Owen (1990) that is extremely helpful. More on it below, but we can indeed select time spans from 'BC 9998-03-20 to AD 9999-12-31' for the earth ephemeris using the HORIZONS interface.\n", "\n", "##On HORIZONS Precession Model\n", "[HORIZONS](http://ssd.jpl.nasa.gov/horizons.cgi) uses two different precession models. One between 1799-Jan-1 to 2202-Jan-1, and the other for outside those bounds.\n", "\n", "> [PRECESSION MODEL](http://ssd.jpl.nasa.gov/?horizons_doc#longterm):\n", "\n", "> For the time-span of 1799-Jan-1 to 2202-Jan-1, the official IAU precession model [16] of Lieske is used. As published, this model is valid for only ~200 years on either side of the J2000.0 epoch. This is due to round-off error in the published coefficients and truncation to a 3rd order polynomial in the expressions for the Euler rotation angles. Therefore, outside this interval, the long-term precession and obliquity model [17] of Owen is used to maintain accuracy in the calculation of apparent (\"of-date\") quantities.\n", "\n", "> This model is a rigorous numerical integration of the equations of motion of the celestial pole using Kinoshita's model for the speed of luni-solar precession.\n", "\n", "Precession (IAU) from 1799-Jan-1 to 2202-Jan-1:\n", "\n", "[16]: Lieske, J., \"[Precession Matrix Based on IAU (1976) System of Astronomical Constants](http://adsabs.harvard.edu/full/1979A%26A....73..282L)\", *Astron. Astrophys.* 73, 282-284, 1979. ([pdf](http://articles.adsabs.harvard.edu/cgi-bin/nph-iarticle_query?1979A%26A....73..282L&defaultprint=YES&filetype=.pdf))\n", "\n", "Precession (long-term) before 1799-Jan-1 and after 2202-Jan-1:\n", "\n", "[17]: Owen, William M., Jr., (JPL) *[A Theory of the Earth's Precession Relative to the Invariable Plane of the Solar System](https://archive.org/details/theoryofearthspr00owen)*, Ph.D. Dissertation, University of Florida, 1990. ([pdf](https://ia700609.us.archive.org/21/items/theoryofearthspr00owen/theoryofearthspr00owen.pdf))\n", "\n", "##On *A Theory of the Earth's Precession Relative to the Invariable Plane of the Solar System*\n", "Chapter 4 by Owen, pg 84, is dedicated to the long-term theory. \n", "\n", "Pg 85, on Kinoshita's luni-solar precession used to integrate the motion of the celestial pole:\n", "\n", "> ... sufficed to give numerical accuracies on the order of a nanoarcsecond after 500,000 years. \n", "\n", "A Laskar file of eccentricity and inclination variables spanning a million years (500,000 to either side of J2000) was used for orbital elements and numerical integrations.\n", "\n", "Pg 242:\n", "\n", "> Chebyshev polynomials to degree nine were fit over 8000-year intervals to the angles obtained from a million-year numerical integration based on the work of Laskar (1990) for the motion of the ecliptic and Kinoshita (1977) for the motion of the equator. \n", "\n", "Fig 4-1, pg 105-110, shows the change in precession angles over these 10000 centuries.\n", "\n", "Ch 4.7, *Sources of Error in the Long-Term Theory*, pg 116:\n", "\n", "> The dominant uncertainty in the long-term theory is in the initial speed of general precession in longitude. The currently-accepted value of 5029\".0966/century may in error by a significant fraction of an arcsecond per century; in other words, by one part in 10\u2074. And since this uncertainty is in a rate, the error in the accumulated angle p\u2090 will grow roughly linearly with time. After 500,000 years, an error of 0\".36/century in p\u2081 would cause an error of half a degree in p\u2090. An error of this magnitude easily swamps all other effects.\n", "\n", "> For this reason it is fruitless at the moment to extend the integration much past the million-year interval covered here.\n", "\n", "Pg253:\n", "\n", "![precession](https://raw.githubusercontent.com/digitalvapor/asterisms/master/notebooks/images/1990_owen_fig5-1.png 'precession over a million years')\n", "\n", "> Figure 5-1: J2000 Equatorial Coordinates of the Ecliptic Pole and of the Pole of the Invariable Plane. The positions of the ecliptic pole at T = -5000, T = 0, and T = +5000 are marked with an open square, open circle, and filled circle, respectively. The position of the pole of the invariable plane is marked with an asterisk.\n", "\n", "Where T is in centuries.\n", "\n", "##On HORIZONS Earth Ephemeris DE431MX\n", "The above-mentioned [time spans](http://ssd.jpl.nasa.gov/eph_spans.cgi?id=A) that HORIZONS limits us to for the Earth stem from Ephemeris File [DE431MX](ftp://ssd.jpl.nasa.gov/pub/eph/planets/ascii/de431/).\n", "\n", "```\n", "ID Name Available Time Span Ephemeris File \n", "399 Earth B.C. 9998-03-20 to A.D. 9999-12-31 DE431MX\n", "```\n", "\n", "From [JPL Planetary and Lunar Ephemerides](http://ssd.jpl.nasa.gov/?planet_eph_export):\n", "\n", "> DE431 :
\n", "> Created April 2013; includes librations and 1980 nutation.
\n", "> Referred to the International Celestial Reference Frame version 2.0.
\n", "> Covers JED -0.3100015.5, (-13200 AUG 15) to JED 8000016.5, (17191 MAR 15).
\n", "> \n", "> DE430 and DE431 are documented in the following document:
\n", "> http://ipnpr.jpl.nasa.gov/progress_report/42-196/196C.pdf\n", "\n", "I haven't read *[The Planetary and Lunar Ephemerides \n", "DE430 and DE431](http://ipnpr.jpl.nasa.gov/progress_report/42-196/196C.pdf)* (2014) as thoroughly as I did *A Theory of the Earth's Precession Relative to the Invariable Plane of the Solar System*. However, I think we can ignore the upper and lower limits of B.C. 9998-03-20 to A.D. 9999-12-31 without further proof that this ephemeris is tied into the precession theory of Owen's.\n", "\n", "Of interest though is pg 10.\n", "\n", "> F. Orientation of the Earth\n", "> \n", "> Only the long-term change of the Earth's orientation is modeled in the ephemeris integration. The Earth orientation model used for the DE430 and DE431 integration is based on the International Astronomical Union (IAU) 1976 precession model [21,22] with an estimated linear correction and on a modified IAU 1980 nutation model [23] including only terms with a period of 18.6 years.\n", "\n", "[21]: J. H. Lieske, T. Lederle, W. Fricke, and B. Morando, \u201cExpression for the Precession\n", "Quantities Based on the IAU (1976) System of Astronomical Constants,\u201d *Astronomy and\n", "Astrophysics*, vol. 58,pp. 1\u201316, 1977.\n", "\n", "[22]: J. H. Lieske, \u201cPrecession Matrix Based on IAU (1976) System of Astronomical Con-\n", "stants,\u201d *Astronomy and Astrophysics*, vol. 73, pp. 282\u2013284, 1979.\n", "\n", "[23]: P. K. Seidelmann, \u201c1980 IAU Theory of Nutation: The Final Report of the IAU Working\n", "Group on Nutation,\u201d *Celestial Mechanics*, vol. 27, pp. 79\u2013106, 1982.\n", "\n", "All of these pre-date Owen's theory. At the time of writing, I am unsure whether I can garner information on purely the precession model by using HORIZONS. I corresponded with Jon D. Giorgini, who cleared up my misconceptions on how DE431, the short-term, and long-term models were implemented by HORIZONS.\n", "\n", "> the long-term precession model is used only to compute a rotation matrix to find apparent positions with respect to the \"true-of-date\" coordinate system. It does not affect the dynamical model (an exception being any asteroid or comet that comes within 0.01 au of Earth, for which Earth oblateness acceleration is included in the numerical integration, requiring knowledge of pole direction). The long-term precession model is not available separately in the sense of being able to get a table of Euler angles from Horizons as a function of time.\n", "\n", "I don't think the bounds of B.C. 9998-03-20 to A.D. 9999-12-31 are limiting. They are due to the \"baked-in\" aspects of the IAU76 precession model with DE431.\n", "\n", "> Note that Horizons is not integrating planetary bodies (only asteroids and comets); it interpolates (reads) planetary positions and velocities from DE431, which was developed and integrated by W. Folkner using the standard IAU76 precession model. The IAU76 precession model, in the sense of physics/dynamics, was \"baked into\" DE431 via the lunar tidal model.\n", "\n", "Giorgini worked with Jay Lieske, who developed the IAU76 precession model. They included a higher precision form of it than was published and adopted by the IAU (or used in DE431).\n", "\n", "> Then, outside the 1799-2200 interval, Owens long-term numerical integration model is used, as indicated in the documentation. Lieske was Owens supervisor for the development of the long-term model.\n", "\n", "> The Chebyshev polynomials published by Owens were converted and stored in a binary direct-access form Horizons accesses to construct a transformation matrix at any instant that transforms from the planetary ephemeris coordinates to a true-equator and equinox of date.\n", "\n", "##How do I use this information?\n", "\n", "I'm not totally sure, but I'd like to try to duplicate the above figure 5-1 from Owen's dissertation. If I can duplicate this for the date range of -500,000 to 500,000 about J2000 as Owen did, by using data from HORIZONS or other, then I think we can conclude that we can use such a model.\n", "\n", "[py-NASA-horizons](https://github.com/tpltnt/py-NASA-horizons) and [HorizonJPL](https://github.com/mihok/horizon-jpl) ([docs](https://docs.google.com/document/d/1g9q3ln9LVAATOZ15986HLTCaqcAj_Jd8e_jOGS3YWrE/pub)) provide nice Python wrappers to the telnet interface. `pip install horizonjpl`." ] }, { "cell_type": "code", "collapsed": false, "input": [ "from horizon.interface import Interface\n", "i = Interface()\n", "i.get_version()" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 1, "text": [ "'3.95'" ] } ], "prompt_number": 1 }, { "cell_type": "code", "collapsed": false, "input": [ "earth = i.get(399)\n", "print(earth['Name'])" ], "language": "python", "metadata": {}, "outputs": [ { "output_type": "stream", "stream": "stdout", "text": [ "Earth\n" ] } ], "prompt_number": 2 }, { "cell_type": "code", "collapsed": false, "input": [ "i.query(399)" ], "language": "python", "metadata": {}, "outputs": [ { "metadata": {}, "output_type": "pyout", "prompt_number": 3, "text": [ "{'Atm. pressure': '1.0 bar',\n", " 'Atmos': '5.1 x 10^18 kg',\n", " 'Escape velocity': '11.186 km s^-1',\n", " 'Fluid core rad': '3480 km',\n", " 'Geometric albedo': '0.367',\n", " 'ID': '399',\n", " 'Inner core rad': '1215 km',\n", " 'Love no., k2': '0.299',\n", " 'Magnetic moment': '0.61 gauss Rp^3',\n", " 'Mass, 10^24 kg': '5.97219+-0.0006',\n", " 'Name': 'Earth',\n", " 'Sidereal period': '365.25636 days',\n", " 'Vis. mag. V(1,0)': '-3.86'}" ] } ], "prompt_number": 3 }, { "cell_type": "markdown", "metadata": {}, "source": [ "After Giorgini cleared up my misconceptions, I realize that I can *not* use HORIZONS in the original way I wanted to. However, I found a second long-term precession model that seemed easier for an astrophysics layman. This follow-up, [vondrak.ipynb](https://github.com/digitalvapor/asterisms/blob/master/notebooks/vondrak.ipynb) ([nbviewer](http://nbviewer.ipython.org/github/digitalvapor/asterisms/blob/master/notebooks/vondrak.ipynb)), implements a long-term precession matrix from beginning to end.\n", "\n", "A helpful read is [The IAU Resolutions on Astronomical Reference Systems, Time Scales, and Earth Rotation Models](http://aa.usno.navy.mil/publications/docs/Circular_179.pdf) (2005)." ] } ], "metadata": {} } ] }