Sky math in PL/SQL

From catalog numbers to "where is it in my sky right now". Seven lessons and five packages: time, coordinates, the Sun and the Moon, the Vedic calendar, and the planets. You finish with SkySutra's first live answers on its home page, and 48 tests that prove the math against published values.

How Part 3 works

Part 2 loaded where each star is on the sky's own map, RA and Dec. A person standing in a backyard wants something else: how high, and which way. Getting from one to the other takes the time, the place, and about 400 lines of PL/SQL, all from one book: Jean Meeus, Astronomical Algorithms (2nd edition), the standard reference for this kind of code.

You won't type 650 lines: you run the finished script, then each lesson opens one piece of it and explains it. That's how the math is checked too: piece by piece, against published answers.

Before you start

You need Part 2 done: the catalog loaded (88 constellations, 9,096 stars). Download this part's files:

Run it now. SQL Workshop › SQL Scripts › Upload skysutra_math.sql, then Run. It creates eight objects (four package specs, four bodies) with no errors. (skysutra_planets.sql comes in Lesson 3.6.) Check:

select object_name, object_type, status
  from user_objects
 where object_name in ('SS_TIME', 'SS_COORD', 'SS_SUN_MOON', 'SS_VEDIC')
 order by object_name, object_type;

Eight rows, all VALID.

Lesson 3.1 — Why the math lives in the database

Outcome: know why SkySutra's astronomy is PL/SQL, and what that buys.

The question. A star chart is drawn in the browser. Why not do the math in JavaScript there?

The decision. One copy of the math, in the database, next to the data it works on.

That last row is the real reason. Because ss_coord.star_altitude is a function, it works inside a query, so the database does the filtering instead of shipping 9,096 stars to a browser to throw most away.

The packages follow Turnout's rules (Turnout Lesson 3.2): a spec that says what, a body that says how, one job per package, names that say what they return. Two more rules come from the subject:

• Degrees everywhere outside, radians only inside. Oracle's sin and cos take radians; people, catalogs and the Naval Observatory use degrees. ss_time.rad and ss_time.deg convert at the last moment.
• deterministic where it's true. A function that always gives the same answer for the same input can say so, and Oracle may reuse an answer instead of computing it again. ss_vedic.ayanamsa reads a table, so it isn't marked: change the row and the answer changes.

Lesson 3.2 — Julian date and sidereal time

Outcome: turn any instant into the two numbers the sky runs on.

The question. Astronomy formulas want "how long since a fixed moment". How do we get that from a TIMESTAMP WITH TIME ZONE?

The decision. SS_TIME: a Julian date for counting, sidereal time for the turning sky.

Julian date

A Julian date is a single count of days since noon on January 1, 4713 BC. Subtract two and you have days, with no months, leap years or time zones in the way (Lesson 1.5). Astronomers pick one modern moment as a reference: J2000.0, noon UTC on January 1, 2000, which is Julian date 2451545.0.

function julian_date (p_at timestamp with time zone) return number deterministic is
    l_diff interval day(9) to second(9) :=
        p_at - timestamp '2000-01-01 12:00:00 UTC';
begin
    return 2451545
         + extract(day from l_diff)
         + (extract(hour from l_diff) * 3600
            + extract(minute from l_diff) * 60
            + extract(second from l_diff)) / 86400;
end;

Subtracting two timestamps with time zones gives an interval, already corrected for zones. Its days, plus its hours, minutes and seconds as a fraction of a day, added to 2451545, is the Julian date. centuries divides the days since J2000 by 36,525, the unit most of Meeus's formulas use.

Sidereal time

The Earth turns once in 23 hours 56 minutes relative to the stars, not 24. Sidereal time is the stars' clock: the right ascension that's due south, right now, at a place. Greenwich's comes from a formula (Meeus 12.4); a place's own is Greenwich's plus its longitude:

return norm360(280.46061837 + 360.98564736629 * l_d
               + 0.000387933 * l_t * l_t - l_t * l_t * l_t / 38710000);

l_d is days since J2000. The 360.9856… is the point: the sky turns a little more than once (360°) per day. norm360 folds any angle back into 0 to 360.

Try it.

select ss_time.julian_date(systimestamp)                          as julian_date,
       round(ss_time.lst(systimestamp, -82.3279) / 15, 3)          as sidereal_hours_wesley_chapel
  from dual;

Dividing by 15 turns degrees into hours. Any star whose RA (in hours) equals that number is due south of Wesley Chapel right now.

Lesson 3.3 — From RA/Dec to altitude and azimuth

Outcome: answer "is it up, and where do I look?" for any star, any place, any time.

