diff options
Diffstat (limited to 'src/solar.c')
| -rw-r--r-- | src/solar.c | 209 |
1 files changed, 0 insertions, 209 deletions
diff --git a/src/solar.c b/src/solar.c deleted file mode 100644 index d8ece88..0000000 --- a/src/solar.c +++ /dev/null @@ -1,209 +0,0 @@ -/* solar.c -- Solar position source - * This file is part of redshift-ng. - * - * redshift-ng is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * redshift-ng 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 General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with redshift-ng. If not, see <http://www.gnu.org/licenses/>. - * - * Copyright (c) 2010 Jon Lund Steffensen <jonlst@gmail.com> -*/ - -/* Ported from javascript code by U.S. Department of Commerce, - * National Oceanic & Atmospheric Administration: - * http://www.srrb.noaa.gov/highlights/sunrise/calcdetails.html - * It is based on equations from "Astronomical Algorithms" by - * Jean Meeus. */ - -#include "common.h" -#include "solar.h" - -#define RAD(x) ((x)*(M_PI/180)) -#define DEG(x) ((x)*(180/M_PI)) - - - -/* Julian day from unix epoch */ -static double -jd_from_epoch(double t) -{ - return (t / 86400.0) + 2440587.5; -} - -/* Julian centuries since J2000.0 from Julian day */ -static double -jcent_from_jd(double jd) -{ - return (jd - 2451545.0) / 36525.0; -} - -/* Julian day from Julian centuries since J2000.0 */ -static double -jd_from_jcent(double t) -{ - return 36525.0*t + 2451545.0; -} - -/* Geometric mean longitude of the sun. - t: Julian centuries since J2000.0 - Return: Geometric mean logitude in radians. */ -static double -sun_geom_mean_lon(double t) -{ - /* FIXME returned value should always be positive */ - return RAD(fmod(280.46646 + t*(36000.76983 + t*0.0003032), 360)); -} - -/* Geometric mean anomaly of the sun. - t: Julian centuries since J2000.0 - Return: Geometric mean anomaly in radians. */ -static double -sun_geom_mean_anomaly(double t) -{ - return RAD(357.52911 + t*(35999.05029 - t*0.0001537)); -} - -/* Eccentricity of earth orbit. - t: Julian centuries since J2000.0 - Return: Eccentricity (unitless). */ -static double -earth_orbit_eccentricity(double t) -{ - return 0.016708634 - t*(0.000042037 + t*0.0000001267); -} - -/* Equation of center of the sun. - t: Julian centuries since J2000.0 - Return: Center(?) in radians */ -static double -sun_equation_of_center(double t) -{ - /* Use the first three terms of the equation. */ - double m = sun_geom_mean_anomaly(t); - double c = sin(m)*(1.914602 - t*(0.004817 + 0.000014*t)) + - sin(2*m)*(0.019993 - 0.000101*t) + - sin(3*m)*0.000289; - return RAD(c); -} - -/* True longitude of the sun. - t: Julian centuries since J2000.0 - Return: True longitude in radians */ -static double -sun_true_lon(double t) -{ - double l_0 = sun_geom_mean_lon(t); - double c = sun_equation_of_center(t); - return l_0 + c; -} - -/* Apparent longitude of the sun. (Right ascension). - t: Julian centuries since J2000.0 - Return: Apparent longitude in radians */ -static double -sun_apparent_lon(double t) -{ - double o = sun_true_lon(t); - return RAD(DEG(o) - 0.00569 - 0.00478*sin(RAD(125.04 - 1934.136*t))); -} - -/* Mean obliquity of the ecliptic - t: Julian centuries since J2000.0 - Return: Mean obliquity in radians */ -static double -mean_ecliptic_obliquity(double t) -{ - double sec = 21.448 - t*(46.815 + t*(0.00059 - t*0.001813)); - return RAD(23.0 + (26.0 + (sec/60.0))/60.0); -} - -/* Corrected obliquity of the ecliptic. - t: Julian centuries since J2000.0 - Return: Currected obliquity in radians */ -static double -obliquity_corr(double t) -{ - double e_0 = mean_ecliptic_obliquity(t); - double omega = 125.04 - t*1934.136; - return RAD(DEG(e_0) + 0.00256*cos(RAD(omega))); -} - -/* Declination of the sun. - t: Julian centuries since J2000.0 - Return: Declination in radians */ -static double -solar_declination(double t) -{ - double e = obliquity_corr(t); - double lambda = sun_apparent_lon(t); - return asin(sin(e)*sin(lambda)); -} - -/* Difference between true solar time and mean solar time. - t: Julian centuries since J2000.0 - Return: Difference in minutes */ -static double -equation_of_time(double t) -{ - double epsilon = obliquity_corr(t); - double l_0 = sun_geom_mean_lon(t); - double e = earth_orbit_eccentricity(t); - double m = sun_geom_mean_anomaly(t); - double y = pow(tan(epsilon/2.0), 2.0); - - double eq_time = y*sin(2*l_0) - 2*e*sin(m) + - 4*e*y*sin(m)*cos(2*l_0) - - 0.5*y*y*sin(4*l_0) - - 1.25*e*e*sin(2*m); - return 4*DEG(eq_time); -} - -/* Angular elevation at the location for the given hour angle. - lat: Latitude of location in degrees - decl: Declination in radians - ha: Hour angle in radians - Return: Angular elevation in radians */ -static double -elevation_from_hour_angle(double lat, double decl, double ha) -{ - return asin(cos(ha)*cos(RAD(lat))*cos(decl) + - sin(RAD(lat))*sin(decl)); -} - -/* Solar angular elevation at the given location and time. - t: Julian centuries since J2000.0 - lat: Latitude of location - lon: Longitude of location - Return: Solar angular elevation in radians */ -static double -solar_elevation_from_time(double t, double lat, double lon) -{ - /* Minutes from midnight */ - double jd = jd_from_jcent(t); - double offset = (jd - round(jd) - 0.5)*1440.0; - - double eq_time = equation_of_time(t); - double ha = RAD((720 - offset - eq_time)/4 - lon); - double decl = solar_declination(t); - return elevation_from_hour_angle(lat, decl, ha); -} - -/* Solar angular elevation at the given location and time. - date: Seconds since unix epoch - lat: Latitude of location - lon: Longitude of location - Return: Solar angular elevation in degrees */ -double -solar_elevation(double date, double lat, double lon) -{ - double jd = jd_from_epoch(date); - return DEG(solar_elevation_from_time(jcent_from_jd(jd), lat, lon)); -} |