<p>In <a href="https://rhodesmill.org/skyfield/" rel="nofollow noreferrer"><code>skyfield</code></a>, to find the times when something happens, you build a function that takes a <code>Time</code> object, outputs an integer value, and has a <code>step_days</code> attribute, and pass that function to the <a href="https://rhodesmill.org/skyfield/searches.html" rel="nofollow noreferrer"><code>find_discrete</code></a> method to find the times when your function transitions from one integer or boolean value to the next during its search period.</p>
<p>In the below script, I create a function <code>ecliptic_longitude_exceeds</code> that returns functions that can be used in the <code>find_discrete</code> method.  It's a simple function, that should work for the Sun or the Moon, because we can trust their ecliptic longitudes to be constantly increasing from 0° to 360°. Not as easy to work with planets (they'll backtrack, requiring more code to detect and ignore the 360° to 0° transition)  We'll use the year 2023 as our time window.</p>
<pre><code>from typing import Callable

import pytz
import skyfield.searchlib
from skyfield import api, timelib
from skyfield.jpllib import SpiceKernel, ChebyshevPosition
from skyfield.vectorlib import VectorSum
from skyfield.framelib import ecliptic_frame


def ecliptic_longitude_exceeds(
        longitude: float,
        target: str|VectorSum|ChebyshevPosition,
        ephemeris: SpiceKernel) -&gt; Callable[[timelib.Time], int]:
    &quot;&quot;&quot;
    Creates functions that check whether target ecliptic longitude exceeds
    value at a specified time
    :param longitude: Ecliptic Longitude in decimal degrees
    :param target: Object to be checked
    :param ephemeris: Ephemeris to be used.
    :return: find_discrete-compatible function
    &quot;&quot;&quot;

    earth = ephemeris['earth']
    target = ephemeris[target] if isinstance(target, str) else target

    def function(time: timelib.Time) -&gt; bool:
        &quot;&quot;&quot;
        :param time: Time of Observation
        :return: True if ecliptic longitude &gt;= longitude False otherwise
        &quot;&quot;&quot;
        observation = earth.at(time).observe(target).apparent()
        _, observed_longitude, _ = observation.frame_latlon(ecliptic_frame)
        return observed_longitude.degrees &gt;= longitude

    function.step_days = 60

    return function


def main():

    ephemeris = api.load('de421.bsp')
    ts = api.load.timescale()
    utc = pytz.timezone(&quot;UTC&quot;)
    start, stop = ts.utc(2023), ts.utc(2024)

    sun_exceeds = ecliptic_longitude_exceeds(
        longitude=56.1074, target=&quot;sun&quot;, ephemeris=ephemeris  
    )

    times, states = skyfield.searchlib.find_discrete(start, stop, sun_exceeds)

    longitude_times = list(time for time, in_state
                           in zip(times.astimezone(utc), states)
                           if in_state)

    print('\n'.join(str(lt) for lt in longitude_times))


if __name__ == '__main__':
    main()

</code></pre>
<p>This script returns:</p>
<pre><code>2023-05-17 06:11:22.543808+00:00
</code></pre>
<p>Quickly checking in my desktop copy of Stellarium puts the Sun's Ecliptic longitude at that time as 56°06'28.8&quot;, which matches up.</p>
<hr />
<p>That said, there was an issue with planets and other solar system bodies. Unlike the Sun. and the Moon, their ecliptic longitudes do not monotonically increase when observed from Earth, so sometimes they'll cross the specified ecliptic longitude in a decreasing direction, and the above script will miss those.</p>
<p>So let's try a different tack, using <code>skyfield</code>'s <code>find_minima</code> function.</p>
<p>In the below script, I create an <code>ecliptic_longitude_angle_difference</code> function that returns functions that can be used with the <code>find_minima</code> function. For the observation times it is sent, the function returns the interior angle between the specified ecliptic longitude and the object's ecliptic longitude.</p>
<p>We'll also switch the target to Mercury (to get that retrograde motion), and widen up the search window to 2020-2024.</p>
<pre><code>from typing import Callable

import numpy as np
import pytz
import skyfield.searchlib
from skyfield import api, timelib
from skyfield.jpllib import SpiceKernel, ChebyshevPosition
from skyfield.vectorlib import VectorSum
from skyfield.framelib import ecliptic_frame


