skyfield Phases of the Moon

[JeffRossMT]

Hello!

I have something odd that could very well be my failure to understand how to get the current phase of the moon via skyfield.

I'm using https://svs.gsfc.nasa.gov/4955 as my authoritative source, using Month: 3, Day: 13, UT Hour: 21 as the input for the Moon Phase and Libration, 2022 page.

Time	Sunday, March 13, 2022, 21:00 UT
Phase	80.4% (11d 3h 25m)
Diameter	1795.6 arcseconds
Distance	399150 km (31.33 Earth diameters)
J2000 Right Ascension, Declination	8h 14m 57s, 24° 35' 3"
Subsolar Longitude, Latitude	48.809°, -1.411°
Sub-Earth Longitude, Latitude	-3.650°, -6.063°
Position Angle	13.478°

Here's the short python3 script I'm using:

#!/usr/bin/env python3

from skyfield import almanac
from skyfield.api import N, W, wgs84, load
from skyfield.framelib import ecliptic_frame
import subprocess
import ephem

m = ephem.Moon()
current_time = subprocess.check_output("psql -X -A -t -q -d postgres -c \"select date_trunc('hour',now() at time zone 'UTC')\"",stderr=subprocess.STDOUT, shell=True)[:-1].decode('utf-8')
print("Current UTC Time: %s" % (current_time))

print("pyephem Previous New Moon at UTC: %s" % (ephem.previous_new_moon(current_time)))
m.compute(current_time)
print("pyepehm Moon Percent Illuminated at Current UTC Time: %s" % (str(round(m.phase,2)) + '%'))

ts = load.timescale()
t = ts.now()
eph = load('de421.bsp')
sun, moon, earth = eph['sun'], eph['moon'], eph['earth']
phase = almanac.moon_phase(eph,t)
print("skyfield Phase of the Moon in degrees: %s" % (str(phase).replace('deg',u'\u00B0').replace(' ','')))
print('skyfield Moon phase: {:.1f} degrees'.format(phase.degrees))

e = earth.at(t)
_, slon, _ = e.observe(sun).apparent().frame_latlon(ecliptic_frame)
_, mlon, _ = e.observe(moon).apparent().frame_latlon(ecliptic_frame)
phase = (mlon.degrees - slon.degrees) % 360.0
print("ecliptic frame Moon phase %s" % (phase))
print("skyfield formatted Moon phase %s" % ('{0:.1f}'.format(phase)))

Running that I get:

jross@aurora-cam:~$ python3 get_moon_position.py
Current UTC Time: 2022-03-13 21:00:00
pyephem Previous New Moon at UTC: 2022/3/2 17:34:45
pyepehm Moon Percent Illuminated at Current UTC Time: 80.41%
skyfield Phase of the Moon in degrees: 127°35'49.2"
skyfield Moon phase: 127.6 degrees
ecliptic frame Moon phase 127.59701283771142
skyfield formatted Moon phase 127.6

pypehem's percent illuminated matches NASA's but I can't see any way to make skyfield's Moon phase result come anywhere near the Waxing Gibbous phase that the Moon is currently at.

What am I missing?

Thanks,

Jeff Ross

[brandon-rhodes]

Those are two different numbers that mean different things. Where in the Skyfield documentation could the description of the Moon phase be improved to clarify that it's the Moon's angle around the sky, not its percent illumination?

To learn the fraction of its surface that's illuminated, try the fraction_illuminated() Skyfield method:

https://rhodesmill.org/skyfield/api-position.html#skyfield.positionlib.ICRF.fraction_illuminated

See if that gives you a result that matches more closely!

[JeffRossMT]

Thank you Brandon! I did a poor job of reading and understanding the documentation. Phase of the Moon in my mind equates to New, Waxing Crescent and so on and that's how pyephem uses it.

Perhaps add a line like this:

jross@MacMini:/Users/jross $ diff -U3 examples.html.orig examples.html
--- examples.html.orig  2022-03-02 22:28:06.000000000 -0700
+++ examples.html   2022-03-18 12:27:59.000000000 -0600
@@ -138,6 +138,7 @@
 90° at the First Quarter,
 180° at the Full Moon,
 and 270° at the Last Quarter.</p>
+<p>To find what portion of the Moon is illuminated, see the <a href = "https://rhodesmill.org/skyfield/api-position.html#skyfield.positionlib.ICRF.fraction_illuminated">fraction illuminated</a> routine.</p>
 <div class="highlight-python notranslate"><div class="highlight"><pre><span></span><span class="kn">from</span> <span class="nn">skyfield.api</span> <span class="kn">import</span> <span class="n">load</span>
 <span class="kn">from</span> <span class="nn">skyfield.framelib</span> <span class="kn">import</span> <span class="n">ecliptic_frame</span>

@@ -456,4 +457,4 @@
 </div>

   </body>
-</html>
\ No newline at end of file
+</html>

What module do I import to use fraction_illuminated()?

Jeff

[brandon-rhodes]

What module do I import to use fraction_illuminated()?

Happily, it's a method on the position you already have; no need to import anything (assuming you have the most recent 1.42 version of Skyfield).

[hanetzer]

Exactly how would one go about using this, in context of the above script? What 'object' or variable is it called from? (fraction_illuminated)

[brandon-rhodes]

Exactly how would one go about using this, in context of the above script?

The script would need to keep the Apparent position object, instead of throwing it away after the frame_latlon() call. It would look something like:

m = e.observe(moon).apparent()
fraction = m.fraction_illuminated(sun)

[hanetzer]

Got it. I ended up going with just using the new moon/conjunction date/times spat out from the almanac library.

[brandon-rhodes]

@JeffRossMT @hanetzer — I've added the fraction of the moon that's illuminated to that example in the documentation, as you can see from the commit just above this comment. The new docs should appear with the next Skyfield release. Thanks for the idea!