/*
Copyright (c) 2014, Vladimir Agafonkin
Copyright (c) 2024, cidoku
All rights reserved.

Redistribution and use in source and binary forms, with or without
modification, are permitted provided that the following conditions are met:

   1. Redistributions of source code must retain the above copyright notice,
   this list of conditions and the following disclaimer.

   2. Redistributions in binary form must reproduce the above copyright notice,
   this list of conditions and the following disclaimer in the documentation
   and/or other materials provided with the distribution.

THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
*/

#ifndef SUNCALC_H
#define SUNCALC_H

// date/time constants and conversions
#define RAD		(M_PI / 180)
#define DAY_SEC	86400
#define MON_SEC (29.53 * DAY_SEC)
#define J1970	2440588
#define J2000	2451545
#define JULIAN(date)	(date) / DAY_SEC - 0.5 + J1970
#define FROM_JULIAN(j)	(time_t)(((j) + 0.5 - J1970) * DAY_SEC)
#define TO_DAYS(date)	JULIAN(date) - J2000

// general calculations for position
#define OBLIQUITY				RAD * 23.4397 // obliquity of the Earth
#define RIGHT_ASCENSION(l, b)	atan2(sin(l) * cos(OBLIQUITY) - tan(b) * sin(OBLIQUITY), cos(l))
#define DECLINATION(l, b)		asin(sin(b) * cos(OBLIQUITY) + cos(b) * sin(OBLIQUITY) * sin(l))
#define AZIMUTH(H, phi, dec)	atan2(sin(H), cos(H) * sin(phi) - tan(dec) * cos(phi))
#define ALTITUDE(H, phi, dec)	asin(sin(phi) * sin(dec) + cos(phi) * cos(dec) * cos(H))
#define SIDEREAL_TIME(d, lw)	RAD * (280.16 + 360.9856235 * d) - lw

// general sun calculations
// sun calculations are based on http://aa.quae.nl/en/reken/zonpositie.html formulas
#define SOLAR_MEAN_ANOMALY(d)	RAD * (357.5291 + 0.98560028 * d)
// M + pi + perihelion of the Earth + equation of center
#define ECLIPTIC_LONGITUDE(M)	M + M_PI + (RAD * (1.9148 * sin(M) + 0.02 * sin(2 * M) + 0.0003 * sin(3 * M))) + (RAD * 102.9372)

// calculations for sun times
#define J0	0.0009
#define JULIAN_CYCLE(d, lw)			round(d - J0 - lw / (2 * M_PI))
#define APPROX_TRANSIT(Ht, lw, n)	J0 + (Ht + lw) / (2 * M_PI) + n
#define SOLAR_TRANSIT_J(ds, M, L)	J2000 + ds + 0.0053 * sin(M) - 0.0069 * sin(2 * L)
#define HOUR_ANGLE(h, phi, d)		acos((sin(h) - sin(phi) * sin(d)) / (cos(phi) * cos(d)))
#define OBSERVER_ANGLE(height)		-2.076 * sqrt(height) / 60
// returns set time for the given sun altitude
#define GET_SET_J(h, lw, phi, dec, n, M, L) SOLAR_TRANSIT_J(APPROX_TRANSIT(HOUR_ANGLE(h, phi, dec), lw, n), M, L)

// sun times configuration: sunrise, sunrise_end, dawn, nautical_dawn,
// night_end, golden_hour_end
#define SUN_TIMES_LEN 6
extern const double sun_times[SUN_TIMES_LEN];
struct s_data {
	time_t solar_noon;
	time_t nadir;
	time_t sunrise;
	time_t sunset;
	time_t sunrise_end;
	time_t sunset_start;
	time_t dawn;
	time_t dusk;
	time_t nautical_dawn;
	time_t nautical_dusk;
	time_t night_end;
	time_t night;
	time_t golden_hour_end;
	time_t golden_hour;
	int length;
};

// calculates sun times for a given date, latitude/longitude, and, optionally,
// the observer height (in meters) relative to the horizon
void get_sun_times(time_t date, double lat, double lng, int height, struct s_data* times);

// moon coords calculation based on http://aa.quae.nl/en/reken/hemelpositie.html formulas
// geocentric ecliptic coordinates of the moon
void moon_coords(double d, double* dec, double* ra, double* dist);

// calculations for illumination parameters of the moon,
// based on http://idlastro.gsfc.nasa.gov/ftp/pro/astro/mphase.pro formulas and
// Chapter 48 of "Astronomical Algorithms" 2nd edition by Jean Meeus (Willmann-Bell, Richmond) 1998.
struct m_light {
	double fraction;
	double phase;
	double angle;
};

void get_moon_illumination(int date, struct m_light* m_light);

#endif // SUNCALC_H