Skip to content

Sun, Moon and Ephemerides

satkit has two sources of solar-system positions: the JPL Development Ephemerides (jplephem), accurate to meters but read from a large file, and analytic low-precision models (lpephem), which need no data and are good to tens of arcseconds (Sun) and a few arcminutes (Moon). The Planetary Ephemerides tutorial compares them.

Low-precision Sun and Moon: lpephem

  • lpephem::sun::pos_gcrf(&t) and lpephem::moon::pos_gcrf(&t): geocentric GCRF positions in meters (pos_mod for the mean-of-date frame).
  • sun::riseset(&t, &coord, sigma): sunrise and sunset on the UTC calendar date of t at an ITRFCoord, as UTC Instants. sigma is the Sun's zenith angle at the event in degrees, None for the standard 90° 50′; 96°, 102° and 108° give civil, nautical and astronomical twilight. Where the Sun stays up or down all day it returns lpephem::Error::NoSunriseOrSunset.
  • moon::phase(&t), moon::illumination(&t) and moon::phase_name(&t): the Sun–Moon elongation, the illuminated fraction, and the named phase.
  • sun::shadowfunc(&sun_pos, &sat_pos): the fraction of the Sun's disc visible from a satellite (1 in sunlight, 0 in the umbra, in between in the penumbra).
  • lpephem::heliocentric_pos(body, &t): approximate heliocentric planet positions from Keplerian elements.
//! Low-precision (analytic) Sun and Moon: positions, sunrise and sunset,
//! Moon phase, and the Earth-shadow function.
//!
//! Runs offline with no data files.
//!
//! ```text
//! cargo run --example sun_moon
//! ```

use satkit::lpephem::{self, moon, sun};
use satkit::{ITRFCoord, Instant};

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let t = Instant::from_datetime(2024, 6, 20, 12, 0, 0.0)?;

    // Geocentric positions in GCRF, meters
    let psun = sun::pos_gcrf(&t);
    let pmoon = moon::pos_gcrf(&t);
    println!("Sun  distance = {:.6} AU", psun.norm() / satkit::consts::AU);
    println!("Moon distance = {:.0} km", pmoon.norm() / 1e3);

    // Moon phase: Sun-Moon elongation in ecliptic longitude
    println!(
        "Moon phase {:.1} deg ({}), {:.1}% illuminated",
        moon::phase(&t).to_degrees(),
        moon::phase_name(&t).name(),
        100.0 * moon::illumination(&t)
    );

    // Sunrise and sunset (UTC) on the input's UTC calendar date. `None`
    // is the standard 90 deg 50' zenith angle; 96, 102 and 108 deg give
    // civil, nautical and astronomical twilight.
    let boston = ITRFCoord::from_geodetic_deg(42.3601, -71.0589, 0.0);
    let (rise, set) = sun::riseset(&t, &boston, None)?;
    println!(
        "Boston sunrise {} UTC, sunset {} UTC",
        rise.strftime("%H:%M")?,
        set.strftime("%H:%M")?
    );
    let (dawn, dusk) = sun::riseset(&t, &boston, Some(96.0))?;
    println!(
        "civil twilight {} to {} UTC",
        dawn.strftime("%H:%M")?,
        dusk.strftime("%H:%M")?
    );

    // Midsummer above the Arctic Circle: no sunset, reported as an error
    let tromso = ITRFCoord::from_geodetic_deg(69.65, 18.96, 0.0);
    match sun::riseset(&t, &tromso, None) {
        Ok((rise, set)) => println!("Tromso: sunrise {rise}, sunset {set}"),
        Err(lpephem::Error::NoSunriseOrSunset) => println!("Tromso: midnight Sun, no sunset"),
        Err(e) => return Err(e.into()),
    }

    // Fraction of the Sun's disc visible from a satellite (0 = umbra)
    let r = 6_878_137.0;
    let sunward = psun * (r / psun.norm());
    let behind = -sunward;
    println!(
        "shadow function: sunward {:.1}, behind the Earth {:.1}",
        sun::shadowfunc(&psun, &sunward),
        sun::shadowfunc(&psun, &behind)
    );

    Ok(())
}
$ cargo run --example sun_moon
Sun  distance = 1.016219 AU
Moon distance = 386478 km
Moon phase 160.9 deg (Full Moon), 97.3% illuminated
Boston sunrise 09:07 UTC, sunset 00:24 UTC
civil twilight 08:32 to 00:59 UTC
Tromso: midnight Sun, no sunset
shadow function: sunward 1.0, behind the Earth 0.0

