How to calculate sunrise and sunset times?
Here is a complete routine to that calculates the sunrise (optionally the sunset, see the code) in C++. It only requires latitude and longitude as input, and is accurate to within seconds for the civil sunrise time. This code is running in a UWP C++ app created with Visual Studio 2017. It includes a subroutine to fill the output minutes and seconds with zeros in the case of single digit results.
String^ toplatformstring(bool fill, long ll) {
// convert ll to platform string
Platform::String^ p_string;
std::string doit;
if (fill == false) {
doit = std::to_string(ll); // convert, don't fill with zeros
}
else {
//convert ll to std doit and fill with zeros
std::stringstream ss;
ss << std::setw(2) << std::setfill('0') << ll;
doit = ss.str();
}
//convert doit to platform string
char const *pchar = doit.c_str();
std::string s_str = std::string(pchar);
std::wstring wid_str = std::wstring(s_str.begin(), s_str.end());
const wchar_t* w_char = wid_str.c_str();
p_string = ref new Platform::String(w_char);
return p_string;
}
//double la = 39.299236; // baltimore
//double lo = -76.609383;
//double la = 37.0; // SF California
//double lo = -122.0;
Platform::String^ MainPage::sunrise(double la, double lo) {
/*double la = 39.300213;
double lo = -76.610516;*/
Platform::String^ sunrisetime;
//// get year, month, day integers
time_t rawtime;
struct tm timeinfo; // get date and time info
time(&rawtime);
localtime_s(&timeinfo, &rawtime);
double xday = timeinfo.tm_mday;
double xmonth = timeinfo.tm_mon;
xmonth = xmonth + 1; // correct for origin 0
//textblockc->Text = xmonth.ToString(); // for debugging
double xyear = timeinfo.tm_year;
//double dayofyear = timeinfo.tm_yday; // day of year also
// calculate the day of the year
// N1 = floor(275 * month / 9)
double xxN1 = floor(275 * xmonth / 9);
// N2 = floor((month + 9) / 12)
double xxN2 = floor((xmonth + 9) / 12);
// N3 = (1 + floor((year - 4 * floor(year / 4) + 2) / 3))
double xxN3 = (1 + floor((xyear - 4 * floor(xyear / 4) + 2) / 3));
// N = N1 - (N2 * N3) + day - 30
double day = xxN1 - (xxN2 * xxN3) + xday - 30;
double zenith = 90.83333333333333;
double D2R = M_PI / 180;
double R2D = 180 / M_PI;
// convert the longitude to hour value and calculate an approximate time
double lnHour = lo / 15;
double t;
//if (sunrise) {
t = day + ((6 - lnHour) / 24);
//} else {
//t = day + ((18 - lnHour) / 24);
//};
//calculate the Sun's mean anomaly
double M = (0.9856 * t) - 3.289;
//calculate the Sun's true longitude
double L = M + (1.916 * sin(M * D2R)) + (0.020 * sin(2 * M * D2R)) + 282.634;
if (L > 360) {
L = L - 360;
}
else if (L < 0) {
L = L + 360;
};
//calculate the Sun's right ascension
double RA = R2D * atan(0.91764 * tan(L * D2R));
if (RA > 360) {
RA = RA - 360;
}
else if (RA < 0) {
RA = RA + 360;
};
//right ascension value needs to be in the same qua
double Lquadrant = (floor(L / (90))) * 90;
double RAquadrant = (floor(RA / 90)) * 90;
RA = RA + (Lquadrant - RAquadrant);
//right ascension value needs to be converted into hours
RA = RA / 15;
//calculate the Sun's declination
double sinDec = 0.39782 * sin(L * D2R);
double cosDec = cos(asin(sinDec));
//calculate the Sun's local hour angle
double cosH = (cos(zenith * D2R) - (sinDec * sin(la * D2R))) / (cosDec * cos(la * D2R));
double H;
//if (sunrise) {
H = 360 - R2D * acos(cosH);
//} else {
//H = R2D * Math.acos(cosH)
//};
H = H / 15;
//calculate local mean time of rising/setting
double T = H + RA - (0.06571 * t) - 6.622;
//adjust back to UTC
double UT = T - lnHour;
if (UT > 24) {
UT = UT - 24;
}
else if (UT < 0) {
UT = UT + 24;
}
//convert UT value to local time zone of latitude/longitude
int offset = (int)(lo / 15); // estimate utc correction
double localT = UT + offset; // -5 for baltimore
//convert to seconds
int seconds = (int)(localT * 3600);
long sec = seconds % 60;
long minutes = seconds % 3600 / 60;
long hours = seconds % 86400 / 3600;
hours = hours % 12;
Platform::String^ ssec = toplatformstring(true, sec);
Platform::String^ mminutes = toplatformstring(true, minutes);
Platform::String^ hhours = toplatformstring(false, hours);
sunrisetime = hhours + ":" + mminutes + ":" + ssec;
return sunrisetime;
}