#!/usr/bin/env python3

# Plot all the planets' orbits, as viewed from a point that
# floats above the Earth's north ecliptic pole and moves with
# the Earth, to demonstrate phenomena like epicycles and the
# Venus pentagram. Idea from Galen Gisler's planetarium show.
#
# Copyright 2020 by Akkana Peck: Share and enjoy under the GPLv2 or later.

# Weird GLib bug: GLib.timeout_add takes an integer number
# of milliseconds, but if you pass it a float, it sets the timeout
# but then doesn't call draw() automatically after configure(),
# resulting in a transparent window since the black background
# doesn't get drawn. I guess it's some mysterious GTK/Cairo bug.
# Of course I could require an integer timeout when parsing arguments,
# but it's such an amusingly weird bug that I've left it as a float.

import ephem

from tkinter import Tk, Canvas, PhotoImage, mainloop, LEFT

import math
import argparse
import sys, os
from datetime import datetime, timezone

ICONDIR = os.path.expanduser("~/Docs/Preso/mars/pix/")

earth = { "name": "Earth", "obj": ephem.Sun(), "color": "#08f",
          "path": [], "xypath": [],
          "line": None, "disk": None,
          "imgname": os.path.join(ICONDIR, "tinyearth.png") }

mars = { "name": "Mars", "obj": ephem.Mars(), "color": "#f80",
         "path": [], "xypath": [], "oppositions": [],
         "line": None, "disk": None,
         "imgname": os.path.join(ICONDIR, "tinymars.png") }

# oppositions will include date, earth hlon, earth dist, mars hlon, mars dist
oppositions = []

table_header = "%-20s %10s %10s" % ("Date", "Distance", "Size")
table_format = "%-20s %10.3f %10.2f"


def find_next_opposition(start_time):
    """Find oppsition and time of closest approach for the given time range.
       Input is the start time, either in ephem.Date or float julian.
       Output is two ephem.Dates: opposition, closest approach
    """
    t = start_time
    timedelta = ephem.hour * 6
    mars = ephem.Mars()
    sun = ephem.Sun()
    min_dist = 20
    oppy_date = None
    closest_date = None
    last_dlon = None

    # Loop til we've found opposition, plus 15 days.
    # Opposition is when dlon changes sign and is very small.
    while not oppy_date or t - oppy_date < 15:
        mars.compute(t)
        sun.compute(t)
        dlon = mars.hlon - sun.hlon

        # Does dlon have the opposite sign from last_dlon?
        if last_dlon and abs(dlon) < .1 and \
           (dlon == 0 or (dlon < 0) != (last_dlon < 0)):
            oppy_date = t
        if mars.earth_distance < min_dist:
            closest_date = t
            min_dist = mars.earth_distance

        if oppy_date and closest_date:
            return ephem.Date(oppy_date), ephem.Date(closest_date)

        last_dlon = dlon
        t += timedelta