def ecliptic_longitude_angle_difference(
        longitude: float,
        target: str | VectorSum | ChebyshevPosition,
        ephemeris: SpiceKernel,
        step_days: int = 25) -&gt; Callable[[timelib.Time], np.ndarray]:
    &quot;&quot;&quot;
    Produces functions for calculating minimal Ecliptic Longitude Distance.
    :param longitude: Ecliptic longitude value to examine.
    :param target: Name or Position of object to be observed from Earth.
    :param ephemeris: Ephemeris used in calculation.
    :param step_days: Initial step size for calculation, default 25 days.
    :return: find_minima-compatible function for determining ecliptic longitude crossings.
    &quot;&quot;&quot;
    earth = ephemeris['earth']
    target = ephemeris[target] if isinstance(target, str) else target

    def function(time: timelib.Time) -&gt; np.ndarray:
        &quot;&quot;&quot;
        :param time: timelib.Time of Observation or observations
        :return: Ecliptic longitude differences at times specified by the timelib.Time object
        &quot;&quot;&quot;
        observation = earth.at(time).observe(target).apparent()
        _, observed_longitude, _ = observation.frame_latlon(ecliptic_frame)

        angle = abs(observed_longitude.degrees - longitude)
        return np.where(angle &lt;= 180, angle, 360 - angle)

    function.step_days = step_days

    return function


def main():
    ephemeris = api.load('de421.bsp')
    ts = api.load.timescale()
    utc = pytz.timezone(&quot;UTC&quot;)
    start, stop = ts.utc(2020), ts.utc(2024)

    target_difference = ecliptic_longitude_angle_difference(
        longitude=56.1074, target=&quot;mercury&quot;, ephemeris=ephemeris
    )

    times, differences = skyfield.searchlib.find_minima(start, stop, target_difference)

    longitude_times = list(time for time, difference in zip(times.astimezone(utc), differences) 
                           if difference &lt; 0.0001)
    print('\n'.join(str(lt) for lt in longitude_times))


if __name__ == '__main__':
    main()
</code></pre>
<p>This produces the following output for Mercury crossing 56.1074° ecliptic longitude:</p>
<pre><code>2020-05-10 02:13:01.258755+00:00
2021-05-02 00:31:27.196698+00:00
2022-04-26 07:10:46.338841+00:00
2022-06-02 14:15:34.412501+00:00
2022-06-04 01:41:10.667429+00:00
2023-06-08 22:02:27.347301+00:00
</code></pre>
<p>Which we can confirm with Stellarium, looking specifically at Mercury's position on 2022-06-02, and a daily image of its progression on the celestial sphere during spring 2022.</p>
<p><a href="https://i.sstatic.net/yMza5.png" rel="nofollow noreferrer"><img src="https://i.sstatic.net/yMza5.png" alt="Mercury's path on the Equatorial Grid in spring 2022." /></a></p>

CC BY-SA 4.0
2023-09-06T23:11:05.367Z
you are amazing it works like a charm <3, i installed Stellarium on my laptop and i will use it as reference ..
dimitri33
CC BY-SA 4.0
2023-09-28T17:58:54.137Z
hey , i tried this for all the planets and as long as they are not in "retrograde" this works perfect .. now the question is can one calculate the same for the "true node"  or nodes in general ? i tried moon_nodes and got "unknown SPICE target
dimitri33
CC BY-SA 4.0
2023-09-06T02:43:48.690Z
2023-12-26T15:20:09.907Z
2023-12-26T15:20:09.907Z
notovny
<p>You had the right idea with <a href="https://astrorigin.com/pyswisseph/sphinx/programmers_manual/crossings_of_planets_over_positions.html" rel="nofollow noreferrer"><code>solcross</code> and <code>mooncross</code></a> from <code>swisseph</code>.
The functions <a href="https://astrorigin.com/pyswisseph/sphinx/programmers_manual/date_and_time/utc_and_julian_day.html" rel="nofollow noreferrer"><code>utc_to_jd</code> and <code>jdet_to_utc</code></a> can help convert from a date+time to a Julian date and back.</p>
<p>This code:</p>
<pre><code>import swisseph as swe
from datetime import datetime, timezone