The question. Aldebaran is RA 68.98°, Dec +16.51°. Where is it in Wesley Chapel's sky at 10 pm?

The decision. One procedure, ss_coord.to_horizon, and SQL-friendly functions around it.

Three steps

1. Precess the catalog position to today. The catalog gives positions for the year 2000. The Earth's axis wobbles (precession, the same wobble that makes the ayanamsa grow, Lesson 1.7), so the whole grid shifts about 50″ a year: 0.37° since 2000. That's most of the Moon's width, so ss_coord.precess moves each star to the current date (Meeus 21).
2. Find the hour angle: how far the star is past due south. Local sidereal time minus the star's RA.
3. Turn the sky on its axis to the place's latitude. Two spherical-trigonometry formulas:

l_h := rad(lst(p_at, p_longitude) - p_ra);   -- hour angle
p_altitude := deg(asin(sin(l_dec) * sin(l_lat) + cos(l_dec) * cos(l_lat) * cos(l_h)));
p_azimuth  := norm360(deg(atan2(-cos(l_dec) * sin(l_h),
                  sin(l_dec) * cos(l_lat) - cos(l_dec) * sin(l_lat) * cos(l_h))));

Azimuth counts from north through east (east is 90°), the way a compass does. atan2 takes two numbers, not their ratio, so it knows which quarter the answer is in; plain atan would confuse northeast with southwest.

star_altitude and star_azimuth do all three steps in one call, and compass turns an azimuth into one of 16 points (N, NNE, NE …).

What it leaves out, on purpose: refraction (the air lifts objects near the horizon by up to half a degree; Part 5 adds it where it matters, at rising and setting), proper motion (Arcturus, the fastest bright star, has moved about 1′ since 2000) and nutation (a nod of the axis under 20″).

Build: the first live answer

Make SkySutra's home page say something true about tonight's sky. In Page Designer, page 1, add a Classic Report region, Bright stars up now, after the SkySutra region, with this SQL:

select s.proper_name as star,
       c.name        as constellation,
       s.magnitude,
       round(ss_coord.star_altitude(s.ra_deg, s.dec_deg, systimestamp, p.latitude, p.longitude)) as altitude,
       ss_coord.compass(ss_coord.star_azimuth(s.ra_deg, s.dec_deg, systimestamp, p.latitude, p.longitude)) as direction
  from ss_stars s
  join ss_constellations c on c.id = s.constellation_id
  join ss_users u on u.id = ss_auth.current_user_id
  join ss_places p on p.id = u.home_place_id
 where s.magnitude < 2
   and s.proper_name is not null
   and ss_coord.star_altitude(s.ra_deg, s.dec_deg, systimestamp, p.latitude, p.longitude) > 0
 order by s.magnitude

Run the app. You see the brightest named stars above your home's horizon, with how high (degrees) and which way. ss_auth.current_user_id (Lesson 1.3) finds your home. It's daytime-blind for now: by day it lists stars hidden by sunlight. Part 5's Tonight page knows when it's dark.

Filter cheaply first

That query is fast because magnitude < 2 leaves 48 stars before any trigonometry runs. Ask the same question of the whole catalog:

select count(*)
  from ss_stars
 where ss_coord.star_altitude(ra_deg, dec_deg, systimestamp, 28.2397, -82.3279) > 0;

That takes about 6 seconds on apex.oracle.com: 9,096 stars, each precessed and turned, in Oracle's NUMBER, which is exact decimal arithmetic done in software. The rule: filter on plain columns first, compute last. (SS_ISS, which steps a satellite through thousands of positions, uses BINARY_DOUBLE, the hardware's floating point, for the same reason.)

Lesson 3.4 — The Sun and the Moon

Outcome: place the Sun and the Moon to a fraction of a degree, and know why that's enough.

The question. Professional ephemerides give planet positions to a thousandth of an arcsecond. Do we need that?

The decision. Meeus's low-precision formulas: good to about 0.01° for the Sun and 0.05° for the Moon. The Moon is 0.5° wide. No one looking up can see an error a tenth of its width.

The Sun

ss_sun_moon.sun_longitude follows Meeus chapter 25: the Sun's mean longitude and mean anomaly grow steadily with time; an "equation of the centre" adds the effect of the Earth's elliptical orbit (up to 1.9°); two small terms adjust for where we see it from. The answer is a longitude on the ecliptic, the Sun's yearly path. ss_coord.ecliptic_to_equatorial turns that into RA and Dec, and Lesson 3.3's to_horizon does the rest.

The Moon

The Moon's orbit is pulled about by the Sun, so its position is a sum of wobbles. Meeus's full tables have 60 terms for longitude; moon_ecliptic uses the 24 largest, plus 11 for latitude and 13 for distance:

