From 12ab2221611336f6d209d643574fdc2d2b1fac16 Mon Sep 17 00:00:00 2001
From: "David A. Madore" <david+git@madore.org>
Date: Tue, 17 Apr 2012 21:09:36 +0200
Subject: Various small fixes + add a simple test (computing solar
coordinates).
---
org/madore/ephem/Ephem.java | 4 +++-
org/madore/ephem/Frames.java | 4 ++--
org/madore/ephem/Test.java | 44 ++++++++++++++++++++++++++++++++++
org/madore/ephem/Time.java | 56 ++++++++++++++++++++++++++++++++++++++++++++
org/madore/ephem/VSOP87.java | 4 ++++
5 files changed, 109 insertions(+), 3 deletions(-)
create mode 100644 org/madore/ephem/Test.java
create mode 100644 org/madore/ephem/Time.java
(limited to 'org/madore')
diff --git a/org/madore/ephem/Ephem.java b/org/madore/ephem/Ephem.java
index 380835c..4afbc09 100644
--- a/org/madore/ephem/Ephem.java
+++ b/org/madore/ephem/Ephem.java
@@ -2,9 +2,11 @@ package org.madore.ephem;
public final class Ephem {
- public static double fromJD(double jd) {
+ public static double fromJd(double jd) {
// Convert Julian date to Julian centuries from J2000
return (jd-2451545)/36525;
}
+ public static final double auSeconds = 499.0047838362; // TDB seconds!
+
}
diff --git a/org/madore/ephem/Frames.java b/org/madore/ephem/Frames.java
index a8a4667..3671b3c 100644
--- a/org/madore/ephem/Frames.java
+++ b/org/madore/ephem/Frames.java
@@ -16,7 +16,7 @@ public final class Frames {
throw new IllegalArgumentException("Vector constructor expects 3 coordinates");
this.v = v;
}
- public Vector(Double... v) {
+ public Vector(Double[] v) {
this(unboxDoubleArray(v));
}
public double dotprod(Vector w) {
@@ -40,7 +40,7 @@ public final class Frames {
throw new IllegalArgumentException("Rotation constructor expects 4 coordinates");
this.q = q;
}
- public Rotation(Double... q) {
+ public Rotation(Double[] q) {
this(unboxDoubleArray(q));
}
public Vector apply(double[] v) {
diff --git a/org/madore/ephem/Test.java b/org/madore/ephem/Test.java
new file mode 100644
index 0000000..70b6663
--- /dev/null
+++ b/org/madore/ephem/Test.java
@@ -0,0 +1,44 @@
+package org.madore.ephem;
+
+import java.util.Date;
+
+public final class Test {
+
+ public static void main(String[] args) {
+ Date now = new Date();
+ long nowTime = now.getTime();
+ int nowUtcMjd = (int)(nowTime/86400000) + 40587;
+ double nowUtcSeconds = (nowTime%86400000)/1000.;
+ double nowTtJd = Time.utcToTt(nowUtcMjd, nowUtcSeconds);
+ System.out.format("UTC = %dms after Unix Epoch\nJD(TT) = %.10f\n", nowTime, nowTtJd);
+ double earthLat = VSOP87.getVal(VSOP87.Planet.EARTH, VSOP87.Variable.LAT, nowTtJd);
+ double earthLong = VSOP87.getVal(VSOP87.Planet.EARTH, VSOP87.Variable.LONG, nowTtJd);
+ double earthDist = VSOP87.getVal(VSOP87.Planet.EARTH, VSOP87.Variable.DIST, nowTtJd);
+ System.out.format("VSOP87 gives: ♁ ecliptic(J2000) lat=%.6f° long=%.6f° dist=%.6fAU\n", earthLat/Comput.degree, earthLong/Comput.degree, earthDist);
+ double delayedJd = nowTtJd - earthDist*Ephem.auSeconds/86400;
+ double sunLat = -VSOP87.getVal(VSOP87.Planet.EARTH, VSOP87.Variable.LAT, delayedJd);
+ double sunLong = Math.PI + VSOP87.getVal(VSOP87.Planet.EARTH, VSOP87.Variable.LONG, delayedJd);
+ double sunDist = VSOP87.getVal(VSOP87.Planet.EARTH, VSOP87.Variable.DIST, delayedJd);
+ System.out.format("⇒ ☉ ecliptic(J2000) lat=%.6f° long=%.6f° dist=%.6fAU\n", sunLat/Comput.degree, sunLong/Comput.degree, sunDist);
+ Frames.Vector sun0 = new Frames.Vector(sunDist*Math.cos(sunLat)*Math.cos(sunLong),
+ sunDist*Math.cos(sunLat)*Math.sin(sunLong),
+ sunDist*Math.sin(sunLat));
+ System.out.format("Check: ☉ ecliptic(J2000) lat=%.6f° long=%.6f°\n",
+ Math.atan(sun0.v[2]/Math.sqrt(sun0.v[0]*sun0.v[0]+sun0.v[1]*sun0.v[1]))/Comput.degree,
+ Math.atan2(sun0.v[1], sun0.v[0])/Comput.degree);
+ Frames.Vector sun;
+ sun = Frames.Rotation.rotx(-Precession.epsilon0*Comput.arcsecond).apply(sun0);
+ System.out.format("☉ equatorial(J2000) lat=%.6f° long=%.6f°\n",
+ Math.atan(sun.v[2]/Math.sqrt(sun.v[0]*sun.v[0]+sun.v[1]*sun.v[1]))/Comput.degree,
+ Math.atan2(sun.v[1], sun.v[0])/Comput.degree);
+ sun = Precession.matrix(Ephem.fromJd(nowTtJd)).apply(sun0);
+ System.out.format("☉ equatorial(mean) lat=%.6f° long=%.6f°\n",
+ Math.atan(sun.v[2]/Math.sqrt(sun.v[0]*sun.v[0]+sun.v[1]*sun.v[1]))/Comput.degree,
+ Math.atan2(sun.v[1], sun.v[0])/Comput.degree);
+ sun = Nutation.matrix(Ephem.fromJd(nowTtJd)).apply(sun);
+ System.out.format("☉ equatorial(true) lat=%.6f° long=%.6f°\n",
+ Math.atan(sun.v[2]/Math.sqrt(sun.v[0]*sun.v[0]+sun.v[1]*sun.v[1]))/Comput.degree,
+ Math.atan2(sun.v[1], sun.v[0])/Comput.degree);
+ }
+
+}
diff --git a/org/madore/ephem/Time.java b/org/madore/ephem/Time.java
new file mode 100644
index 0000000..3d2cb3e
--- /dev/null
+++ b/org/madore/ephem/Time.java
@@ -0,0 +1,56 @@
+package org.madore.ephem;
+
+public final class Time {
+
+ public static final class LeapSecond {
+ public int mjd;
+ public int offset;
+ public LeapSecond(int mjd, int offset) {
+ this.mjd = mjd;
+ this.offset = offset;
+ }
+ }
+
+ public static final LeapSecond[] leapSecondsTable = new LeapSecond[] {
+ new LeapSecond(41317, 10),
+ new LeapSecond(41499, 11),
+ new LeapSecond(41683, 12),
+ new LeapSecond(42048, 13),
+ new LeapSecond(42413, 14),
+ new LeapSecond(42778, 15),
+ new LeapSecond(43144, 16),
+ new LeapSecond(43509, 17),
+ new LeapSecond(43874, 18),
+ new LeapSecond(44239, 19),
+ new LeapSecond(44786, 20),
+ new LeapSecond(45151, 21),
+ new LeapSecond(45516, 22),
+ new LeapSecond(46247, 23),
+ new LeapSecond(47161, 24),
+ new LeapSecond(47892, 25),
+ new LeapSecond(48257, 26),
+ new LeapSecond(48804, 27),
+ new LeapSecond(49169, 28),
+ new LeapSecond(49534, 29),
+ new LeapSecond(50083, 30),
+ new LeapSecond(50630, 31),
+ new LeapSecond(51179, 32),
+ new LeapSecond(53736, 33),
+ new LeapSecond(54832, 34),
+ new LeapSecond(56109, 35),
+ };
+
+ public static double utcToTt(int utcMjd, double utcSeconds) {
+ // Returns JD in TT
+ int i0 = 0; int i1 = leapSecondsTable.length;
+ while ( i1-i0 > 1 ) {
+ int i = (i0+i1)/2;
+ if ( utcMjd >= leapSecondsTable[i].mjd )
+ i0 = i;
+ else
+ i1 = i;
+ }
+ return ( utcSeconds + leapSecondsTable[i0].offset + 32.184 )/86400. + utcMjd + 2400000.5;
+ }
+
+}
diff --git a/org/madore/ephem/VSOP87.java b/org/madore/ephem/VSOP87.java
index 3f06cf8..623401c 100644
--- a/org/madore/ephem/VSOP87.java
+++ b/org/madore/ephem/VSOP87.java
@@ -72,4 +72,8 @@ public final class VSOP87 {
return data.get(pl).get(v);
}
+ public static double getVal(Planet pl, Variable v, double tdbJd) {
+ return getFunc(pl, v).v(Ephem.fromJd(tdbJd)/10);
+ }
+
}
--
cgit v1.2.3