def jd_from_dt(dt_ut):
    &quot;&quot;&quot;Converts a datetime in UTC to a Julian date in TT.&quot;&quot;&quot;
    (jd_tt, jd_ut) = swe.utc_to_jd(dt_ut.year, dt_ut.month, dt_ut.day,
            dt_ut.hour, dt_ut.minute, dt_ut.second)
    return jd_tt

def dt_from_jd(jd_tt):
    &quot;&quot;&quot;Converts a Julian date in TT to a datetime in UTC.&quot;&quot;&quot;
    (yr, mo, day, h, m, s) = swe.jdet_to_utc(jd_tt)
    dt_ut = datetime(yr, mo, day, h, m, round(s), tzinfo=timezone.utc)
    return dt_ut

# examples from question

post_jd = jd_from_dt(datetime(2023, 9, 5, 23, tzinfo=timezone.utc))

sun1_lon = 150 + 13 + 7/60 + 9/3600
sun1a_jd = swe.solcross(sun1_lon, post_jd - 30)
sun1b_jd = swe.solcross(sun1_lon, post_jd + 335)
print('Sun at lon={0:.4f}:\n  {1}\n  {2}'.format(
        sun1_lon, dt_from_jd(sun1a_jd), dt_from_jd(sun1b_jd)))

moon1_lon = 60 + 1 + 33/60 + 23/3600
moon1a_jd = swe.mooncross(moon1_lon, post_jd - 3)
moon1b_jd = swe.mooncross(moon1_lon, post_jd + 24)
print('Moon at lon={0:.4f}:\n  {1}\n  {2}'.format(
        moon1_lon, dt_from_jd(moon1a_jd), dt_from_jd(moon1b_jd)))

moon2_lon = 60 + 5 + 45/60
moon2a_jd = swe.mooncross(moon2_lon, post_jd - 3)
moon2b_jd = swe.mooncross(moon2_lon, post_jd + 24)
print('Moon at lon={0:.4f}:\n  {1}\n  {2}'.format(
        moon2_lon, dt_from_jd(moon2a_jd), dt_from_jd(moon2b_jd)))

# example from notovny's answer

jan1_jd = jd_from_dt(datetime(2023, 1, 1, 0, tzinfo=timezone.utc))

sun2_lon = 56.1074
sun2a_jd = swe.solcross(sun2_lon, jan1_jd)
sun2b_jd = swe.solcross(sun2_lon, jan1_jd + 365)
print('Sun at lon={0:.4f}:\n  {1}\n  {2}'.format(
        sun2_lon, dt_from_jd(sun2a_jd), dt_from_jd(sun2b_jd)))
</code></pre>
<p>produces this output:</p>
<pre><code>Sun at lon=163.1192:
  2023-09-05 22:54:26+00:00
  2024-09-05 04:39:52+00:00
Moon at lon=61.5564:
  2023-09-05 22:58:23+00:00
  2023-10-03 07:50:07+00:00
Moon at lon=65.7500:
  2023-09-06 06:44:16+00:00
  2023-10-03 15:23:09+00:00
Sun at lon=56.1074:
  2023-05-17 06:11:23+00:00
  2024-05-16 11:57:53+00:00
</code></pre>

CC BY-SA 4.0
2023-09-07T15:52:14.790Z
thank you Mike, this changes things a little bit <3
dimitri33
CC BY-SA 4.0
2023-09-06T22:52:44.763Z
2023-09-07T15:09:53.600Z
2023-09-07T15:09:53.600Z
Mike G
<p>the Sun is currently in Virgo ♍, 13°7'9&quot; or longitude 56.0762 degrees
how does one calculate when exactly the Sun will reach this position again. Since the moon is doing a spin in 29.x days lets use that as example also, &quot;she&quot; is in Gemini ♊, 1°33'23&quot; or  56.1074 degrees
,  can we calculate when will it reach Gemini 05° 45' 00&quot; or 56.1263 degrees ( a random coordination)
Any math or programming language will be just fine even tho I messed with swisseph ( solcross, mooncross) and skyfield and no result ...</p>