class OrbitViewWindow():
    def __init__(self, auscale, timestep, time_increment=1,
                 start_time=None, stopped=False):
        """time_increment is in days.
           start_time is anything that can be turned into a ephem.Date object.
        """
        self.auscale = auscale
        self.timestep = timestep
        if stopped:
            self.stepping = False
        else:
            self.stepping = True

        if start_time:
            self.time = ephem.Date(start_time)
        else:
            # Default to starting 30 days before present
            self.time = ephem.Date(datetime.now(tz=timezone.utc)) \
                - ephem.hour * 24 * 30
        print("Start time:", ephem.Date(self.time))

        self.opp_date, self.closest_date = find_next_opposition(self.time)
        # print("Next opposition:", self.opp_date)
        # print("Next closest:", self.closest_date)

        self.time_increment = ephem.hour * time_increment * 24

        self.linewidth = 3

        self.width = 1024
        self.height = 768

        self.halfwidth = self.width/2.
        self.halfheight = self.height/2.
        self.dist_scale = self.halfheight / self.auscale

        tkmaster =  Tk()
        tkmaster.title("Mars Oppositions")
        self.canvas = Canvas(tkmaster, bg="black",
                             width=self.width, height=self.height)
        # Start with just the Sun
        try:
            self.sunimg = PhotoImage(file=os.path.join(ICONDIR, "tinysun.png"))
            self.canvas.create_image(self.width/2, self.height/2,
                                     image=self.sunimg)
        except:
            sunrad = 20
            self.canvas.create_oval(self.width/2 - sunrad,
                                    self.height/2 - sunrad,
                                    self.width/2 + sunrad,
                                    self.height/2 + sunrad,
                                    fill="yellow")

        self.canvas.pack()

        # Canvas requires that the app save a reference to PhotoImages:
        # the canvas doesn't keep the references.
        try:
            earth["tinyimg"] = PhotoImage(file=earth["imgname"])
        except:
            earth["tinyimg"] = None
        try:
            mars["tinyimg"] = PhotoImage(file=mars["imgname"])
        except:
            mars["tinyimg"] = None

        tkmaster.bind("<KeyPress-q>", sys.exit)
        tkmaster.bind("<KeyPress-space>", self.toggle_stepping)

        print(table_header)

        # Schedule the first draw
        self.step_draw()

    def toggle_stepping(self, key):
        self.stepping = not self.stepping

    def step_draw(self):
        """Calculate and draw the next position of each planet.
        """
        # If we don't call step_draw at all, we'll never get further key events
        # that could restart the animation. So just step at a much slower pace.
        if not self.stepping:
            self.canvas.after(500, self.step_draw)
            return

        # Adding a float to ephem.Date turns it into a float.
        # You can get back an ephem.Date with: ephem.Date(self.time).
        self.time += self.time_increment

        for p in (earth, mars):
            p["obj"].compute(self.time)

            # ephem treats Earth specially, what a hassle!
            # There is no ephem.Earth body; ephem.Sun gives the Earth's
            # hlon as hlon, but I guess we need to use earth_distance.
            oppy = False
            if p["name"] == "Earth":
                hlon = p["obj"].hlon
                sundist = p["obj"].earth_distance
                earthdist = 0
                size = 0
            else:
                hlon = p["obj"].hlon
                sundist = p["obj"].sun_distance
                earthdist = p["obj"].earth_distance
                size = p["obj"].size

                if abs(self.time - self.opp_date) <= .5:
                    oppy = True
                    if self.opp_date < self.closest_date:
                        print(table_format % (self.opp_date, earthdist, size),
                              "Opposition")
                        print(table_format % (self.closest_date,
                                              earthdist, size),
                              "Closest approach")
                    else:
                        print(table_format % (self.closest_date,
                                              earthdist, size),
                              "Closest approach")
                        print(table_format % (self.opp_date, earthdist, size),
                              "Opposition")

            xn, yn = self.planet_x_y(hlon, sundist)
            radius = 10

            if oppy:
                # Create outline circles for Mars and Earth at opposition.
                # xn, yn should be Mars since Earth was done first.
                # Create the open circle at the bottom of the stacking order.
                # There may be a way to do this by passing in tags,
                # but I can't find any documentation on tags.
                self.canvas.tag_lower(
                    self.canvas.create_oval(xn-radius, yn-radius,
                                            xn+radius, yn+radius,
                                            outline=p["color"], width=3)
                )

                earthx = earth["xypath"][-2]
                earthy = earth["xypath"][-1]
                self.canvas.tag_lower(
                    self.canvas.create_oval(earthx-radius, earthy-radius,
                                            earthx+radius, earthy+radius,
                                            outline=earth["color"], width=3)
                )

                localtz = datetime.now().astimezone().tzinfo
                oppdate = ephem.to_timezone(self.opp_date, localtz)
                opp_str = oppdate.strftime("%Y-%m-%d") + \
                    '\n%.3f AU\n%.1f"' % (earthdist, size)
                if xn < self.width/2:
                    if yn < self.height / 2:
                        anchor = "se"
                    else:
                        anchor = "ne"
                    xtxt = xn - radius
                else:
                    if yn < self.height / 2:
                        anchor = "sw"
                    else:
                        anchor = "nw"
                    xtxt = xn + radius
                ytxt = yn

                txtobj = self.canvas.create_text(xtxt, ytxt,
                                                 fill="white", justify=LEFT,
                                                 font=('sans', 14, 'bold'),
                                                 anchor=anchor,
                                                 text=opp_str)
                # Make sure it's not offscreen
                xt1, yt1, xt2, yt2 = self.canvas.bbox(txtobj)
                if xt1 < 0:
                    xtxt -= xt1
                elif xt2 > self.width:
                    xtxt -= (xt2 - self.width)
                if yt1 < 0:
                    ytxt -= yt1
                elif yt2 > self.height:
                     ytxt -= yt2 - self.height
                self.canvas.coords(txtobj, xtxt, ytxt)

                # Done with this opposition: find the next one.
                self.opp_date, self.closest_date \
                    = find_next_opposition(self.time + 500)

            p["xypath"].append(int(xn))
            p["xypath"].append(int(yn))
            if p["line"]:
                self.canvas.coords(p["line"], p["xypath"])
                if p["tinyimg"]:
                    self.canvas.coords(p["disk"], xn, yn)
                else:
                    self.canvas.coords(p["disk"], xn-radius, yn-radius,
                                       xn+radius, yn+radius)

            else:
                p["line"] = self.canvas.create_line(xn, yn, xn, yn,
                                                    width=self.linewidth,
                                                    fill=p["color"])
                # Use images if there are any, else circles
                if p["tinyimg"]:
                    p["disk"] = self.canvas.create_image(xn-radius, yn-radius,
                                                         image=p["tinyimg"])
                else:
                    p["disk"] = self.canvas.create_oval(xn-radius, yn-radius,
                                                        xn+radius, yn+radius,
                                                        fill=p["color"])

            p["path"].append((hlon, sundist, earthdist, size))

        if self.stepping:
            self.canvas.after(self.timestep, self.step_draw)

    def planet_x_y(self, hlon, dist):
        return (dist * self.dist_scale * math.cos(hlon) + self.halfwidth,
                dist * self.dist_scale * math.sin(hlon) + self.halfheight)


