Turn Earth by the Earth Rotation Angle, so its lit face stays Horizons' at AD 1 as it is today
Earth was turned by its IAU W, taken at the clock's UT plus today's 69.184 s. That W is a straight
line fitted to the present: 360.9856235 degrees a day, which, once its pole's -0.641 degrees a
century in right ascension is counted, runs 6.3e-6 degrees a day slow of Earth's real turning. The
followsUt comment said the clock's date "already says how far it has turned"; it did not. Against
Horizons (observer quantity 14 from the Sun, TIME_TYPE=UT, Earth one light-time back) the drawn
sub-solar point was 2.3 degrees off at AD 1000 and 4.5 at AD 1.
Earth is now turned by the IERS Earth Rotation Angle (IERS Conventions 2010, eq. 5.15) at the
clock's date, counted from the node the IAU's W starts at, 90 degrees past the pole's right
ascension. The pole is unchanged. Drawn minus Horizons, in degrees:
date before after
2025-06-01 12:00 (unit) +0.06 -0.003
AD 1000, JD 2086455 (unit) -2.3 -0.001
AD 1, JD 1721600 (unit) -4.5 +0.051
live app, :4301, same probe as the review's
JD 2460900.25 +0.089 +0.005
JD 2086300.5 -2.281 +0.010
JD 1800000 -4.049 +0.056
JD 1721450.75 -4.530 +0.072
At noon UTC on 1 June 2025 the Sun now stands over 0.52 W on the drawn sphere, where the equation
of time puts it at 0.53 W (0.43 W before).
TT_MINUS_UTC_DAYS had no other use and is removed; its comment also counted 37 leap seconds where
UTC has taken 27 on top of the 10 s it started from in 1972.
Tests: body-orientation.spec 'lights Earth's face where Horizons does at the far end of the clock
too: AD 1000 and AD 1' (within 0.15 degrees), and the renderer's AD 1000 Earth test now checks the
drawn face against Horizons instead of against the IAU W the old code used. Control: turning Earth
by its IAU W at UT + 69.184 s again fails both named tests (2 failed, 834 passed). The README says
which model turns Earth.
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
This commit is contained in:
@@ -1,7 +1,7 @@
|
||||
import * as THREE from 'three/webgpu';
|
||||
import { describe, expect, it, vi } from 'vitest';
|
||||
|
||||
import { DEFAULT_EPOCH_JD, GM_SUN_AU3_PER_DAY2, TT_MINUS_UTC_DAYS, ttMinusUtSeconds } from '../../shared/astro/constants';
|
||||
import { DEFAULT_EPOCH_JD, GM_SUN_AU3_PER_DAY2, ttMinusUtSeconds } from '../../shared/astro/constants';
|
||||
import { keplerRates } from '../../shared/astro/kepler';
|
||||
import { eclipticToEquatorial, laplacePlaneToEquatorial, OBLIQUITY_J2000_DEG } from '../../shared/astro/coordinates';
|
||||
import { orientationAt } from '../../shared/astro/rotational-elements';
|
||||
@@ -734,15 +734,17 @@ describe('solar-system bodies against Horizons', () => {
|
||||
expect((drawnPrimeMeridian('jupiter').angleTo(iauPrimeMeridian('jupiter', jdUt + 1574.1 / 86400)) * 180) / Math.PI).toBeLessThan(0.01);
|
||||
});
|
||||
|
||||
it('turns Earth by the UT the clock names, which is its turning: at AD 1000 its W is not moved on by ΔT', () => {
|
||||
const jdUt = 2086307.5;
|
||||
renderer.update(jdUt);
|
||||
expect((drawnPrimeMeridian('earth').angleTo(iauPrimeMeridian('earth', jdUt + TT_MINUS_UTC_DAYS)) * 180) / Math.PI).toBeLessThan(0.01);
|
||||
it('turns Earth by the UT the clock names, which is its turning: at AD 1000 the Sun stands over Horizons’ point', () => {
|
||||
// Horizons' sub-solar longitude from the Sun (observer quantity 14, TIME_TYPE=UT) on JD 2086455,
|
||||
// 1.0510 E, is Earth as it was 8.454 minutes before. Turned by the IAU's W at UT the drawn face
|
||||
// was 2.3 degrees off; taken at TDB, 6.6.
|
||||
renderer.update(2086455 - 8.45437443 / 1440);
|
||||
expect(apart(facing('earth', new THREE.Vector3()).eastDeg, 1.05101)).toBeLessThan(0.15);
|
||||
});
|
||||
|
||||
it('lights Earth where the Sun really stands: within 4 degrees of Greenwich at noon UTC', () => {
|
||||
// The equation of time is all that separates them: on 1 June 2025 it puts the Sun over 0.53 W,
|
||||
// and the drawn sphere has it over 0.43 W.
|
||||
// and the drawn sphere has it over 0.52 W.
|
||||
renderer.update(JUNE_1_2025_NOON_UTC);
|
||||
expect(Math.abs(facing('earth', new THREE.Vector3()).eastDeg)).toBeLessThan(4);
|
||||
});
|
||||
|
||||
@@ -19,13 +19,6 @@ export const DEFAULT_EPOCH_JD = 2451545.0;
|
||||
*/
|
||||
export const GM_SUN_AU3_PER_DAY2 = 0.01720209895 * 0.01720209895;
|
||||
|
||||
/**
|
||||
* TT - UTC today, in days: 32.184 s plus the 37 leap seconds UTC has taken since 1972, the last at
|
||||
* the end of 2016. TDB, which ephemerides run on, stays within 2 ms of TT. See {@link ttMinusUtSeconds}
|
||||
* for other dates.
|
||||
*/
|
||||
export const TT_MINUS_UTC_DAYS = 69.184 / 86400;
|
||||
|
||||
/** The first day of each month UTC took a leap second at the start of, from its 10 s of 1972. */
|
||||
const LEAP_SECONDS_FROM = [
|
||||
[1972, 7], [1973, 1], [1974, 1], [1975, 1], [1976, 1], [1977, 1], [1978, 1], [1979, 1], [1980, 1], [1981, 7],
|
||||
|
||||
@@ -47,6 +47,15 @@ describe('bodyPageView', () => {
|
||||
expect(Math.abs(geodetic - 22.260426)).toBeLessThan(0.05);
|
||||
});
|
||||
|
||||
it('lights Earth’s face where Horizons does at the far end of the clock too: AD 1000 and AD 1', () => {
|
||||
// Horizons' sub-solar longitude from the Sun (observer quantity 14, TIME_TYPE=UT), Earth taken
|
||||
// one light-time back: 1.0510 E on JD 2086455 and 1.5606 E on JD 1721600. The IAU's W, taken at
|
||||
// UT, drew them 2.3 and 4.5 degrees west of that.
|
||||
for (const [jdUt, lightMinutes, eastDeg] of [[2086455, 8.45437443, 1.05101], [1721600, 8.45020842, 1.560644]]) {
|
||||
expect(Math.abs(subSolarPoint(EARTH, jdUt - lightMinutes / 1440).eastDeg - eastDeg)).toBeLessThan(0.15);
|
||||
}
|
||||
});
|
||||
|
||||
it('takes a moon’s Sun from where it and its planet are: the Moon’s sub-solar point is Horizons’', () => {
|
||||
const moon = subSolarPoint(MOON, JUNE_1_2025_NOON_UTC);
|
||||
// 116.2859 E and 1.5030 N, seen from Earth's centre.
|
||||
|
||||
@@ -1,6 +1,6 @@
|
||||
import * as THREE from 'three/webgpu';
|
||||
|
||||
import { tdbFromUtc, TT_MINUS_UTC_DAYS } from '../astro/constants';
|
||||
import { tdbFromUtc } from '../astro/constants';
|
||||
import { CartesianCoordinates, eclipticToEquatorial, laplacePlaneToEquatorial } from '../astro/coordinates';
|
||||
import { meanElementsAt, positionAtEpoch } from '../astro/kepler';
|
||||
import { orientationAt } from '../astro/rotational-elements';
|
||||
@@ -54,14 +54,18 @@ export function poleFrame(pole: { raDeg: number; decDeg: number }, target = new
|
||||
*
|
||||
* The clock is UT and the IAU's elements run on TDB, 69.184 s ahead today and 1 574 s at AD 1000;
|
||||
* in 69 s Earth turns 0.29 degrees, Jupiter 0.70 and Phobos 0.90, so the date is taken to TDB here
|
||||
* (see `tdbFromUtc`). Earth, `followsUt`, is the one exception: its turning is what UT counts,
|
||||
* so the clock's date already says how far it has turned, and its W, fitted to today, is taken at
|
||||
* that date plus today's TT - UTC. Taken at TDB, it would turn ΔT further: 44 degrees at AD 1.
|
||||
* (see `tdbFromUtc`). Earth, `followsUt`, is the one exception: its turning is what UT counts, so
|
||||
* it is turned by the IERS Earth Rotation Angle at the clock's date (IERS Conventions 2010, eq.
|
||||
* 5.15), counted from the node its W starts at, 90 degrees past its pole's right ascension. The
|
||||
* IAU's W for Earth, fitted to today, runs 6.3e-6 degrees a day slow of that once its pole's drift
|
||||
* is counted: taken at UT, it left Earth's lit face 2.3 degrees off Horizons at AD 1000 and 4.5 at
|
||||
* AD 1. Taken at TDB, it would have turned ΔT further, 44 degrees at AD 1.
|
||||
*/
|
||||
export function bodyOrientation(elements: RotationalElements, jdUtc: number, target = new THREE.Quaternion(), followsUt = false): THREE.Quaternion {
|
||||
const { poleRaDeg, poleDecDeg, primeMeridianDeg } = orientationAt(elements, followsUt ? jdUtc + TT_MINUS_UTC_DAYS : tdbFromUtc(jdUtc));
|
||||
const { poleRaDeg, poleDecDeg, primeMeridianDeg } = orientationAt(elements, tdbFromUtc(jdUtc));
|
||||
const turnDeg = followsUt ? 360 * (0.779057273264 + 1.00273781191135448 * (jdUtc - 2451545)) - 90 - poleRaDeg : primeMeridianDeg;
|
||||
return poleFrame({ raDeg: poleRaDeg, decDeg: poleDecDeg }, target)
|
||||
.multiply(scratchTurn.setFromAxisAngle(Z_AXIS, primeMeridianDeg * DEG_TO_RAD))
|
||||
.multiply(scratchTurn.setFromAxisAngle(Z_AXIS, turnDeg * DEG_TO_RAD))
|
||||
.multiply(MAP_TO_BODY);
|
||||
}
|
||||
|
||||
|
||||
Reference in New Issue
Block a user