fix for #3400
replace low_accuracy math functions (sint_t, cos_t, atan_t, ...) with standard libm functions that have higher accuracy.
This commit is contained in:
parent
5eadbe7ecd
commit
166316e0c5
@ -412,8 +412,8 @@ int getSunriseUTC(int year, int month, int day, float lat, float lon, bool sunse
|
|||||||
//1. first calculate the day of the year
|
//1. first calculate the day of the year
|
||||||
float N1 = 275 * month / 9;
|
float N1 = 275 * month / 9;
|
||||||
float N2 = (month + 9) / 12;
|
float N2 = (month + 9) / 12;
|
||||||
float N3 = (1 + floor_t((year - 4 * floor_t(year / 4) + 2) / 3));
|
float N3 = (1.0f + floorf((year - 4 * floorf(year / 4) + 2.0f) / 3.0f));
|
||||||
float N = N1 - (N2 * N3) + day - 30;
|
float N = N1 - (N2 * N3) + day - 30.0f;
|
||||||
|
|
||||||
//2. convert the longitude to hour value and calculate an approximate time
|
//2. convert the longitude to hour value and calculate an approximate time
|
||||||
float lngHour = lon / 15.0f;
|
float lngHour = lon / 15.0f;
|
||||||
@ -423,37 +423,37 @@ int getSunriseUTC(int year, int month, int day, float lat, float lon, bool sunse
|
|||||||
float M = (0.9856f * t) - 3.289f;
|
float M = (0.9856f * t) - 3.289f;
|
||||||
|
|
||||||
//4. calculate the Sun's true longitude
|
//4. calculate the Sun's true longitude
|
||||||
float L = fmod_t(M + (1.916f * sin_t(DEG_TO_RAD*M)) + (0.02f * sin_t(2*DEG_TO_RAD*M)) + 282.634f, 360.0f);
|
float L = fmodf(M + (1.916f * sinf(DEG_TO_RAD*M)) + (0.02f * sinf(2*DEG_TO_RAD*M)) + 282.634f, 360.0f);
|
||||||
|
|
||||||
//5a. calculate the Sun's right ascension
|
//5a. calculate the Sun's right ascension
|
||||||
float RA = fmod_t(RAD_TO_DEG*atan_t(0.91764f * tan_t(DEG_TO_RAD*L)), 360.0f);
|
float RA = fmodf(RAD_TO_DEG*atan(0.91764f * tan(DEG_TO_RAD*L)), 360.0f);
|
||||||
|
|
||||||
//5b. right ascension value needs to be in the same quadrant as L
|
//5b. right ascension value needs to be in the same quadrant as L
|
||||||
float Lquadrant = floor_t( L/90) * 90;
|
float Lquadrant = floorf( L/90) * 90;
|
||||||
float RAquadrant = floor_t(RA/90) * 90;
|
float RAquadrant = floorf(RA/90) * 90;
|
||||||
RA = RA + (Lquadrant - RAquadrant);
|
RA = RA + (Lquadrant - RAquadrant);
|
||||||
|
|
||||||
//5c. right ascension value needs to be converted into hours
|
//5c. right ascension value needs to be converted into hours
|
||||||
RA /= 15.0f;
|
RA /= 15.0f;
|
||||||
|
|
||||||
//6. calculate the Sun's declination
|
//6. calculate the Sun's declination
|
||||||
float sinDec = 0.39782f * sin_t(DEG_TO_RAD*L);
|
float sinDec = 0.39782f * sinf(DEG_TO_RAD*L);
|
||||||
float cosDec = cos_t(asin_t(sinDec));
|
float cosDec = cosf(asinf(sinDec));
|
||||||
|
|
||||||
//7a. calculate the Sun's local hour angle
|
//7a. calculate the Sun's local hour angle
|
||||||
float cosH = (sin_t(DEG_TO_RAD*ZENITH) - (sinDec * sin_t(DEG_TO_RAD*lat))) / (cosDec * cos_t(DEG_TO_RAD*lat));
|
float cosH = (sinf(DEG_TO_RAD*ZENITH) - (sinDec * sinf(DEG_TO_RAD*lat))) / (cosDec * cosf(DEG_TO_RAD*lat));
|
||||||
if (cosH > 1 && !sunset) return 0; // the sun never rises on this location (on the specified date)
|
if ((cosH > 1.0f) && !sunset) return 0; // the sun never rises on this location (on the specified date)
|
||||||
if (cosH < -1 && sunset) return 0; // the sun never sets on this location (on the specified date)
|
if ((cosH < -1.0f) && sunset) return 0; // the sun never sets on this location (on the specified date)
|
||||||
|
|
||||||
//7b. finish calculating H and convert into hours
|
//7b. finish calculating H and convert into hours
|
||||||
float H = sunset ? RAD_TO_DEG*acos_t(cosH) : 360 - RAD_TO_DEG*acos_t(cosH);
|
float H = sunset ? RAD_TO_DEG*acosf(cosH) : 360 - RAD_TO_DEG*acosf(cosH);
|
||||||
H /= 15.0f;
|
H /= 15.0f;
|
||||||
|
|
||||||
//8. calculate local mean time of rising/setting
|
//8. calculate local mean time of rising/setting
|
||||||
float T = H + RA - (0.06571f * t) - 6.622f;
|
float T = H + RA - (0.06571f * t) - 6.622f;
|
||||||
|
|
||||||
//9. adjust back to UTC
|
//9. adjust back to UTC
|
||||||
float UT = fmod_t(T - lngHour, 24.0f);
|
float UT = fmodf(T - lngHour, 24.0f);
|
||||||
|
|
||||||
// return in minutes from midnight
|
// return in minutes from midnight
|
||||||
return UT*60;
|
return UT*60;
|
||||||
|
Loading…
Reference in New Issue
Block a user