Skip to content

Commit 7abaa16

Browse files
cliniqueclinique
authored
[astro] Factor standard Meeus calculations (#20138)
* Factor standard meeus calculations Signed-off-by: clinique <gael@lhopital.org>
1 parent 2cf478f commit 7abaa16

7 files changed

Lines changed: 213 additions & 146 deletions

File tree

‎bundles/org.openhab.binding.astro/src/main/java/org/openhab/binding/astro/internal/calc/AstroCalc.java‎

Lines changed: 0 additions & 54 deletions
This file was deleted.

‎bundles/org.openhab.binding.astro/src/main/java/org/openhab/binding/astro/internal/calc/EclipseCalc.java‎

Lines changed: 24 additions & 25 deletions
Original file line numberDiff line numberDiff line change
@@ -19,6 +19,7 @@
1919
import java.util.Set;
2020

2121
import org.eclipse.jdt.annotation.NonNullByDefault;
22+
import org.openhab.binding.astro.internal.calc.moon.LunationArguments;
2223
import org.openhab.binding.astro.internal.model.EclipseKind;
2324
import org.openhab.binding.astro.internal.model.Position;
2425
import org.openhab.binding.astro.internal.util.DateTimeUtils;
@@ -30,7 +31,7 @@
3031
* @author Gaël L'hopital - Extracted from MoonCalc
3132
*/
3233
@NonNullByDefault
33-
public abstract class EclipseCalc extends AstroCalc {
34+
public abstract class EclipseCalc {
3435
public record Eclipse(EclipseKind kind, double when) {
3536
public LocalizedEclipse withPosition(Position position) {
3637
return new LocalizedEclipse(this, position.getElevationAsDouble());
@@ -59,9 +60,9 @@ public double calculate(double midnightJd, EclipseKind eclipse) {
5960
double tz = 0;
6061
double eclipseJd = 0;
6162
do {
62-
double k = varK(midnightJd, tz);
63+
double k = LunationArguments.varK(midnightJd, tz);
6364
tz += 1;
64-
eclipseJd = getEclipse(Math.floor(k) + getJDAjust(), eclipse);
65+
eclipseJd = getEclipse(k, getJDAjust(), eclipse);
6566
} while (eclipseJd <= midnightJd);
6667
return eclipseJd;
6768
}
@@ -75,32 +76,30 @@ public double calculate(double midnightJd, EclipseKind eclipse) {
7576
/**
7677
* Calculates the eclipse.
7778
*/
78-
protected double getEclipse(double kMod, EclipseKind eclipse) {
79-
double t = kMod / 1236.85;
80-
double f = varF(kMod, t);
79+
protected double getEclipse(double k, double adjust, EclipseKind eclipse) {
80+
LunationArguments la = new LunationArguments(k, adjust);
8181

82-
if (sinDeg(Math.abs(f)) > .36) {
82+
if (sinDeg(Math.abs(la.f)) > .36) {
8383
return 0;
8484
}
8585

86-
double o = varO(kMod, t);
87-
double f1 = f - .02665 * sinDeg(o);
88-
double a1 = 299.77 + .107408 * kMod - .009173 * t * t;
89-
double e = varE(t);
90-
double m = varM(kMod, t);
91-
double m1 = varM1(kMod, t);
92-
double p = .207 * e * sinDeg(m) + .0024 * e * sinDeg(2 * m) - .0392 * sinDeg(m1) + .0116 * sinDeg(2 * m1)
93-
- .0073 * e * sinDeg(m1 + m) + .0067 * e * sinDeg(m1 - m) + .0118 * sinDeg(2 * f1);
94-
double q = 5.2207 - .0048 * e * cosDeg(m) + .002 * e * cosDeg(2 * m) - .3299 * cosDeg(m1)
95-
- .006 * e * cosDeg(m1 + m) + .0041 * e * cosDeg(m1 - m);
86+
double f1 = la.f - .02665 * sinDeg(la.o);
87+
double a1 = 299.77 + .107408 * la.kMod - .009173 * la.t2;
88+
double p = .207 * la.e * sinDeg(la.m) + .0024 * la.e * sinDeg(2 * la.m) - .0392 * sinDeg(la.m1)
89+
+ .0116 * sinDeg(2 * la.m1) - .0073 * la.e * sinDeg(la.m1 + la.m) + .0067 * la.e * sinDeg(la.m1 - la.m)
90+
+ .0118 * sinDeg(2 * f1);
91+
double q = 5.2207 - .0048 * la.e * cosDeg(la.m) + .002 * la.e * cosDeg(2 * la.m) - .3299 * cosDeg(la.m1)
92+
- .006 * la.e * cosDeg(la.m1 + la.m) + .0041 * la.e * cosDeg(la.m1 - la.m);
9693
double g = (p * cosDeg(f1) + q * sinDeg(f1)) * (1 - .0048 * cosDeg(Math.abs(f1)));
97-
double u = .0059 + .0046 * e * cosDeg(m) - .0182 * cosDeg(m1) + .0004 * cosDeg(2 * m1) - .0005 * cosDeg(m + m1);
98-
double jd = varJde(kMod, t);
99-
jd += .0161 * sinDeg(2 * m1) - .0097 * sinDeg(2 * f1) + .0073 * e * sinDeg(m1 - m) - .005 * e * sinDeg(m1 + m)
100-
- .0023 * sinDeg(m1 - 2 * f1) + .0021 * e * sinDeg(2 * m);
101-
jd += .0012 * sinDeg(m1 + 2 * f1) + .0006 * e * sinDeg(2 * m1 + m) - .0004 * sinDeg(3 * m1)
102-
- .0003 * e * sinDeg(m + 2 * f1) + .0003 * sinDeg(a1) - .0002 * e * sinDeg(m - 2 * f1)
103-
- .0002 * e * sinDeg(2 * m1 - m) - .0002 * sinDeg(o);
104-
return astroAdjust(eclipse, e, m, m1, g, u, jd);
94+
double u = .0059 + .0046 * la.e * cosDeg(la.m) - .0182 * cosDeg(la.m1) + .0004 * cosDeg(2 * la.m1)
95+
- .0005 * cosDeg(la.m + la.m1);
96+
97+
double jd = la.jde;
98+
jd += .0161 * sinDeg(2 * la.m1) - .0097 * sinDeg(2 * f1) + .0073 * la.e * sinDeg(la.m1 - la.m)
99+
- .005 * la.e * sinDeg(la.m1 + la.m) - .0023 * sinDeg(la.m1 - 2 * f1) + .0021 * la.e * sinDeg(2 * la.m);
100+
jd += .0012 * sinDeg(la.m1 + 2 * f1) + .0006 * la.e * sinDeg(2 * la.m1 + la.m) - .0004 * sinDeg(3 * la.m1)
101+
- .0003 * la.e * sinDeg(la.m + 2 * f1) + .0003 * sinDeg(a1) - .0002 * la.e * sinDeg(la.m - 2 * f1)
102+
- .0002 * la.e * sinDeg(2 * la.m1 - la.m) - .0002 * sinDeg(la.o);
103+
return astroAdjust(eclipse, la.e, la.m, la.m1, g, u, jd);
105104
}
106105
}

‎bundles/org.openhab.binding.astro/src/main/java/org/openhab/binding/astro/internal/calc/MoonCalc.java‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -42,7 +42,7 @@
4242
* zodiac based on http://lexikon.astronomie.info/java/sunmoon/
4343
*/
4444
@NonNullByDefault
45-
public class MoonCalc extends AstroCalc {
45+
public class MoonCalc {
4646
private static final double FL = 1.0 - AstroConstants.WGS84_EARTH_FLATTENING;
4747
private static final EclipseCalc ECLIPSE_CALC = new MoonEclipseCalc();
4848

‎bundles/org.openhab.binding.astro/src/main/java/org/openhab/binding/astro/internal/calc/MoonDistanceCalc.java‎

Lines changed: 3 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -15,10 +15,9 @@
1515
import static java.lang.Math.*;
1616

1717
import org.eclipse.jdt.annotation.NonNullByDefault;
18+
import org.openhab.binding.astro.internal.calc.moon.LunarArguments;
1819
import org.openhab.binding.astro.internal.model.DistanceType;
1920
import org.openhab.binding.astro.internal.model.MoonDistance;
20-
import org.openhab.binding.astro.internal.util.AstroConstants;
21-
import org.openhab.binding.astro.internal.util.DateTimeUtils;
2221

2322
/**
2423
* Moon Distance Calculator
@@ -50,15 +49,8 @@ public class MoonDistanceCalc {
5049
* Calculates the distance from the moon to earth in metres
5150
*/
5251
public static MoonDistance calculate(double jd) {
53-
double t = DateTimeUtils.toJulianCenturies(jd);
54-
double t2 = t * t;
55-
double t3 = t2 * t;
56-
double t4 = t3 * t;
57-
double d = toRadians(297.8502042 + 445267.11151686 * t - .00163 * t2 + t3 / 545868 - t4 / 113065000);
58-
double m = toRadians(AstroConstants.E05_0 + 35999.0502909 * t - .0001536 * t2 + t3 / 24490000);
59-
double m1 = toRadians(134.9634114 + 477198.8676313 * t + .008997 * t2 + t3 / 69699 - t4 / 14712000);
60-
double f = toRadians(93.2720993 + 483202.0175273 * t - .0034029 * t2 - t3 / 3526000 + t4 / 863310000);
61-
return new MoonDistance(jd, 385000560 + getCoefficient(d, m, m1, f));
52+
LunarArguments la = new LunarArguments(jd);
53+
return new MoonDistance(jd, 385000560 + getCoefficient(la.d, la.m, la.m1, la.f));
6254
}
6355

6456
public static MoonDistance get(DistanceType type, double julianDate) {

‎bundles/org.openhab.binding.astro/src/main/java/org/openhab/binding/astro/internal/calc/MoonPhaseCalc.java‎

Lines changed: 54 additions & 55 deletions
Original file line numberDiff line numberDiff line change
@@ -20,18 +20,18 @@
2020
import java.util.stream.Collectors;
2121

2222
import org.eclipse.jdt.annotation.NonNullByDefault;
23+
import org.openhab.binding.astro.internal.calc.moon.LunarArguments;
24+
import org.openhab.binding.astro.internal.calc.moon.LunationArguments;
2325
import org.openhab.binding.astro.internal.model.MoonPhase;
2426
import org.openhab.binding.astro.internal.model.MoonPhaseSet;
25-
import org.openhab.binding.astro.internal.util.AstroConstants;
26-
import org.openhab.binding.astro.internal.util.DateTimeUtils;
2727

2828
/**
2929
* Moon Phase Calculator
3030
*
3131
* @author Gaël L'hopital - Initial contribution
3232
*/
3333
@NonNullByDefault
34-
public class MoonPhaseCalc extends AstroCalc {
34+
public class MoonPhaseCalc {
3535
public static MoonPhaseSet calculate(InstantSource instantSource, double julianDate, MoonPhaseSet previousMP,
3636
ZoneId zone) {
3737
final MoonPhaseSet result;
@@ -58,16 +58,12 @@ public static MoonPhaseSet calculate(InstantSource instantSource, double julianD
5858
* Calculates the illumination.
5959
*/
6060
private static double getIllumination(double jd) {
61-
double t = DateTimeUtils.toJulianCenturies(jd);
62-
double t2 = t * t;
63-
double t3 = t2 * t;
64-
double t4 = t3 * t;
65-
double d = 297.8502042 + 445267.11151686 * t - .00163 * t2 + t3 / 545868 - t4 / 113065000;
66-
double m = AstroConstants.E05_0 + 35999.0502909 * t - .0001536 * t2 + t3 / 24490000;
67-
double m1 = 134.9634114 + 477198.8676313 * t + .008997 * t2 + t3 / 69699 - t4 / 14712000;
68-
double i = 180 - d - 6.289 * sinDeg(m1) + 2.1 * sinDeg(m) - 1.274 * sinDeg(2 * d - m1) - .658 * sinDeg(2 * d)
69-
- .241 * sinDeg(2 * m1) - .110 * sinDeg(d);
70-
return (1 + cosDeg(i)) / 2;
61+
LunarArguments la = new LunarArguments(jd);
62+
double i = Math.PI - la.d - Math.toRadians(6.289) * Math.sin(la.m1) + Math.toRadians(2.1) * Math.sin(la.m)
63+
- Math.toRadians(1.274) * Math.sin(2 * la.d - la.m1) - Math.toRadians(0.658) * Math.sin(2 * la.d)
64+
- Math.toRadians(0.241) * Math.sin(2 * la.m1) - Math.toRadians(0.110) * Math.sin(la.d);
65+
66+
return (1 + Math.cos(i)) / 2;
7167
}
7268

7369
/**
@@ -77,7 +73,7 @@ private static double getPhase(double jd, MoonPhase phase, boolean forward) {
7773
double tz = 0;
7874
double phaseJd = 0;
7975
do {
80-
double k = varK(jd, tz);
76+
double k = LunationArguments.varK(jd, tz);
8177
tz += forward ? 1 : -1;
8278
phaseJd = calcMoonPhase(k, phase);
8379
} while (forward ? phaseJd <= jd : phaseJd > jd);
@@ -88,55 +84,58 @@ private static double getPhase(double jd, MoonPhase phase, boolean forward) {
8884
* Calculates the moon phase.
8985
*/
9086
private static double calcMoonPhase(double k, MoonPhase phase) {
91-
double kMod = Math.floor(k) + phase.cycleProgress;
92-
double t = kMod / 1236.85;
93-
double e = varE(t);
94-
double m = varM(kMod, t);
95-
double m1 = varM1(kMod, t);
96-
double f = varF(kMod, t);
97-
double o = varO(kMod, t);
98-
double jd = varJde(kMod, t);
87+
LunationArguments la = new LunationArguments(k, phase.cycleProgress);
88+
double jd = la.jde;
9989
switch (phase) {
10090
case NEW:
101-
jd += -.4072 * sinDeg(m1) + .17241 * e * sinDeg(m) + .01608 * sinDeg(2 * m1) + .01039 * sinDeg(2 * f)
102-
+ .00739 * e * sinDeg(m1 - m) - .00514 * e * sinDeg(m1 + m) + .00208 * e * e * sinDeg(2 * m)
103-
- .00111 * sinDeg(m1 - 2 * f) - .00057 * sinDeg(m1 + 2 * f);
104-
jd += .00056 * e * sinDeg(2 * m1 + m) - .00042 * sinDeg(3 * m1) + .00042 * e * sinDeg(m + 2 * f)
105-
+ .00038 * e * sinDeg(m - 2 * f) - .00024 * e * sinDeg(2 * m1 - m) - .00017 * sinDeg(o)
106-
- .00007 * sinDeg(m1 + 2 * m) + .00004 * sinDeg(2 * m1 - 2 * f);
107-
jd += .00004 * sinDeg(3 * m) + .00003 * sinDeg(m1 + m - 2 * f) + .00003 * sinDeg(2 * m1 + 2 * f)
108-
- .00003 * sinDeg(m1 + m + 2 * f) + .00003 * sinDeg(m1 - m + 2 * f)
109-
- .00002 * sinDeg(m1 - m - 2 * f) - .00002 * sinDeg(3 * m1 + m);
110-
jd += .00002 * sinDeg(4 * m1);
91+
jd += -.4072 * sinDeg(la.m1) + .17241 * la.e * sinDeg(la.m) + .01608 * sinDeg(2 * la.m1)
92+
+ .01039 * sinDeg(2 * la.f) + .00739 * la.e * sinDeg(la.m1 - la.m)
93+
- .00514 * la.e * sinDeg(la.m1 + la.m) + .00208 * la.e * la.e * sinDeg(2 * la.m)
94+
- .00111 * sinDeg(la.m1 - 2 * la.f) - .00057 * sinDeg(la.m1 + 2 * la.f);
95+
jd += .00056 * la.e * sinDeg(2 * la.m1 + la.m) - .00042 * sinDeg(3 * la.m1)
96+
+ .00042 * la.e * sinDeg(la.m + 2 * la.f) + .00038 * la.e * sinDeg(la.m - 2 * la.f)
97+
- .00024 * la.e * sinDeg(2 * la.m1 - la.m) - .00017 * sinDeg(la.o)
98+
- .00007 * sinDeg(la.m1 + 2 * la.m) + .00004 * sinDeg(2 * la.m1 - 2 * la.f);
99+
jd += .00004 * sinDeg(3 * la.m) + .00003 * sinDeg(la.m1 + la.m - 2 * la.f)
100+
+ .00003 * sinDeg(2 * la.m1 + 2 * la.f) - .00003 * sinDeg(la.m1 + la.m + 2 * la.f)
101+
+ .00003 * sinDeg(la.m1 - la.m + 2 * la.f) - .00002 * sinDeg(la.m1 - la.m - 2 * la.f)
102+
- .00002 * sinDeg(3 * la.m1 + la.m);
103+
jd += .00002 * sinDeg(4 * la.m1);
111104
break;
112105
case FULL:
113-
jd += -.40614 * sinDeg(m1) + .17302 * e * sinDeg(m) + .01614 * sinDeg(2 * m1) + .01043 * sinDeg(2 * f)
114-
+ .00734 * e * sinDeg(m1 - m) - .00515 * e * sinDeg(m1 + m) + .00209 * e * e * sinDeg(2 * m)
115-
- .00111 * sinDeg(m1 - 2 * f) - .00057 * sinDeg(m1 + 2 * f);
116-
jd += .00056 * e * sinDeg(2 * m1 + m) - .00042 * sinDeg(3 * m1) + .00042 * e * sinDeg(m + 2 * f)
117-
+ .00038 * e * sinDeg(m - 2 * f) - .00024 * e * sinDeg(2 * m1 - m) - .00017 * sinDeg(o)
118-
- .00007 * sinDeg(m1 + 2 * m) + .00004 * sinDeg(2 * m1 - 2 * f);
119-
jd += .00004 * sinDeg(3 * m) + .00003 * sinDeg(m1 + m - 2 * f) + .00003 * sinDeg(2 * m1 + 2 * f)
120-
- .00003 * sinDeg(m1 + m + 2 * f) + .00003 * sinDeg(m1 - m + 2 * f)
121-
- .00002 * sinDeg(m1 - m - 2 * f) - .00002 * sinDeg(3 * m1 + m);
122-
jd += .00002 * sinDeg(4 * m1);
106+
jd += -.40614 * sinDeg(la.m1) + .17302 * la.e * sinDeg(la.m) + .01614 * sinDeg(2 * la.m1)
107+
+ .01043 * sinDeg(2 * la.f) + .00734 * la.e * sinDeg(la.m1 - la.m)
108+
- .00515 * la.e * sinDeg(la.m1 + la.m) + .00209 * la.e * la.e * sinDeg(2 * la.m)
109+
- .00111 * sinDeg(la.m1 - 2 * la.f) - .00057 * sinDeg(la.m1 + 2 * la.f);
110+
jd += .00056 * la.e * sinDeg(2 * la.m1 + la.m) - .00042 * sinDeg(3 * la.m1)
111+
+ .00042 * la.e * sinDeg(la.m + 2 * la.f) + .00038 * la.e * sinDeg(la.m - 2 * la.f)
112+
- .00024 * la.e * sinDeg(2 * la.m1 - la.m) - .00017 * sinDeg(la.o)
113+
- .00007 * sinDeg(la.m1 + 2 * la.m) + .00004 * sinDeg(2 * la.m1 - 2 * la.f);
114+
jd += .00004 * sinDeg(3 * la.m) + .00003 * sinDeg(la.m1 + la.m - 2 * la.f)
115+
+ .00003 * sinDeg(2 * la.m1 + 2 * la.f) - .00003 * sinDeg(la.m1 + la.m + 2 * la.f)
116+
+ .00003 * sinDeg(la.m1 - la.m + 2 * la.f) - .00002 * sinDeg(la.m1 - la.m - 2 * la.f)
117+
- .00002 * sinDeg(3 * la.m1 + la.m);
118+
jd += .00002 * sinDeg(4 * la.m1);
123119
break;
124120
default:
125-
jd += -.62801 * sinDeg(m1) + .17172 * e * sinDeg(m) - .01183 * e * sinDeg(m1 + m)
126-
+ .00862 * sinDeg(2 * m1) + .00804 * sinDeg(2 * f) + .00454 * e * sinDeg(m1 - m)
127-
+ .00204 * e * e * sinDeg(2 * m) - .0018 * sinDeg(m1 - 2 * f) - .0007 * sinDeg(m1 + 2 * f);
128-
jd += -.0004 * sinDeg(3 * m1) - .00034 * e * sinDeg(2 * m1 - m) + .00032 * e * sinDeg(m + 2 * f)
129-
+ .00032 * e * sinDeg(m - 2 * f) - .00028 * e * e * sinDeg(m1 + 2 * m)
130-
+ .00027 * e * sinDeg(2 * m1 + m) - .00017 * sinDeg(o);
131-
jd += -.00005 * sinDeg(m1 - m - 2 * f) + .00004 * sinDeg(2 * m1 + 2 * f)
132-
- .00004 * sinDeg(m1 + m + 2 * f) + .00004 * sinDeg(m1 - 2 * m)
133-
+ .00003 * sinDeg(m1 + m - 2 * f) + .00003 * sinDeg(3 * m) + .00002 * sinDeg(2 * m1 - 2 * f);
134-
jd += .00002 * sinDeg(m1 - m + 2 * f) - .00002 * sinDeg(3 * m1 + m);
135-
double w = .00306 - .00038 * e * cosDeg(m) + .00026 * cosDeg(m1) - .00002 * cosDeg(m1 - m)
136-
+ .00002 * cosDeg(m1 + m) + .00002 * cosDeg(2 * f);
121+
jd += -.62801 * sinDeg(la.m1) + .17172 * la.e * sinDeg(la.m) - .01183 * la.e * sinDeg(la.m1 + la.m)
122+
+ .00862 * sinDeg(2 * la.m1) + .00804 * sinDeg(2 * la.f) + .00454 * la.e * sinDeg(la.m1 - la.m)
123+
+ .00204 * la.e * la.e * sinDeg(2 * la.m) - .0018 * sinDeg(la.m1 - 2 * la.f)
124+
- .0007 * sinDeg(la.m1 + 2 * la.f);
125+
jd += -.0004 * sinDeg(3 * la.m1) - .00034 * la.e * sinDeg(2 * la.m1 - la.m)
126+
+ .00032 * la.e * sinDeg(la.m + 2 * la.f) + .00032 * la.e * sinDeg(la.m - 2 * la.f)
127+
- .00028 * la.e * la.e * sinDeg(la.m1 + 2 * la.m) + .00027 * la.e * sinDeg(2 * la.m1 + la.m)
128+
- .00017 * sinDeg(la.o);
129+
jd += -.00005 * sinDeg(la.m1 - la.m - 2 * la.f) + .00004 * sinDeg(2 * la.m1 + 2 * la.f)
130+
- .00004 * sinDeg(la.m1 + la.m + 2 * la.f) + .00004 * sinDeg(la.m1 - 2 * la.m)
131+
+ .00003 * sinDeg(la.m1 + la.m - 2 * la.f) + .00003 * sinDeg(3 * la.m)
132+
+ .00002 * sinDeg(2 * la.m1 - 2 * la.f);
133+
jd += .00002 * sinDeg(la.m1 - la.m + 2 * la.f) - .00002 * sinDeg(3 * la.m1 + la.m);
134+
double w = .00306 - .00038 * la.e * cosDeg(la.m) + .00026 * cosDeg(la.m1)
135+
- .00002 * cosDeg(la.m1 - la.m) + .00002 * cosDeg(la.m1 + la.m) + .00002 * cosDeg(2 * la.f);
137136
jd += MoonPhase.FIRST_QUARTER.equals(phase) ? w : -w;
138137
}
139-
return moonCorrection(jd, t, kMod);
138+
return moonCorrection(jd, la.t, la.kMod);
140139
}
141140

142141
private static double moonCorrection(double jd, double t, double k) {

0 commit comments

Comments
 (0)
⚡