def print_table():
    """Super quickie hack to print out a table for the current opposition.
       Unrelated to any of the other code in this script.
       Ought to be generalized to take start and stop times, etc.
    """
    start_date = ephem.Date(datetime.now()) - ephem.hour*24*10
    end_date = start_date + ephem.hour*24*40
    opp_date, closest_date = find_next_opposition(start_date)
    print("Opposition:", opp_date)
    print("Closest approach:", closest_date)

    # Define "opposition season"
    d = opp_date - ephem.hour * 24 * 15
    end_date = opp_date + ephem.hour * 24 * 20

    mars = ephem.Mars()
    print(table_header)
    while d < start_date + 60:
        mars.compute(d)
        d += ephem.hour * 24
        if abs(d - opp_date) <= .5:
            print(table_format % (ephem.Date(d), mars.earth_distance,
                                  mars.size), "** OPPOSITION")
        elif abs(d - closest_date) <= .5:
            print(table_format % (ephem.Date(d), mars.earth_distance,
                                  mars.size), "** CLOSEST APPROACH")
        else:
            print(table_format % (ephem.Date(d), mars.earth_distance,
                                  mars.size))


if __name__ == '__main__':
    parser = argparse.ArgumentParser(
        description="""Draw planet orbits from the north ecliptic pole.

Key bindings:
  space   Start/stop animation
  q       quit""",
        formatter_class=argparse.RawTextHelpFormatter)
    parser.add_argument('-a', "--au", dest="auscale", type=float, default=1.7,
                        action="store",
                        help="""Scale of the window in astronomical units.""")

    parser.add_argument("-s", "--start", dest="start", default=None,
                        help="Start date, YYYY-MM-DD, "
                             "default: 30 days beforetoday")

    parser.add_argument('-t', "--timestep", dest="timestep",
                        type=int, default=30,
                        help="""Time step in milliseconds (default 30).
Controls how fast the orbits are drawn.""")

    parser.add_argument('-S', "--stopped", dest="stopped", action="store_true",
                        help="Bring up the window but don't immediately "
                              "start animating: wait for the spacebar.")

    parser.add_argument('-T', "--table", dest="table", action="store_true",
                        help="Forget all that graphic stuff and just "
                             "print a table of sizes around opposition")

    args = parser.parse_args(sys.argv[1:])
    if args.start:
        args.start = datetime.strptime(args.start, "%Y-%m-%d")

    if args.table:
        print_table()
        sys.exit(0)

    win = OrbitViewWindow(auscale=args.auscale, start_time=args.start,
                          timestep=args.timestep, stopped=args.stopped)

    mainloop()

