/*
|
Zeitrechnung - a class library to determine calendar events
|
Copyright (c) 1984-2023 Ulrich Hilger, http://uhilger.de
|
|
This program is free software: you can redistribute it and/or modify
|
it under the terms of the GNU Affero General Public License as published by
|
the Free Software Foundation, either version 3 of the License, or
|
(at your option) any later version.
|
|
This program is distributed in the hope that it will be useful,
|
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
GNU Affero General Public License for more details.
|
|
You should have received a copy of the GNU Affero General Public License
|
along with this program. If not, see <http://www.gnu.org/licenses/>.
|
*/
|
package de.uhilger.zeitrechnung.kalender;
|
|
import de.uhilger.zeitrechnung.Definition;
|
import de.uhilger.zeitrechnung.Ort;
|
import de.uhilger.zeitrechnung.Zeit;
|
|
/**
|
* Abstrakte Basisklasse fuer Klassen, die ein Kalendersystem implementieren.
|
*
|
* Hier sind neben allerlei relevanten Rechenmethoden die grundlegenden
|
* astronomischen Algorithmen für die Zeit- und Kalenderrechnung implementiert.
|
*
|
* @author Ulrich Hilger
|
*/
|
public abstract class BasisKalender implements Zeitrechnung {
|
|
/* Implementierung der Schnittstelle Zeitrechnung */
|
|
@Override
|
public long ganzzahlQuotient(double x, double y) {
|
return (long) Math.floor(x / y);
|
}
|
|
@Override
|
public long modulo(long x, long y) {
|
return (long) (x - y * Math.floor((double) x / y));
|
}
|
|
public double modulo(double x, double y) {
|
return x - y * Math.floor(x / y);
|
}
|
|
@Override
|
public long tagNach(long datum, int t) {
|
return tagAmOderNach(datum + 7, t);
|
}
|
|
@Override
|
public long tagAmOderNach(long datum, int t) {
|
return datum - wochentagVonGenerisch(datum - t);
|
}
|
|
@Override
|
public long wochentagVonGenerisch(long datum) {
|
return modulo(datum, 7);
|
}
|
|
@Override
|
public long nterTag(int n, int t, long datum) {
|
return n > 0 ?
|
tagVor(datum, t) + 7 * n :
|
tagNach(datum, t) + 7 * n;
|
}
|
|
@Override
|
public long tagVor(long datum, int t) {
|
return tagAmOderVor(datum - 1, t);
|
}
|
|
@Override
|
public long tagAmOderVor(long datum, int t) {
|
return datum - wochentagVonGenerisch(datum - t);
|
}
|
|
@Override
|
public long letzterTag(int t, long datum) {
|
return nterTag(-1, t, datum);
|
}
|
|
/* ----- */
|
|
public double moduloAngepasst(double x, double y) {
|
return y + modulo(x, -y);
|
}
|
|
/* ---- Zeit ----- */
|
|
public double zuMoment(int stunde, int minute, double sekunde) {
|
return stunde / 24d + minute / (24d * 60) + sekunde / (24d * 60 * 60);
|
}
|
|
public double zuMoment(Zeit z) {
|
return BasisKalender.this.zuMoment(z.getStunde(), z.getMinute(), z.getSekunde());
|
}
|
|
public Zeit vonMoment(double t) {
|
Zeit z = new Zeit();
|
z.setStunde((int)Math.floor(modulo(t * 24, 24)));
|
z.setMinute((int)Math.floor(modulo(t * 24 * 60, 60)));
|
z.setSekunde(modulo(t * 24 * 60 * 60, 60));
|
return z;
|
}
|
|
/* ----------- Mond ----------- */
|
|
/** durchschnittliche Dauer eines Mondphasenzyklus (synodischer Monat) in Tagen */
|
public static final double MITTLERER_SYNODISCHER_MONAT = 29.530588853;
|
|
public double mondphase(double t) {
|
return modulo(mondLaenge(t) - solareLaenge(t), 360);
|
}
|
|
public double mondHoehe(double t, Ort ort) {
|
double phi = ort.getBreite();
|
double psi = ort.getLaenge();
|
double varepsilon = schiefstand(t);
|
double lambda = mondLaenge(t);
|
double beta = mondBreite(t);
|
double alpha = arcTanGrad(
|
(sinGrad(lambda) * kosGrad(varepsilon) - tanGrad(beta) * sinGrad(varepsilon)) /
|
kosGrad(lambda),
|
(int)ganzzahlQuotient(lambda, (double) (90)) + 1
|
);
|
double delta = arcSinGrad(sinGrad(beta) * kosGrad(varepsilon) +
|
kosGrad(beta) * sinGrad(varepsilon) * sinGrad(lambda));
|
double theta0 = siderischVonMoment(t);
|
double capH = modulo(theta0 + psi - alpha, 360);
|
double altitude = arcSinGrad(sinGrad(phi) * sinGrad(delta) + kosGrad(phi) * kosGrad(delta) * kosGrad(capH));
|
return modulo(altitude + (double) (180), 360) - (double) (180);
|
}
|
|
public double mondLaenge(double t) {
|
double c = julJahrhunderte(t);
|
double meanMoon = grad(poly(c, llon.coeffMeanMoon));
|
double elongation = grad(poly(c, llon.coeffElongation));
|
double solarAnomaly = grad(poly(c, llon.coeffSolarAnomaly));
|
double lunarAnomaly = grad(poly(c, llon.coeffLunarAnomaly));
|
double moonNode = grad(poly(c, llon.coeffMoonNode));
|
double capE = poly(c, llon.coeffCapE);
|
double sigma = 0;
|
for(int i = 0; i < llon.argsLunarElongation.length; ++i) {
|
double x = llon.argsSolarAnomaly[i];
|
sigma += llon.sineCoefficients[i] *
|
Math.pow(capE, Math.abs(x)) *
|
sinGrad( llon.argsLunarElongation[i] * elongation +
|
x * solarAnomaly +
|
llon.argsLunarAnomaly[i] * lunarAnomaly +
|
llon.argsMoonFromNode[i] * moonNode);
|
}
|
double correction = ((double) (1) / 1000000) * sigma;
|
double venus = ((double) (3958) / 1000000) * sinGrad(119.75 + c * 131.849);
|
double jupiter = ((double) (318) / 1000000) * sinGrad(53.09 + c * 479264.29);
|
double flatEarth = ((double) (1962) / 1000000) * sinGrad(meanMoon - moonNode);
|
return modulo(meanMoon + correction + venus + jupiter + flatEarth + nutation(t), 360);
|
}
|
private static class llon {
|
private static final double[] coeffMeanMoon = new double[] {218.3164591, 481267.88134236, -0.0013268, 1d/538841, -1d/65194000};
|
private static final double[] coeffElongation = new double[] {297.8502042, 445267.1115168, -0.00163, 1d/545868, -1d/113065000};
|
private static final double[] coeffSolarAnomaly = new double[] {357.5291092, 35999.0502909, -0.0001536, 1d/24490000};
|
private static final double[] coeffLunarAnomaly = new double[] {134.9634114, 477198.8676313, 0.008997, 1d/69699, -1d/14712000};
|
private static final double[] coeffMoonNode = new double[] {93.2720993, 483202.0175273, -0.0034029, -1d/3526000, 1d/863310000};
|
private static final double[] coeffCapE = new double[] {1, -0.002516, -0.0000074};
|
private static final byte[] argsLunarElongation = new byte[] {
|
0, 2, 2, 0, 0, 0, 2, 2, 2, 2, 0, 1, 0, 2, 0, 0, 4, 0, 4, 2, 2, 1,
|
1, 2, 2, 4, 2, 0, 2, 2, 1, 2, 0, 0, 2, 2, 2, 4, 0, 3, 2, 4, 0, 2,
|
2, 2, 4, 0, 4, 1, 2, 0, 1, 3, 4, 2, 0, 1, 2
|
};
|
private static final byte[] argsSolarAnomaly = new byte[] {
|
0, 0, 0, 0, 1, 0, 0, -1, 0, -1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1,
|
0, 1, -1, 0, 0, 0, 1, 0, -1, 0, -2, 1, 2, -2, 0, 0, -1, 0, 0, 1,
|
-1, 2, 2, 1, -1, 0, 0, -1, 0, 1, 0, 1, 0, 0, -1, 2, 1, 0
|
};
|
private static final byte[] argsLunarAnomaly = new byte[] {
|
1, -1, 0, 2, 0, 0, -2, -1, 1, 0, -1, 0, 1, 0, 1, 1, -1, 3, -2,
|
-1, 0, -1, 0, 1, 2, 0, -3, -2, -1, -2, 1, 0, 2, 0, -1, 1, 0,
|
-1, 2, -1, 1, -2, -1, -1, -2, 0, 1, 4, 0, -2, 0, 2, 1, -2, -3,
|
2, 1, -1, 3
|
};
|
private static final byte[] argsMoonFromNode = new byte[] {
|
0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, -2, 2, -2, 0, 0, 0, 0, 0,
|
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, -2, 2, 0, 2, 0, 0, 0, 0,
|
0, 0, -2, 0, 0, 0, 0, -2, -2, 0, 0, 0, 0, 0, 0, 0
|
};
|
private static final int[] sineCoefficients = new int[] {
|
6288774, 1274027, 658314, 213618, -185116, -114332,
|
58793, 57066, 53322, 45758, -40923, -34720, -30383,
|
15327, -12528, 10980, 10675, 10034, 8548, -7888,
|
-6766, -5163, 4987, 4036, 3994, 3861, 3665, -2689,
|
-2602, 2390, -2348, 2236, -2120, -2069, 2048, -1773,
|
-1595, 1215, -1110, -892, -810, 759, -713, -700, 691,
|
596, 549, 537, 520, -487, -399, -381, 351, -340, 330,
|
327, -323, 299, 294
|
};
|
}
|
|
// aus "Astronomical Algorithms" von Jean Meeus,
|
// Willmann-Bell, Inc., 1998.
|
|
public double mondBreite(double t) {
|
double c = julJahrhunderte(t);
|
double longitude = grad(poly(c, llat.coeffLongitude));
|
double elongation = grad(poly(c, llat.coeffElongation));
|
double solarAnomaly = grad(poly(c, llat.coeffSolarAnomaly));
|
double lunarAnomaly = grad(poly(c, llat.coeffLunarAnomaly));
|
double moonNode = grad(poly(c, llat.coeffMoonNode));
|
double capE = poly(c, llat.coeffCapE);
|
double latitude = 0;
|
for(int i = 0; i < llat.argsLunarElongation.length; ++i) {
|
double x = llat.argsSolarAnomaly[i];
|
latitude += llat.sineCoefficients[i] *
|
Math.pow(capE, Math.abs(x)) *
|
sinGrad( llat.argsLunarElongation[i] * elongation +
|
x * solarAnomaly +
|
llat.argsLunarAnomaly[i] * lunarAnomaly +
|
llat.argsMoonNode[i] * moonNode);
|
}
|
latitude *= (double) (1) / 1000000;
|
double venus = ((double) (175) / 1000000) * (sinGrad((double) (119.75) + c * 131.849 + moonNode) + sinGrad((double) (119.75) + c * 131.849 - moonNode));
|
double flatEarth = ((double) (-2235) / 1000000) * sinGrad(longitude) +
|
((double) (127) / 1000000) * sinGrad(longitude - lunarAnomaly) +
|
((double) (-115) / 1000000) * sinGrad(longitude + lunarAnomaly);
|
double extra = ((double) (382) / 1000000) * sinGrad((double) (313.45) + c * (double) (481266.484));
|
return modulo(latitude + venus + flatEarth + extra, 360);
|
}
|
private static class llat {
|
private static final double[] coeffLongitude = new double[] {218.3164591, 481267.88134236, -0.0013268, 1d/538841, -1d/65194000};
|
private static final double[] coeffElongation = new double[] {297.8502042, 445267.1115168, -0.00163, 1d/545868, -1d/113065000};
|
private static final double[] coeffSolarAnomaly = new double[] {357.5291092, 35999.0502909, -0.0001536, 1d/24490000};
|
private static final double[] coeffLunarAnomaly = new double[] {134.9634114, 477198.8676313, 0.008997, 1d/69699, -1d/14712000};
|
private static final double[] coeffMoonNode = new double[] {93.2720993, 483202.0175273, -0.0034029, -1d/3526000, 1d/863310000};
|
private static final double[] coeffCapE = new double[] {1, -0.002516, -0.0000074};
|
private static final byte[] argsLunarElongation = new byte[] {
|
0, 0, 0, 2, 2, 2, 2, 0, 2, 0, 2, 2, 2, 2, 2, 2, 2, 0, 4, 0, 0, 0,
|
1, 0, 0, 0, 1, 0, 4, 4, 0, 4, 2, 2, 2, 2, 0, 2, 2, 2, 2, 4, 2, 2,
|
0, 2, 1, 1, 0, 2, 1, 2, 0, 4, 4, 1, 4, 1, 4, 2};
|
private static final byte[] argsSolarAnomaly = new byte[] {
|
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 1, -1, -1, -1, 1, 0, 1,
|
0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 1, 1,
|
0, -1, -2, 0, 1, 1, 1, 1, 1, 0, -1, 1, 0, -1, 0, 0, 0, -1, -2};
|
private static final byte[] argsLunarAnomaly = new byte[] {
|
0, 1, 1, 0, -1, -1, 0, 2, 1, 2, 0, -2, 1, 0, -1, 0, -1, -1, -1,
|
0, 0, -1, 0, 1, 1, 0, 0, 3, 0, -1, 1, -2, 0, 2, 1, -2, 3, 2, -3,
|
-1, 0, 0, 1, 0, 1, 1, 0, 0, -2, -1, 1, -2, 2, -2, -1, 1, 1, -2,
|
0, 0};
|
private static final byte[] argsMoonNode = new byte[] {
|
1, 1, -1, -1, 1, -1, 1, 1, -1, -1, -1, -1, 1, -1, 1, 1, -1, -1,
|
-1, 1, 3, 1, 1, 1, -1, -1, -1, 1, -1, 1, -3, 1, -3, -1, -1, 1,
|
-1, 1, -1, 1, 1, 1, 1, -1, 3, -1, -1, 1, -1, -1, 1, -1, 1, -1,
|
-1, -1, -1, -1, -1, 1};
|
private static final int[] sineCoefficients = new int[] {
|
5128122, 280602, 277693, 173237, 55413, 46271, 32573,
|
17198, 9266, 8822, 8216, 4324, 4200, -3359, 2463, 2211,
|
2065, -1870, 1828, -1794, -1749, -1565, -1491, -1475,
|
-1410, -1344, -1335, 1107, 1021, 833, 777, 671, 607,
|
596, 491, -451, 439, 422, 421, -366, -351, 331, 315,
|
302, -283, -229, 223, 223, -220, -220, -185, 181,
|
-177, 176, 166, -164, 132, -119, 115, 107};
|
}
|
|
public double arcTanGrad(double x, int quad) {
|
double alpha = bogenmassZuGrad(Math.atan(x));
|
return modulo(quad == 1 || quad == 4 ? alpha : alpha + (double) (180), 360);
|
}
|
|
public double siderischVonMoment(double t) {
|
double c = (t - j2000()) / 36525;
|
return modulo(poly(c, sfm.siderealCoeff), 360);
|
}
|
private static class sfm {
|
private static final double[] siderealCoeff = new double[] {280.46061837, 36525 * 360.98564736629, 0.000387933, 1d/38710000};
|
}
|
|
public double neumondNach(double tee) {
|
return nterNeumond(1 + Math.round(tee / MITTLERER_SYNODISCHER_MONAT - mondphase(tee) / (double) (360)));
|
}
|
|
public double neumondVor(double tee) {
|
return nterNeumond(Math.round(tee / MITTLERER_SYNODISCHER_MONAT - mondphase(tee) / (double) (360)));
|
}
|
|
public double nterNeumond(long n) {
|
double k = n - 24724;
|
double c = k / 1236.85;
|
double approx = poly(c, nm.coeffApprox);
|
double capE = poly(c, nm.coeffCapE);
|
double solarAnomaly = poly(c, nm.coeffSolarAnomaly);
|
double lunarAnomaly = poly(c, nm.coeffLunarAnomaly);
|
double moonArgument = poly(c, nm.coeffMoonArgument);
|
double capOmega = poly(c, nm.coeffCapOmega);
|
double correction = -0.00017 * sinGrad(capOmega);
|
for(int ix = 0; ix < nm.sineCoeff.length; ++ix) {
|
correction += nm.sineCoeff[ix] * Math.pow(capE, nm.EFactor[ix]) *
|
sinGrad(nm.solarCoeff[ix] * solarAnomaly +
|
nm.lunarCoeff[ix] * lunarAnomaly +
|
nm.moonCoeff[ix] * moonArgument);
|
}
|
double additional = 0;
|
for(int ix = 0; ix < nm.addConst.length; ++ix) {
|
additional += nm.addFactor[ix] *
|
sinGrad(nm.addConst[ix] + nm.addCoeff[ix] * k);
|
}
|
double extra = 0.000325 * sinGrad(poly(c, nm.extra));
|
return universalVonDynamisch(approx + correction + extra + additional);
|
}
|
private static class nm {
|
private static final double[] coeffApprox = new double[] {730125.59765, MITTLERER_SYNODISCHER_MONAT * 1236.85, 0.0001337, -0.000000150, 0.00000000073};
|
private static final double[] coeffCapE = new double[] {1, -0.002516, -0.0000074};
|
private static final double[] coeffSolarAnomaly = new double[] {2.5534, 29.10535669 * 1236.85, -0.0000218, -0.00000011};
|
private static final double[] coeffLunarAnomaly = new double[] {201.5643, 385.81693528 * 1236.85, 0.0107438, 0.00001239, -0.000000058};
|
private static final double[] coeffMoonArgument = new double[] {160.7108, 390.67050274 * 1236.85, -0.0016341, -0.00000227, 0.000000011};
|
private static final double[] coeffCapOmega = new double[] {124.7746, -1.56375580 * 1236.85, 0.0020691, 0.00000215};
|
private static final byte[] EFactor = new byte[] {0, 1, 0, 0, 1, 1, 2, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
|
private static final byte[] solarCoeff = new byte[] {0, 1, 0, 0, -1, 1, 2, 0, 0, 1, 0, 1, 1, -1, 2, 0, 3, 1, 0, 1, -1, -1, 1, 0};
|
private static final byte[] lunarCoeff = new byte[] {1, 0, 2, 0, 1, 1, 0, 1, 1, 2, 3, 0, 0, 2, 1, 2, 0, 1, 2, 1, 1, 1, 3, 4};
|
private static final byte[] moonCoeff = new byte[] {0, 0, 0, 2, 0, 0, 0, -2, 2, 0, 0, 2, -2, 0, 0, -2, 0, -2, 2, 2, 2, -2, 0, 0};
|
private static final double[] sineCoeff = new double[] {
|
-0.40720, 0.17241, 0.01608, 0.01039, 0.00739, -0.00514, 0.00208,
|
-0.00111, -0.00057, 0.00056, -0.00042, 0.00042, 0.00038, -0.00024,
|
-0.00007, 0.00004, 0.00004, 0.00003, 0.00003, -0.00003, 0.00003,
|
-0.00002, -0.00002, 0.00002
|
};
|
private static final double[] addConst = new double[] {
|
251.88, 251.83, 349.42, 84.66, 141.74, 207.14, 154.84, 34.52, 207.19,
|
291.34, 161.72, 239.56, 331.55
|
};
|
private static final double[] addCoeff = new double[] {
|
0.016321, 26.641886, 36.412478, 18.206239, 53.303771, 2.453732,
|
7.306860, 27.261239, 0.121824, 1.844379, 24.198154, 25.513099, 3.592518
|
};
|
private static final double[] addFactor = new double[] {
|
0.000165, 0.000164, 0.000126, 0.000110, 0.000062, 0.000060, 0.000056,
|
0.000047, 0.000042, 0.000040, 0.000037, 0.000035, 0.000023
|
};
|
private static final double[] extra = new double[] {
|
299.77, 132.8475848, -0.009173
|
};
|
}
|
|
public double universalVonDynamisch(double tee) {
|
return tee - ephemeridenKorrektur(tee);
|
}
|
|
public double universalVonStandard(double teeS, Ort locale) {
|
return teeS - locale.getZeitzone() / 24;
|
}
|
|
/* ----------- Sonnenauf- und -untergang ----------- */
|
|
public double sonnenaufgang(long date, Ort ort) throws Exception {
|
return morgen(date, ort, alpha(ort));
|
}
|
|
public double sonnenuntergang(long date, Ort ort)
|
throws Exception
|
{
|
return abend(date, ort, alpha(ort));
|
}
|
|
public double alpha(Ort ort) {
|
double h = Math.max(0, ort.getHoehe());
|
final double capR = (double) 6.372E+6;
|
double dip = arcKosGrad(capR / (capR + h));
|
return winkel(0, 50, 0) + dip;
|
}
|
|
public double morgen(long date, Ort ort, double alpha) throws Exception {
|
double approx;
|
try {
|
approx = zeitVonHorizont(date + .25, ort, alpha);
|
} catch(Exception ex) {
|
approx = date;
|
}
|
double result = zeitVonHorizont(approx, ort, alpha);
|
return standardVonLokal(result, ort);
|
}
|
|
public double abend(long date, Ort ort, double alpha)
|
throws Exception
|
{
|
double approx;
|
try {
|
approx = zeitVonHorizont(date + .75, ort, alpha);
|
} catch(Exception ex) {
|
approx = date + .99d;
|
}
|
double result = zeitVonHorizont(approx, ort, alpha);
|
return standardVonLokal(result, ort);
|
}
|
|
public double zeitVonHorizont(double approx, Ort ort, double alpha) throws Exception {
|
double phi = ort.getBreite();
|
double t = universalVonLokal(approx, ort);
|
double delta = arcSinGrad(sinGrad(schiefstand(t)) * sinGrad(solareLaenge(t)));
|
boolean morgen = modulo(approx, 1) < 0.5;
|
double sinusAbstand = tanGrad(phi) * tanGrad(delta) +
|
sinGrad(alpha) / (kosGrad(delta) * kosGrad(phi));
|
double offset = modulo(0.5 + arcSinGrad(sinusAbstand) / (double) 360, 1) - 0.5;
|
if(Math.abs(sinusAbstand) > 1) {
|
throw new Exception();
|
}
|
return lokalVonScheinbar(Math.floor(approx) + (morgen ? .25 - offset : .75 + offset));
|
}
|
|
public double universalVonLokal(double tl, Ort ort) {
|
return tl - ort.getLaenge() / (double) 360;
|
}
|
|
public double standardVonLokal(double tl, Ort ort) {
|
return standardVonUniversal(universalVonLokal(tl, ort), ort);
|
}
|
|
public double lokalVonScheinbar(double t) {
|
return t - zeitgleichung(t);
|
}
|
|
public double zeitgleichung(double t) {
|
double c = julJahrhunderte(t);
|
double laenge = poly(c, et.koeffLaenge);
|
double anomalie = poly(c, et.koeffAnomalie);
|
double exzentrizitaet = poly(c, et.koeffExzentrizitaet);
|
double varepsilon = schiefstand(t);
|
double y = quadrat(tanGrad(varepsilon / 2));
|
double equation = (1d / 2d / Math.PI) * (y * sinGrad(2 * laenge) +
|
-2 * exzentrizitaet * sinGrad(anomalie) +
|
4 * exzentrizitaet * y * sinGrad(anomalie) * kosGrad(2 * laenge) +
|
-0.5 * y * y * sinGrad(4 * laenge) +
|
-1.25 * exzentrizitaet * exzentrizitaet * sinGrad(2 * anomalie));
|
return vorzeichen(equation) * Math.min(Math.abs(equation), stunde(12));
|
}
|
private static class et {
|
private static final double[] koeffLaenge = new double[] {280.46645, 36000.76983, 0.0003032};
|
private static final double[] koeffAnomalie = new double[] {357.52910, 35999.05030, -0.0001559, -0.00000048};
|
private static final double[] koeffExzentrizitaet = new double[] {0.016708617, -0.000042037, -0.0000001236};
|
}
|
|
public double schiefstand(double t) {
|
double c = julJahrhunderte(t);
|
return winkel(23, 26, 21.448) + poly(c, coeffObliquity);
|
}
|
private final double[] coeffObliquity = new double[] {0, winkel(0, 0, -46.8150), winkel(0, 0, -0.00059), winkel(0, 0, 0.001813)};
|
|
public int vorzeichen(double x) {
|
if(x < 0)
|
return -1;
|
else if(x > 0)
|
return 1;
|
else
|
return 0;
|
}
|
|
public double quadrat(double x) {
|
return x * x;
|
}
|
|
public double kosGrad(double theta) {
|
return Math.cos(gradZuBogenmass(theta));
|
}
|
|
public double arcSinGrad(double x) {
|
return bogenmassZuGrad(Math.asin(x));
|
}
|
|
public double tanGrad(double theta) {
|
return Math.tan(gradZuBogenmass(theta));
|
}
|
|
public double arcKosGrad(double x) {
|
return bogenmassZuGrad(Math.acos(x));
|
}
|
|
public double bogenmassZuGrad(double theta) {
|
return grad(theta / Math.PI * 180);
|
}
|
|
public double winkel(double d, double m, double s) {
|
return d + (m + s / 60) / 60;
|
}
|
|
/* ---------------- Jahreszeiten ----- */
|
|
/** durchschnittliche Dauer eines Umlaufs der Erde um die Sonne in Tagen */
|
public static final double MITTLERES_TROPISCHES_JAHR = 365.242189;
|
|
public double standardVonUniversal(double t, Ort ort) {
|
return t + ort.getZeitzone() / 24;
|
}
|
|
public double solareLaengeNach(double t, double phi) {
|
double varepsilon = 0.00001;
|
double rate = MITTLERES_TROPISCHES_JAHR / (double) 360;
|
double tau = t + rate * modulo(phi - solareLaenge(t), 360);
|
double l = Math.max(t, tau - 5);
|
double u = tau + 5;
|
|
double lo = l, hi = u, x = (hi + lo) / 2;
|
while(hi - lo > varepsilon) {
|
if(modulo(solareLaenge(x) - phi, 360) < (double) 180)
|
hi = x;
|
else
|
lo = x;
|
|
x = (hi + lo) / 2;
|
}
|
return x;
|
}
|
|
public double solareLaenge(double t) {
|
double c = julJahrhunderte(t);
|
double sigma = 0;
|
for(int i = 0; i < sLaenge.koeffizienten.length; ++i) {
|
sigma += sLaenge.koeffizienten[i] * sinGrad(sLaenge.multiplikatoren[i] * c + sLaenge.summanden[i]);
|
}
|
double laenge = (double) 282.7771834 +
|
36000.76953744 * c +
|
0.000005729577951308232 * sigma;
|
return modulo(laenge + aberration(t) + nutation(t), 360);
|
}
|
|
public double geschaetzteSolareLaengeVor(double tee, double phi) {
|
double rate = MITTLERES_TROPISCHES_JAHR / (double) (360);
|
double tau = tee - rate * modulo(solareLaenge(tee) - phi, 360);
|
double capDelta = modulo(solareLaenge(tau) - phi + (double) (180), 360) - (double) (180);
|
return Math.min(tee, tau - rate * capDelta);
|
}
|
|
public double julJahrhunderte(double t) {
|
return (dynamischVonUniversal(t) - j2000()) / 36525;
|
}
|
|
public double dynamischVonUniversal(double tee) {
|
return tee + ephemeridenKorrektur(tee);
|
}
|
|
public double ephemeridenKorrektur(double t) {
|
double[] koeffizient17tes = new double[] {196.58333, -4.0675, 0.0219167};
|
double[] koeffizient19tes = new double[] {-0.00002, 0.000297, 0.025184, -0.181133, 0.553040, -0.861938, 0.677066, -0.212591};
|
double[] koeffizient18tes = new double[] {-0.000009, 0.003844, 0.083563, 0.865736, 4.867575, 15.845535, 31.332267, 38.291999, 28.316289, 11.636204, 2.043794};
|
ISOKalender w = new ISOKalender();
|
long jahr = w.jahrVonTagen((long)Math.floor(t));
|
double c = differenz(w.zuTagen(1900, Definition.JANUAR, 1),
|
w.zuTagen(jahr, Definition.JULI, 1)) / 36525d;
|
double result;
|
if(1988 <= jahr && jahr <= 2019) {
|
result = (jahr - 1933) / (24d * 60 * 60);
|
} else if (1900 <= jahr && jahr <= 1987) {
|
result = poly(c, koeffizient19tes);
|
} else if (1800 <= jahr && jahr <= 1899) {
|
result = poly(c, koeffizient18tes);
|
} else if (1620 <= jahr && jahr <= 1799) {
|
result = poly(jahr - 1600, koeffizient17tes) / (24 * 60 * 60);
|
} else {
|
double x = stunde(12) + differenz(w.zuTagen(1810, Definition.JANUAR, 1),
|
w.zuTagen(jahr, Definition.JANUAR, 1));
|
return (x * x / 41048480 - 15) / (24 * 60 * 60);
|
}
|
return result;
|
}
|
|
public double nutation(double t) {
|
double[] koeffa = new double[] {124.90, -1934.134, 0.002063};
|
double[] koeffb = new double[] {201.11, 72001.5377, 0.00057};
|
double c = julJahrhunderte(t);
|
double capA = poly(c, koeffa);
|
double capB = poly(c, koeffb);
|
return (double) -0.004778 * sinGrad(capA) +
|
(double) -0.0003667 * sinGrad(capB);
|
}
|
|
public static long differenz(long datum1, long datum2) {
|
return datum2 - datum1;
|
}
|
|
public double stunde(double x) {
|
return x / 24;
|
}
|
|
public double j2000() {
|
ISOKalender w = new ISOKalender();
|
return stunde(12) + w.zuTagen(2000, Definition.JANUAR, 1);
|
}
|
|
public double sinGrad(double theta) {
|
return Math.sin(gradZuBogenmass(theta));
|
}
|
|
public double gradZuBogenmass(double theta) {
|
return grad(theta) * Math.PI / 180;
|
}
|
|
public double grad(double theta) {
|
return modulo(theta, 360);
|
}
|
|
public double aberration(double t) {
|
double c = julJahrhunderte(t);
|
return (double) 0.0000974 * gradKosinus((double) 177.63 + (double) 35999.01848 * c) - (double) 0.005575;
|
}
|
|
public double gradKosinus(double theta) {
|
return Math.cos(gradZuBogenmass(theta));
|
}
|
|
public double poly(double x, double[] a) {
|
double ergebnis = a[0];
|
for(int i = 1; i < a.length; ++i) {
|
ergebnis += a[i] * Math.pow(x, i);
|
}
|
return ergebnis;
|
}
|
|
private static class sLaenge {
|
private static final int[] koeffizienten = new int[] {
|
403406, 195207, 119433, 112392, 3891, 2819, 1721,
|
660, 350, 334, 314, 268, 242, 234, 158, 132, 129, 114,
|
99, 93, 86, 78, 72, 68, 64, 46, 38, 37, 32, 29, 28, 27, 27,
|
25, 24, 21, 21, 20, 18, 17, 14, 13, 13, 13, 12, 10, 10, 10,
|
10
|
};
|
private static final double[] multiplikatoren = new double[] {
|
0.9287892, 35999.1376958, 35999.4089666,
|
35998.7287385, 71998.20261, 71998.4403,
|
36000.35726, 71997.4812, 32964.4678,
|
-19.4410, 445267.1117, 45036.8840, 3.1008,
|
22518.4434, -19.9739, 65928.9345,
|
9038.0293, 3034.7684, 33718.148, 3034.448,
|
-2280.773, 29929.992, 31556.493, 149.588,
|
9037.750, 107997.405, -4444.176, 151.771,
|
67555.316, 31556.080, -4561.540,
|
107996.706, 1221.655, 62894.167,
|
31437.369, 14578.298, -31931.757,
|
34777.243, 1221.999, 62894.511,
|
-4442.039, 107997.909, 119.066, 16859.071,
|
-4.578, 26895.292, -39.127, 12297.536,
|
90073.778
|
};
|
private static final double[] summanden = new double[] {
|
270.54861, 340.19128, 63.91854, 331.26220,
|
317.843, 86.631, 240.052, 310.26, 247.23,
|
260.87, 297.82, 343.14, 166.79, 81.53,
|
3.50, 132.75, 182.95, 162.03, 29.8,
|
266.4, 249.2, 157.6, 257.8, 185.1, 69.9,
|
8.0, 197.1, 250.4, 65.3, 162.7, 341.5,
|
291.6, 98.5, 146.7, 110.0, 5.2, 342.6,
|
230.9, 256.1, 45.3, 242.9, 115.2, 151.8,
|
285.3, 53.3, 126.6, 205.7, 85.9,
|
146.1
|
};
|
}
|
}
|
