#![forbid(
missing_docs,
anonymous_parameters,
unused_extern_crates,
unused_import_braces,
missing_copy_implementations,
trivial_casts,
variant_size_differences,
missing_debug_implementations,
trivial_numeric_casts
)]
#![deny(unused_qualifications, unsafe_code)]
#![deny(clippy::all)]
#![warn(clippy::pedantic)]
#![allow(
clippy::many_single_char_names,
clippy::unreadable_literal,
clippy::excessive_precision,
clippy::must_use_candidate
)]
#![cfg_attr(all(test, feature = "no_std"), allow(unused_imports))]
#![cfg_attr(feature = "no_std", no_std)]
#![warn(box_pointers, unused_results)]
pub mod vsop87a;
pub mod vsop87b;
pub mod vsop87c;
pub mod vsop87d;
pub mod vsop87e;
mod earth_moon;
mod jupiter;
mod mars;
mod mercury;
mod neptune;
mod saturn;
mod uranus;
mod venus;
#[cfg(feature = "no_std")]
use core::f64::consts::PI;
#[cfg(feature = "no_std")]
use libm::{acos, asin, atan, cos, sin, sincos, sqrt};
#[cfg(not(feature = "no_std"))]
use std::f64::consts::PI;
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct KeplerianElements {
ecc: f64,
sma: f64,
incl: f64,
lan: f64,
lper: f64,
l0: f64,
}
impl KeplerianElements {
pub fn eccentricity(&self) -> f64 {
self.ecc
}
pub fn semimajor_axis(&self) -> f64 {
self.sma
}
pub fn inclination(&self) -> f64 {
self.incl
}
pub fn ascending_node(&self) -> f64 {
self.lan
}
pub fn periapsis(&self) -> f64 {
self.lper
}
pub fn mean_anomaly(&self) -> f64 {
self.l0
}
}
impl From<VSOP87Elements> for KeplerianElements {
fn from(elts: VSOP87Elements) -> Self {
#[cfg(feature = "no_std")]
{
let ecc = sqrt(elts.h * elts.h + elts.k * elts.k);
let i = acos(1_f64 - 2_f64 * (elts.p * elts.p + elts.q * elts.q));
let lan = atan(elts.p / elts.q);
let lper = asin(elts.h / ecc);
Self {
ecc,
sma: elts.a,
incl: i,
lan,
lper,
l0: elts.l,
}
}
#[cfg(not(feature = "no_std"))]
{
let ecc = (elts.h * elts.h + elts.k * elts.k).sqrt();
let i = (1_f64 - 2_f64 * (elts.p * elts.p + elts.q * elts.q)).acos();
let lan = (elts.p / elts.q).atan();
let lper = (elts.h / ecc).asin();
Self {
ecc,
sma: elts.a,
incl: i,
lan,
lper,
l0: elts.l,
}
}
}
}
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct RectangularCoordinates {
pub x: f64,
pub y: f64,
pub z: f64,
}
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct SphericalCoordinates {
lon: f64,
lat: f64,
dist: f64,
}
impl SphericalCoordinates {
pub fn longitude(&self) -> f64 {
self.lon
}
pub fn latitude(&self) -> f64 {
self.lat
}
pub fn distance(&self) -> f64 {
self.dist
}
}
#[inline]
fn calculate_t(jde: f64) -> f64 {
(jde - 2_451_545_f64) / 365_250_f64
}
#[inline]
fn calculate_var(t: f64, a: &[f64], b: &[f64], c: &[f64]) -> f64 {
#[cfg(all(
any(target_arch = "x86", target_arch = "x86_64"),
feature = "simd",
not(feature = "no_std")
))]
#[allow(unsafe_code)]
{
if is_x86_feature_detected!("avx") {
unsafe { calculate_var_avx(t, a, b, c) }
} else {
calculate_var_fallback(t, a, b, c)
}
}
#[cfg(any(
all(not(target_arch = "x86"), not(target_arch = "x86_64")),
not(feature = "simd"),
feature = "no_std"
))]
{
calculate_var_fallback(t, a, b, c)
}
}
#[inline]
fn calculate_var_fallback(t: f64, a: &[f64], b: &[f64], c: &[f64]) -> f64 {
#[cfg(not(feature = "no_std"))]
{
a.iter()
.zip(b)
.zip(c)
.fold(0_f64, |term, ((a, b), c)| term + a * (b + c * t).cos())
}
#[cfg(feature = "no_std")]
{
a.iter()
.zip(b)
.zip(c)
.fold(0_f64, |term, ((a, b), c)| term + a * cos(b + c * t))
}
}
#[target_feature(enable = "avx")]
#[cfg(all(
any(target_arch = "x86", target_arch = "x86_64"),
feature = "simd",
not(feature = "no_std")
))]
#[allow(unsafe_code)]
unsafe fn calculate_var_avx(t: f64, a: &[f64], b: &[f64], c: &[f64]) -> f64 {
#[cfg(feature = "no_std")]
use core::{f64, mem};
#[cfg(not(feature = "no_std"))]
use std::{f64, mem};
#[cfg(all(feature = "no_std", target_arch = "x86_64"))]
use core::arch::x86_64::{_mm256_add_pd, _mm256_mul_pd, _mm256_set1_pd, _mm256_set_pd};
#[cfg(all(not(feature = "no_std"), target_arch = "x86_64"))]
use std::arch::x86_64::{_mm256_add_pd, _mm256_mul_pd, _mm256_set1_pd, _mm256_set_pd};
#[cfg(all(feature = "no_std", target_arch = "x86"))]
use core::arch::x86::{_mm256_add_pd, _mm256_mul_pd, _mm256_set1_pd, _mm256_set_pd};
#[cfg(all(not(feature = "no_std"), target_arch = "x86"))]
use std::arch::x86::{_mm256_add_pd, _mm256_mul_pd, _mm256_set1_pd, _mm256_set_pd};
unsafe fn vector_term(
(a1, b1, c1): (f64, f64, f64),
(a2, b2, c2): (f64, f64, f64),
(a3, b3, c3): (f64, f64, f64),
(a4, b4, c4): (f64, f64, f64),
t: f64,
) -> (f64, f64, f64, f64) {
let a = _mm256_set_pd(a1, a2, a3, a4);
let b = _mm256_set_pd(b1, b2, b3, b4);
let c = _mm256_set_pd(c1, c2, c3, c4);
let t = _mm256_set1_pd(t);
let ct = _mm256_mul_pd(c, t);
let bct = _mm256_add_pd(b, ct);
let bct_unpacked: (f64, f64, f64, f64) = mem::transmute(bct);
let bct = _mm256_set_pd(
bct_unpacked.3.cos(),
bct_unpacked.2.cos(),
bct_unpacked.1.cos(),
bct_unpacked.0.cos(),
);
let term = _mm256_mul_pd(a, bct);
let term_unpacked: (f64, f64, f64, f64) = mem::transmute(term);
term_unpacked
}
let iter1 = a.chunks_exact(4);
let iter2 = b.chunks_exact(4);
let iter3 = c.chunks_exact(4);
let remainder = match (iter1.remainder(), iter2.remainder(), iter3.remainder()) {
(&[a1, a2, a3], &[b1, b2, b3], &[c1, c2, c3]) => {
let (_term4, term3, term2, term1) = vector_term(
(a1, b1, c1),
(a2, b2, c2),
(a3, b3, c3),
(f64::NAN, f64::NAN, f64::NAN),
t,
);
term1 + term2 + term3
}
(&[a1, a2], &[b1, b2], &[c1, c2]) => a1 * (b1 + c1 * t).cos() + a2 * (b2 + c2 * t).cos(),
(&[a], &[b], &[c]) => a * (b + c * t).cos(),
(&[], &[], &[]) => 0_f64,
_ => unreachable!(),
};
let res = iter1
.zip(iter2)
.zip(iter3)
.map(|vars| match vars {
((&[a1, a2, a3, a4], &[b1, b2, b3, b4]), &[c1, c2, c3, c4]) => {
let (term4, term3, term2, term1) =
vector_term((a1, b1, c1), (a2, b2, c2), (a3, b3, c3), (a4, b4, c4), t);
term1 + term2 + term3 + term4
}
_ => unreachable!(),
})
.sum::<f64>();
res + remainder
}
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct VSOP87Elements {
pub a: f64,
pub l: f64,
pub k: f64,
pub h: f64,
pub q: f64,
pub p: f64,
}
impl From<KeplerianElements> for VSOP87Elements {
fn from(elts: KeplerianElements) -> Self {
#[cfg(feature = "no_std")]
{
let (lper_sin, lper_cos) = sincos(elts.lper);
let (lan_sin, lan_cos) = sincos(elts.lan);
let incl_sin = sin(elts.incl / 2.0);
Self {
a: elts.sma,
l: elts.l0,
k: elts.ecc * lper_cos,
h: elts.ecc * lper_sin,
q: incl_sin * lan_cos,
p: incl_sin * lan_sin,
}
}
#[cfg(not(feature = "no_std"))]
{
let (lper_sin, lper_cos) = elts.lper.sin_cos();
let (lan_sin, lan_cos) = elts.lan.sin_cos();
let incl_sin = (elts.incl / 2.0).sin();
Self {
a: elts.sma,
l: elts.l0,
k: elts.ecc * lper_cos,
h: elts.ecc * lper_sin,
q: incl_sin * lan_cos,
p: incl_sin * lan_sin,
}
}
}
}
pub fn mercury(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &mercury::A0[0], &mercury::A0[1], &mercury::A0[2]);
let a1 = calculate_var(t, &mercury::A1[0], &mercury::A1[1], &mercury::A1[2]);
let a2 = calculate_var(t, &mercury::A2[0], &mercury::A2[1], &mercury::A2[2]);
let l0 = calculate_var(t, &mercury::L0[0], &mercury::L0[1], &mercury::L0[2]);
let l1 = calculate_var(t, &mercury::L1[0], &mercury::L1[1], &mercury::L1[2]);
let l2 = calculate_var(t, &mercury::L2[0], &mercury::L2[1], &mercury::L2[2]);
let l3 = calculate_var(t, &mercury::L3[0], &mercury::L3[1], &mercury::L3[2]);
let k0 = calculate_var(t, &mercury::K0[0], &mercury::K0[1], &mercury::K0[2]);
let k1 = calculate_var(t, &mercury::K1[0], &mercury::K1[1], &mercury::K1[2]);
let k2 = calculate_var(t, &mercury::K2[0], &mercury::K2[1], &mercury::K2[2]);
let k3 = calculate_var(t, &mercury::K3[0], &mercury::K3[1], &mercury::K3[2]);
let k4 = calculate_var(t, &mercury::K4[0], &mercury::K4[1], &mercury::K4[2]);
let k5 = calculate_var(t, &mercury::K5[0], &mercury::K5[1], &mercury::K5[2]);
let h0 = calculate_var(t, &mercury::H0[0], &mercury::H0[1], &mercury::H0[2]);
let h1 = calculate_var(t, &mercury::H1[0], &mercury::H1[1], &mercury::H1[2]);
let h2 = calculate_var(t, &mercury::H2[0], &mercury::H2[1], &mercury::H2[2]);
let h3 = calculate_var(t, &mercury::H3[0], &mercury::H3[1], &mercury::H3[2]);
let h4 = calculate_var(t, &mercury::H4[0], &mercury::H4[1], &mercury::H4[2]);
let h5 = calculate_var(t, &mercury::H5[0], &mercury::H5[1], &mercury::H5[2]);
let q0 = calculate_var(t, &mercury::Q0[0], &mercury::Q0[1], &mercury::Q0[2]);
let q1 = calculate_var(t, &mercury::Q1[0], &mercury::Q1[1], &mercury::Q1[2]);
let q2 = calculate_var(t, &mercury::Q2[0], &mercury::Q2[1], &mercury::Q2[2]);
let q3 = calculate_var(t, &mercury::Q3[0], &mercury::Q3[1], &mercury::Q3[2]);
let q4 = calculate_var(t, &mercury::Q4[0], &mercury::Q4[1], &mercury::Q4[2]);
let q5 = calculate_var(t, &mercury::Q5[0], &mercury::Q5[1], &mercury::Q5[2]);
let p0 = calculate_var(t, &mercury::P0[0], &mercury::P0[1], &mercury::P0[2]);
let p1 = calculate_var(t, &mercury::P1[0], &mercury::P1[1], &mercury::P1[2]);
let p2 = calculate_var(t, &mercury::P2[0], &mercury::P2[1], &mercury::P2[2]);
let p3 = calculate_var(t, &mercury::P3[0], &mercury::P3[1], &mercury::P3[2]);
let p4 = calculate_var(t, &mercury::P4[0], &mercury::P4[1], &mercury::P4[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4 + k5 * t5;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4 + h5 * t5;
let q = q0 + q1 * t + q2 * t2 + q3 * t3 + q4 * t4 + q5 * t5;
let p = p0 + p1 * t + p2 * t2 + p3 * t3 + p4 * t4;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
pub fn venus(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &venus::A0[0], &venus::A0[1], &venus::A0[2]);
let a1 = calculate_var(t, &venus::A1[0], &venus::A1[1], &venus::A1[2]);
let a2 = calculate_var(t, &venus::A2[0], &venus::A2[1], &venus::A2[2]);
let l0 = calculate_var(t, &venus::L0[0], &venus::L0[1], &venus::L0[2]);
let l1 = calculate_var(t, &venus::L1[0], &venus::L1[1], &venus::L1[2]);
let l2 = calculate_var(t, &venus::L2[0], &venus::L2[1], &venus::L2[2]);
let l3 = calculate_var(t, &venus::L3[0], &venus::L3[1], &venus::L3[2]);
let k0 = calculate_var(t, &venus::K0[0], &venus::K0[1], &venus::K0[2]);
let k1 = calculate_var(t, &venus::K1[0], &venus::K1[1], &venus::K1[2]);
let k2 = calculate_var(t, &venus::K2[0], &venus::K2[1], &venus::K2[2]);
let k3 = calculate_var(t, &venus::K3[0], &venus::K3[1], &venus::K3[2]);
let k4 = calculate_var(t, &venus::K4[0], &venus::K4[1], &venus::K4[2]);
let k5 = calculate_var(t, &venus::K5[0], &venus::K5[1], &venus::K5[2]);
let h0 = calculate_var(t, &venus::H0[0], &venus::H0[1], &venus::H0[2]);
let h1 = calculate_var(t, &venus::H1[0], &venus::H1[1], &venus::H1[2]);
let h2 = calculate_var(t, &venus::H2[0], &venus::H2[1], &venus::H2[2]);
let h3 = calculate_var(t, &venus::H3[0], &venus::H3[1], &venus::H3[2]);
let h4 = calculate_var(t, &venus::H4[0], &venus::H4[1], &venus::H4[2]);
let h5 = calculate_var(t, &venus::H5[0], &venus::H5[1], &venus::H5[2]);
let q0 = calculate_var(t, &venus::Q0[0], &venus::Q0[1], &venus::Q0[2]);
let q1 = calculate_var(t, &venus::Q1[0], &venus::Q1[1], &venus::Q1[2]);
let q2 = calculate_var(t, &venus::Q2[0], &venus::Q2[1], &venus::Q2[2]);
let q3 = calculate_var(t, &venus::Q3[0], &venus::Q3[1], &venus::Q3[2]);
let q4 = calculate_var(t, &venus::Q4[0], &venus::Q4[1], &venus::Q4[2]);
let q5 = calculate_var(t, &venus::Q5[0], &venus::Q5[1], &venus::Q5[2]);
let p0 = calculate_var(t, &venus::P0[0], &venus::P0[1], &venus::P0[2]);
let p1 = calculate_var(t, &venus::P1[0], &venus::P1[1], &venus::P1[2]);
let p2 = calculate_var(t, &venus::P2[0], &venus::P2[1], &venus::P2[2]);
let p3 = calculate_var(t, &venus::P3[0], &venus::P3[1], &venus::P3[2]);
let p4 = calculate_var(t, &venus::P4[0], &venus::P4[1], &venus::P4[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4 + k5 * t5;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4 + h5 * t5;
let q = q0 + q1 * t + q2 * t2 + q3 * t3 + q4 * t4 + q5 * t5;
let p = p0 + p1 * t + p2 * t2 + p3 * t3 + p4 * t4;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
#[allow(clippy::too_many_lines)]
pub fn earth_moon(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(
t,
&earth_moon::A0[0],
&earth_moon::A0[1],
&earth_moon::A0[2],
);
let a1 = calculate_var(
t,
&earth_moon::A1[0],
&earth_moon::A1[1],
&earth_moon::A1[2],
);
let a2 = calculate_var(
t,
&earth_moon::A2[0],
&earth_moon::A2[1],
&earth_moon::A2[2],
);
let l0 = calculate_var(
t,
&earth_moon::L0[0],
&earth_moon::L0[1],
&earth_moon::L0[2],
);
let l1 = calculate_var(
t,
&earth_moon::L1[0],
&earth_moon::L1[1],
&earth_moon::L1[2],
);
let l2 = calculate_var(
t,
&earth_moon::L2[0],
&earth_moon::L2[1],
&earth_moon::L2[2],
);
let l3 = calculate_var(
t,
&earth_moon::L3[0],
&earth_moon::L3[1],
&earth_moon::L3[2],
);
let l4 = calculate_var(
t,
&earth_moon::L4[0],
&earth_moon::L4[1],
&earth_moon::L4[2],
);
let l5 = calculate_var(
t,
&earth_moon::L5[0],
&earth_moon::L5[1],
&earth_moon::L5[2],
);
let k0 = calculate_var(
t,
&earth_moon::K0[0],
&earth_moon::K0[1],
&earth_moon::K0[2],
);
let k1 = calculate_var(
t,
&earth_moon::K1[0],
&earth_moon::K1[1],
&earth_moon::K1[2],
);
let k2 = calculate_var(
t,
&earth_moon::K2[0],
&earth_moon::K2[1],
&earth_moon::K2[2],
);
let k3 = calculate_var(
t,
&earth_moon::K3[0],
&earth_moon::K3[1],
&earth_moon::K3[2],
);
let k4 = calculate_var(
t,
&earth_moon::K4[0],
&earth_moon::K4[1],
&earth_moon::K4[2],
);
let k5 = calculate_var(
t,
&earth_moon::K5[0],
&earth_moon::K5[1],
&earth_moon::K5[2],
);
let h0 = calculate_var(
t,
&earth_moon::H0[0],
&earth_moon::H0[1],
&earth_moon::H0[2],
);
let h1 = calculate_var(
t,
&earth_moon::H1[0],
&earth_moon::H1[1],
&earth_moon::H1[2],
);
let h2 = calculate_var(
t,
&earth_moon::H2[0],
&earth_moon::H2[1],
&earth_moon::H2[2],
);
let h3 = calculate_var(
t,
&earth_moon::H3[0],
&earth_moon::H3[1],
&earth_moon::H3[2],
);
let h4 = calculate_var(
t,
&earth_moon::H4[0],
&earth_moon::H4[1],
&earth_moon::H4[2],
);
let h5 = calculate_var(
t,
&earth_moon::H5[0],
&earth_moon::H5[1],
&earth_moon::H5[2],
);
let q0 = calculate_var(
t,
&earth_moon::Q0[0],
&earth_moon::Q0[1],
&earth_moon::Q0[2],
);
let q1 = calculate_var(
t,
&earth_moon::Q1[0],
&earth_moon::Q1[1],
&earth_moon::Q1[2],
);
let q2 = calculate_var(
t,
&earth_moon::Q2[0],
&earth_moon::Q2[1],
&earth_moon::Q2[2],
);
let q3 = calculate_var(
t,
&earth_moon::Q3[0],
&earth_moon::Q3[1],
&earth_moon::Q3[2],
);
let q4 = calculate_var(
t,
&earth_moon::Q4[0],
&earth_moon::Q4[1],
&earth_moon::Q4[2],
);
let q5 = calculate_var(
t,
&earth_moon::Q5[0],
&earth_moon::Q5[1],
&earth_moon::Q5[2],
);
let p0 = calculate_var(
t,
&earth_moon::P0[0],
&earth_moon::P0[1],
&earth_moon::P0[2],
);
let p1 = calculate_var(
t,
&earth_moon::P1[0],
&earth_moon::P1[1],
&earth_moon::P1[2],
);
let p2 = calculate_var(
t,
&earth_moon::P2[0],
&earth_moon::P2[1],
&earth_moon::P2[2],
);
let p3 = calculate_var(
t,
&earth_moon::P3[0],
&earth_moon::P3[1],
&earth_moon::P3[2],
);
let p4 = calculate_var(
t,
&earth_moon::P4[0],
&earth_moon::P4[1],
&earth_moon::P4[2],
);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3 + l4 * t4 + l5 * t5) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4 + k5 * t5;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4 + h5 * t5;
let q = q0 + q1 * t + q2 * t2 + q3 * t3 + q4 * t4 + q5 * t5;
let p = p0 + p1 * t + p2 * t2 + p3 * t3 + p4 * t4;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
pub fn mars(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &mars::A0[0], &mars::A0[1], &mars::A0[2]);
let a1 = calculate_var(t, &mars::A1[0], &mars::A1[1], &mars::A1[2]);
let a2 = calculate_var(t, &mars::A2[0], &mars::A2[1], &mars::A2[2]);
let l0 = calculate_var(t, &mars::L0[0], &mars::L0[1], &mars::L0[2]);
let l1 = calculate_var(t, &mars::L1[0], &mars::L1[1], &mars::L1[2]);
let l2 = calculate_var(t, &mars::L2[0], &mars::L2[1], &mars::L2[2]);
let l3 = calculate_var(t, &mars::L3[0], &mars::L3[1], &mars::L3[2]);
let l4 = calculate_var(t, &mars::L4[0], &mars::L4[1], &mars::L4[2]);
let l5 = calculate_var(t, &mars::L5[0], &mars::L5[1], &mars::L5[2]);
let k0 = calculate_var(t, &mars::K0[0], &mars::K0[1], &mars::K0[2]);
let k1 = calculate_var(t, &mars::K1[0], &mars::K1[1], &mars::K1[2]);
let k2 = calculate_var(t, &mars::K2[0], &mars::K2[1], &mars::K2[2]);
let k3 = calculate_var(t, &mars::K3[0], &mars::K3[1], &mars::K3[2]);
let k4 = calculate_var(t, &mars::K4[0], &mars::K4[1], &mars::K4[2]);
let k5 = calculate_var(t, &mars::K5[0], &mars::K5[1], &mars::K5[2]);
let h0 = calculate_var(t, &mars::H0[0], &mars::H0[1], &mars::H0[2]);
let h1 = calculate_var(t, &mars::H1[0], &mars::H1[1], &mars::H1[2]);
let h2 = calculate_var(t, &mars::H2[0], &mars::H2[1], &mars::H2[2]);
let h3 = calculate_var(t, &mars::H3[0], &mars::H3[1], &mars::H3[2]);
let h4 = calculate_var(t, &mars::H4[0], &mars::H4[1], &mars::H4[2]);
let h5 = calculate_var(t, &mars::H5[0], &mars::H5[1], &mars::H5[2]);
let q0 = calculate_var(t, &mars::Q0[0], &mars::Q0[1], &mars::Q0[2]);
let q1 = calculate_var(t, &mars::Q1[0], &mars::Q1[1], &mars::Q1[2]);
let q2 = calculate_var(t, &mars::Q2[0], &mars::Q2[1], &mars::Q2[2]);
let q3 = calculate_var(t, &mars::Q3[0], &mars::Q3[1], &mars::Q3[2]);
let q4 = calculate_var(t, &mars::Q4[0], &mars::Q4[1], &mars::Q4[2]);
let q5 = calculate_var(t, &mars::Q5[0], &mars::Q5[1], &mars::Q5[2]);
let p0 = calculate_var(t, &mars::P0[0], &mars::P0[1], &mars::P0[2]);
let p1 = calculate_var(t, &mars::P1[0], &mars::P1[1], &mars::P1[2]);
let p2 = calculate_var(t, &mars::P2[0], &mars::P2[1], &mars::P2[2]);
let p3 = calculate_var(t, &mars::P3[0], &mars::P3[1], &mars::P3[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3 + l4 * t4 + l5 * t5) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4 + k5 * t5;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4 + h5 * t5;
let q = q0 + q1 * t + q2 * t2 + q3 * t3 + q4 * t4 + q5 * t5;
let p = p0 + p1 * t + p2 * t2 + p3 * t3;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
pub fn jupiter(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &jupiter::A0[0], &jupiter::A0[1], &jupiter::A0[2]);
let a1 = calculate_var(t, &jupiter::A1[0], &jupiter::A1[1], &jupiter::A1[2]);
let a2 = calculate_var(t, &jupiter::A2[0], &jupiter::A2[1], &jupiter::A2[2]);
let a3 = calculate_var(t, &jupiter::A3[0], &jupiter::A3[1], &jupiter::A3[2]);
let a4 = calculate_var(t, &jupiter::A4[0], &jupiter::A4[1], &jupiter::A4[2]);
let a5 = calculate_var(t, &jupiter::A5[0], &jupiter::A5[1], &jupiter::A5[2]);
let l0 = calculate_var(t, &jupiter::L0[0], &jupiter::L0[1], &jupiter::L0[2]);
let l1 = calculate_var(t, &jupiter::L1[0], &jupiter::L1[1], &jupiter::L1[2]);
let l2 = calculate_var(t, &jupiter::L2[0], &jupiter::L2[1], &jupiter::L2[2]);
let l3 = calculate_var(t, &jupiter::L3[0], &jupiter::L3[1], &jupiter::L3[2]);
let l4 = calculate_var(t, &jupiter::L4[0], &jupiter::L4[1], &jupiter::L4[2]);
let l5 = calculate_var(t, &jupiter::L5[0], &jupiter::L5[1], &jupiter::L5[2]);
let k0 = calculate_var(t, &jupiter::K0[0], &jupiter::K0[1], &jupiter::K0[2]);
let k1 = calculate_var(t, &jupiter::K1[0], &jupiter::K1[1], &jupiter::K1[2]);
let k2 = calculate_var(t, &jupiter::K2[0], &jupiter::K2[1], &jupiter::K2[2]);
let k3 = calculate_var(t, &jupiter::K3[0], &jupiter::K3[1], &jupiter::K3[2]);
let k4 = calculate_var(t, &jupiter::K4[0], &jupiter::K4[1], &jupiter::K4[2]);
let h0 = calculate_var(t, &jupiter::H0[0], &jupiter::H0[1], &jupiter::H0[2]);
let h1 = calculate_var(t, &jupiter::H1[0], &jupiter::H1[1], &jupiter::H1[2]);
let h2 = calculate_var(t, &jupiter::H2[0], &jupiter::H2[1], &jupiter::H2[2]);
let h3 = calculate_var(t, &jupiter::H3[0], &jupiter::H3[1], &jupiter::H3[2]);
let h4 = calculate_var(t, &jupiter::H4[0], &jupiter::H4[1], &jupiter::H4[2]);
let q0 = calculate_var(t, &jupiter::Q0[0], &jupiter::Q0[1], &jupiter::Q0[2]);
let q1 = calculate_var(t, &jupiter::Q1[0], &jupiter::Q1[1], &jupiter::Q1[2]);
let q2 = calculate_var(t, &jupiter::Q2[0], &jupiter::Q2[1], &jupiter::Q2[2]);
let q3 = calculate_var(t, &jupiter::Q3[0], &jupiter::Q3[1], &jupiter::Q3[2]);
let p0 = calculate_var(t, &jupiter::P0[0], &jupiter::P0[1], &jupiter::P0[2]);
let p1 = calculate_var(t, &jupiter::P1[0], &jupiter::P1[1], &jupiter::P1[2]);
let p2 = calculate_var(t, &jupiter::P2[0], &jupiter::P2[1], &jupiter::P2[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2 + a3 * t3 + a4 * t4 + a5 * t5;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3 + l4 * t4 + l5 * t5) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4;
let q = q0 + q1 * t + q2 * t2 + q3 * t3;
let p = p0 + p1 * t + p2 * t2;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
pub fn saturn(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &saturn::A0[0], &saturn::A0[1], &saturn::A0[2]);
let a1 = calculate_var(t, &saturn::A1[0], &saturn::A1[1], &saturn::A1[2]);
let a2 = calculate_var(t, &saturn::A2[0], &saturn::A2[1], &saturn::A2[2]);
let a3 = calculate_var(t, &saturn::A3[0], &saturn::A3[1], &saturn::A3[2]);
let a4 = calculate_var(t, &saturn::A4[0], &saturn::A4[1], &saturn::A4[2]);
let a5 = calculate_var(t, &saturn::A5[0], &saturn::A5[1], &saturn::A5[2]);
let l0 = calculate_var(t, &saturn::L0[0], &saturn::L0[1], &saturn::L0[2]);
let l1 = calculate_var(t, &saturn::L1[0], &saturn::L1[1], &saturn::L1[2]);
let l2 = calculate_var(t, &saturn::L2[0], &saturn::L2[1], &saturn::L2[2]);
let l3 = calculate_var(t, &saturn::L3[0], &saturn::L3[1], &saturn::L3[2]);
let l4 = calculate_var(t, &saturn::L4[0], &saturn::L4[1], &saturn::L4[2]);
let l5 = calculate_var(t, &saturn::L5[0], &saturn::L5[1], &saturn::L5[2]);
let k0 = calculate_var(t, &saturn::K0[0], &saturn::K0[1], &saturn::K0[2]);
let k1 = calculate_var(t, &saturn::K1[0], &saturn::K1[1], &saturn::K1[2]);
let k2 = calculate_var(t, &saturn::K2[0], &saturn::K2[1], &saturn::K2[2]);
let k3 = calculate_var(t, &saturn::K3[0], &saturn::K3[1], &saturn::K3[2]);
let k4 = calculate_var(t, &saturn::K4[0], &saturn::K4[1], &saturn::K4[2]);
let k5 = calculate_var(t, &saturn::K5[0], &saturn::K5[1], &saturn::K5[2]);
let h0 = calculate_var(t, &saturn::H0[0], &saturn::H0[1], &saturn::H0[2]);
let h1 = calculate_var(t, &saturn::H1[0], &saturn::H1[1], &saturn::H1[2]);
let h2 = calculate_var(t, &saturn::H2[0], &saturn::H2[1], &saturn::H2[2]);
let h3 = calculate_var(t, &saturn::H3[0], &saturn::H3[1], &saturn::H3[2]);
let h4 = calculate_var(t, &saturn::H4[0], &saturn::H4[1], &saturn::H4[2]);
let h5 = calculate_var(t, &saturn::H5[0], &saturn::H5[1], &saturn::H5[2]);
let q0 = calculate_var(t, &saturn::Q0[0], &saturn::Q0[1], &saturn::Q0[2]);
let q1 = calculate_var(t, &saturn::Q1[0], &saturn::Q1[1], &saturn::Q1[2]);
let q2 = calculate_var(t, &saturn::Q2[0], &saturn::Q2[1], &saturn::Q2[2]);
let q3 = calculate_var(t, &saturn::Q3[0], &saturn::Q3[1], &saturn::Q3[2]);
let q4 = calculate_var(t, &saturn::Q4[0], &saturn::Q4[1], &saturn::Q4[2]);
let p0 = calculate_var(t, &saturn::P0[0], &saturn::P0[1], &saturn::P0[2]);
let p1 = calculate_var(t, &saturn::P1[0], &saturn::P1[1], &saturn::P1[2]);
let p2 = calculate_var(t, &saturn::P2[0], &saturn::P2[1], &saturn::P2[2]);
let p3 = calculate_var(t, &saturn::P3[0], &saturn::P3[1], &saturn::P3[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2 + a3 * t3 + a4 * t4 + a5 * t5;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3 + l4 * t4 + l5 * t5) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4 + k5 * t5;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4 + h5 * t5;
let q = q0 + q1 * t + q2 * t2 + q3 * t3 + q4 * t4;
let p = p0 + p1 * t + p2 * t2 + p3 * t3;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
pub fn uranus(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &uranus::A0[0], &uranus::A0[1], &uranus::A0[2]);
let a1 = calculate_var(t, &uranus::A1[0], &uranus::A1[1], &uranus::A1[2]);
let a2 = calculate_var(t, &uranus::A2[0], &uranus::A2[1], &uranus::A2[2]);
let a3 = calculate_var(t, &uranus::A3[0], &uranus::A3[1], &uranus::A3[2]);
let a4 = calculate_var(t, &uranus::A4[0], &uranus::A4[1], &uranus::A4[2]);
let a5 = calculate_var(t, &uranus::A5[0], &uranus::A5[1], &uranus::A5[2]);
let l0 = calculate_var(t, &uranus::L0[0], &uranus::L0[1], &uranus::L0[2]);
let l1 = calculate_var(t, &uranus::L1[0], &uranus::L1[1], &uranus::L1[2]);
let l2 = calculate_var(t, &uranus::L2[0], &uranus::L2[1], &uranus::L2[2]);
let l3 = calculate_var(t, &uranus::L3[0], &uranus::L3[1], &uranus::L3[2]);
let l4 = calculate_var(t, &uranus::L4[0], &uranus::L4[1], &uranus::L4[2]);
let l5 = calculate_var(t, &uranus::L5[0], &uranus::L5[1], &uranus::L5[2]);
let k0 = calculate_var(t, &uranus::K0[0], &uranus::K0[1], &uranus::K0[2]);
let k1 = calculate_var(t, &uranus::K1[0], &uranus::K1[1], &uranus::K1[2]);
let k2 = calculate_var(t, &uranus::K2[0], &uranus::K2[1], &uranus::K2[2]);
let k3 = calculate_var(t, &uranus::K3[0], &uranus::K3[1], &uranus::K3[2]);
let k4 = calculate_var(t, &uranus::K4[0], &uranus::K4[1], &uranus::K4[2]);
let h0 = calculate_var(t, &uranus::H0[0], &uranus::H0[1], &uranus::H0[2]);
let h1 = calculate_var(t, &uranus::H1[0], &uranus::H1[1], &uranus::H1[2]);
let h2 = calculate_var(t, &uranus::H2[0], &uranus::H2[1], &uranus::H2[2]);
let h3 = calculate_var(t, &uranus::H3[0], &uranus::H3[1], &uranus::H3[2]);
let h4 = calculate_var(t, &uranus::H4[0], &uranus::H4[1], &uranus::H4[2]);
let q0 = calculate_var(t, &uranus::Q0[0], &uranus::Q0[1], &uranus::Q0[2]);
let q1 = calculate_var(t, &uranus::Q1[0], &uranus::Q1[1], &uranus::Q1[2]);
let q2 = calculate_var(t, &uranus::Q2[0], &uranus::Q2[1], &uranus::Q2[2]);
let q3 = calculate_var(t, &uranus::Q3[0], &uranus::Q3[1], &uranus::Q3[2]);
let p0 = calculate_var(t, &uranus::P0[0], &uranus::P0[1], &uranus::P0[2]);
let p1 = calculate_var(t, &uranus::P1[0], &uranus::P1[1], &uranus::P1[2]);
let p2 = calculate_var(t, &uranus::P2[0], &uranus::P2[1], &uranus::P2[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2 + a3 * t3 + a4 * t4 + a5 * t5;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3 + l4 * t4 + l5 * t5) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4;
let q = q0 + q1 * t + q2 * t2 + q3 * t3;
let p = p0 + p1 * t + p2 * t2;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}
pub fn neptune(jde: f64) -> VSOP87Elements {
let t = calculate_t(jde);
let a0 = calculate_var(t, &neptune::A0[0], &neptune::A0[1], &neptune::A0[2]);
let a1 = calculate_var(t, &neptune::A1[0], &neptune::A1[1], &neptune::A1[2]);
let a2 = calculate_var(t, &neptune::A2[0], &neptune::A2[1], &neptune::A2[2]);
let a3 = calculate_var(t, &neptune::A3[0], &neptune::A3[1], &neptune::A3[2]);
let a4 = calculate_var(t, &neptune::A4[0], &neptune::A4[1], &neptune::A4[2]);
let a5 = calculate_var(t, &neptune::A5[0], &neptune::A5[1], &neptune::A5[2]);
let l0 = calculate_var(t, &neptune::L0[0], &neptune::L0[1], &neptune::L0[2]);
let l1 = calculate_var(t, &neptune::L1[0], &neptune::L1[1], &neptune::L1[2]);
let l2 = calculate_var(t, &neptune::L2[0], &neptune::L2[1], &neptune::L2[2]);
let l3 = calculate_var(t, &neptune::L3[0], &neptune::L3[1], &neptune::L3[2]);
let l4 = calculate_var(t, &neptune::L4[0], &neptune::L4[1], &neptune::L4[2]);
let l5 = calculate_var(t, &neptune::L5[0], &neptune::L5[1], &neptune::L5[2]);
let k0 = calculate_var(t, &neptune::K0[0], &neptune::K0[1], &neptune::K0[2]);
let k1 = calculate_var(t, &neptune::K1[0], &neptune::K1[1], &neptune::K1[2]);
let k2 = calculate_var(t, &neptune::K2[0], &neptune::K2[1], &neptune::K2[2]);
let k3 = calculate_var(t, &neptune::K3[0], &neptune::K3[1], &neptune::K3[2]);
let k4 = calculate_var(t, &neptune::K4[0], &neptune::K4[1], &neptune::K4[2]);
let k5 = calculate_var(t, &neptune::K5[0], &neptune::K5[1], &neptune::K5[2]);
let h0 = calculate_var(t, &neptune::H0[0], &neptune::H0[1], &neptune::H0[2]);
let h1 = calculate_var(t, &neptune::H1[0], &neptune::H1[1], &neptune::H1[2]);
let h2 = calculate_var(t, &neptune::H2[0], &neptune::H2[1], &neptune::H2[2]);
let h3 = calculate_var(t, &neptune::H3[0], &neptune::H3[1], &neptune::H3[2]);
let h4 = calculate_var(t, &neptune::H4[0], &neptune::H4[1], &neptune::H4[2]);
let h5 = calculate_var(t, &neptune::H5[0], &neptune::H5[1], &neptune::H5[2]);
let q0 = calculate_var(t, &neptune::Q0[0], &neptune::Q0[1], &neptune::Q0[2]);
let q1 = calculate_var(t, &neptune::Q1[0], &neptune::Q1[1], &neptune::Q1[2]);
let q2 = calculate_var(t, &neptune::Q2[0], &neptune::Q2[1], &neptune::Q2[2]);
let q3 = calculate_var(t, &neptune::Q3[0], &neptune::Q3[1], &neptune::Q3[2]);
let p0 = calculate_var(t, &neptune::P0[0], &neptune::P0[1], &neptune::P0[2]);
let p1 = calculate_var(t, &neptune::P1[0], &neptune::P1[1], &neptune::P1[2]);
let p2 = calculate_var(t, &neptune::P2[0], &neptune::P2[1], &neptune::P2[2]);
let t2 = t * t;
let t3 = t2 * t;
let t4 = t2 * t2;
let t5 = t2 * t3;
let a = a0 + a1 * t + a2 * t2 + a3 * t3 + a4 * t4 + a5 * t5;
let l = (l0 + l1 * t + l2 * t2 + l3 * t3 + l4 * t4 + l5 * t5) % (2_f64 * PI);
let k = k0 + k1 * t + k2 * t2 + k3 * t3 + k4 * t4 + k5 * t5;
let h = h0 + h1 * t + h2 * t2 + h3 * t3 + h4 * t4 + h5 * t5;
let q = q0 + q1 * t + q2 * t2 + q3 * t3;
let p = p0 + p1 * t + p2 * t2;
VSOP87Elements {
a,
l: if l > 0_f64 { l } else { 2_f64 * PI + l },
k,
h,
q,
p,
}
}