p_lon := l_lp
    + 6.288774 * s(l_mp)              + 1.274027 * s(2*l_d - l_mp)
    + 0.658314 * s(2*l_d)             + 0.213618 * s(2*l_mp)
    ...

The first term alone is the elliptical orbit (6.3°); the second, the Sun's pull (1.3°). The ones left out are each under 0.004°.

Parallax: the Moon is close

Stars are so far away that it doesn't matter where on Earth you stand. The Moon is about 384,000 km away, only 60 Earth radii, so from your backyard it sits up to 1° lower than from the Earth's centre. moon_horizon corrects for it:

p_altitude := p_altitude - deg(asin(6378.14 / l_km)) * cos(rad(p_altitude));

At 10 pm on October 22 the Moon was 57.81° up from the Earth's centre and 57.30° from Wesley Chapel. Half a degree: one Moon-width. That's why the Naval Observatory's navigation data lists parallax separately.

Lesson 3.5 — Phase, tithi and nakshatra

Outcome: compute "Today in the sky" from two longitudes.

The question. A Moon-phase app, a Hindu calendar and a nakshatra table look like three different subjects. What do they share?

The decision. One angle: the elongation, the Moon's longitude minus the Sun's.

function elongation (p_at timestamp with time zone) return number deterministic is
begin
    return ss_time.norm360(moon_longitude(p_at) - sun_longitude(p_at));
end;

The nakshatra needs one more step: the Moon's sidereal longitude (Lesson 0.4), which is its longitude minus the ayanamsa read from SS_AYANAMSA (Lesson 1.7). Each nakshatra is 13°20′, so:

return floor(sidereal(ss_sun_moon.moon_longitude(p_at), p_at) / (360 / 27)) + 1;

Two small extras: rahu_longitude gives the Moon's ascending node (Rahu; Ketu is opposite), and ss_vedic.vara names the weekday (Ravivar for Sunday through Shanivar for Saturday).

Try it. Today in the sky, for right now:

select ss_sun_moon.moon_phase_name(systimestamp)                as phase,
       round(ss_sun_moon.moon_illumination(systimestamp) * 100) as lit_pct,
       ss_vedic.tithi_name(ss_vedic.tithi(systimestamp))        as tithi,
       n.name || ' (' || n.english_name || ')'                  as moon_nakshatra,
       ss_vedic.vara(trunc(cast(systimestamp at time zone 'America/New_York' as date))) as vara
  from ss_nakshatras n
 where n.seq = ss_vedic.moon_nakshatra(systimestamp);

On the morning of October 24, 2026 that gave: Waxing gibbous, 96% lit, Shukla Chaturdashi, Uttara Bhadrapada, Shanivar. Check yours against any published panchang.

Lesson 3.6 — The planets

Outcome: place Mercury to Neptune to within 2 arcminutes, say how bright each is, and say when it rises and sets.

The question. The plan was to ask the Naval Observatory for planet positions. On October 24, 2026, apex.oracle.com stopped reaching it. Now what?

The decision. Compute them, like the Sun and the Moon. A planet follows an ellipse around the Sun that changes only slowly; six numbers describe it, and Paul Schlyter's well-known method ("How to compute planetary positions") gives those six numbers for any date, plus the few corrections that matter for the eye.

Run it. SQL Scripts › Upload skysutra_planets.sql, Run. Then upload and run the updated skysutra_rise_set.sql, which teaches SS_ALMANAC the planets (more below).

Five steps from orbit to sky

After step 5 the planet has an ecliptic longitude and latitude, like the Moon in Lesson 3.4, so the same ecliptic_to_equatorial and to_horizon finish the job (ss_planets.horizon).

-- planet from the Sun, plus the Sun from the Earth, is the planet from the Earth
l_x := l_r * c(l_hlon) * c(l_hlat) + l_sr * c(l_slon);
l_y := l_r * s(l_hlon) * c(l_hlat) + l_sr * s(l_slon);
l_z := l_r * s(l_hlat);

How good is it? Against Astropy, for dates from 2020 to 2030, every planet is within 2 arcminutes: a fifteenth of the Moon's width.

Brightness and glare

Two more things a stargazer asks:

• How bright? ss_planets.magnitude uses each planet's distance from the Sun and from us and how much of its lit side faces us (Saturn adds its rings, brighter the more they're tilted towards us). It's approximate, within about 0.2 magnitudes: enough to say Venus outshines everything but the Moon.
• Can I see it? A planet near the Sun in the sky is lost in its glare. ss_planets.sun_distance gives the angle between them; under about 15° the planet rises and sets with the Sun.

