diff options
Diffstat (limited to 'src/solar.c')
| -rw-r--r-- | src/solar.c | 69 |
1 files changed, 40 insertions, 29 deletions
diff --git a/src/solar.c b/src/solar.c index cbee27d..8a57507 100644 --- a/src/solar.c +++ b/src/solar.c @@ -20,6 +20,12 @@ +#if !defined(CLOCK_REALTIME_COARSE) +# define CLOCK_REALTIME_COARSE CLOCK_REALTIME +#endif + + + /** * Get current Julian Centuries time (100 Julian days since J2000.) * @@ -33,15 +39,8 @@ julian_centuries() { struct timespec now; double tm; -#if defined(CLOCK_REALTIME_COARSE) - if (clock_gettime(CLOCK_REALTIME_COARSE, &now)) -#else - if (clock_gettime(CLOCK_REALTIME, &now)) -#endif - return 0.0; - tm = (double)(now.tv_nsec); - tm /= 1000000000.0; - tm += (double)(now.tv_sec); + if (clock_gettime(CLOCK_REALTIME_COARSE, &now)) return 0.0; + tm = (double)(now.tv_nsec) / 1000000000.0 + (double)(now.tv_sec); tm = (tm / 86400.0 + 2440587.5 - 2451545.0) / 36525.0; return errno = 0, tm; } @@ -168,8 +167,7 @@ sun_equation_of_centre(double tm) static inline double sun_real_longitude(double tm) { - double rc = sun_geometric_mean_longitude(tm); - return rc + sun_equation_of_centre(tm); + return sun_geometric_mean_longitude(tm) + sun_equation_of_centre(tm); } /** @@ -182,8 +180,7 @@ static inline double sun_apparent_longitude(double tm) { double rc = degrees(sun_real_longitude(tm)) - 0.00569; - rc -= 0.00478 * sin(radians(-1934.136 * tm + 125.04)); - return radians(rc); + return radians(rc - 0.00478 * sin(radians(-1934.136 * tm + 125.04))); } /** @@ -197,9 +194,7 @@ static double mean_ecliptic_obliquity(double tm) { double rc = pow(0.001813 * tm, 3.0) - pow(0.00059 * tm, 2.0) - 46.815 * tm + 21.448; - rc = 26 + rc / 60; - rc = 23 + rc / 60; - return radians(rc); + return radians(23.0 + (26.0 + rc / 60.0) / 60.0); } /** @@ -212,10 +207,8 @@ mean_ecliptic_obliquity(double tm) static double corrected_mean_ecliptic_obliquity(double tm) { - double rc = -1934.136 * tm + 125.04; - rc = 0.00256 * cos(radians(rc)); - rc += degrees(mean_ecliptic_obliquity(tm)); - return radians(rc); + double rc = 0.00256 * cos(radians(-1934.136 * tm + 125.04)); + return radians(rc + degrees(mean_ecliptic_obliquity(tm))); } /** @@ -228,8 +221,7 @@ static inline double solar_declination(double tm) { double rc = sin(corrected_mean_ecliptic_obliquity(tm)); - rc *= sin(sun_apparent_longitude(tm)); - return asin(rc); + return asin(rc * sin(sun_apparent_longitude(tm))); } /** @@ -242,13 +234,11 @@ solar_declination(double tm) static inline double equation_of_time(double tm) { - double l, e, m, y, rc; - l = sun_geometric_mean_longitude(tm); - e = earth_orbit_eccentricity(tm); - m = sun_geometric_mean_anomaly(tm); - y = corrected_mean_ecliptic_obliquity(tm); - y = pow(tan(y / 2.0), 2.0); - rc = y * sin(2.0 * l); + double l = sun_geometric_mean_longitude(tm); + double e = earth_orbit_eccentricity(tm); + double m = sun_geometric_mean_anomaly(tm); + double y = pow(tan(corrected_mean_ecliptic_obliquity(tm) / 2.0), 2.0); + double rc = y * sin(2.0 * l); rc += (4.0 * y * cos(2.0 * l) - 2.0) * e * sin(m); rc -= pow(0.5 * y, 2.0) * sin(4.0 * l); rc -= pow(1.25 * e, 2.0) * sin(2.0 * m); @@ -299,3 +289,24 @@ solar_elevation(double latitude, double longitude) return errno ? -1 : degrees(solar_elevation_from_time(tm, latitude, longitude)); } + +/** + * Exit if time the is before year 0 in J2000. + * + * @return 0 on success, -1 on error. + */ +#if defined(TIMETRAVELLER) +int +check_timetravel(void) +{ + struct timespec now; + if (clock_gettime(CLOCK_REALTIME, &now)) return -1; + if (now.tv_nsec < (time_t)946728000L) + fprintf(stderr, "We have detected that you are a time-traveller" + "(or your clock is not configured correctly.)" + "Please recompile with -DTIMETRAVELLER" + "(or correct your clock.)"), exit(1); + return 0; +} +#endif + |