1 /*
2 SPDX-FileCopyrightText: 2009 Prakash Mohan <prakash.mohan@kdemail.net>
3
4 SPDX-License-Identifier: GPL-2.0-or-later
5 */
6
7 #include "ksalmanac.h"
8
9 #include "geolocation.h"
10 #include "ksnumbers.h"
11 #include "kstarsdata.h"
12
KSAlmanac()13 KSAlmanac::KSAlmanac()
14 {
15 KStarsData *data = KStarsData::Instance();
16
17 /*dt = KStarsDateTime::currentDateTime();
18 geo = data->geo();
19 dt.setTime(QTime());
20 dt = geo->LTtoUT(dt);*/
21
22 // Jasem 2015-08-24 Do NOT use KStarsDataTime for local time, it is only for UTC
23 QDateTime midnight = QDateTime(data->lt().date(), QTime());
24 geo = data->geo();
25 dt = geo->LTtoUT(KStarsDateTime(midnight));
26 update();
27 }
28
KSAlmanac(const KStarsDateTime & midnight,const GeoLocation * g)29 KSAlmanac::KSAlmanac(const KStarsDateTime &midnight, const GeoLocation *g)
30 {
31 geo = g ? g : KStarsData::Instance()->geo();
32
33 dt = midnight.isValid() ?
34 midnight.timeSpec() == Qt::LocalTime ?
35 geo->LTtoUT(midnight) :
36 midnight :
37 KStarsData::Instance()->ut();
38
39 update();
40 }
41
update()42 void KSAlmanac::update()
43 {
44 RiseSetTime(&m_Sun, &SunRise, &SunSet, &SunRiseT, &SunSetT);
45 RiseSetTime(&m_Moon, &MoonRise, &MoonSet, &MoonRiseT, &MoonSetT);
46 // qDebug() << "Sun rise: " << SunRiseT.toString() << " Sun set: " << SunSetT.toString() << " Moon rise: " << MoonRiseT.toString() << " Moon set: " << MoonSetT.toString();
47 findDawnDusk();
48 findMoonPhase();
49 }
50
RiseSetTime(SkyObject * o,double * riseTime,double * setTime,QTime * RiseTime,QTime * SetTime)51 void KSAlmanac::RiseSetTime(SkyObject *o, double *riseTime, double *setTime, QTime *RiseTime, QTime *SetTime)
52 {
53 // Compute object rise and set times
54 const KStarsDateTime today = dt;
55 const GeoLocation *_geo = geo;
56 *RiseTime = o->riseSetTime(
57 today, _geo,
58 true); // FIXME: Should we add a day here so that we report future rise time? Not doing so produces the right results for the moon. Not sure about the sun.
59 *SetTime = o->riseSetTime(today, _geo, false);
60 *riseTime = -1.0 * RiseTime->secsTo(QTime(0, 0, 0, 0)) / 86400.0;
61 *setTime = -1.0 * SetTime->secsTo(QTime(0, 0, 0, 0)) / 86400.0;
62
63 // Check to see if the object is circumpolar
64 // NOTE: Since we are working on a local copy of the Sun / Moon,
65 // we freely change the geolocation / time without setting
66 // them back.
67
68 KSNumbers num(dt.djd());
69 CachingDms LST = geo->GSTtoLST(dt.gst());
70 o->updateCoords(&num, true, geo->lat(), &LST, true);
71 if (o->checkCircumpolar(geo->lat()))
72 {
73 if (o->alt().Degrees() > 0.0)
74 {
75 //Circumpolar, signal it this way:
76 *riseTime = 0.0;
77 *setTime = 1.0;
78 }
79 else
80 {
81 //never rises, signal it this way:
82 *riseTime = 0.0;
83 *setTime = -1.0;
84 }
85 }
86 }
87
findDawnDusk(double altitude)88 void KSAlmanac::findDawnDusk(double altitude)
89 {
90 KStarsDateTime today = dt;
91 KSNumbers num(today.djd());
92 CachingDms LST = geo->GSTtoLST(today.gst());
93
94 // Relocate our local Sun to this almanac time - local midnight
95 m_Sun.updateCoords(&num, true, geo->lat(), &LST, true);
96
97 // Granularity
98 int const h_inc = 5;
99
100 // Compute the altitude of the Sun twelve hours before this almanac time
101 int const start_h = -1200, end_h = +1200;
102 double last_alt = findAltitude(&m_Sun, start_h/100.0);
103
104 int dawn = -1300, dusk = -1300, min_alt_time = -1300;
105 double max_alt = -100.0, min_alt = +100.0;
106
107 // Compute the dawn and dusk times in a [-12,+12] hours around the day midnight of this almanac as well as min and max altitude
108 // See the header comment about dawn and dusk positions
109 for (int h = start_h + h_inc; h <= end_h; h += h_inc)
110 {
111 // Compute the Sun's altitude in an increasing hour interval
112 double const alt = findAltitude(&m_Sun, h/100.0);
113
114 // Deduce whether the Sun is rising or setting
115 bool const rising = alt - last_alt > 0;
116
117 // Extend min/max altitude interval, push minimum time down
118 if (max_alt < alt)
119 {
120 max_alt = alt;
121 }
122 else if (alt < min_alt)
123 {
124 min_alt = alt;
125 min_alt_time = h;
126 }
127
128 // Dawn is when the Sun is rising and crosses an altitude of -18 degrees
129 if (dawn < 0)
130 if (rising)
131 if (last_alt <= altitude && altitude <= alt)
132 dawn = h - h_inc * (alt == last_alt ? 0 : (alt - altitude)/(alt - last_alt));
133
134 // Dusk is when the Sun is setting and crosses an altitude of -18 degrees
135 if (dusk < 0)
136 if (!rising)
137 if (last_alt >= altitude && altitude >= alt)
138 dusk = h - h_inc * (alt == last_alt ? 0 : (alt - altitude)/(alt - last_alt));
139
140 last_alt = alt;
141 }
142
143 // If the Sun did not cross the astronomical twilight altitude, use the minimal altitude
144 if (dawn < start_h || dusk < start_h)
145 {
146 DawnAstronomicalTwilight = static_cast <double> (min_alt_time) / 2400.0;
147 DuskAstronomicalTwilight = static_cast <double> (min_alt_time) / 2400.0;
148 }
149 // If the Sun did cross the astronomical twilight, use the computed time offsets
150 else
151 {
152 DawnAstronomicalTwilight = static_cast <double> (dawn) / 2400.0;
153 DuskAstronomicalTwilight = static_cast <double> (dusk) / 2400.0;
154 }
155
156 SunMaxAlt = max_alt;
157 SunMinAlt = min_alt;
158 }
159
findMoonPhase()160 void KSAlmanac::findMoonPhase()
161 {
162 const KStarsDateTime today = dt;
163 KSNumbers num(today.djd());
164 CachingDms LST = geo->GSTtoLST(today.gst());
165
166 m_Sun.updateCoords(&num, true, geo->lat(), &LST, true); // We can abuse our own copy of the sun and/or moon
167 m_Moon.updateCoords(&num, true, geo->lat(), &LST, true);
168 m_Moon.findPhase(&m_Sun);
169 MoonPhase = m_Moon.phase().Degrees();
170 }
171
setDate(const KStarsDateTime & utc_midnight)172 void KSAlmanac::setDate(const KStarsDateTime &utc_midnight)
173 {
174 dt = utc_midnight;
175 update();
176 }
177
setLocation(const GeoLocation * geo_)178 void KSAlmanac::setLocation(const GeoLocation *geo_)
179 {
180 geo = geo_;
181 update();
182 }
183
sunZenithAngleToTime(double z) const184 double KSAlmanac::sunZenithAngleToTime(double z) const
185 {
186 // TODO: Correct for movement of the sun
187 double HA = acos((cos(z * dms::DegToRad) - m_Sun.dec().sin() * geo->lat()->sin()) /
188 (m_Sun.dec().cos() * geo->lat()->cos()));
189 double HASunset = acos((-m_Sun.dec().sin() * geo->lat()->sin()) / (m_Sun.dec().cos() * geo->lat()->cos()));
190 return SunSet + (HA - HASunset) / 24.0;
191 }
192
findAltitude(const SkyPoint * p,double hour)193 double KSAlmanac::findAltitude(const SkyPoint *p, double hour)
194 {
195 return SkyPoint::findAltitude(p, dt, geo, hour).Degrees();
196 }
197