JPL ephemerides: jplephem

jplephem::geocentric_pos(body, &t), geocentric_state, barycentric_pos and barycentric_state return positions (m) and velocities (m/s) of a SolarSystem body in the ICRF/GCRF axes. The first query loads the DE440 file linux_p1550p2650.440 (1550–2650), downloading it (102 MB, SHA-256 verified) if it is not in a data directory; SATKIT_JPLEPHEM_FILE selects another file, such as the 14 MB DE421 (Selecting a JPL ephemeris file). Offline with no copy on disk, the queries return jplephem::Error.

//! JPL planetary ephemerides, compared with the low-precision models.
//!
//! Needs the JPL DE440 ephemeris `linux_p1550p2650.440` (102 MB), which is
//! downloaded and SHA-256 verified on first use, then cached in the data
//! directory. With `SATKIT_OFFLINE=1` and no copy on disk the queries fail
//! with `jplephem::Error`.
//!
//! ```text
//! cargo run --example ephemeris
//! ```

use satkit::{jplephem, lpephem, Instant, SolarSystem};

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let t = Instant::from_datetime(2024, 6, 20, 12, 0, 0.0)?;

    // Geocentric position (m) and velocity (m/s) of the Moon, in the GCRF
    let (p, v) = jplephem::geocentric_state(SolarSystem::Moon, &t)?;
    println!(
        "Moon: {:.0} km away, moving at {:.3} km/s",
        p.norm() / 1e3,
        v.norm() / 1e3
    );

    // Positions relative to the solar-system barycenter
    for body in [SolarSystem::Venus, SolarSystem::Mars, SolarSystem::Jupiter] {
        let p = jplephem::barycentric_pos(body, &t)?;
        println!(
            "{body:?}: {:.4} AU from the barycenter",
            p.norm() / satkit::consts::AU
        );
    }

    // The low-precision Sun and Moon are much faster; how far off are they?
    let angle = |a: satkit::Vector3, b: satkit::Vector3| {
        (a.dot(&b) / (a.norm() * b.norm()))
            .clamp(-1.0, 1.0)
            .acos()
            .to_degrees()
            * 3600.0
    };
    let sun_jpl = jplephem::geocentric_pos(SolarSystem::Sun, &t)?;
    let sun_lp = lpephem::sun::pos_gcrf(&t);
    println!(
        "low-precision Sun:  {:.0} arcsec from DE440",
        angle(sun_jpl, sun_lp)
    );
    let moon_lp = lpephem::moon::pos_gcrf(&t);
    println!(
        "low-precision Moon: {:.0} arcsec from DE440",
        angle(p, moon_lp)
    );

    // Constants from the ephemeris file header
    if let Some(au_km) = jplephem::consts("AU") {
        println!("AU in the DE440 header: {au_km} km");
    }

    Ok(())
}
$ cargo run --example ephemeris
Moon: 386605 km away, moving at 1.021 km/s
Venus: 0.7165 AU from the barycenter
Mars: 1.3874 AU from the barycenter
Jupiter: 5.0174 AU from the barycenter
low-precision Sun:  17 arcsec from DE440
low-precision Moon: 306 arcsec from DE440
AU in the DE440 header: 149597870.7 km

See also