Kstars

modcalcvizequinox.cpp
1/*
2 SPDX-FileCopyrightText: 2007 Jason Harris <kstars@30doradus.org>
3
4 SPDX-License-Identifier: GPL-2.0-or-later
5*/
6
7#include "modcalcvizequinox.h"
8
9#include "dms.h"
10#include "kstarsdata.h"
11#include "ksnotification.h"
12#include "skyobjects/kssun.h"
13#include "widgets/dmsbox.h"
14
15#include <KLineEdit>
16#include <KPlotAxis>
17#include <KPlotObject>
18#include <KPlotPoint>
19
20#include <cmath>
21#include <QtMath>
22
23modCalcEquinox::modCalcEquinox(QWidget *parentSplit) : QFrame(parentSplit), dSpring(), dSummer(), dAutumn(), dWinter()
24{
25 setupUi(this);
26
27 connect(Year, SIGNAL(valueChanged(int)), this, SLOT(slotCompute()));
28 connect(InputFileBatch, SIGNAL(urlSelected(QUrl)), this, SLOT(slotCheckFiles()));
29 connect(OutputFileBatch, SIGNAL(urlSelected(QUrl)), this, SLOT(slotCheckFiles()));
30 connect(RunButtonBatch, SIGNAL(clicked()), this, SLOT(slotRunBatch()));
31 connect(ViewButtonBatch, SIGNAL(clicked()), this, SLOT(slotViewBatch()));
32
33 Plot->axis(KPlotWidget::LeftAxis)->setLabel(i18n("Sun's Declination"));
34 Plot->setTopPadding(40);
35 //Don't draw Top & Bottom axes; custom axes drawn as plot objects
36 Plot->axis(KPlotWidget::BottomAxis)->setVisible(false);
37 Plot->axis(KPlotWidget::TopAxis)->setVisible(false);
38
39 //This will call slotCompute():
40 Year->setValue(KStarsData::Instance()->lt().date().year());
41
42 RunButtonBatch->setEnabled(false);
43 ViewButtonBatch->setEnabled(false);
44
45 show();
46}
47
48double modCalcEquinox::dmonth(int i)
49{
50 Q_ASSERT(i >= 0 && i < 12 && "Month must be in 0 .. 11 range");
51 return DMonth[i];
52}
53
54void modCalcEquinox::slotCheckFiles()
55{
56 RunButtonBatch->setEnabled(!InputFileBatch->lineEdit()->text().isEmpty() &&
57 !OutputFileBatch->lineEdit()->text().isEmpty());
58}
59
60void modCalcEquinox::slotRunBatch()
61{
62 QString inputFileName = InputFileBatch->url().toLocalFile();
63
64 if (QFile::exists(inputFileName))
65 {
66 QFile f(inputFileName);
67 if (!f.open(QIODevice::ReadOnly))
68 {
69 KSNotification::sorry(i18n("Could not open file %1.", f.fileName()), i18n("Could Not Open File"));
70 inputFileName.clear();
71 return;
72 }
73
74 QTextStream istream(&f);
75 processLines(istream);
76
77 ViewButtonBatch->setEnabled(true);
78
79 f.close();
80 }
81 else
82 {
83 KSNotification::sorry(i18n("Invalid file: %1", inputFileName), i18n("Invalid file"));
84 inputFileName.clear();
85 return;
86 }
87}
88
89void modCalcEquinox::processLines(QTextStream &istream)
90{
91 QFile fOut(OutputFileBatch->url().toLocalFile());
92 fOut.open(QIODevice::WriteOnly);
93 QTextStream ostream(&fOut);
94 int originalYear = Year->value();
95
96 //Write header to output file
97 ostream << i18n("# Timing of Equinoxes and Solstices\n") << i18n("# computed by KStars\n#\n")
98 << i18n("# Vernal Equinox\t\tSummer Solstice\t\t\tAutumnal Equinox\t\tWinter Solstice\n#\n");
99
100 while (!istream.atEnd())
101 {
102 QString line = istream.readLine();
103 bool ok = false;
104 int year = line.toInt(&ok);
105
106 //for now I will simply change the value of the Year widget to trigger
107 //computation of the Equinoxes and Solstices.
108 if (ok)
109 {
110 //triggers slotCompute(), which sets values of dSpring et al.:
111 Year->setValue(year);
112
113 //Write to output file
114 ostream << QLocale().toString(dSpring.date(), QLocale::LongFormat) << "\t"
115 << QLocale().toString(dSummer.date(), QLocale::LongFormat) << "\t"
116 << QLocale().toString(dAutumn.date(), QLocale::LongFormat) << "\t"
117 << QLocale().toString(dWinter.date(), QLocale::LongFormat) << '\n';
118 }
119 }
120
121 if (Year->value() != originalYear)
122 Year->setValue(originalYear);
123}
124
125void modCalcEquinox::slotViewBatch()
126{
127 QFile fOut(OutputFileBatch->url().toLocalFile());
128 fOut.open(QIODevice::ReadOnly);
129 QTextStream istream(&fOut);
130 QStringList text;
131
132 while (!istream.atEnd())
133 text.append(istream.readLine());
134
135 fOut.close();
136
137 KMessageBox::informationList(nullptr, i18n("Results of Sidereal time calculation"), text,
138 OutputFileBatch->url().toLocalFile());
139}
140
141void modCalcEquinox::slotCompute()
142{
143 KStarsData *data = KStarsData::Instance();
144 KSSun Sun;
145 int year0 = Year->value();
146
147 KStarsDateTime dt(QDate(year0, 1, 1), QTime(0, 0, 0));
148 long double jd0 = dt.djd(); //save JD on Jan 1st
149 for (int imonth = 0; imonth < 12; imonth++)
150 {
151 KStarsDateTime kdt(QDate(year0, imonth + 1, 1), QTime(0, 0, 0));
152 DMonth[imonth] = kdt.djd() - jd0;
153 }
154
155 Plot->removeAllPlotObjects();
156
157 //Add the celestial equator, just a single line bisecting the plot horizontally
158 KPlotObject *ce = new KPlotObject(data->colorScheme()->colorNamed("EqColor"), KPlotObject::Lines, 2.0);
159 ce->addPoint(0.0, 0.0);
160 ce->addPoint(366.0, 0.0);
161 Plot->addPlotObject(ce);
162
163 // Tropic of cancer
164 KPlotObject *tcan = new KPlotObject(data->colorScheme()->colorNamed("EqColor"), KPlotObject::Lines, 2.0);
165 tcan->addPoint(0.0, 23.5);
166 tcan->addPoint(366.0, 23.5);
167 Plot->addPlotObject(tcan);
168
169 // Tropic of capricorn
170 KPlotObject *tcap = new KPlotObject(data->colorScheme()->colorNamed("EqColor"), KPlotObject::Lines, 2.0);
171 tcap->addPoint(0.0, -23.5);
172 tcap->addPoint(366.0, -23.5);
173 Plot->addPlotObject(tcap);
174
175 //Add Ecliptic. This is more complicated than simply incrementing the
176 //ecliptic longitude, because we want the x-axis to be time, not RA.
177 //For each day in the year, compute the Sun's position.
178 KPlotObject *ecl = new KPlotObject(data->colorScheme()->colorNamed("EclColor"), KPlotObject::Lines, 2);
179 ecl->setLinePen(QPen(ecl->pen().color(), 4));
180
181 Plot->setLimits(1.0, double(dt.date().daysInYear()), -30.0, 30.0);
182
183 //Add top and bottom axis lines, and custom tickmarks at each month
184 addDateAxes();
185
186 for (int i = 1; i <= dt.date().daysInYear(); i++)
187 {
188 KSNumbers num(dt.djd());
189 Sun.findPosition(&num);
190 ecl->addPoint(double(i), Sun.dec().Degrees());
191
192 dt = dt.addDays(1);
193 }
194 Plot->addPlotObject(ecl);
195
196 // Calculate dates for all solstices and equinoxes
197 findSolsticeAndEquinox(Year->value());
198
199 //Display the Date/Time of each event in the text fields
200 VEquinox->setText(QLocale().toString(dSpring, QLocale::LongFormat));
201 SSolstice->setText(QLocale().toString(dSummer, QLocale::LongFormat));
202 AEquinox->setText(QLocale().toString(dAutumn, QLocale::LongFormat));
203 WSolstice->setText(QLocale().toString(dWinter, QLocale::LongFormat));
204
205 //Add vertical dotted lines at times of the equinoxes and solstices
206 KPlotObject *poSpring = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
207 poSpring->setLinePen(QPen(Qt::white, 1.0, Qt::DotLine));
208 poSpring->addPoint(dSpring.djd() - jd0, Plot->dataRect().top());
209 poSpring->addPoint(dSpring.djd() - jd0, Plot->dataRect().bottom());
210 Plot->addPlotObject(poSpring);
211 KPlotObject *poSummer = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
212 poSummer->setLinePen(QPen(Qt::white, 1.0, Qt::DotLine));
213 poSummer->addPoint(dSummer.djd() - jd0, Plot->dataRect().top());
214 poSummer->addPoint(dSummer.djd() - jd0, Plot->dataRect().bottom());
215 Plot->addPlotObject(poSummer);
216 KPlotObject *poAutumn = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
217 poAutumn->setLinePen(QPen(Qt::white, 1.0, Qt::DotLine));
218 poAutumn->addPoint(dAutumn.djd() - jd0, Plot->dataRect().top());
219 poAutumn->addPoint(dAutumn.djd() - jd0, Plot->dataRect().bottom());
220 Plot->addPlotObject(poAutumn);
221 KPlotObject *poWinter = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
222 poWinter->setLinePen(QPen(Qt::white, 1.0, Qt::DotLine));
223 poWinter->addPoint(dWinter.djd() - jd0, Plot->dataRect().top());
224 poWinter->addPoint(dWinter.djd() - jd0, Plot->dataRect().bottom());
225 Plot->addPlotObject(poWinter);
226}
227
228//Add custom top/bottom axes with tickmarks for each month
229void modCalcEquinox::addDateAxes()
230{
231 KPlotObject *poTopAxis = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
232 poTopAxis->addPoint(0.0, Plot->dataRect().bottom()); //y-axis is reversed!
233 poTopAxis->addPoint(366.0, Plot->dataRect().bottom());
234 Plot->addPlotObject(poTopAxis);
235
236 KPlotObject *poBottomAxis = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
237 poBottomAxis->addPoint(0.0, Plot->dataRect().top() + 0.02);
238 poBottomAxis->addPoint(366.0, Plot->dataRect().top() + 0.02);
239 Plot->addPlotObject(poBottomAxis);
240
241 //Tick mark for each month
242 for (int imonth = 0; imonth < 12; imonth++)
243 {
244 KPlotObject *poMonth = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
245 poMonth->addPoint(dmonth(imonth), Plot->dataRect().top());
246 poMonth->addPoint(dmonth(imonth), Plot->dataRect().top() + 1.4);
247 Plot->addPlotObject(poMonth);
248 poMonth = new KPlotObject(Qt::white, KPlotObject::Lines, 1);
249 poMonth->addPoint(dmonth(imonth), Plot->dataRect().bottom());
250 poMonth->addPoint(dmonth(imonth), Plot->dataRect().bottom() - 1.4);
251 Plot->addPlotObject(poMonth);
252 }
253}
254
255// Calculate and store dates for all solstices and equinoxes
256void modCalcEquinox::findSolsticeAndEquinox(uint32_t year)
257{
258 // All the magic numbers are taken from Meeus - 1991 chapter 27
259 qreal m, jdme[4];
260 if(year > 1000)
261 {
262 m = (year-2000.0) / 1000.0;
263 jdme[0] = (-0.00057 * qPow(m, 4)) + (-0.00411 * qPow(m, 3)) + (0.05169 * qPow(m, 2)) + (365242.37404 * m) + 2451623.80984;
264 jdme[1] = (-0.0003 * qPow(m, 4)) + (0.00888 * qPow(m, 3)) + (0.00325 * qPow(m, 2)) + (365241.62603 * m) + 2451716.56767;
265 jdme[2] = (0.00078 * qPow(m, 4)) + (0.00337 * qPow(m, 3)) + (-0.11575 * qPow(m, 2)) + (365242.01767 * m) + 2451810.21715;
266 jdme[3] = (0.00032 * qPow(m, 4)) + (-0.00823 * qPow(m, 3)) + (-0.06223 * qPow(m, 2)) + (365242.74049 * m) + 2451900.05952;
267 }
268 else
269 {
270 m = year / 1000.0;
271 jdme[0] = (-0.00071 * qPow(m, 4)) + (0.00111 * qPow(m, 3)) + (0.06134 * qPow(m, 2)) + (365242.13740 * m) + 1721139.29189;
272 jdme[1] = (0.00025 * qPow(m, 4)) + (0.00907 * qPow(m, 3)) + (-0.05323 * qPow(m, 2)) + (365241.72562 * m) + 1721233.25401;
273 jdme[2] = (0.00074 * qPow(m, 4)) + (-0.00297 * qPow(m, 3)) + (-0.11677 * qPow(m, 2)) + (365242.49558 * m) + 1721325.70455;
274 jdme[3] = (-0.00006 * qPow(m, 4)) + (-0.00933 * qPow(m, 3)) + (-0.00769 * qPow(m, 2)) + (365242.88257 * m) + 1721414.39987;
275 }
276
277 qreal t[4];
278 for(int i = 0; i < 4; i++)
279 {
280 t[i] = (jdme[i] - 2451545)/36525;
281 }
282
283 qreal a[4][24] = {{485, 203, 199, 182, 156, 136, 77, 74, 70, 58, 52, 50, 45, 44, 29, 18, 17, 16, 14, 12, 12, 12, 9, 8}, {485, 203, 199, 182, 156, 136, 77, 74, 70, 58, 52, 50, 45, 44, 29, 18, 17, 16, 14, 12, 12, 12, 9, 8}, {485, 203, 199, 182, 156, 136, 77, 74, 70, 58, 52, 50, 45, 44, 29, 18, 17, 16, 14, 12, 12, 12, 9, 8}, {485, 203, 199, 182, 156, 136, 77, 74, 70, 58, 52, 50, 45, 44, 29, 18, 17, 16, 14, 12, 12, 12, 9, 8}};
284
285 qreal b[4][24] = {{324.96, 337.23, 342.08, 27.85, 73.14, 171.52, 222.54, 296.72, 243.58, 119.81, 297.17, 21.02, 247.54, 325.15, 60.93, 155.12, 288.79, 198.04, 199.76, 95.39, 287.11, 320.81, 227.73, 15.45}, {324.96, 337.23, 342.08, 27.85, 73.14, 171.52, 222.54, 296.72, 243.58, 119.81, 297.17, 21.02, 247.54, 325.15, 60.93, 155.12, 288.79, 198.04, 199.76, 95.39, 287.11, 320.81, 227.73, 15.45}, {324.96, 337.23, 342.08, 27.85, 73.14, 171.52, 222.54, 296.72, 243.58, 119.81, 297.17, 21.02, 247.54, 325.15, 60.93, 155.12, 288.79, 198.04, 199.76, 95.39, 287.11, 320.81, 227.73, 15.45}, {324.96, 337.23, 342.08, 27.85, 73.14, 171.52, 222.54, 296.72, 243.58, 119.81, 297.17, 21.02, 247.54, 325.15, 60.93, 155.12, 288.79, 198.04, 199.76, 95.39, 287.11, 320.81, 227.73, 15.45}};
286
287 qreal c[] = {1934.136, 32964.467, 20.186, 445267.112, 45036.886, 22518.443, 65928.934, 3034.906, 9037.513, 33718.147, 150.678, 2281.226, 29929.562, 31555.956, 4443.417, 67555.328, 4562.452, 62894.029, 31436.921, 14577.848, 31931.756, 34777.259, 1222.114, 16859.074};
288
289 qreal d[4][24];
290
291 for (int i = 0; i < 4; i++)
292 {
293 for (int j = 0; j < 24; j++)
294 {
295 d[i][j] = c[j] * t[i];
296 }
297 }
298
299 for (int i = 0; i < 4; i++)
300 {
301 for (int j = 0; j < 24; j++)
302 {
303 d[i][j] = qCos(qDegreesToRadians(b[i][j] + d[i][j]));
304 }
305 }
306
307 for (int i = 0; i < 4; i++)
308 {
309 for (int j = 0; j < 24; j++)
310 {
311 d[i][j] = d[i][j] * a[i][j];
312 }
313 }
314
315 qreal e[4], f[4], g[4], jd[4];
316
317 for (int i = 0; i < 4; i++)
318 {
319 e[i] = 0;
320 for (int j = 0; j < 24; j++)
321 {
322 e[i] += d[i][j];
323 }
324 }
325
326 for (int i = 0; i < 4; i++)
327 {
328 f[i] = (35999.373*t[i]-2.47);
329 }
330
331 for (int i = 0; i < 4; i++)
332 {
333 g[i] = 1 + (0.0334 * qCos(qDegreesToRadians(f[i]))) + (0.0007 * qCos(qDegreesToRadians(2*f[i])));
334 }
335
336 for (int i = 0; i < 4; i++)
337 {
338 jd[i] = jdme[i] + (0.00001*e[i]/g[i]);
339 }
340
341 // Get correction
342 qreal correction = FindCorrection(year);
343
344 for (int i = 0; i < 4; i++)
345 {
346 KStarsDateTime dt(jd[i]);
347
348 // Apply correction
349 dt = dt.addSecs(correction);
350
351 switch(i)
352 {
353 case 0 :
354 dSpring = dt;
355 break;
356 case 1 :
357 dSummer = dt;
358 break;
359 case 2 :
360 dAutumn = dt;
361 break;
362 case 3 :
363 dWinter = dt;
364 break;
365 }
366 }
367}
368
369qreal modCalcEquinox::FindCorrection(uint32_t year)
370{
371 uint32_t tblFirst = 1620, tblLast = 2002;
372
373 // Corrections taken from Meeus -1991 chapter 10
374 qreal tbl[] = {/*1620*/ 121,112,103, 95, 88, 82, 77, 72, 68, 63, 60, 56, 53, 51, 48, 46, 44, 42, 40, 38,
375 /*1660*/ 35, 33, 31, 29, 26, 24, 22, 20, 18, 16, 14, 12, 11, 10, 9, 8, 7, 7, 7, 7,
376 /*1700*/ 7, 7, 8, 8, 9, 9, 9, 9, 9, 10, 10, 10, 10, 10, 10, 10, 10, 11, 11, 11,
377 /*1740*/ 11, 11, 12, 12, 12, 12, 13, 13, 13, 14, 14, 14, 14, 15, 15, 15, 15, 15, 16, 16,
378 /*1780*/ 16, 16, 16, 16, 16, 16, 15, 15, 14, 13,
379 /*1800*/ 13.1, 12.5, 12.2, 12.0, 12.0, 12.0, 12.0, 12.0, 12.0, 11.9, 11.6, 11.0, 10.2, 9.2, 8.2,
380 /*1830*/ 7.1, 6.2, 5.6, 5.4, 5.3, 5.4, 5.6, 5.9, 6.2, 6.5, 6.8, 7.1, 7.3, 7.5, 7.6,
381 /*1860*/ 7.7, 7.3, 6.2, 5.2, 2.7, 1.4, -1.2, -2.8, -3.8, -4.8, -5.5, -5.3, -5.6, -5.7, -5.9,
382 /*1890*/ -6.0, -6.3, -6.5, -6.2, -4.7, -2.8, -0.1, 2.6, 5.3, 7.7, 10.4, 13.3, 16.0, 18.2, 20.2,
383 /*1920*/ 21.1, 22.4, 23.5, 23.8, 24.3, 24.0, 23.9, 23.9, 23.7, 24.0, 24.3, 25.3, 26.2, 27.3, 28.2,
384 /*1950*/ 29.1, 30.0, 30.7, 31.4, 32.2, 33.1, 34.0, 35.0, 36.5, 38.3, 40.2, 42.2, 44.5, 46.5, 48.5,
385 /*1980*/ 50.5, 52.5, 53.8, 54.9, 55.8, 56.9, 58.3, 60.0, 61.6, 63.0, 63.8, 64.3};
386
387 qreal deltaT = 0;
388 qreal t = (year - 2000.0)/100.0;
389
390 if(year >= tblFirst && year <= tblLast)
391 {
392 if(year % 2)
393 {
394 deltaT = (tbl[(year-tblFirst-1)/2] + tbl[(year-tblFirst+1)/2]) / 2;
395 }
396 else
397 {
398 deltaT = tbl[(year-tblFirst)/2];
399 }
400 }
401 else if(year < 948)
402 {
403 deltaT = 2177 + 497*t + 44.1*qPow(t, 2);
404 }
405 else if(year >= 948)
406 {
407 deltaT = 102 + 102*t + 25.3*qPow(t, 2);
408 if (year>=2000 && year <=2100)
409 {
410 // Special correction to avoid discontinuity in 2000
411 deltaT += 0.37 * ( year - 2100.0 );
412 }
413 }
414 return -deltaT;
415}
ColorScheme::colorNamed
QColor colorNamed(const QString &name) const
Retrieve a color by name.
Definition colorscheme.cpp:89
KPlotObject::addPoint
void addPoint(const QPointF &p, const QString &label=QString(), double barWidth=0.0)
KPlotObject::setLinePen
void setLinePen(const QPen &p)
KPlotObject::pen
const QPen & pen() const
KSNumbers
There are several time-dependent values used in position calculations, that are not specific to an ob...
Definition ksnumbers.h:43
KSPlanetBase::findPosition
void findPosition(const KSNumbers *num, const CachingDms *lat=nullptr, const CachingDms *LST=nullptr, const KSPlanetBase *Earth=nullptr)
Find position, including correction for Figure-of-the-Earth.
Definition ksplanetbase.cpp:109
KSSun
Child class of KSPlanetBase; encapsulates information about the Sun.
Definition kssun.h:24
KStarsData
KStarsData is the backbone of KStars.
Definition kstarsdata.h:72
KStarsData::colorScheme
ColorScheme * colorScheme()
Definition kstarsdata.h:172
KStarsDateTime
Extension of QDateTime for KStars KStarsDateTime can represent the date/time as a Julian Day,...
Definition kstarsdatetime.h:36
KStarsDateTime::djd
long double djd() const
Definition kstarsdatetime.h:163
SkyPoint::dec
const CachingDms & dec() const
Definition skypoint.h:269
dms::Degrees
const double & Degrees() const
Definition dms.h:141
i18n
QString i18n(const char *text, const TYPE &arg...)
Baloo::toString
char * toString(const EngineQuery &query)
KMessageBox::informationList
void informationList(QWidget *parent, const QString &text, const QStringList &strlist, const QString &title=QString(), const QString &dontShowAgainName=QString(), Options options=Notify)
QDateTime::date
QDate date() const const
QFile::exists
bool exists() const const
QList::append
void append(QList< T > &&value)
QLocale::toString
QString toString(QDate date, FormatType format) const const
QPen::color
QColor color() const const
QString::clear
void clear()
QString::toInt
int toInt(bool *ok, int base) const const
QTextStream::atEnd
bool atEnd() const const
QTextStream::readLine
QString readLine(qint64 maxlen)
QtFuture::connect
QFuture< ArgsType< Signal > > connect(Sender *sender, Signal signal)
This file is part of the KDE documentation.
Documentation copyright © 1996-2024 The KDE developers.
Generated on Mon Nov 4 2024 16:38:44 by doxygen 1.12.0 written by Dimitri van Heesch, © 1997-2006

KDE's Doxygen guidelines are available online.