Changeset View
Changeset View
Standalone View
Standalone View
src/astroseasons.cpp
1 | /* | 1 | /* | ||
---|---|---|---|---|---|
2 | This file is part of the kholidays library. | 2 | This file is part of the kholidays library. | ||
3 | 3 | | |||
4 | Copyright (c) 2004,2006-2007 Allen Winter <winter@kde.org> | 4 | Copyright (c) 2004,2006-2007 Allen Winter <winter@kde.org> | ||
5 | Copyright (c) 2018 Daniel Vrátil <dvratil@kde.org> | ||||
5 | 6 | | |||
6 | This library is free software; you can redistribute it and/or | 7 | This library is free software; you can redistribute it and/or | ||
7 | modify it under the terms of the GNU Library General Public | 8 | modify it under the terms of the GNU Library General Public | ||
8 | License as published by the Free Software Foundation; either | 9 | License as published by the Free Software Foundation; either | ||
9 | version 2 of the License, or (at your option) any later version. | 10 | version 2 of the License, or (at your option) any later version. | ||
10 | 11 | | |||
11 | This library is distributed in the hope that it will be useful, | 12 | This library is distributed in the hope that it will be useful, | ||
12 | but WITHOUT ANY WARRANTY; without even the implied warranty of | 13 | but WITHOUT ANY WARRANTY; without even the implied warranty of | ||
13 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | 14 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | ||
14 | Library General Public License for more details. | 15 | Library General Public License for more details. | ||
15 | 16 | | |||
16 | You should have received a copy of the GNU Library General Public License | 17 | You should have received a copy of the GNU Library General Public License | ||
17 | along with this library; see the file COPYING.LIB. If not, write to | 18 | along with this library; see the file COPYING.LIB. If not, write to | ||
18 | the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, | 19 | the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, | ||
19 | Boston, MA 02110-1301, USA. | 20 | Boston, MA 02110-1301, USA. | ||
20 | */ | 21 | */ | ||
21 | 22 | | |||
22 | #include "astroseasons.h" | 23 | #include "astroseasons.h" | ||
23 | 24 | | |||
25 | #include <cmath> | ||||
26 | #include <numeric> | ||||
27 | | ||||
24 | #include <QDate> | 28 | #include <QDate> | ||
25 | #include <QCoreApplication> | 29 | #include <QCoreApplication> | ||
26 | 30 | | |||
27 | using namespace KHolidays; | 31 | using namespace KHolidays; | ||
28 | 32 | | |||
33 | /* Correction tables and calculation algorithms below are based on the book | ||||
34 | * Astronomical Algorithms by Jean Meeus, (c) 1991 by Willman-Bell, Inc., | ||||
35 | * | ||||
36 | * The correction tables are precise for years -1000 to +3000 but according | ||||
37 | * to the book they can be used for several centuries before and after that | ||||
38 | * period with the error still being quite small. As we only care about date | ||||
39 | * and not the time, the precision of the algorithm is good enough even for | ||||
40 | * greater time span, therefor no checks for the given year are performed | ||||
41 | * anywhere in the code. | ||||
42 | */ | ||||
43 | namespace { | ||||
44 | | ||||
45 | double meanJDE(AstroSeasons::Season season, int year) | ||||
46 | { | ||||
47 | if (year <= 1000) { | ||||
48 | // Astronomical Algorithms, Jean Meeus, chapter 26, table 26.A | ||||
49 | // mean season Julian dates for years -1000 to 1000 | ||||
50 | const double y = year / 1000.0; | ||||
51 | switch (season) { | ||||
52 | case AstroSeasons::MarchEquinox: | ||||
53 | return 1721139.29189 + (365242.13740 * y) + (0.06134 * std::pow(y, 2)) + (0.00111 * std::pow(y, 3)) - (0.00071 * std::pow(y, 4)); | ||||
54 | case AstroSeasons::JuneSolstice: | ||||
55 | return 1721233.25401 + (365241.72562 * y) - (0.05323 * std::pow(y, 2)) + (0.00907 * std::pow(y, 3)) + (0.00025 * std::pow(y, 4)); | ||||
56 | case AstroSeasons::SeptemberEquinox: | ||||
57 | return 1721325.70455 + (365242.49558 * y) - (0.11677 * std::pow(y, 2)) - (0.00297 * std::pow(y, 3)) + (0.00074 * std::pow(y, 4)); | ||||
58 | case AstroSeasons::DecemberSolstice: | ||||
59 | return 1721414.39987 + (365242.88257 * y) - (0.00769 * std::pow(y, 2)) - (0.00933 * std::pow(y, 3)) - (0.00006 * std::pow(y, 4)); | ||||
60 | case AstroSeasons::None: | ||||
61 | Q_ASSERT(false); | ||||
62 | return 0; | ||||
63 | } | ||||
64 | } else { | ||||
65 | // Astronomical Algorithms, Jean Meeus, chapter 26, table 26.B | ||||
66 | // mean season Julian dates for years 1000 to 3000 | ||||
67 | const double y = (year - 2000) / 1000.0; | ||||
68 | switch (season) { | ||||
69 | case AstroSeasons::MarchEquinox: | ||||
70 | return 2451623.80984 + (365242.37404 * y) + (0.05169 * std::pow(y, 2)) - (0.00411 * std::pow(y, 3)) - (0.00057 * std::pow(y, 4)); | ||||
71 | case AstroSeasons::JuneSolstice: | ||||
72 | return 2451716.56767 + (365241.62603 * y) + (0.00325 * std::pow(y, 2)) + (0.00888 * std::pow(y, 3)) - (0.00030 * std::pow(y, 4)); | ||||
73 | case AstroSeasons::SeptemberEquinox: | ||||
74 | return 2451810.21715 + (365242.01767 * y) - (0.11575 * std::pow(y, 2)) + (0.00337 * std::pow(y, 3)) + (0.00078 * std::pow(y, 4)); | ||||
75 | case AstroSeasons::DecemberSolstice: | ||||
76 | return 2451900.05952 + (365242.74049 * y) - (0.06223 * std::pow(y, 2)) - (0.00823 * std::pow(y, 3)) + (0.00032 * std::pow(y, 4)); | ||||
77 | case AstroSeasons::None: | ||||
78 | Q_ASSERT(false); | ||||
79 | return 0; | ||||
80 | } | ||||
81 | } | ||||
82 | | ||||
83 | return 0; | ||||
84 | } | ||||
85 | | ||||
86 | double periodicTerms(double t) | ||||
87 | { | ||||
88 | // Astronomical Algorithms, Jean Meeus, chapter 26, table 26.C | ||||
89 | // The table gives the periodic terms in degrees, but the values are converted to radians | ||||
90 | // at compile time so that they can be passed to std::cos() | ||||
91 | struct Periodic { | ||||
92 | constexpr Periodic(int a, double b_deg, double c_deg) | ||||
93 | : a(a), b_rad(b_deg * (M_PI / 180.0)), c_rad(c_deg * (M_PI / 180.0)) | ||||
94 | {} | ||||
95 | | ||||
96 | int a; | ||||
97 | double b_rad; | ||||
98 | double c_rad; | ||||
99 | } periodic[] = { | ||||
100 | {485, 324.96, 1934.136}, {203, 337.23, 32964.467}, {199, 342.08, 20.186}, {182, 27.85, 445267.112}, | ||||
101 | {156, 73.14, 45036.886}, {136, 171.52, 22518.443}, { 77, 222.54, 65928.934}, { 74, 296.72, 3034.906}, | ||||
102 | { 70, 243.58, 9037.513}, { 58, 119.81, 33718.147}, { 52, 297.17, 150.678}, { 50, 21.02, 2281.226}, | ||||
103 | { 45, 247.54, 29929.562}, { 44, 325.15, 31555.956}, { 29, 60.93, 4443.417}, { 18, 155.12, 67555.328}, | ||||
104 | { 17, 288.79, 4562.452}, { 16, 198.04, 62894.029}, { 14, 199.76, 31436.921}, { 12, 95.39, 14577.848}, | ||||
105 | { 12, 287.11, 31931.756}, { 12, 320.81, 34777.259}, { 9, 227.73, 1222.114}, { 8, 15.45, 16859.074} | ||||
106 | }; | ||||
107 | | ||||
108 | return std::accumulate(std::begin(periodic), std::end(periodic), 0.0, | ||||
109 | [t](double s, const Periodic &p) { return s + p.a * std::cos(p.b_rad + p.c_rad * t); }); | ||||
110 | } | ||||
111 | | ||||
112 | // Returns julian date of given season in given year | ||||
113 | double seasonJD(AstroSeasons::Season season, int year) | ||||
114 | { | ||||
115 | // Astronimical Algorithms, Jean Meeus, chapter 26 | ||||
116 | const auto jde0 = meanJDE(season, year); | ||||
117 | const auto T = (jde0 - 2451545.0) / 36525; | ||||
118 | const auto W_deg = 35999.373 * T + 2.47; | ||||
119 | const auto W_rad = W_deg * (M_PI / 180.0); | ||||
120 | const auto dLambda = 1 + (0.0334 * std::cos(W_rad)) + (0.0007 * std::cos(2 * W_rad)); | ||||
121 | const auto S = periodicTerms(T); | ||||
122 | return jde0 + (0.00001 * S) / dLambda; | ||||
123 | } | ||||
124 | | ||||
125 | } | ||||
126 | | ||||
127 | QDate AstroSeasons::seasonDate(Season season, int year) | ||||
128 | { | ||||
129 | if (season == None) { | ||||
130 | return {}; | ||||
131 | } | ||||
132 | const qint64 jd = round(seasonJD(season, year)); | ||||
133 | return QDate::fromJulianDay(jd); | ||||
134 | } | ||||
135 | | ||||
29 | QString AstroSeasons::seasonNameAtDate(const QDate &date) | 136 | QString AstroSeasons::seasonNameAtDate(const QDate &date) | ||
30 | { | 137 | { | ||
31 | return seasonName(seasonAtDate(date)); | 138 | return seasonName(seasonAtDate(date)); | ||
32 | } | 139 | } | ||
33 | 140 | | |||
34 | QString AstroSeasons::seasonName(AstroSeasons::Season season) | 141 | QString AstroSeasons::seasonName(AstroSeasons::Season season) | ||
35 | { | 142 | { | ||
36 | switch (season) { | 143 | switch (season) { | ||
37 | case JuneSolstice: | 144 | case JuneSolstice: | ||
38 | return QCoreApplication::translate("AstroSeasons", "June Solstice"); | 145 | return QCoreApplication::translate("AstroSeasons", "June Solstice"); | ||
39 | case DecemberSolstice: | 146 | case DecemberSolstice: | ||
40 | return QCoreApplication::translate("AstroSeasons", "December Solstice"); | 147 | return QCoreApplication::translate("AstroSeasons", "December Solstice"); | ||
41 | case MarchEquinox: | 148 | case MarchEquinox: | ||
42 | return QCoreApplication::translate("AstroSeasons", "March Equinox"); | 149 | return QCoreApplication::translate("AstroSeasons", "March Equinox"); | ||
43 | case SeptemberEquinox: | 150 | case SeptemberEquinox: | ||
44 | return QCoreApplication::translate("AstroSeasons", "September Equinox"); | 151 | return QCoreApplication::translate("AstroSeasons", "September Equinox"); | ||
45 | case None: | 152 | case None: | ||
46 | return QString(); | 153 | return QString(); | ||
47 | } | 154 | } | ||
48 | return QString(); | 155 | return QString(); | ||
49 | } | 156 | } | ||
50 | 157 | | |||
51 | AstroSeasons::Season AstroSeasons::seasonAtDate(const QDate &date) | 158 | AstroSeasons::Season AstroSeasons::seasonAtDate(const QDate &date) | ||
52 | { | 159 | { | ||
53 | // see http://www.hermetic.ch/cal_sw/ve/ve.php | 160 | for (auto season : { JuneSolstice, DecemberSolstice, MarchEquinox, SeptemberEquinox }) { | ||
54 | Season retSeason = None; | 161 | if (seasonDate(season, date.year()) == date) { | ||
55 | 162 | return season; | |||
56 | const int year = date.year(); | | |||
57 | //Use dumb method for now | | |||
58 | if (date == QDate(year, 6, 22)) { | | |||
59 | return JuneSolstice; | | |||
60 | } | 163 | } | ||
61 | if (date == QDate(year, 12, 22)) { | | |||
62 | return DecemberSolstice; | | |||
63 | } | 164 | } | ||
64 | if (date == QDate(year, 3, 22)) { | 165 | return None; | ||
65 | return MarchEquinox; | | |||
66 | } | | |||
67 | if (date == QDate(year, 9, 22)) { | | |||
68 | return SeptemberEquinox; | | |||
69 | } | | |||
70 | | ||||
71 | return retSeason; | | |||
72 | } | 166 | } |