Rise and set

SS_ALMANAC (the package that fills SS_RISE_SET) used to stop at the Sun and the Moon. The updated script adds the planets with the same method: step through the day ten minutes at a time, and close in on each horizon crossing by halving the gap. A planet is a point of light, so its horizon is just refraction, 0.5667° below the true horizon. Against Astropy, planet rise and set times agree within 5 seconds: run the first block of skysutra_rise_set_tests.sql in SQL Commands and it ends 21 passed, 0 failed.

The same package now finds the Moon's four main phases, the moments the elongation (Lesson 3.5) passes 0°, 90°, 180° and 270°, and stores them in SS_MOON_PHASES:

begin
    ss_almanac.refresh_moon_phases(trunc(sysdate), 90);
    commit;
end;

Filling SS_RISE_SET for a place now takes about 18 seconds a week (nine bodies, seven days), so it's a job for the nightly automation in Part 10, one place at a time. SQL Commands stops a request after about a minute.

Build: planets on the home page

Page 1, a second Classic Report, Planets, after Bright stars up now:

select b.name as planet,
       b.vedic_name,
       ss_planets.magnitude(b.code, systimestamp) as magnitude,
       round(ss_planets.altitude(b.code, systimestamp, p.latitude, p.longitude)) as altitude,
       ss_coord.compass(ss_planets.azimuth(b.code, systimestamp, p.latitude, p.longitude)) as direction,
       case when ss_planets.sun_distance(b.code, systimestamp) < 15 then 'Too close to the Sun' end as note
  from ss_bodies b
  join ss_users u on u.id = ss_auth.current_user_id
  join ss_places p on p.id = u.home_place_id
 where b.kind = 'PLANET'
 order by 3

At 9 am on October 24, 2026 in Wesley Chapel that listed Venus (Shukra, −3.5, too close to the Sun), Jupiter (Guru, −1.7, high in the south) and Mars (Mangala, +1.0) first, all in daylight, which Part 5 will take into account.

Lesson 3.7 — Testing against published answers

Outcome: run 48 tests, read them, and know where each expected answer came from.

The question. The numbers look plausible. How do we know they're right, and that they stay right after a change?

The decision. Fixed dates, published answers, and a tolerance for each, in one script you can run any time (Turnout Lesson 3.7's pattern).

near('Sun altitude, 1 pm EDT Oct 23 (USNO)', l_alt, 50.046, 0.02);

near passes if our number is within the tolerance of the published one, and prints both either way.

Run it. Open skysutra_math_tests.sql, copy the first block (from declare to the /), paste it into SQL Commands, Run. It ends:

48 passed, 0 failed.

Read a few lines and notice the size of the differences: stars within 0.005°, the Sun within 0.002°, the Moon within 0.05°, the planets within 0.03°. The ayanamsa is 0.006° (22″) below the published Lahiri value, because we leave out nutation; that moves a nakshatra boundary by under a minute of time.

The live check. The second part of the file asks the Naval Observatory directly, with apex_web_service.make_rest_request, and compares twelve stars. It needs the database to reach aa.usno.navy.mil. On October 23 it did; on October 24 apex.oracle.com couldn't (ORA-29259: end-of-input reached in UTL_HTTP call), while other services worked. That's why the tests above don't call it: a test that depends on someone else's server tells you about their server, not your code. Part 4 deals with services that go away.

Checkpoint

SkySutra can place every star, the Sun, the Moon and the planets in anyone's sky, and knows the Moon's phase, tithi and nakshatra. Part 4 brings in what can't be computed: the ISS's latest orbit, space weather, the picture of the day.

Test it

1. The eight package objects are VALID.
2. skysutra_math_tests.sql: 48 passed, 0 failed.
3. The home page shows Bright stars up now and Planets for your home, and the altitudes look right: compare one star and one planet with a sky app on your phone.
4. skysutra_rise_set_tests.sql, first block: 21 passed, 0 failed.

Files

Exercises

• Add a test of to_horizon: Meeus example 13.b gives Venus's altitude and azimuth over Washington from its RA and Dec. Careful: Meeus measures azimuth from the south, SkySutra from the north.
• Add the Moon to the home page: a second region with moon_horizon for your home, its altitude, compass direction and moon_phase_name.
• Write ss_coord.angular_distance(ra1, dec1, ra2, dec2) and use it to find the star closest to the Moon right now.
• Stretch: copy star_altitude into a test function that uses BINARY_DOUBLE throughout, and time the full-catalog count from Lesson 3.3 with each.

Review questions

Try each one before you open it. The graded quiz is at the end of the page.