I have a software tool which calculate the synodic period with degree inputs with two planets. Synodic Period is the temporal interval that it takes for an object to reappear at the same point in relation to two other objects (linear nodes).

But I'm not sure how the calculation done behind the scene.

Here's the settings:

Planets: Sun and Moon (Geocentric Ecliptic)

Synodic Degree: 180

Begin: 12:00AM, 4 January 2022 GMT-5

End result: 1:23PM, 19 January 2022 GMT-5

How to calculate end result datetime with synodic degree between Sun and Moon?

  • $\begingroup$ It's not easy to accurately calculate the positions of Solar system bodies. We can't guess which method your software tool uses. Maybe it uses VSOP, or the SPICE Toolkit from JPL / NASA. $\endgroup$
    – PM 2Ring
    Jul 2, 2022 at 15:32
  • $\begingroup$ @PM2Ring Can it be done using astro.com/swisseph/swephprg.htm Swisspeh ? Accuracy is not my concern. If you know overall theory how to calculate, that's would be great. I would like to know more. $\endgroup$
    – Arcana
    Jul 3, 2022 at 16:15
  • $\begingroup$ Yes I had the stuffs I need practically. But I don't know how to compute the equation theoretically. $\endgroup$
    – Arcana
    Jul 4, 2022 at 2:11
  • $\begingroup$ I show how to compute the ecliptic longitude (and right ascension) of the Sun to reasonable accuracy in this answer. Computing the Moon's position is much harder. As I said there, lunar theory is complicated! To do what you want, you can make a rough estimate using the mean synodic month length (~29.53059 days), and then use software to do a more accurate calculation. $\endgroup$
    – PM 2Ring
    Jul 4, 2022 at 3:16
  • $\begingroup$ Let's assume the accuracy for Sun and Moon are correct, what's the formula to calculate the end time, or theoretically find the end time with brute force on ephemeris data? $\endgroup$
    – Arcana
    Jul 4, 2022 at 8:17

1 Answer 1


I'm pretty sure that you want the difference in ecliptic longitude of the Moon and Sun. That's what's used in the definition of the phases of the Moon. That is, New Moon (or Dark Moon) is at $0°$, First Quarter is at $90°$, Full Moon is at $180°$, and Last Quarter is at $270°$. This Moon phase angle is also known as the elongation of the Moon.

In Chat you said:

I write my own programming to run every day/hour/minute to get every degree angles from moon or sun. Then compute the formula. I'm not sure this is the way.

That's extreme brute force! Fortunately, there is an easier way.

Here's a plot of the Moon's elongation for the first synodic month of 2022, produced using JPL Horizons. It uses a time step of six hours. The numbers on the $X$ axis are Julian days.

Moon elongation, Jan 2022

As you can see, the graph is almost a straight line. If you have the $(x, y)$ coordinates of two points on that curve you can estimate the $x$ value for any $y$ value using linear interpolation.

We just need to rearrange the standard equation of a line in two-point form:

$$x = x_0 + \left(\frac{x_1-x_0}{y_1-y_0}\right)(y-y_0)$$

If $(x_0, y_0)$ and $(x_1, y_1)$ are two points on the curve, then $(x, y)$ lies on the line through those points, so it will be very close to the curve. So we can look up the $y$ value in the ephemeris corresponding to that new $x$ value.

We can then repeat the process, replacing whichever of $(x_0, y_0)$ and $(x_1, y_1)$ is furthest from the curve by our new $(x, y)$. We stop when $y$ is close enough to the $y$ value that we want.

There's a slight complication with using this method because elongation is an angle, so it wraps around at $360°=0°$. A simple way to deal with this is to use angles in the range $-180°$ to $180°$, and shift the elongation angle we're searching for to $0°$. Then we just need to use estimates $(x_0, y_0)$ and $(x_1, y_1)$ that are within two weeks of the correct time.

