diff --git a/js/ephem.js b/js/ephem.js new file mode 100644 index 0000000..ff5c03a --- /dev/null +++ b/js/ephem.js @@ -0,0 +1,171 @@ +// SOLARGOD ephemeris engine — JPL SSD "Approximate Positions of the Planets" +// (Standish). Tables extracted 2026-07-15 from +// https://ssd.jpl.nasa.gov/planets/approx_pos.html — authoritative, verbatim. +// Angles in degrees, rates per Julian century; `e` is dimensionless (the page's +// "rad" header for e is its own typo). No n-body integration — this IS the spec. + +import { centuriesSinceJ2000, normDegPM180, DEG } from './lib.js'; + +// Row layout: el = [a(au), e, I(deg), L(deg), ϖ long.peri(deg), Ω long.node(deg)] +// rate = same order, per century. +const TABLE1 = { + mercury: { el: [0.38709927, 0.20563593, 7.00497902, 252.25032350, 77.45779628, 48.33076593], + rate: [0.00000037, 0.00001906, -0.00594749, 149472.67411175, 0.16047689, -0.12534081] }, + venus: { el: [0.72333566, 0.00677672, 3.39467605, 181.97909950, 131.60246718, 76.67984255], + rate: [0.00000390, -0.00004107, -0.00078890, 58517.81538729, 0.00268329, -0.27769418] }, + earth: { el: [1.00000261, 0.01671123, -0.00001531, 100.46457166, 102.93768193, 0.0], + rate: [0.00000562, -0.00004392, -0.01294668, 35999.37244981, 0.32327364, 0.0] }, + mars: { el: [1.52371034, 0.09339410, 1.84969142, -4.55343205, -23.94362959, 49.55953891], + rate: [0.00001847, 0.00007882, -0.00813131, 19140.30268499, 0.44441088, -0.29257343] }, + jupiter: { el: [5.20288700, 0.04838624, 1.30439695, 34.39644051, 14.72847983, 100.47390909], + rate: [-0.00011607, -0.00013253, -0.00183714, 3034.74612775, 0.21252668, 0.20469106] }, + saturn: { el: [9.53667594, 0.05386179, 2.48599187, 49.95424423, 92.59887831, 113.66242448], + rate: [-0.00125060, -0.00050991, 0.00193609, 1222.49362201, -0.41897216, -0.28867794] }, + uranus: { el: [19.18916464, 0.04725744, 0.77263783, 313.23810451, 170.95427630, 74.01692503], + rate: [-0.00196176, -0.00004397, -0.00242939, 428.48202785, 0.40805281, 0.04240589] }, + neptune: { el: [30.06992276, 0.00859048, 1.77004347, -55.12002969, 44.96476227, 131.78422574], + rate: [0.00026291, 0.00005105, 0.00035372, 218.45945325, -0.32241464, -0.00508664] }, +}; + +// Table 2a — valid 3000 BC–3000 AD. Shipped for the extended-range config swap +// (CONFIG.useExtendedRange); unused by v1 UI. +const TABLE2A = { + mercury: { el: [0.38709843, 0.20563661, 7.00559432, 252.25166724, 77.45771895, 48.33961819], + rate: [0.00000000, 0.00002123, -0.00590158, 149472.67486623, 0.15940013, -0.12214182] }, + venus: { el: [0.72332102, 0.00676399, 3.39777545, 181.97970850, 131.76755713, 76.67261496], + rate: [-0.00000026, -0.00005107, 0.00043494, 58517.81560260, 0.05679648, -0.27274174] }, + earth: { el: [1.00000018, 0.01673163, -0.00054346, 100.46691572, 102.93005885, -5.11260389], + rate: [-0.00000003, -0.00003661, -0.01337178, 35999.37306329, 0.31795260, -0.24123856] }, + mars: { el: [1.52371243, 0.09336511, 1.85181869, -4.56813164, -23.91744784, 49.71320984], + rate: [0.00000097, 0.00009149, -0.00724757, 19140.29934243, 0.45223625, -0.26852431] }, + jupiter: { el: [5.20248019, 0.04853590, 1.29861416, 34.33479152, 14.27495244, 100.29282654], + rate: [-0.00002864, 0.00018026, -0.00322699, 3034.90371757, 0.18199196, 0.13024619] }, + saturn: { el: [9.54149883, 0.05550825, 2.49424102, 50.07571329, 92.86136063, 113.63998702], + rate: [-0.00003065, -0.00032044, 0.00451969, 1222.11494724, 0.54179478, -0.25015002] }, + uranus: { el: [19.18797948, 0.04685740, 0.77298127, 314.20276625, 172.43404441, 73.96250215], + rate: [-0.00020455, -0.00001550, -0.00180155, 428.49512595, 0.09266985, 0.05739699] }, + neptune: { el: [30.06952752, 0.00895439, 1.77005520, 304.22289287, 46.68158724, 131.78635853], + rate: [0.00006447, 0.00000818, 0.00022400, 218.46515314, 0.01009938, -0.00606302] }, +}; + +// Table 2b — extra terms added to M (deg) for Jupiter–Neptune, Table 2a only: +// M += b·T² + c·cos(f·T) + s·sin(f·T). [b, c, s, f] +const TABLE2B = { + jupiter: [-0.00012452, 0.06064060, -0.35635438, 38.35125000], + saturn: [0.00025899, -0.13434469, 0.87320147, 38.35125000], + uranus: [0.00058331, -0.97731848, 0.17689245, 7.67025000], + neptune: [-0.00041348, 0.68346318, -0.10162547, 7.67025000], +}; + +export const PLANET_IDS = ['mercury', 'venus', 'earth', 'mars', 'jupiter', 'saturn', 'uranus', 'neptune']; + +let useExtended = false; +export function setExtendedRange(on) { useExtended = !!on; } + +// Solve Kepler's equation E − e·sinE = M (radians) by Newton–Raphson. +// Seed E₀ = M + e·sinM; converges in a few iters even for e→0.97 (verify gate). +export function solveKepler(M, e, tol = 1e-9, maxIter = 20) { + let E = M + e * Math.sin(M); + for (let i = 0; i < maxIter; i++) { + const dE = (E - e * Math.sin(E) - M) / (1 - e * Math.cos(E)); + E -= dE; + if (Math.abs(dE) < tol) break; + } + return E; +} + +// The reusable core: resolved classical elements → heliocentric J2000 ecliptic +// {x,y,z} in AU. a[au], e, and iDeg/OmDeg(Ω)/wDeg(ω arg.peri)/Mdeg all in degrees. +// Shared by planets (element+rate) and small bodies (epoch elements) — brief §3. +export function eclFromElements(a, e, iDeg, OmDeg, wDeg, Mdeg, out = {}) { + const M = normDegPM180(Mdeg) * DEG; + const E = solveKepler(M, e); + // Orbital-plane coords (perifocal), AU. + const xp = a * (Math.cos(E) - e); + const yp = a * Math.sqrt(1 - e * e) * Math.sin(E); + const w = wDeg * DEG, Om = OmDeg * DEG, I = iDeg * DEG; + const cw = Math.cos(w), sw = Math.sin(w); + const cO = Math.cos(Om), sO = Math.sin(Om); + const cI = Math.cos(I), sI = Math.sin(I); + // Perifocal → ecliptic rotation (brief §3 step 5 = standard R_z(Ω)R_x(I)R_z(ω)). + out.x = (cw * cO - sw * sO * cI) * xp + (-sw * cO - cw * sO * cI) * yp; + out.y = (cw * sO + sw * cO * cI) * xp + (-sw * sO + cw * cO * cI) * yp; + out.z = (sw * sI) * xp + (cw * sI) * yp; + return out; +} + +// Resolve a planet's elements at jd (Table 1, 1800–2050) → {a,e,I,Om,w,M} deg. +function resolvePlanet(bodyId, jd) { + const T = centuriesSinceJ2000(jd); + const row = TABLE1[bodyId]; + if (!row) throw new Error(`ephem: unknown planet ${bodyId}`); + const a = row.el[0] + row.rate[0] * T; + const e = row.el[1] + row.rate[1] * T; + const I = row.el[2] + row.rate[2] * T; + const L = row.el[3] + row.rate[3] * T; + const wb = row.el[4] + row.rate[4] * T; // ϖ longitude of perihelion + const Om = row.el[5] + row.rate[5] * T; // Ω longitude of ascending node + const w = wb - Om; // ω argument of perihelion + const M = L - wb; // mean anomaly (normalize AFTER subtract) + return { a, e, I, Om, w, M }; +} +// Extended-range resolver (Table 2a + the 2b M-correction for Jupiter–Neptune, +// valid 3000 BC–3000 AD). f·T uses T in centuries; f is deg-per-century so f·T +// is degrees → convert to rad for the trig. Only used when useExtended is set. +function resolvePlanetExt(bodyId, jd) { + const T = centuriesSinceJ2000(jd); + const row = TABLE2A[bodyId]; + const a = row.el[0] + row.rate[0] * T; + const e = row.el[1] + row.rate[1] * T; + const I = row.el[2] + row.rate[2] * T; + const L = row.el[3] + row.rate[3] * T; + const wb = row.el[4] + row.rate[4] * T; + const Om = row.el[5] + row.rate[5] * T; + const w = wb - Om; + let M = L - wb; + const b2 = TABLE2B[bodyId]; + if (b2) { + const [b, c, s, f] = b2; + M += b * T * T + c * Math.cos(f * T * DEG) + s * Math.sin(f * T * DEG); + } + return { a, e, I, Om, w, M }; +} + +// Public: heliocentric J2000 ecliptic position of a planet at jd, AU. +export function helioEcl(bodyId, jd, out = {}) { + const p = useExtended ? resolvePlanetExt(bodyId, jd) : resolvePlanet(bodyId, jd); + return eclFromElements(p.a, p.e, p.I, p.Om, p.w, p.M, out); +} + +// Mean motion [deg/day] from semi-major axis via Kepler's third law (Gaussian). +export function meanMotionDegPerDay(a) { return 0.9856076686 / Math.pow(a, 1.5); } + +// Small-body position (asteroids/comets/Pluto): epoch elements propagated by +// M = M0 + n·(jd − epoch). el = {a,e,i,om(Ω),w(ω),ma(M0),epoch}. +export function smallBodyEcl(el, jd, out = {}) { + const n = meanMotionDegPerDay(el.a); + const M = el.ma + n * (jd - el.epoch); + return eclFromElements(el.a, el.e, el.i, el.om, el.w, M, out); +} + +// One full orbit sampled in mean anomaly, as a Float64Array of xyz triples +// (nSamples+1 points; last == first so a LineLoop/LineSegments closes cleanly — +// verify gate "orbitPath endpoints join"). Elements frozen at jd (they drift +// slowly; caller rebuilds only when |ΔT| > 0.1 Cy, brief §11). +export function orbitPath(bodyId, jd, nSamples = 256) { + const p = useExtended ? resolvePlanetExt(bodyId, jd) : resolvePlanet(bodyId, jd); + return samplePath(p.a, p.e, p.I, p.Om, p.w, nSamples); +} +export function orbitPathFromElements(el, nSamples = 256) { + return samplePath(el.a, el.e, el.i, el.om, el.w, nSamples); +} +function samplePath(a, e, I, Om, w, nSamples) { + const out = new Float64Array((nSamples + 1) * 3); + const tmp = {}; + for (let k = 0; k <= nSamples; k++) { + const M = (360 * k) / nSamples; // 0..360; k=nSamples reproduces k=0 + eclFromElements(a, e, I, Om, w, M, tmp); + out[k * 3] = tmp.x; out[k * 3 + 1] = tmp.y; out[k * 3 + 2] = tmp.z; + } + return out; +} diff --git a/js/layers/planets.js b/js/layers/planets.js new file mode 100644 index 0000000..cdce17d --- /dev/null +++ b/js/layers/planets.js @@ -0,0 +1,93 @@ +// SOLARGOD planets layer (Stage 1) — the eight planets + Pluto as textured, +// exaggerated spheres on compressed orbit paths, with labels. Positions come from +// ctx.bodyWorld (computed centrally in main.js); this layer only renders. Orbit +// geometry is rebuilt only when the MEGA/TRUE exponent is animating or elements +// drift past CONFIG.orbitRebuildCy (brief §11). + +export default function create(ctx) { + const { THREE, scene, CONFIG, lib, ui, scale, ephem, bodies, worldGroup, + clock, focus, toLocal, effectiveE, makeLabel, registerPick, loadTextureInto, PLANET_ORDER } = ctx; + + const ROSTER = [...PLANET_ORDER, 'pluto']; + const orbitMat = new THREE.LineBasicMaterial({ color: CONFIG.colors.orbit, transparent: true, opacity: 0.6 }); + const N = CONFIG.orbitSamples; + const entries = []; + const tmp = new THREE.Vector3(); + const _v = {}; + + for (const id of ROSTER) { + const b = bodies[id]; + // tilt group holds the spinning sphere so axial tilt + rotation compose. + const group = new THREE.Group(); + const mat = new THREE.MeshStandardMaterial({ color: new THREE.Color(b.color), roughness: 1, metalness: 0 }); + const mesh = new THREE.Mesh(new THREE.SphereGeometry(1, 48, 32), mat); + group.rotation.z = (b.axialTiltDeg || 0) * lib.DEG; + group.add(mesh); + scene.add(group); + loadTextureInto(mat, b.texture); + registerPick(mesh, id); + const lbl = makeLabel(group, b.name, 'body-label'); + + // orbit path (absolute view-world, lives under worldGroup) + const geo = new THREE.BufferGeometry(); + geo.setAttribute('position', new THREE.BufferAttribute(new Float32Array((N + 1) * 3), 3)); + const line = new THREE.LineLoop(geo, orbitMat); + line.frustumCulled = false; + worldGroup.add(line); + + const e = { id, b, group, mesh, mat, geo, line, lbl, auPath: null, lastT: NaN }; + buildAuPath(e); + fillView(e); + entries.push(e); + } + + function buildAuPath(e) { + e.auPath = e.b.ephemId + ? ephem.orbitPath(e.b.ephemId, clock.jd, N) + : ephem.orbitPathFromElements(e.b.elements, N); + e.lastT = lib.centuriesSinceJ2000(clock.jd); + } + function fillView(e) { + const pos = e.geo.attributes.position.array; + const au = e.auPath; + for (let k = 0; k <= N; k++) { + _v.x = au[k * 3]; _v.y = au[k * 3 + 1]; _v.z = au[k * 3 + 2]; + scale.viewFromEcl(_v, _v); + pos[k * 3] = _v.x; pos[k * 3 + 1] = _v.y; pos[k * 3 + 2] = _v.z; + } + e.geo.attributes.position.needsUpdate = true; + e.geo.computeBoundingSphere(); + } + + let showOrbits = true, showLabels = true; + ui.addLayer('orbits', 'Orbit paths', true, (on) => { showOrbits = on; for (const e of entries) e.line.visible = on; }); + ui.addLayer('labels', 'Body labels', true, (on) => { showLabels = on; for (const e of entries) e.lbl.div.style.display = on ? '' : 'none'; }); + ui.setStatus('orbits', `${entries.length} orbits`, 'ok'); + ui.setStatus('labels', 'names shown', 'ok'); + + return { + id: 'planets', + onClockTick(simMs, jd, isLive) { + const transitioning = scale.isTransitioning(); + for (const e of entries) { + const abs = ctx.bodyWorld[e.id]; + if (!abs) continue; + // position (floating-origin JS subtraction) + exaggerated radius + toLocal(abs, e.group.position); + const r = scale.drawRadius(e.b.radiusKm, effectiveE(e.id)); + e.mesh.scale.setScalar(r); + // spin: fraction of a rotation since J2000 (retrograde via sign) + if (e.b.rotationHours) { + const rot = ((jd - lib.J2000_JD) * 24) / e.b.rotationHours; + e.mesh.rotation.y = rot * 2 * Math.PI; + } + // focus highlight + e.lbl.div.classList.toggle('focused', e.id === focus.id); + // orbit rebuild only when needed + const needAu = Math.abs(lib.centuriesSinceJ2000(jd) - e.lastT) > CONFIG.orbitRebuildCy; + if (needAu) buildAuPath(e); + if (needAu || transitioning) fillView(e); + } + }, + }; +} diff --git a/verify.html b/verify.html new file mode 100644 index 0000000..dc5d2e7 --- /dev/null +++ b/verify.html @@ -0,0 +1,145 @@ + + + + + SOLARGOD ▸ verify + + + + +

SOLARGOD ▸ VERIFY

+
lib.js + ephem.js + scale.js assertions · live Horizons via serve.py proxy. Serve with python3 serve.py then open /verify.html.
+
#gateresultdetail
+
running…
+ + + +