Last active
April 12, 2018 12:32
-
-
Save thure/8d167eb16598a1c96069 to your computer and use it in GitHub Desktop.
Heliocentric coordinates from Kepler orbit and time
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
var moment = require('moment'); | |
module.exports = function(/* orbit, [Date, time String], [format String], [etc] */){ | |
// Takes an orbit definition of the kind used in GEOF stellar-objects and | |
// remaining arguments as moment arguments, returning heliocentric | |
// rectangular coordinates in α's unit (probably AU). | |
// TODO: Does not yet take into account changes in orbit over time. | |
// TODO: Converts between degrees and radians too often, I think. | |
// TODO: May not actually need to keep longitudes between 0 and 2π||360. | |
// TODO: Might be useful to get spherical coordinates (λ,β,Δ?) too. | |
// TODO: Eventually: Rotation, precession, and axial tilt. | |
// TODO: Eventually: Can this do satellites? | |
// Technique based upon http://www.stargazing.net/kepler/ellipse.html | |
var args = Array.prototype.slice.call(arguments) | |
, orbit = args.shift(); | |
var π = Math.PI | |
, rad_deg = π / 180 // rad per deg | |
, deg_rad = 180 / π // deg per rad | |
, ms_d = 24 * 60 * 60 * 1000 // ms per day | |
, cos = function(deg){return Math.cos(rad_deg * deg);} | |
, sin = function(deg){return Math.sin(rad_deg * deg);} | |
, pow = Math.pow | |
, t = moment.apply(this, args).diff(moment(orbit._0.epoch), 'milliseconds') / ms_d | |
, α = orbit.α | |
, e = orbit.e | |
, i = orbit.i | |
, Ω = orbit.Ω | |
, ϖ = orbit.ϖ | |
, P = orbit.orbital_period | |
, M_0 = orbit._0.L - ϖ | |
; | |
// find the mean motion | |
var n = 360 / P; | |
// find the mean anomaly | |
var M = (M_0 + n * t) % 360; | |
M = M < 0 ? M + 360 : M; | |
// find the true anomaly | |
// this is an approximation using equation of the center | |
var ν = M + deg_rad * [ (2 * e - pow(e,3)/4) * sin(M) | |
+ 5/4 * pow(e,2) * sin(2*M) | |
+ 13/12 * pow(e,3) * sin(3*M) ]; | |
// find the radius vector | |
var r = α * (1 - pow(e,2)) / [1 + e * cos(ν)]; | |
var λ = (ν + ϖ - Ω) % 360; | |
λ = λ < 0 ? λ + 360 : λ; | |
var x = r * [cos(Ω) * cos(λ) - sin(Ω) * sin(λ) * cos(i)] | |
, y = r * [sin(Ω) * cos(λ) + cos(Ω) * sin(λ) * cos(i)] | |
, z = r * [sin(λ) * sin(i)] | |
; | |
// var accuracy = (Math.sqrt(Math.pow(x,2) + Math.pow(y,2) + Math.pow(z,2)) - r); | |
return {x: x, y: y, z: z}; | |
}; |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment