diff --git a/examples/pyerfa/README.md b/examples/pyerfa/README.md new file mode 100644 index 0000000..7d18dac --- /dev/null +++ b/examples/pyerfa/README.md @@ -0,0 +1,18 @@ +# pyerfa Examples + +Each sub-directory contains a self-contained example. The order in +which the examples are to appear is specified in `order.json` (an +array of directory names in the expected order). + +In each example directory you'll find: + +* `config.toml` - must conform to the specification outlined here: + https://docs.pyscript.net/latest/user-guide/configuration/ This is + parsed and ultimately turned into a JSON representation as part of + the package's API object. +* `setup.py` - Python code for contextual and environmental setup, + NOT SEEN BY THE END USER, but is run before the `code.py` code is + evaluated. Allows us to create useful (IPython) shims, avoid + repeating boilerplate and whatnot. +* `code.py` - the actual code added to the editor which forms the + practical example of using the package. diff --git a/examples/pyerfa/angles_and_separations/code.py b/examples/pyerfa/angles_and_separations/code.py new file mode 100644 index 0000000..31d8010 --- /dev/null +++ b/examples/pyerfa/angles_and_separations/code.py @@ -0,0 +1,71 @@ +# --------------------------------------------------------------------- +# Sexagesimal parsing, angular separations, and great-circle distances. +# --------------------------------------------------------------------- +# +# ERFA includes a rich set of angle utilities. Here we'll parse +# right-ascension and declination strings, then use erfa.seps to +# compute the angular separation between bright stars. +import numpy as np +import erfa + + +heading("Parsing RA/Dec strings into radians") +note( + "erfa.tf2a converts a sign + hours/minutes/seconds tuple to " + "radians (hours form), and erfa.af2a does the same for " + "degrees/arcminutes/arcseconds." +) + +stars = [ + # (name, RA h/m/s, Dec sign/d/m/s) + ("Sirius", ( 6, 45, 8.92), ("-", 16, 42, 58.0)), + ("Canopus", ( 6, 23, 57.11), ("-", 52, 41, 44.4)), + ("Vega", (18, 36, 56.34), ("+", 38, 47, 1.3)), + ("Betelgeuse", (5, 55, 10.31), ("+", 7, 24, 25.4)), + ("Rigel", ( 5, 14, 32.27), ("-", 8, 12, 6.0)), +] + +# Build arrays of RA and Dec in radians. +names = [s[0] for s in stars] +ra_rad = np.array([erfa.tf2a("+", h, m, s) for _, (h, m, s), _ in stars]) +dec_rad = np.array([erfa.af2a(sign, d, m, s) for _, _, (sign, d, m, s) in stars]) + +note("Star catalogue (RA, Dec converted to degrees for display):") +rows = ["
{text}
"), append=True) + diff --git a/examples/pyerfa/calendar_and_constants/code.py b/examples/pyerfa/calendar_and_constants/code.py new file mode 100644 index 0000000..ef8cee4 --- /dev/null +++ b/examples/pyerfa/calendar_and_constants/code.py @@ -0,0 +1,59 @@ +""" +A first look at PyERFA: the Python wrapper for the ERFA C library +of fundamental astronomy routines (the open-source counterpart of +the IAU's SOFA library). + +Every ERFA function is exposed as a NumPy universal function, so +they accept scalars or arrays interchangeably. See the docs at +https://pyerfa.readthedocs.io/ for the full catalogue. +""" +from IPython.core.display import display, HTML + +import numpy as np +import erfa + + +# ERFA exposes a number of useful astronomical constants from erfam.h. +heading("Some ERFA constants") +note( + "These come straight from the ERFA C library and are handy " + "when working with time and angles." +) +constants = { + "DAYSEC (seconds in a day)": erfa.DAYSEC, + "DJY (days in a Julian year)": erfa.DJY, + "DAU (astronomical unit, metres)": erfa.DAU, + "CMPS (speed of light, m/s)": erfa.CMPS, + "DR2D (radians to degrees)": erfa.DR2D, +} +rows = "".join(f"{text}
"), append=True) diff --git a/examples/pyerfa/order.json b/examples/pyerfa/order.json new file mode 100644 index 0000000..2d70c9f --- /dev/null +++ b/examples/pyerfa/order.json @@ -0,0 +1,5 @@ +[ + "calendar_and_constants", + "planet_positions", + "angles_and_separations" +] diff --git a/examples/pyerfa/planet_positions/code.py b/examples/pyerfa/planet_positions/code.py new file mode 100644 index 0000000..b0b48bb --- /dev/null +++ b/examples/pyerfa/planet_positions/code.py @@ -0,0 +1,61 @@ +# --------------------------------------------------------------------- +# Plotting the inner planets' orbits using erfa.plan94. +# --------------------------------------------------------------------- +# +# erfa.plan94 gives heliocentric position and velocity (in AU and +# AU/day) for planets 1..8 (Mercury through Neptune) using an +# analytical approximation accurate enough for many purposes. + +import numpy as np +import matplotlib.pyplot as plt +import erfa + + +heading("Tracing the inner planets over one Earth-year") +note( + "We sample plan94 every few days for a full year starting " + "on JD 2460000.5 and project the heliocentric positions onto " + "the ecliptic plane (X, Y in AU)." +) + +planets = { + 1: ("Mercury", "tab:gray"), + 2: ("Venus", "tab:orange"), + 3: ("Earth", "tab:blue"), + 4: ("Mars", "tab:red"), +} + +jd_base = 2460000.5 +days = np.arange(0, 700, 2.0) # enough to close Mars' orbit too + +fig, ax = plt.subplots(figsize=(6, 6)) +ax.plot(0, 0, marker="*", color="gold", markersize=18, label="Sun") + +for planet_id, (name, color) in planets.items(): + pv = erfa.plan94(jd_base, days, planet_id) + # pv is a structured array with fields 'p' (position) and 'v' (velocity). + x = pv["p"][:, 0] + y = pv["p"][:, 1] + ax.plot(x, y, color=color, linewidth=1.2, label=name) + # Mark the starting position. + ax.plot(x[0], y[0], "o", color=color, markersize=5) + +ax.set_aspect("equal") +ax.set_xlabel("X (AU)") +ax.set_ylabel("Y (AU)") +ax.set_title("Inner planets, heliocentric ecliptic plane") +ax.legend(loc="upper right", fontsize=9) +ax.grid(alpha=0.3) +fig.tight_layout() +display(fig, append=True) + +# How far is each planet from the Sun on day zero? +heading("Heliocentric distances on the start date") +pv0 = erfa.plan94(jd_base, 0.0, 3) # Earth as a check: should be ~1 AU +note(f"Earth's distance from the Sun on JD {jd_base}: " + f"{np.linalg.norm(pv0['p']):.4f} AU") + +for planet_id, (name, _) in planets.items(): + pv = erfa.plan94(jd_base, 0.0, planet_id) + r = np.linalg.norm(pv["p"]) + note(f"{name}: {r:.4f} AU") diff --git a/examples/pyerfa/planet_positions/config.toml b/examples/pyerfa/planet_positions/config.toml new file mode 100644 index 0000000..1cad054 --- /dev/null +++ b/examples/pyerfa/planet_positions/config.toml @@ -0,0 +1 @@ +packages = ["pyerfa", "numpy", "matplotlib"] diff --git a/examples/pyerfa/planet_positions/setup.py b/examples/pyerfa/planet_positions/setup.py new file mode 100644 index 0000000..4b05c08 --- /dev/null +++ b/examples/pyerfa/planet_positions/setup.py @@ -0,0 +1,19 @@ +"""Lightweight setup for example 2 (no IPython shim needed).""" +import js +from pyscript import window, HTML, display as _display + +js.alert = window.alert + + +def display(*args, **kwargs): + return _display( + *args, **kwargs, target=__pyscript_display_target__, + ) + + +def heading(text, level=2): + display(HTML(f"{text}
"), append=True)