Here's some Python code which performs this process. It gets ephemeris data from Horizons, so it's quite accurate. However, it uses ecliptic longitudes for an observer at the centre of the Earth, and ignores atmospheric refraction. (Horizons can actually calculate longitudes for an observer on the surface of the Earth, if you give it the observer's longitude, latitude and altitude. It can also include an estimate of the effects of refraction).

""" Fetch Sun & Moon ecliptic longitude data from Horizons
    to compute the Moon's elongation (phase) for a given UTC date & time.
    Also, find the date & time when the elongation has changed by a given delta.

    Written by PM 2Ring 2022.07.07
    Run as Python!

import re, requests
from functools import lru_cache
from datetime import timedelta

url = "https://ssd.jpl.nasa.gov/api/horizons_file.api"
api_version = "1.0"

base_cmd = """

def fetch_data(target, day, verbose=False):
    ttype = "JD" if isinstance(day, float) else "CAL"
    cmd = f"""

    if verbose:
    req = requests.post(url, data={'format': 'text'}, files={'input': ('cmd', cmd)})
    m = re.search(r"API VERSION:\s*(\S*)", req.text)
    if m is None:
        print("Malformed command file, aborting")
        return None

    version = m.group(1)
    if version != api_version:
        print(f"Warning: API version is {version}, but this script is for {api_version}")

    m = re.search(r"(?s)\\\$\\\$SOE(.*)\\\$\\\$EOE", req.text)
    if m is None:
        print("NO EPHEMERIS")
        return None

    if verbose:
        print(' %' * 30)

    data = m.group(1)[1:]
    row = data.split(', ')
    return row[0], row[1], float(row[4])

# Horizons body IDs
# See https://ssd.jpl.nasa.gov/api/horizons.api?format=text&OBJ_DATA=YES&MAKE_EPHEM=NO&COMMAND=MB
# Sun: 10
# Moon: 301

def horizons_elong(day):
    """ Fetch Sun & Moon ecliptic longitude from Horizons
        and thence compute the Moon's elongation
    cal, jd, sun_pos = fetch_data(10, day)
    moon_pos = fetch_data(301, day)[2]
    elong = (moon_pos - sun_pos) % 360
    print(cal, jd)
    print("Sun  ", sun_pos)
    print("Moon ", moon_pos)
    print("elong", elong, "\n")
    return elong, jd

# Mean synodic month length
synodic_mean = 29.5305889

# Find x given y such that (x, y) is on the line
# through (x0, y0) & (x1, y1)
def interpolate(y, x0, y0, x1, y1):
    dy = y1 - y0
    if dy == 0:
        return x0
    return x0 + (y - y0) * (x1 - x0) / dy

def main(start='2000-Jan-06 18:13:38.114', delta=0, auto_update=False):
    # If start is purely numeric, treat it as a JD number
        start = float(start)
    except ValueError:
        start = start.strip()

    start_elong, start_jd = horizons_elong(start)

    if delta == 0:

    start_jd = float(start_jd)

    delta = float(delta)
    print("delta", delta)

    goal = (start_elong + delta) % 360
    print("goal ", goal, "\n")

    # Shift coordinates to put goal at 0
    def func(jd):
        return (180 + horizons_elong(jd)[0] - goal) % 360 - 180

    # Initial estimate
    x0 = start_jd + (delta / 360) * synodic_mean

    # Step forward an hour
    x1 = x0 + 1 / 24
    y0, y1 = func(x0), func(x1)

    # Refine the estimate
    for i in range(9):
        print("\nRound", i)

        x = interpolate(0, x0, y0, x1, y1)
        y = func(x)

        if abs(y) < 5e-7:

        # Replace the point which is furthest from the X axis
        if abs(y1) > abs(y0):
            x1, y1 = x, y
            x0, y0 = x, y

    d = x - start_jd
    print(f"\n{d} days")

Dates may be given in most of the time formats that Horizons accepts, except for Modified Julian Day number. However, they must be given in UTC, you cannot specify a different timezone or time scale. Note that Horizons uses the Julian calendar for dates prior to 1582-Oct-15.

If you give a standard Julian day number, do not give the JD prefix, just enter the plain number.

In theory, this program should work for any dates in the range 9999 BC to 9999 AD, but there seems to be a minor Horizons bug affecting dates earlier than JD $-999999$ (BC 7451-Feb-25).

delta can be any numeric value. It can also be a numeric expression, eg 5*360 + 90. If delta is zero, it just prints the Sun & Moon longitudes and Moon elongation for that date.

This program uses the current value of the mean synodic month length, $29.5305889$. The mean length is actually increasing, from ~$29.53056$ in 10,000 BC to ~$29.53061$ in 10,000 AD.

Here's the output of the program for the data given in the question.

 2022-Jan-04 05:00:00.000 2459583.708333333
Sun   283.8004537
Moon  303.8864265
elong 20.085972800000036

delta 180.0
goal  200.08597280000004

 2022-Jan-18 23:22:01.440 2459598.473627778
Sun   298.8430904
Moon  129.8813202
elong 191.0382298

 2022-Jan-19 00:22:01.440 2459598.515294444
Sun   298.8854946
Moon  130.3958492
elong 191.51035459999997

Round 0
 2022-Jan-19 18:31:51.408 2459599.272122778
Sun   299.6556721
Moon  139.7942933
elong 200.1386212


Round 1
 2022-Jan-19 18:25:12.408 2459599.267504722
Sun   299.6509729
Moon  139.7366397
elong 200.0856668


Round 2
 2022-Jan-19 18:25:14.714 2459599.267531412
Sun   299.6510001
Moon  139.7369729
elong 200.08597280000004


15.559198077768087 days
15 days, 13:25:14.713919

And here's a link to a live version running on the SageMathCell server. The program itself is plain Python. It uses some Sage features for getting the input data, but they can easily be replaced with standard Python functions if you want to run the program on the command line. It uses the popular requests library to fetch data from Horizons.

  • $\begingroup$ Thanks looking great. I have questions regarding to other celestial bodies. I realized planets before Earth (Mercury and Venus), has different maximum delta value depends on date. And planets after Earth (Mars above) doesn't has limit. Are those two types of calculation differ from Moon itself? Would like to learn more if you can provide website references. Thank you. $\endgroup$
    – Arcana
    Sep 16, 2022 at 6:10
  • $\begingroup$ @Arcana Yes, it's a bit different calculating elongation of planets because of apparent retrograde motion. $\endgroup$
    – PM 2Ring
    Sep 16, 2022 at 8:53
  • $\begingroup$ Is it very complex to calculate other planets? $\endgroup$
    – Arcana
    Sep 17, 2022 at 11:09
  • $\begingroup$ @Arcana No, it's just a bit messy. My technique assumes the elongation is approximately linear. I recommend using a graph to find the stationary points. $\endgroup$
    – PM 2Ring
    Sep 17, 2022 at 18:58

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .