Planet_Position_2021_teacher.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "Planet Position 2021 teacher",
      "provenance": [],
      "authorship_tag": "ABX9TyMew4NyL6+gpyg68LDvDRPn",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/tracieschroeder/Astronomy/blob/main/Planet_Position_2021_teacher.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "id": "LWfkGcyqpEFD"
      },
      "outputs": [],
      "source": [
        "import numpy as np\n",
        "import pandas as pd\n",
        "import matplotlib.pyplot as plt"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "![link text](https://rhodesmill.org/skyfield/_static/logo.png)\n",
        "\n",
        "Skyfield is a pure-Python astronomy package that is compatible with both Python 2 and 3 and makes it easy to generate high precision research-grade positions for planets and Earth satellites.\n",
        "\n",
        "Documentation can be found here: https://rhodesmill.org/skyfield/"
      ],
      "metadata": {
        "id": "BVP6Lbvjrj2s"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "pip install skyfield"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "qFNH_Xs5pIqh",
        "outputId": "4d73299c-e745-4ec1-c005-e482a66f3642"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Collecting skyfield\n",
            "  Downloading skyfield-1.42-py3-none-any.whl (436 kB)\n",
            "\u001b[K     |████████████████████████████████| 436 kB 7.8 MB/s \n",
            "\u001b[?25hCollecting sgp4>=2.2\n",
            "  Downloading sgp4-2.21-cp37-cp37m-manylinux_2_5_x86_64.manylinux1_x86_64.manylinux_2_17_x86_64.manylinux2014_x86_64.whl (228 kB)\n",
            "\u001b[K     |████████████████████████████████| 228 kB 39.4 MB/s \n",
            "\u001b[?25hRequirement already satisfied: numpy in /usr/local/lib/python3.7/dist-packages (from skyfield) (1.21.6)\n",
            "Collecting jplephem>=2.13\n",
            "  Downloading jplephem-2.17-py3-none-any.whl (46 kB)\n",
            "\u001b[K     |████████████████████████████████| 46 kB 2.3 MB/s \n",
            "\u001b[?25hRequirement already satisfied: certifi>=2017.4.17 in /usr/local/lib/python3.7/dist-packages (from skyfield) (2021.10.8)\n",
            "Installing collected packages: sgp4, jplephem, skyfield\n",
            "Successfully installed jplephem-2.17 sgp4-2.21 skyfield-1.42\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "For this plot, we are going to import two of the skyfield libraries. \n",
        "\n",
        "almanac contains a ton of information about solar system bodies at given times. For example, we are going to ask the almanac to compare the exact times for the sunrise and sunset each day as well as the planet position at that time. https://rhodesmill.org/skyfield/almanac.html\n",
        "\n",
        "load, wgs84 saves the pulled data into a working directory, while wgs84 refers to the World Geodetic System 1984. This is the reference coordinate system used by GPS to find locations on Earth. https://rhodesmill.org/skyfield/api-topos.html#skyfield.toposlib.wgs84 "
      ],
      "metadata": {
        "id": "zMLKZ4aKs3tC"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "from skyfield import almanac\n",
        "from skyfield.api import load, wgs84"
      ],
      "metadata": {
        "id": "KVLzjORlpSsp"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Now we need to choose an ephemeris. An ephemeris is a book with tables that gives the trajectory of astronomical objects in the sky. It is a lot of data compiled in one place. \n",
        "\n",
        "de421.bsp is a fairly short ephemeris created in 2008 by NASA's JPL that contains data from 1900 to 2050.\n",
        "\n",
        "If you would like to look in the the ephemeris world further: https://rhodesmill.org/skyfield/planets.html"
      ],
      "metadata": {
        "id": "co4ZCKDoxlhl"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "planets = load('de421.bsp')"
      ],
      "metadata": {
        "id": "_kWv0Mh_xb_Q",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "2e28ed36-908a-4b6c-ab26-2861212289e2"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stderr",
          "text": [
            "[#################################] 100% de421.bsp\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "You can check to see if a certain ephemeris contains data on a specific planet by using the print command."
      ],
      "metadata": {
        "id": "7hFXzZzazP7H"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "print(planets)"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "9nUGv95FzDzt",
        "outputId": "a89262b7-a0be-4c94-fa0d-e720d583b3ed"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "SPICE kernel file 'de421.bsp' has 15 segments\n",
            "  JD 2414864.50 - JD 2471184.50  (1899-07-28 through 2053-10-08)\n",
            "      0 -> 1    SOLAR SYSTEM BARYCENTER -> MERCURY BARYCENTER\n",
            "      0 -> 2    SOLAR SYSTEM BARYCENTER -> VENUS BARYCENTER\n",
            "      0 -> 3    SOLAR SYSTEM BARYCENTER -> EARTH BARYCENTER\n",
            "      0 -> 4    SOLAR SYSTEM BARYCENTER -> MARS BARYCENTER\n",
            "      0 -> 5    SOLAR SYSTEM BARYCENTER -> JUPITER BARYCENTER\n",
            "      0 -> 6    SOLAR SYSTEM BARYCENTER -> SATURN BARYCENTER\n",
            "      0 -> 7    SOLAR SYSTEM BARYCENTER -> URANUS BARYCENTER\n",
            "      0 -> 8    SOLAR SYSTEM BARYCENTER -> NEPTUNE BARYCENTER\n",
            "      0 -> 9    SOLAR SYSTEM BARYCENTER -> PLUTO BARYCENTER\n",
            "      0 -> 10   SOLAR SYSTEM BARYCENTER -> SUN\n",
            "      3 -> 301  EARTH BARYCENTER -> MOON\n",
            "      3 -> 399  EARTH BARYCENTER -> EARTH\n",
            "      1 -> 199  MERCURY BARYCENTER -> MERCURY\n",
            "      2 -> 299  VENUS BARYCENTER -> VENUS\n",
            "      4 -> 499  MARS BARYCENTER -> MARS\n"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "source": [
        "earth = planets['earth']\n",
        "mercury = planets['mercury']"
      ],
      "metadata": {
        "id": "8a652BKBspQO"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "Start by defining the time period we want to plot. \n",
        "The timescale() object manages this conversion for us. \n",
        "Then we can set it to look at every day of the year 2021. "
      ],
      "metadata": {
        "id": "EvGIT0xugrmh"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "ts = load.timescale()\n",
        "start, end = ts.utc(2021, 1, 1), ts.utc(2021, 12, 31)\n",
        "\n",
        "ra, dec, distance = mercury.radec()"
      ],
      "metadata": {
        "id": "52WjlEvGgpv3",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 200
        },
        "outputId": "dc5a0989-10ec-41f2-e250-7f88b156fa4a"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "error",
          "ename": "AttributeError",
          "evalue": "ignored",
          "traceback": [
            "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
            "\u001b[0;31mAttributeError\u001b[0m                            Traceback (most recent call last)",
            "\u001b[0;32m<ipython-input-7-6b76c43deec1>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      2\u001b[0m \u001b[0mstart\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mend\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mts\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mutc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m2021\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m1\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mts\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mutc\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m2021\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m12\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m31\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      3\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 4\u001b[0;31m \u001b[0mra\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdec\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mdistance\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mmercury\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mradec\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
            "\u001b[0;31mAttributeError\u001b[0m: 'VectorSum' object has no attribute 'radec'"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "Now, we plot.\n",
        "\n",
        "We want right ascension (ra.hours) plotted on the x axis and declination (dec.degees) on the y axis. "
      ],
      "metadata": {
        "id": "B2Jxp4JMZJQD"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "fig = plt.figure(figsize=(15, 10))\n",
        "plt.scatter(ra.hours, dec.degrees)\n",
        "\n",
        "plt.title(\"Venus Position from Council Grove, Kansas\")\n",
        "plt.xlabel(\"Right Ascension\")\n",
        "plt.ylabel(\"Declination\");"
      ],
      "metadata": {
        "id": "yQg8RqI8lyiB"
      },
      "execution_count": null,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "source": [
        "#Challenge\n",
        "You have plotted the planets position using right ascension and declination. These plots will be the same for any location on Earth. (If you doubt that statement, go back to where we added Council Grove into our code. Change the location and replot.)\n",
        "\n",
        "To see what the planet path would look like as viewed from Council Grove, change the graph to plot altitude and azimuth."
      ],
      "metadata": {
        "id": "EI1swN37r9RZ"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "#Credits\n",
        "\n",
        "This notebook was adapted from Skyfield's [Charting an Apparition of Venus](https://rhodesmill.org/skyfield/example-plots.html) example plot by Tracie Schroeder."
      ],
      "metadata": {
        "id": "pjpHyBuQdS0p"
      }
    }
  ]
}