CC BY-SA 4.0
2023-09-05T23:30:28.617Z
The longitudes in Virgo and Gemini cannot both be 56 degrees. The Sun is approaching the Sep equinox, so the longitude should be close to 180 degrees. (Otherwise, the question is still valid about crossing a specific longitude.)
JohnHoltz
CC BY-SA 4.0
2023-09-06T00:45:29.263Z
thanks john, i was using a python code to do that so i guess its wrong .. :(
dimitri33
CC BY-SA 4.0
2023-09-06T05:59:51.667Z
At the time of posting, the sun was in Leo, at an ecliptic long of +162 degrees, and the moon was in  Taurus at 62 degrees. The sun is only in "virgo" in the astrological sense, which is no longer related to the actual positions of the stars in the sky,
James K
<the-sun><the-moon>
Sun or moon, finding the date and time of crossing a specific longitude
dimitri33
# Sun or moon, finding the date and time of crossing a specific longitude

(Asked by: username_0 on 2023-09-05)

the Sun is currently in Virgo ♍, 13°7'9" or longitude 56.0762 degrees
how does one calculate when exactly the Sun will reach this position again. Since the moon is doing a spin in 29.x days lets use that as example also, "she" is in Gemini ♊, 1°33'23" or 56.1074 degrees
, can we calculate when will it reach Gemini 05° 45' 00" or 56.1263 degrees ( a random coordination)
Any math or programming language will be just fine even tho I messed with swisseph ( solcross, mooncross) and skyfield and no result ...

**Tags:** the-sun, the-moon

> **Comments:**
> **username_1**: The longitudes in Virgo and Gemini cannot both be 56 degrees. The Sun is approaching the Sep equinox, so the longitude should be close to 180 degrees. (Otherwise, the question is still valid about crossing a specific longitude.)  
> **username_0**: thanks john, i was using a python code to do that so i guess its wrong .. :(  
> **username_2**: At the time of posting, the sun was in Leo, at an ecliptic long of +162 degrees, and the moon was in  Taurus at 62 degrees. The sun is only in "virgo" in the astrological sense, which is no longer related to the actual positions of the stars in the sky,  

---

## 2 answers

### Answer by username_3 on 2023-09-06. Score: 3

In [`skyfield`](https://rhodesmill.org/skyfield/), to find the times when something happens, you build a function that takes a `Time` object, outputs an integer value, and has a `step_days` attribute, and pass that function to the [`find_discrete`](https://rhodesmill.org/skyfield/searches.html) method to find the times when your function transitions from one integer or boolean value to the next during its search period.

In the below script, I create a function `ecliptic_longitude_exceeds` that returns functions that can be used in the `find_discrete` method. It's a simple function, that should work for the Sun or the Moon, because we can trust their ecliptic longitudes to be constantly increasing from 0° to 360°. Not as easy to work with planets (they'll backtrack, requiring more code to detect and ignore the 360° to 0° transition) We'll use the year 2023 as our time window.

```
from typing import Callable

import pytz
import skyfield.searchlib
from skyfield import api, timelib
from skyfield.jpllib import SpiceKernel, ChebyshevPosition
from skyfield.vectorlib import VectorSum
from skyfield.framelib import ecliptic_frame

def ecliptic_longitude_exceeds(
        longitude: float,
        target: str|VectorSum|ChebyshevPosition,
        ephemeris: SpiceKernel) -> Callable[[timelib.Time], int]:
    """
    Creates functions that check whether target ecliptic longitude exceeds
    value at a specified time
    :param longitude: Ecliptic Longitude in decimal degrees
    :param target: Object to be checked
    :param ephemeris: Ephemeris to be used.
    :return: find_discrete-compatible function
    """

    earth = ephemeris['earth']
    target = ephemeris[target] if isinstance(target, str) else target

    def function(time: timelib.Time) -> bool:
        """
        :param time: Time of Observation
        :return: True if ecliptic longitude >= longitude False otherwise
        """
        observation = earth.at(time).observe(target).apparent()
        _, observed_longitude, _ = observation.frame_latlon(ecliptic_frame)
        return observed_longitude.degrees >= longitude

    function.step_days = 60

    return function

def main():

    ephemeris = api.load('de421.bsp')
    ts = api.load.timescale()
    utc = pytz.timezone("UTC")
    start, stop = ts.utc(2023), ts.utc(2024)

    sun_exceeds = ecliptic_longitude_exceeds(
        longitude=56.1074, target="sun", ephemeris=ephemeris  
    )

    times, states = skyfield.searchlib.find_discrete(start, stop, sun_exceeds)

    longitude_times = list(time for time, in_state
                           in zip(times.astimezone(utc), states)
                           if in_state)

    print('\n'.join(str(lt) for lt in longitude_times))

if __name__ == '__main__':
    main()

```

This script returns:

```
2023-05-17 06:11:22.543808+00:00

```

Quickly checking in my desktop copy of Stellarium puts the Sun's Ecliptic longitude at that time as 56°06'28.8", which matches up.

---

That said, there was an issue with planets and other solar system bodies. Unlike the Sun. and the Moon, their ecliptic longitudes do not monotonically increase when observed from Earth, so sometimes they'll cross the specified ecliptic longitude in a decreasing direction, and the above script will miss those.

So let's try a different tack, using `skyfield`'s `find_minima` function.

In the below script, I create an `ecliptic_longitude_angle_difference` function that returns functions that can be used with the `find_minima` function. For the observation times it is sent, the function returns the interior angle between the specified ecliptic longitude and the object's ecliptic longitude.

We'll also switch the target to Mercury (to get that retrograde motion), and widen up the search window to 2020-2024.

```
from typing import Callable

import numpy as np
import pytz
import skyfield.searchlib
from skyfield import api, timelib
from skyfield.jpllib import SpiceKernel, ChebyshevPosition
from skyfield.vectorlib import VectorSum
from skyfield.framelib import ecliptic_frame

def ecliptic_longitude_angle_difference(
        longitude: float,
        target: str | VectorSum | ChebyshevPosition,
        ephemeris: SpiceKernel,
        step_days: int = 25) -> Callable[[timelib.Time], np.ndarray]:
    """
    Produces functions for calculating minimal Ecliptic Longitude Distance.
    :param longitude: Ecliptic longitude value to examine.
    :param target: Name or Position of object to be observed from Earth.
    :param ephemeris: Ephemeris used in calculation.
    :param step_days: Initial step size for calculation, default 25 days.
    :return: find_minima-compatible function for determining ecliptic longitude crossings.
    """
    earth = ephemeris['earth']
    target = ephemeris[target] if isinstance(target, str) else target

    def function(time: timelib.Time) -> np.ndarray:
        """
        :param time: timelib.Time of Observation or observations
        :return: Ecliptic longitude differences at times specified by the timelib.Time object
        """
        observation = earth.at(time).observe(target).apparent()
        _, observed_longitude, _ = observation.frame_latlon(ecliptic_frame)

        angle = abs(observed_longitude.degrees - longitude)
        return np.where(angle <= 180, angle, 360 - angle)

    function.step_days = step_days

    return function

def main():
    ephemeris = api.load('de421.bsp')
    ts = api.load.timescale()
    utc = pytz.timezone("UTC")
    start, stop = ts.utc(2020), ts.utc(2024)

    target_difference = ecliptic_longitude_angle_difference(
        longitude=56.1074, target="mercury", ephemeris=ephemeris
    )

    times, differences = skyfield.searchlib.find_minima(start, stop, target_difference)

    longitude_times = list(time for time, difference in zip(times.astimezone(utc), differences) 
                           if difference < 0.0001)
    print('\n'.join(str(lt) for lt in longitude_times))

if __name__ == '__main__':
    main()

```

This produces the following output for Mercury crossing 56.1074° ecliptic longitude:

```
2020-05-10 02:13:01.258755+00:00
2021-05-02 00:31:27.196698+00:00
2022-04-26 07:10:46.338841+00:00
2022-06-02 14:15:34.412501+00:00
2022-06-04 01:41:10.667429+00:00
2023-06-08 22:02:27.347301+00:00

```

Which we can confirm with Stellarium, looking specifically at Mercury's position on 2022-06-02, and a daily image of its progression on the celestial sphere during spring 2022.

[![Mercury's path on the Equatorial Grid in spring 2022.](https://i.sstatic.net/yMza5.png)](https://i.sstatic.net/yMza5.png)

> **Comments:**
> **username_0**: you are amazing it works like a charm <3, i installed Stellarium on my laptop and i will use it as reference ..  
> **username_0**: hey , i tried this for all the planets and as long as they are not in "retrograde" this works perfect .. now the question is can one calculate the same for the "true node"  or nodes in general ? i tried moon_nodes and got "unknown SPICE target  

---

### Answer by username_4 on 2023-09-06. Score: 2

You had the right idea with [`solcross` and `mooncross`](https://astrorigin.com/pyswisseph/sphinx/programmers_manual/crossings_of_planets_over_positions.html) from `swisseph`.
The functions [`utc_to_jd` and `jdet_to_utc`](https://astrorigin.com/pyswisseph/sphinx/programmers_manual/date_and_time/utc_and_julian_day.html) can help convert from a date+time to a Julian date and back.

This code:

```
import swisseph as swe
from datetime import datetime, timezone

def jd_from_dt(dt_ut):
    """Converts a datetime in UTC to a Julian date in TT."""
    (jd_tt, jd_ut) = swe.utc_to_jd(dt_ut.year, dt_ut.month, dt_ut.day,
            dt_ut.hour, dt_ut.minute, dt_ut.second)
    return jd_tt

def dt_from_jd(jd_tt):
    """Converts a Julian date in TT to a datetime in UTC."""
    (yr, mo, day, h, m, s) = swe.jdet_to_utc(jd_tt)
    dt_ut = datetime(yr, mo, day, h, m, round(s), tzinfo=timezone.utc)
    return dt_ut

# examples from question

post_jd = jd_from_dt(datetime(2023, 9, 5, 23, tzinfo=timezone.utc))

sun1_lon = 150 + 13 + 7/60 + 9/3600
sun1a_jd = swe.solcross(sun1_lon, post_jd - 30)
sun1b_jd = swe.solcross(sun1_lon, post_jd + 335)
print('Sun at lon={0:.4f}:\n  {1}\n  {2}'.format(
        sun1_lon, dt_from_jd(sun1a_jd), dt_from_jd(sun1b_jd)))

moon1_lon = 60 + 1 + 33/60 + 23/3600
moon1a_jd = swe.mooncross(moon1_lon, post_jd - 3)
moon1b_jd = swe.mooncross(moon1_lon, post_jd + 24)
print('Moon at lon={0:.4f}:\n  {1}\n  {2}'.format(
        moon1_lon, dt_from_jd(moon1a_jd), dt_from_jd(moon1b_jd)))

moon2_lon = 60 + 5 + 45/60
moon2a_jd = swe.mooncross(moon2_lon, post_jd - 3)
moon2b_jd = swe.mooncross(moon2_lon, post_jd + 24)
print('Moon at lon={0:.4f}:\n  {1}\n  {2}'.format(
        moon2_lon, dt_from_jd(moon2a_jd), dt_from_jd(moon2b_jd)))

# example from username_3's answer

jan1_jd = jd_from_dt(datetime(2023, 1, 1, 0, tzinfo=timezone.utc))

sun2_lon = 56.1074
sun2a_jd = swe.solcross(sun2_lon, jan1_jd)
sun2b_jd = swe.solcross(sun2_lon, jan1_jd + 365)
print('Sun at lon={0:.4f}:\n  {1}\n  {2}'.format(
        sun2_lon, dt_from_jd(sun2a_jd), dt_from_jd(sun2b_jd)))

```

produces this output:

```
Sun at lon=163.1192:
  2023-09-05 22:54:26+00:00
  2024-09-05 04:39:52+00:00
Moon at lon=61.5564:
  2023-09-05 22:58:23+00:00
  2023-10-03 07:50:07+00:00
Moon at lon=65.7500:
  2023-09-06 06:44:16+00:00
  2023-10-03 15:23:09+00:00
Sun at lon=56.1074:
  2023-05-17 06:11:23+00:00
  2024-05-16 11:57:53+00:00

```

> **Comments:**
> **username_0**: thank you username_4, this changes things a little bit <3  

---

