Commit | Line | Data |
---|---|---|
9f34a2a0 JB |
1 | ;;; solar.el --- calendar functions for solar events. |
2 | ||
75af4a4a | 3 | ;; Copyright (C) 1992, 1993, 1995 Free Software Foundation, Inc. |
9f34a2a0 JB |
4 | |
5 | ;; Author: Edward M. Reingold <reingold@cs.uiuc.edu> | |
087c56fa | 6 | ;; Denis B. Roegel <Denis.Roegel@loria.fr> |
68e60225 ER |
7 | ;; Keywords: calendar |
8 | ;; Human-Keywords: sunrise, sunset, equinox, solstice, calendar, diary, | |
9 | ;; holidays | |
9f34a2a0 JB |
10 | |
11 | ;; This file is part of GNU Emacs. | |
12 | ||
59243403 RS |
13 | ;; GNU Emacs is free software; you can redistribute it and/or modify |
14 | ;; it under the terms of the GNU General Public License as published by | |
15 | ;; the Free Software Foundation; either version 2, or (at your option) | |
16 | ;; any later version. | |
17 | ||
9f34a2a0 | 18 | ;; GNU Emacs is distributed in the hope that it will be useful, |
59243403 RS |
19 | ;; but WITHOUT ANY WARRANTY; without even the implied warranty of |
20 | ;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the | |
21 | ;; GNU General Public License for more details. | |
22 | ||
23 | ;; You should have received a copy of the GNU General Public License | |
b578f267 EN |
24 | ;; along with GNU Emacs; see the file COPYING. If not, write to the |
25 | ;; Free Software Foundation, Inc., 59 Temple Place - Suite 330, | |
26 | ;; Boston, MA 02111-1307, USA. | |
9f34a2a0 JB |
27 | |
28 | ;;; Commentary: | |
29 | ||
6ff099c3 ER |
30 | ;; This collection of functions implements the features of calendar.el, |
31 | ;; diary.el, and holiday.el that deal with times of day, sunrise/sunset, and | |
ac138dc0 | 32 | ;; equinoxes/solstices. |
9f34a2a0 JB |
33 | |
34 | ;; Based on the ``Almanac for Computers 1984,'' prepared by the Nautical | |
75af4a4a ER |
35 | ;; Almanac Office, United States Naval Observatory, Washington, 1984, on |
36 | ;; ``Astronomical Formulae for Calculators,'' 3rd ed., by Jean Meeus, | |
087c56fa ER |
37 | ;; Willmann-Bell, Inc., 1985, on ``Astronomical Algorithms'' by Jean Meeus, |
38 | ;; Willmann-Bell, Inc., 1991, and on ``Planetary Programs and Tables from | |
39 | ;; -4000 to +2800'' by Pierre Bretagnon and Jean-Louis Simon, Willmann-Bell, | |
40 | ;; Inc., 1986. | |
75af4a4a | 41 | |
9f34a2a0 | 42 | ;; |
087c56fa ER |
43 | ;; Accuracy: |
44 | ;; 1. Sunrise/sunset times will be accurate to the minute for years | |
45 | ;; 1951--2050. For other years the times will be within +/- 2 minutes. | |
9f34a2a0 | 46 | ;; |
087c56fa ER |
47 | ;; 2. Equinox/solstice times will be accurate to the minute for years |
48 | ;; 1951--2050. For other years the times will be within +/- 1 minute. | |
9f34a2a0 JB |
49 | |
50 | ;; Comments, corrections, and improvements should be sent to | |
51 | ;; Edward M. Reingold Department of Computer Science | |
52 | ;; (217) 333-6733 University of Illinois at Urbana-Champaign | |
53 | ;; reingold@cs.uiuc.edu 1304 West Springfield Avenue | |
54 | ;; Urbana, Illinois 61801 | |
55 | ||
56 | ;;; Code: | |
57 | ||
58 | (if (fboundp 'atan) | |
59 | (require 'lisp-float-type) | |
6ff099c3 | 60 | (error "Solar/lunar calculations impossible since floating point is unavailable.")) |
9f34a2a0 | 61 | |
a92ade89 | 62 | (require 'cal-dst) |
087c56fa | 63 | (require 'cal-julian) |
a92ade89 JB |
64 | |
65 | ;;;###autoload | |
66 | (defvar calendar-time-display-form | |
67 | '(12-hours ":" minutes am-pm | |
68 | (if time-zone " (") time-zone (if time-zone ")")) | |
69 | "*The pseudo-pattern that governs the way a time of day is formatted. | |
70 | ||
71 | A pseudo-pattern is a list of expressions that can involve the keywords | |
72 | `12-hours', `24-hours', and `minutes', all numbers in string form, | |
73 | and `am-pm' and `time-zone', both alphabetic strings. | |
74 | ||
75 | For example, the form | |
76 | ||
77 | '(24-hours \":\" minutes | |
78 | (if time-zone \" (\") time-zone (if time-zone \")\")) | |
79 | ||
80 | would give military-style times like `21:07 (UTC)'.") | |
81 | ||
82 | ;;;###autoload | |
83 | (defvar calendar-latitude nil | |
6ff099c3 ER |
84 | "*Latitude of `calendar-location-name' in degrees. |
85 | ||
86 | The value can be either a decimal fraction (one place of accuracy is | |
87 | sufficient), + north, - south, such as 40.7 for New York City, or the value | |
88 | can be a vector [degrees minutes north/south] such as [40 50 north] for New | |
89 | York City. | |
90 | ||
75af4a4a | 91 | This variable should be set in site-local.el.") |
a92ade89 JB |
92 | |
93 | ;;;###autoload | |
94 | (defvar calendar-longitude nil | |
6ff099c3 ER |
95 | "*Longitude of `calendar-location-name' in degrees. |
96 | ||
97 | The value can be either a decimal fraction (one place of accuracy is | |
98 | sufficient), + east, - west, such as -73.9 for New York City, or the value | |
99 | can be a vector [degrees minutes east/west] such as [73 55 west] for New | |
100 | York City. | |
101 | ||
087c56fa | 102 | This variable should be set in site-start.el.") |
6ff099c3 ER |
103 | |
104 | (defsubst calendar-latitude () | |
105 | "Convert calendar-latitude to a signed decimal fraction, if needed." | |
106 | (if (numberp calendar-latitude) | |
107 | calendar-latitude | |
108 | (let ((lat (+ (aref calendar-latitude 0) | |
109 | (/ (aref calendar-latitude 1) 60.0)))) | |
110 | (if (equal (aref calendar-latitude 2) 'north) | |
111 | lat | |
112 | (- lat))))) | |
113 | ||
114 | (defsubst calendar-longitude () | |
115 | "Convert calendar-longitude to a signed decimal fraction, if needed." | |
116 | (if (numberp calendar-longitude) | |
117 | calendar-longitude | |
118 | (let ((long (+ (aref calendar-longitude 0) | |
119 | (/ (aref calendar-longitude 1) 60.0)))) | |
120 | (if (equal (aref calendar-longitude 2) 'east) | |
121 | long | |
122 | (- long))))) | |
a92ade89 JB |
123 | |
124 | ;;;###autoload | |
125 | (defvar calendar-location-name | |
126 | '(let ((float-output-format "%.1f")) | |
127 | (format "%s%s, %s%s" | |
6ff099c3 ER |
128 | (if (numberp calendar-latitude) |
129 | (abs calendar-latitude) | |
130 | (+ (aref calendar-latitude 0) | |
131 | (/ (aref calendar-latitude 1) 60.0))) | |
132 | (if (numberp calendar-latitude) | |
133 | (if (> calendar-latitude 0) "N" "S") | |
134 | (if (equal (aref calendar-latitude 2) 'north) "N" "S")) | |
135 | (if (numberp calendar-longitude) | |
136 | (abs calendar-longitude) | |
137 | (+ (aref calendar-longitude 0) | |
138 | (/ (aref calendar-longitude 1) 60.0))) | |
139 | (if (numberp calendar-longitude) | |
140 | (if (> calendar-longitude 0) "E" "W") | |
562a94a0 | 141 | (if (equal (aref calendar-longitude 2) 'east) "E" "W")))) |
a92ade89 | 142 | "*Expression evaluating to name of `calendar-longitude', calendar-latitude'. |
6ff099c3 ER |
143 | For example, \"New York City\". Default value is just the latitude, longitude |
144 | pair. | |
145 | ||
087c56fa ER |
146 | This variable should be set in site-start.el.") |
147 | ||
148 | (defvar solar-error 0.5 | |
149 | "*Tolerance (in minutes) for sunrise/sunset calculations. | |
150 | ||
151 | A larger value makes the calculations for sunrise/sunset faster, but less | |
152 | accurate. The default is half a minute (30 seconds), so that sunrise/sunset | |
153 | times will be correct to the minute. | |
154 | ||
155 | It is useless to set the value smaller than 4*delta, where delta is the | |
156 | accuracy in the longitude of the sun (given by the function | |
157 | `solar-ecliptic-coordinates') in degrees since (delta/360) x (86400/60) = 4 x | |
158 | delta. At present, delta = 0.01 degrees, so the value of the variable | |
159 | `solar-error' should be at least 0.04 minutes (about 2.5 seconds).") | |
9f34a2a0 | 160 | |
3d9dece2 | 161 | (defvar solar-n-hemi-seasons |
fc68b552 BF |
162 | '("Vernal Equinox" "Summer Solstice" "Autumnal Equinox" "Winter Solstice") |
163 | "List of season changes for the northern hemisphere.") | |
164 | ||
3d9dece2 | 165 | (defvar solar-s-hemi-seasons |
fc68b552 BF |
166 | '("Autumnal Equinox" "Winter Solstice" "Vernal Equinox" "Summer Solstice") |
167 | "List of season changes for the southern hemisphere.") | |
168 | ||
087c56fa ER |
169 | (defvar solar-sidereal-time-greenwich-midnight |
170 | nil | |
171 | "Sidereal time at Greenwich at midnight (universal time).") | |
172 | ||
173 | (defvar solar-spring-or-summer-season nil | |
174 | "T if spring or summer and nil otherwise. | |
175 | Needed for polar areas, in order to know whether the day lasts 0 or 24 hours.") | |
176 | ||
9f34a2a0 JB |
177 | (defun solar-setup () |
178 | "Prompt user for latitude, longitude, and time zone." | |
179 | (beep) | |
180 | (if (not calendar-longitude) | |
181 | (setq calendar-longitude | |
182 | (solar-get-number | |
183 | "Enter longitude (decimal fraction; + east, - west): "))) | |
184 | (if (not calendar-latitude) | |
185 | (setq calendar-latitude | |
186 | (solar-get-number | |
187 | "Enter latitude (decimal fraction; + north, - south): "))) | |
188 | (if (not calendar-time-zone) | |
189 | (setq calendar-time-zone | |
190 | (solar-get-number | |
e2fe2f52 | 191 | "Enter difference from Coordinated Universal Time (in minutes): ")))) |
9f34a2a0 JB |
192 | |
193 | (defun solar-get-number (prompt) | |
194 | "Return a number from the minibuffer, prompting with PROMPT. | |
195 | Returns nil if nothing was entered." | |
196 | (let ((x (read-string prompt ""))) | |
197 | (if (not (string-equal x "")) | |
198 | (string-to-int x)))) | |
199 | ||
087c56fa ER |
200 | ;; The condition-case stuff is needed to catch bogus arithmetic |
201 | ;; exceptions that occur on some machines (like Sparcs) | |
202 | (defun solar-sin-degrees (x) | |
203 | (condition-case nil | |
204 | (sin (degrees-to-radians (mod x 360.0))) | |
205 | (solar-sin-degrees x))) | |
206 | (defun solar-cosine-degrees (x) | |
207 | (condition-case nil | |
208 | (cos (degrees-to-radians (mod x 360.0))) | |
209 | (solar-cosine-degrees x))) | |
210 | (defun solar-tangent-degrees (x) | |
211 | (condition-case nil | |
212 | (tan (degrees-to-radians (mod x 360.0))) | |
213 | (solar-tangent-degrees x))) | |
214 | ||
9f34a2a0 JB |
215 | (defun solar-xy-to-quadrant (x y) |
216 | "Determines the quadrant of the point X, Y." | |
217 | (if (> x 0) | |
218 | (if (> y 0) 1 4) | |
219 | (if (> y 0) 2 3))) | |
220 | ||
221 | (defun solar-degrees-to-quadrant (angle) | |
222 | "Determines the quadrant of ANGLE." | |
3a2e3ab5 | 223 | (1+ (floor (mod angle 360) 90))) |
9f34a2a0 JB |
224 | |
225 | (defun solar-arctan (x quad) | |
226 | "Arctangent of X in quadrant QUAD." | |
227 | (let ((deg (radians-to-degrees (atan x)))) | |
228 | (cond ((equal quad 2) (+ deg 180)) | |
229 | ((equal quad 3) (+ deg 180)) | |
230 | ((equal quad 4) (+ deg 360)) | |
231 | (t deg)))) | |
232 | ||
087c56fa ER |
233 | (defun solar-atn2 (x y) |
234 | "Arctan of point X, Y." | |
235 | (if (= y 0) | |
236 | (if (> x 0) 90 270) | |
237 | (solar-arctan (/ x y) y))) | |
238 | ||
9f34a2a0 | 239 | (defun solar-arccos (x) |
087c56fa | 240 | "Arcos of X." |
9f34a2a0 JB |
241 | (let ((y (sqrt (- 1 (* x x))))) |
242 | (solar-arctan (/ y x) (solar-xy-to-quadrant x y)))) | |
243 | ||
244 | (defun solar-arcsin (y) | |
087c56fa | 245 | "Arcsin of Y." |
9f34a2a0 JB |
246 | (let ((x (sqrt (- 1 (* y y))))) |
247 | (solar-arctan (/ y x) (solar-xy-to-quadrant x y)))) | |
248 | ||
087c56fa ER |
249 | (defsubst solar-degrees-to-hours (degrees) |
250 | "Convert DEGREES to hours." | |
251 | (/ degrees 15.0)) | |
9f34a2a0 | 252 | |
6ff099c3 | 253 | (defsubst solar-hours-to-days (hour) |
087c56fa | 254 | "Convert HOUR to decimal fraction of a day." |
6ff099c3 | 255 | (/ hour 24.0)) |
9f34a2a0 | 256 | |
087c56fa ER |
257 | (defun solar-right-ascension (longitude obliquity) |
258 | "Right ascension of the sun, in hours, given LONGITUDE and OBLIQUITY. | |
259 | Both arguments are in degrees." | |
9f34a2a0 JB |
260 | (solar-degrees-to-hours |
261 | (solar-arctan | |
087c56fa | 262 | (* (solar-cosine-degrees obliquity) (solar-tangent-degrees longitude)) |
9f34a2a0 JB |
263 | (solar-degrees-to-quadrant longitude)))) |
264 | ||
087c56fa ER |
265 | (defun solar-declination (longitude obliquity) |
266 | "Declination of the sun, in degrees, given LONGITUDE and OBLIQUITY. | |
267 | Both arguments are in degrees." | |
9f34a2a0 | 268 | (solar-arcsin |
087c56fa | 269 | (* (solar-sin-degrees obliquity) |
9f34a2a0 JB |
270 | (solar-sin-degrees longitude)))) |
271 | ||
087c56fa ER |
272 | (defun solar-sunrise-and-sunset (time latitude longitude) |
273 | "Sunrise, sunset and length of day. | |
274 | Parameters are the midday TIME and the LATITUDE, LONGITUDE of the location. | |
275 | ||
276 | TIME is a pair with the first component being the number of Julian centuries | |
277 | elapsed at 0 Universal Time, and the second component being the universal | |
278 | time. For instance, the pair corresponding to November 28, 1995 at 16 UT is | |
279 | (-0.040945 16), -0.040945 being the number of julian centuries elapsed between | |
280 | Jan 1, 2000 at 12 UT and November 28, 1995 at 0 UT. | |
281 | ||
282 | Coordinates are included because this function is called with latitude=10 | |
283 | degrees to find out if polar regions have 24 hours of sun or only night." | |
284 | (let* ((rise-time (solar-moment -1 latitude longitude time)) | |
285 | (set-time (solar-moment 1 latitude longitude time)) | |
286 | (day-length)) | |
287 | (if (not (and rise-time set-time)) | |
288 | (if (or (and (> latitude 0) solar-spring-or-summer-season) | |
289 | (and (< latitude 0) (not solar-spring-or-summer-season))) | |
290 | (setq day-length 24) | |
291 | (setq day-length 0)) | |
292 | (setq day-length (- set-time rise-time))) | |
293 | (list (+ rise-time (/ calendar-time-zone 60.0)) | |
294 | (+ set-time (/ calendar-time-zone 60.0)) day-length))) | |
295 | ||
296 | (defun solar-moment (direction latitude longitude time) | |
297 | "Sunrise/sunset at location. | |
298 | Sunrise if DIRECTION =-1 or sunset if =1 at LATITUDE, LONGITUDE, with midday | |
299 | being TIME. | |
300 | ||
301 | TIME is a pair with the first component being the number of Julian centuries | |
302 | elapsed at 0 Universal Time, and the second component being the universal | |
303 | time. For instance, the pair corresponding to November 28, 1995 at 16 UT is | |
304 | (-0.040945 16), -0.040945 being the number of julian centuries elapsed between | |
305 | Jan 1, 2000 at 12 UT and November 28, 1995 at 0 UT. | |
306 | ||
307 | Uses binary search." | |
308 | (let* ((ut (car (cdr time))) | |
309 | (possible 1) ; we assume that rise or set are possible | |
310 | (utmin (+ ut (* direction 12.0))) | |
311 | (utmax ut) ; the time searched is between utmin and utmax | |
312 | ; utmin and utmax are in hours | |
313 | (utmoment-old 0.0) ; rise or set approximation | |
314 | (utmoment 1.0) ; rise or set approximation | |
315 | (hut 0) ; sun height at utmoment | |
316 | (t0 (car time)) | |
317 | (hmin (car (cdr | |
318 | (solar-horizontal-coordinates (list t0 utmin) | |
319 | latitude longitude t)))) | |
320 | (hmax (car (cdr | |
321 | (solar-horizontal-coordinates (list t0 utmax) | |
322 | latitude longitude t))))) | |
323 | ; -0.61 degrees is the height of the middle of the sun, when it rises | |
324 | ; or sets. | |
325 | (if (< hmin -0.61) | |
326 | (if (> hmax -0.61) | |
327 | (while ;(< i 20) ; we perform a simple dichotomy | |
328 | ; (> (abs (+ hut 0.61)) epsilon) | |
329 | (>= (abs (- utmoment utmoment-old)) | |
330 | (/ solar-error 60)) | |
331 | (setq utmoment-old utmoment) | |
332 | (setq utmoment (/ (+ utmin utmax) 2)) | |
333 | (setq hut (car (cdr | |
334 | (solar-horizontal-coordinates | |
335 | (list t0 utmoment) latitude longitude t)))) | |
336 | (if (< hut -0.61) (setq utmin utmoment)) | |
337 | (if (> hut -0.61) (setq utmax utmoment)) | |
338 | ) | |
339 | (setq possible 0)) the sun never rises | |
340 | (setq possible 0)) ; the sun never sets | |
341 | (if (equal possible 0) nil utmoment))) | |
9f34a2a0 | 342 | |
c68488d2 | 343 | (defun solar-time-string (time time-zone) |
75af4a4a | 344 | "Printable form for decimal fraction TIME in TIME-ZONE. |
c68488d2 ER |
345 | Format used is given by `calendar-time-display-form'." |
346 | (let* ((time (round (* 60 time))) | |
a92ade89 JB |
347 | (24-hours (/ time 60)) |
348 | (minutes (format "%02d" (% time 60))) | |
349 | (12-hours (format "%d" (1+ (% (+ 24-hours 11) 12)))) | |
9f34a2a0 | 350 | (am-pm (if (>= 24-hours 12) "pm" "am")) |
9f34a2a0 JB |
351 | (24-hours (format "%02d" 24-hours))) |
352 | (mapconcat 'eval calendar-time-display-form ""))) | |
353 | ||
087c56fa ER |
354 | |
355 | (defun solar-daylight (time) | |
356 | "Printable form for time expressed in hours." | |
357 | (format "%d:%02d" | |
358 | (floor time) | |
359 | (floor (* 60 (- time (floor time)))))) | |
360 | ||
361 | (defun solar-exact-local-noon (date) | |
362 | "Date and Universal Time of local noon at *local date* date. | |
363 | ||
364 | The date may be different from the one asked for, but it will be the right | |
365 | local date. The second component of date should be an integer." | |
366 | (let* ((nd date) | |
367 | (ut (- 12.0 (/ calendar-longitude 15))) | |
368 | (te (solar-time-equation date ut))) | |
369 | (setq ut (- ut te)) | |
370 | (if (>= ut 24) | |
371 | (progn | |
372 | (setq nd (list (car date) (+ 1 (car (cdr date))) | |
373 | (car (cdr (cdr date))))) | |
374 | (setq ut (- ut 24)))) | |
375 | (if (< ut 0) | |
376 | (progn | |
377 | (setq nd (list (car date) (- (car (cdr date)) 1) | |
378 | (car (cdr (cdr date))))) | |
379 | (setq ut (+ ut 24)))) | |
380 | (setq nd (calendar-gregorian-from-absolute | |
381 | (calendar-absolute-from-gregorian nd))) | |
382 | ; date standardization | |
383 | (list nd ut))) | |
384 | ||
9f34a2a0 | 385 | (defun solar-sunrise-sunset (date) |
087c56fa ER |
386 | "List of *local* times of sunrise, sunset, and daylight on Gregorian DATE. |
387 | ||
388 | Corresponding value is nil if there is no sunrise/sunset." | |
389 | (let* (; first, get the exact moment of local noon. | |
390 | (exact-local-noon (solar-exact-local-noon date)) | |
391 | ; get the the time from the 2000 epoch. | |
392 | (t0 (solar-julian-ut-centuries (car exact-local-noon))) | |
393 | ; store the sidereal time at Greenwich at midnight of UT time. | |
394 | ; find if summer or winter slightly above the equator | |
395 | (equator-rise-set | |
396 | (progn (setq solar-sidereal-time-greenwich-midnight | |
397 | (solar-sidereal-time t0)) | |
398 | (solar-sunrise-and-sunset | |
399 | (list t0 (car (cdr exact-local-noon))) | |
400 | 10.0 | |
401 | calendar-longitude))) | |
402 | ; store the spring/summer information, | |
403 | ; compute sunrise and sunset (two first components of rise-set). | |
404 | ; length of day is the third component (it is only the difference | |
405 | ; between sunset and sunrise when there is a sunset and a sunrise) | |
406 | (rise-set | |
407 | (progn | |
408 | (setq solar-spring-or-summer-season | |
409 | (if (> (car (cdr (cdr equator-rise-set))) 12) 1 0)) | |
410 | (solar-sunrise-and-sunset | |
411 | (list t0 (car (cdr exact-local-noon))) | |
412 | calendar-latitude | |
413 | calendar-longitude))) | |
414 | (rise (car rise-set)) | |
415 | (adj-rise (if rise (dst-adjust-time date rise) nil)) | |
416 | (set (car (cdr rise-set))) | |
417 | (adj-set (if set (dst-adjust-time date set) nil)) | |
418 | (length (car (cdr (cdr rise-set)))) ) | |
419 | (list | |
420 | (and rise (calendar-date-equal date (car adj-rise)) (cdr adj-rise)) | |
421 | (and set (calendar-date-equal date (car adj-set)) (cdr adj-set)) | |
422 | (solar-daylight length)))) | |
423 | ||
424 | (defun solar-sunrise-sunset-string (date) | |
425 | "String of *local* times of sunrise, sunset, and daylight on Gregorian DATE." | |
426 | (let ((l (solar-sunrise-sunset date))) | |
427 | (format | |
428 | "%s, %s at %s (%s hours daylight)" | |
429 | (if (car l) | |
430 | (concat "Sunrise " (apply 'solar-time-string (car l))) | |
431 | "No sunrise") | |
432 | (if (car (cdr l)) | |
433 | (concat "sunset " (apply 'solar-time-string (car (cdr l)))) | |
434 | "no sunset") | |
435 | (eval calendar-location-name) | |
436 | (car (cdr (cdr l)))))) | |
437 | ||
438 | (defun solar-julian-ut-centuries (date) | |
439 | "Number of Julian centuries elapsed since 1 Jan, 2000 at noon U.T. for Gregorian DATE." | |
440 | (/ (- (calendar-absolute-from-gregorian date) | |
441 | (calendar-absolute-from-gregorian '(1 1.5 2000))) | |
442 | 36525.0)) | |
443 | ||
444 | (defun solar-ephemeris-time(time) | |
445 | "Ephemeris Time at moment TIME. | |
446 | ||
447 | TIME is a pair with the first component being the number of Julian centuries | |
448 | elapsed at 0 Universal Time, and the second component being the universal | |
449 | time. For instance, the pair corresponding to November 28, 1995 at 16 UT is | |
450 | (-0.040945 16), -0.040945 being the number of julian centuries elapsed between | |
451 | Jan 1, 2000 at 12 UT and November 28, 1995 at 0 UT. | |
452 | ||
453 | Result is in julian centuries of ephemeris time." | |
454 | (let* ((t0 (car time)) | |
455 | (ut (car (cdr time))) | |
456 | (t1 (+ t0 (/ (/ ut 24.0) 36525))) | |
457 | (y (+ 2000 (* 100 t1))) | |
458 | (dt (* 86400 (solar-ephemeris-correction (floor y))))) | |
459 | (+ t1 (/ (/ dt 86400) 36525)))) | |
9f34a2a0 | 460 | |
75af4a4a ER |
461 | (defun solar-date-next-longitude (d l) |
462 | "First moment on or after Julian day number D when sun's longitude is a | |
463 | multiple of L degrees at calendar-location-name with that location's | |
464 | local time (including any daylight savings rules). | |
465 | ||
466 | L must be an integer divisor of 360. | |
467 | ||
468 | Result is in local time expressed astronomical (Julian) day numbers. | |
469 | ||
470 | The values of calendar-daylight-savings-starts, | |
471 | calendar-daylight-savings-starts-time, calendar-daylight-savings-ends, | |
472 | calendar-daylight-savings-ends-time, calendar-daylight-time-offset, and | |
473 | calendar-time-zone are used to interpret local time." | |
474 | (let* ((long) | |
475 | (start d) | |
476 | (start-long (solar-longitude d)) | |
477 | (next (mod (* l (1+ (floor (/ start-long l)))) 360)) | |
478 | (end (+ d (* (/ l 360.0) 400))) | |
479 | (end-long (solar-longitude end))) | |
480 | (while ;; bisection search for nearest minute | |
481 | (< 0.00001 (- end start)) | |
482 | ;; start <= d < end | |
483 | ;; start-long <= next < end-long when next != 0 | |
484 | ;; when next = 0, we look for the discontinuity (start-long is near 360 | |
485 | ;; and end-long is small (less than l). | |
486 | (setq d (/ (+ start end) 2.0)) | |
487 | (setq long (solar-longitude d)) | |
488 | (if (or (and (/= next 0) (< long next)) | |
489 | (and (= next 0) (< l long))) | |
490 | (progn | |
491 | (setq start d) | |
492 | (setq start-long long)) | |
493 | (setq end d) | |
494 | (setq end-long long))) | |
495 | (/ (+ start end) 2.0))) | |
496 | ||
087c56fa ER |
497 | (defun solar-horizontal-coordinates |
498 | (time latitude longitude for-sunrise-sunset) | |
499 | "Azimuth and height of the sun at TIME, LATITUDE, and LONGITUDE. | |
500 | ||
501 | TIME is a pair with the first component being the number of Julian centuries | |
502 | elapsed at 0 Universal Time, and the second component being the universal | |
503 | time. For instance, the pair corresponding to November 28, 1995 at 16 UT is | |
504 | (-0.040945 16), -0.040945 being the number of julian centuries elapsed between | |
505 | Jan 1, 2000 at 12 UT and November 28, 1995 at 0 UT. | |
506 | ||
507 | The azimuth is given in degrees as well as the height (between -180 and 180)." | |
508 | (let* ((ut (car (cdr time))) | |
509 | (ec (solar-equatorial-coordinates time for-sunrise-sunset)) | |
510 | (st (+ solar-sidereal-time-greenwich-midnight | |
511 | (* ut 1.00273790935))) | |
512 | (ah (- (* st 15) (* 15 (car ec)) (* -1 calendar-longitude))) | |
513 | ; hour angle (in degrees) | |
514 | (de (car (cdr ec))) | |
515 | (azimuth (solar-atn2 (solar-sin-degrees ah) | |
516 | (- (* (solar-cosine-degrees ah) | |
517 | (solar-sin-degrees latitude)) | |
518 | (* (solar-tangent-degrees de) | |
519 | (solar-cosine-degrees latitude))))) | |
520 | (height (solar-arcsin | |
521 | (+ (* (solar-sin-degrees latitude) (solar-sin-degrees de)) | |
522 | (* (solar-cosine-degrees latitude) | |
523 | (solar-cosine-degrees de) | |
524 | (solar-cosine-degrees ah)))))) | |
525 | (if (> height 180) (setq height (- height 360))) | |
526 | (list azimuth height))) | |
527 | ||
528 | (defun solar-equatorial-coordinates (time for-sunrise-sunset) | |
529 | "Right ascension (in hours) and declination (in degrees) of the sun at TIME. | |
530 | ||
531 | TIME is a pair with the first component being the number of Julian centuries | |
532 | elapsed at 0 Universal Time, and the second component being the universal | |
533 | time. For instance, the pair corresponding to November 28, 1995 at 16 UT is | |
534 | (-0.040945 16), -0.040945 being the number of julian centuries elapsed between | |
535 | Jan 1, 2000 at 12 UT and November 28, 1995 at 0 UT." | |
536 | (let* ((tm (solar-ephemeris-time time)) | |
537 | (ec (solar-ecliptic-coordinates tm for-sunrise-sunset))) | |
538 | (list (solar-right-ascension (car ec) (car (cdr ec))) | |
539 | (solar-declination (car ec) (car (cdr ec)))))) | |
540 | ||
541 | (defun solar-ecliptic-coordinates (time for-sunrise-sunset) | |
542 | "Apparent longitude of the sun, ecliptic inclination, (both in degrees) | |
543 | equation of time (in hours) and nutation in longitude (in seconds) | |
544 | at moment `time', expressed in julian centuries of Ephemeris Time | |
545 | since January 1st, 2000, at 12 ET." | |
546 | (let* ((l (+ 280.46645 | |
547 | (* 36000.76983 time) | |
548 | (* 0.0003032 time time))) ; sun mean longitude | |
549 | (ml (+ 218.3165 | |
550 | (* 481267.8813 time))) ; moon mean longitude | |
551 | (m (+ 357.52910 | |
552 | (* 35999.05030 time) | |
553 | (* -0.0001559 time time) | |
554 | (* -0.00000048 time time time))) ; sun mean anomaly | |
555 | (i (+ 23.43929111 (* -0.013004167 time) | |
556 | (* -0.00000016389 time time) | |
557 | (* 0.0000005036 time time time))); mean inclination | |
558 | (c (+ (* (+ 1.914600 | |
559 | (* -0.004817 time) | |
560 | (* -0.000014 time time)) | |
561 | (solar-sin-degrees m)) | |
562 | (* (+ 0.019993 (* -0.000101 time)) | |
563 | (solar-sin-degrees (* 2 m))) | |
564 | (* 0.000290 | |
565 | (solar-sin-degrees (* 3 m))))) ; center equation | |
566 | (L (+ l c)) ; total longitude | |
567 | (omega (+ 125.04 | |
568 | (* -1934.136 time))) ; longitude of moon's ascending node | |
569 | ; on the ecliptic | |
570 | (nut (if (not for-sunrise-sunset) | |
571 | (+ (* -17.20 (solar-sin-degrees omega)) | |
572 | (* -1.32 (solar-sin-degrees (* 2 l))) | |
573 | (* -0.23 (solar-sin-degrees (* 2 ml))) | |
574 | (* 0.21 (solar-sin-degrees (* 2 omega)))) | |
575 | nil)) | |
576 | ; nut = nutation in longitude, measured in seconds of angle. | |
577 | (ecc (if (not for-sunrise-sunset) | |
578 | (+ 0.016708617 | |
579 | (* -0.000042037 time) | |
580 | (* -0.0000001236 time time)) ; eccentricity of earth's orbit | |
581 | nil)) | |
582 | (app (+ L | |
583 | -0.00569 | |
584 | (* -0.00478 | |
585 | (solar-sin-degrees omega)))) ; apparent longitude of sun | |
586 | (y (if (not for-sunrise-sunset) | |
587 | (* (solar-tangent-degrees (/ i 2)) | |
588 | (solar-tangent-degrees (/ i 2))) | |
589 | nil)) | |
590 | (time-eq (if (not for-sunrise-sunset) | |
591 | (/ (* 12 (+ (* y (solar-sin-degrees (* 2 l))) | |
592 | (* -2 ecc (solar-sin-degrees m)) | |
593 | (* 4 ecc y (solar-sin-degrees m) | |
594 | (solar-cosine-degrees (* 2 l))) | |
595 | (* -0.5 y y (solar-sin-degrees (* 4 l))) | |
596 | (* -1.25 ecc ecc (solar-sin-degrees (* 2 m))))) | |
597 | 3.1415926535) | |
598 | nil))) | |
599 | ; equation of time, in hours | |
600 | (list app i time-eq nut))) | |
601 | ||
75af4a4a ER |
602 | (defun solar-longitude (d) |
603 | "Longitude of sun on astronomical (Julian) day number D. | |
087c56fa | 604 | Accurary is about 0.0006 degree (about 365.25*24*60*0.0006/360 = 1 minutes). |
75af4a4a ER |
605 | |
606 | The values of calendar-daylight-savings-starts, | |
607 | calendar-daylight-savings-starts-time, calendar-daylight-savings-ends, | |
608 | calendar-daylight-savings-ends-time, calendar-daylight-time-offset, and | |
609 | calendar-time-zone are used to interpret local time." | |
610 | (let* ((a-d (calendar-absolute-from-astro d)) | |
75af4a4a | 611 | ;; get Universal Time |
087c56fa ER |
612 | (date (calendar-astro-from-absolute |
613 | (- a-d | |
614 | (if (dst-in-effect a-d) | |
615 | (/ calendar-daylight-time-offset 24.0 60.0) 0) | |
616 | (/ calendar-time-zone 60.0 24.0)))) | |
75af4a4a ER |
617 | ;; get Ephemeris Time |
618 | (date (+ date (solar-ephemeris-correction | |
619 | (extract-calendar-year | |
620 | (calendar-gregorian-from-absolute | |
621 | (floor | |
622 | (calendar-absolute-from-astro | |
623 | date))))))) | |
087c56fa ER |
624 | (U (/ (- date 2451545) 3652500)) |
625 | (longitude | |
626 | (+ 4.9353929 | |
627 | (* 62833.1961680 U) | |
628 | (* 0.0000001 | |
629 | (apply '+ | |
630 | (mapcar '(lambda (x) | |
631 | (* (car x) | |
632 | (sin (mod | |
633 | (+ (car (cdr x)) | |
634 | (* (car (cdr (cdr x))) U)) | |
635 | (* 2 pi))))) | |
636 | solar-data-list))))) | |
637 | (aberration | |
638 | (* 0.0000001 (- (* 17 (cos (+ 3.10 (* 62830.14 U)))) 973))) | |
639 | (A1 (mod (+ 2.18 (* U (+ -3375.70 (* 0.36 U)))) (* 2 pi))) | |
640 | (A2 (mod (+ 3.51 (* U (+ 125666.39 (* 0.10 U)))) (* 2 pi))) | |
641 | (nutation (* -0.0000001 (+ (* 834 (sin A1)) (* 64 (sin A2)))))) | |
642 | (mod (radians-to-degrees (+ longitude aberration nutation)) 360.0))) | |
643 | ||
644 | (defconst solar-data-list | |
645 | '((403406 4.721964 1.621043) | |
646 | (195207 5.937458 62830.348067) | |
647 | (119433 1.115589 62830.821524) | |
648 | (112392 5.781616 62829.634302) | |
649 | (3891 5.5474 125660.5691) | |
650 | (2819 1.5120 125660.984) | |
651 | (1721 4.1897 62832.4766) | |
652 | (0 1.163 0.813) | |
653 | (660 5.415 125659.31) | |
654 | (350 4.315 57533.85) | |
655 | (334 4.553 -33.931) | |
656 | (314 5.198 777137.715) | |
657 | (268 5.989 78604.191) | |
658 | (242 2.911 5.412) | |
659 | (234 1.423 39302.098) | |
660 | (158 0.061 -34.861) | |
661 | (132 2.317 115067.698) | |
662 | (129 3.193 15774.337) | |
663 | (114 2.828 5296.670) | |
664 | (99 0.52 58849.27) | |
665 | (93 4.65 5296.11) | |
666 | (86 4.35 -3980.70) | |
667 | (78 2.75 52237.69) | |
668 | (72 4.50 55076.47) | |
669 | (68 3.23 261.08) | |
670 | (64 1.22 15773.85) | |
671 | (46 0.14 188491.03) | |
672 | (38 3.44 -7756.55) | |
673 | (37 4.37 264.89) | |
674 | (32 1.14 117906.27) | |
675 | (29 2.84 55075.75) | |
676 | (28 5.96 -7961.39) | |
677 | (27 5.09 188489.81) | |
678 | (27 1.72 2132.19) | |
679 | (25 2.56 109771.03) | |
680 | (24 1.92 54868.56) | |
681 | (21 0.09 25443.93) | |
682 | (21 5.98 -55731.43) | |
683 | (20 4.03 60697.74) | |
684 | (18 4.47 2132.79) | |
685 | (17 0.79 109771.63) | |
686 | (14 4.24 -7752.82) | |
687 | (13 2.01 188491.91) | |
688 | (13 2.65 207.81) | |
689 | (13 4.98 29424.63) | |
690 | (12 0.93 -7.99) | |
691 | (10 2.21 46941.14) | |
692 | (10 3.59 -68.29) | |
693 | (10 1.50 21463.25) | |
694 | (10 2.55 157208.40))) | |
9f34a2a0 JB |
695 | |
696 | (defun solar-ephemeris-correction (year) | |
087c56fa ER |
697 | "Ephemeris time minus Universal Time during Gregorian year. |
698 | Result is in days. | |
699 | ||
700 | For the years 1800-1987, the maximum error is 1.9 seconds. | |
75af4a4a ER |
701 | For the other years, the maximum error is about 30 seconds." |
702 | (cond ((and (<= 1988 year) (< year 2020)) | |
703 | (/ (+ year -2000 67.0) 60.0 60.0 24.0)) | |
704 | ((and (<= 1900 year) (< year 1988)) | |
705 | (let* ((theta (/ (- (calendar-astro-from-absolute | |
706 | (calendar-absolute-from-gregorian | |
707 | (list 7 1 year))) | |
708 | (calendar-astro-from-absolute | |
709 | (calendar-absolute-from-gregorian | |
710 | '(1 1 1900)))) | |
711 | 36525.0)) | |
712 | (theta2 (* theta theta)) | |
713 | (theta3 (* theta2 theta)) | |
714 | (theta4 (* theta2 theta2)) | |
715 | (theta5 (* theta3 theta2))) | |
716 | (+ -0.00002 | |
717 | (* 0.000297 theta) | |
718 | (* 0.025184 theta2) | |
719 | (* -0.181133 theta3) | |
720 | (* 0.553040 theta4) | |
721 | (* -0.861938 theta5) | |
722 | (* 0.677066 theta3 theta3) | |
723 | (* -0.212591 theta4 theta3)))) | |
724 | ((and (<= 1800 year) (< year 1900)) | |
725 | (let* ((theta (/ (- (calendar-astro-from-absolute | |
726 | (calendar-absolute-from-gregorian | |
727 | (list 7 1 year))) | |
728 | (calendar-astro-from-absolute | |
729 | (calendar-absolute-from-gregorian | |
730 | '(1 1 1900)))) | |
731 | 36525.0)) | |
732 | (theta2 (* theta theta)) | |
733 | (theta3 (* theta2 theta)) | |
734 | (theta4 (* theta2 theta2)) | |
735 | (theta5 (* theta3 theta2))) | |
736 | (+ -0.000009 | |
737 | (* 0.003844 theta) | |
738 | (* 0.083563 theta2) | |
739 | (* 0.865736 theta3) | |
740 | (* 4.867575 theta4) | |
741 | (* 15.845535 theta5) | |
742 | (* 31.332267 theta3 theta3) | |
743 | (* 38.291999 theta4 theta3) | |
744 | (* 28.316289 theta4 theta4) | |
745 | (* 11.636204 theta4 theta5) | |
746 | (* 2.043794 theta5 theta5)))) | |
747 | ((and (<= 1620 year) (< year 1800)) | |
748 | (let ((x (/ (- year 1600) 10.0))) | |
749 | (/ (+ (* 2.19167 x x) (* -40.675 x) 196.58333) 60.0 60.0 24.0))) | |
750 | (t (let* ((tmp (- (calendar-astro-from-absolute | |
751 | (calendar-absolute-from-gregorian | |
752 | (list 1 1 year))) | |
753 | 2382148)) | |
754 | (second (- (/ (* tmp tmp) 41048480.0) 15))) | |
755 | (/ second 60.0 60.0 24.0))))) | |
9f34a2a0 | 756 | |
087c56fa ER |
757 | (defun solar-sidereal-time (t0) |
758 | "Sidereal time (in hours) in Greenwich. | |
759 | ||
760 | At T0=Julian centuries of universal time. | |
761 | T0 must correspond to 0 hours UT." | |
762 | (let* ((mean-sid-time (+ 6.6973746 | |
763 | (* 2400.051337 t0) | |
764 | (* 0.0000258622 t0 t0) | |
765 | (* -0.0000000017222 t0 t0 t0))) | |
766 | (et (solar-ephemeris-time (list t0 0.0))) | |
767 | (nut-i (solar-ecliptic-coordinates et nil)) | |
768 | (nut (car (cdr (cdr (cdr nut-i))))) ; nutation | |
769 | (i (car (cdr nut-i)))) ; inclination | |
770 | (mod (+ (mod (+ mean-sid-time | |
771 | (/ (/ (* nut (solar-cosine-degrees i)) 15) 3600)) 24.0) | |
772 | 24.0) | |
773 | 24.0))) | |
774 | ||
775 | (defun solar-time-equation (date ut) | |
776 | "Equation of time expressed in hours at Gregorian DATE at Universal time UT." | |
777 | (let* ((et (solar-date-to-et date ut)) | |
778 | (ec (solar-ecliptic-coordinates et nil))) | |
779 | (car (cdr (cdr ec))))) | |
780 | ||
781 | (defun solar-date-to-et (date ut) | |
782 | "Ephemeris Time at Gregorian DATE at Universal Time UT (in hours). | |
783 | Expressed in julian centuries of Ephemeris Time." | |
784 | (let ((t0 (solar-julian-ut-centuries date))) | |
785 | (solar-ephemeris-time (list t0 ut)))) | |
786 | ||
9f34a2a0 JB |
787 | ;;;###autoload |
788 | (defun sunrise-sunset (&optional arg) | |
087c56fa | 789 | "Local time of sunrise and sunset for today. Accurate to a few seconds. |
ec4dfb6b | 790 | If called with an optional prefix argument, prompt for date. |
9f34a2a0 | 791 | |
ec4dfb6b RS |
792 | If called with an optional double prefix argument, prompt for longitude, |
793 | latitude, time zone, and date, and always use standard time. | |
9f34a2a0 JB |
794 | |
795 | This function is suitable for execution in a .emacs file." | |
796 | (interactive "p") | |
6a6e6405 | 797 | (or arg (setq arg 1)) |
fbfed6f0 JB |
798 | (if (and (< arg 16) |
799 | (not (and calendar-latitude calendar-longitude calendar-time-zone))) | |
800 | (solar-setup)) | |
9f34a2a0 | 801 | (let* ((calendar-longitude |
fbfed6f0 | 802 | (if (< arg 16) calendar-longitude |
9f34a2a0 JB |
803 | (solar-get-number |
804 | "Enter longitude (decimal fraction; + east, - west): "))) | |
805 | (calendar-latitude | |
fbfed6f0 | 806 | (if (< arg 16) calendar-latitude |
9f34a2a0 JB |
807 | (solar-get-number |
808 | "Enter latitude (decimal fraction; + north, - south): "))) | |
809 | (calendar-time-zone | |
fbfed6f0 | 810 | (if (< arg 16) calendar-time-zone |
9f34a2a0 | 811 | (solar-get-number |
e2fe2f52 | 812 | "Enter difference from Coordinated Universal Time (in minutes): "))) |
9f34a2a0 | 813 | (calendar-location-name |
fbfed6f0 JB |
814 | (if (< arg 16) calendar-location-name |
815 | (let ((float-output-format "%.1f")) | |
816 | (format "%s%s, %s%s" | |
6ff099c3 ER |
817 | (if (numberp calendar-latitude) |
818 | (abs calendar-latitude) | |
819 | (+ (aref calendar-latitude 0) | |
820 | (/ (aref calendar-latitude 1) 60.0))) | |
821 | (if (numberp calendar-latitude) | |
822 | (if (> calendar-latitude 0) "N" "S") | |
823 | (if (equal (aref calendar-latitude 2) 'north) "N" "S")) | |
824 | (if (numberp calendar-longitude) | |
825 | (abs calendar-longitude) | |
826 | (+ (aref calendar-longitude 0) | |
827 | (/ (aref calendar-longitude 1) 60.0))) | |
828 | (if (numberp calendar-longitude) | |
829 | (if (> calendar-longitude 0) "E" "W") | |
562a94a0 | 830 | (if (equal (aref calendar-longitude 2) 'east) |
6ff099c3 | 831 | "E" "W")))))) |
9f34a2a0 | 832 | (calendar-standard-time-zone-name |
fbfed6f0 | 833 | (if (< arg 16) calendar-standard-time-zone-name |
a92ade89 | 834 | (cond ((= calendar-time-zone 0) "UTC") |
fbfed6f0 | 835 | ((< calendar-time-zone 0) |
a92ade89 JB |
836 | (format "UTC%dmin" calendar-time-zone)) |
837 | (t (format "UTC+%dmin" calendar-time-zone))))) | |
ec4dfb6b RS |
838 | (calendar-daylight-savings-starts |
839 | (if (< arg 16) calendar-daylight-savings-starts)) | |
840 | (calendar-daylight-savings-ends | |
841 | (if (< arg 16) calendar-daylight-savings-ends)) | |
fbfed6f0 | 842 | (date (if (< arg 4) (calendar-current-date) (calendar-read-date))) |
9f34a2a0 | 843 | (date-string (calendar-date-string date t)) |
087c56fa | 844 | (time-string (solar-sunrise-sunset-string date)) |
9f34a2a0 JB |
845 | (msg (format "%s: %s" date-string time-string)) |
846 | (one-window (one-window-p t))) | |
80897760 | 847 | (if (<= (length msg) (frame-width)) |
9f34a2a0 JB |
848 | (message msg) |
849 | (with-output-to-temp-buffer "*temp*" | |
850 | (princ (concat date-string "\n" time-string))) | |
851 | (message (substitute-command-keys | |
852 | (if one-window | |
853 | (if pop-up-windows | |
854 | "Type \\[delete-other-windows] to remove temp window." | |
855 | "Type \\[switch-to-buffer] RET to remove temp window.") | |
856 | "Type \\[switch-to-buffer-other-window] RET to restore old contents of temp window.")))))) | |
857 | ||
858 | (defun calendar-sunrise-sunset () | |
859 | "Local time of sunrise and sunset for date under cursor. | |
087c56fa | 860 | Accurate to a few seconds." |
9f34a2a0 JB |
861 | (interactive) |
862 | (if (not (and calendar-latitude calendar-longitude calendar-time-zone)) | |
863 | (solar-setup)) | |
6a6e6405 | 864 | (let ((date (calendar-cursor-to-date t))) |
abd93e66 RS |
865 | (message "%s: %s" |
866 | (calendar-date-string date t t) | |
087c56fa | 867 | (solar-sunrise-sunset-string date)))) |
9f34a2a0 JB |
868 | |
869 | (defun diary-sunrise-sunset () | |
870 | "Local time of sunrise and sunset as a diary entry. | |
087c56fa | 871 | Accurate to a few seconds." |
9f34a2a0 JB |
872 | (if (not (and calendar-latitude calendar-longitude calendar-time-zone)) |
873 | (solar-setup)) | |
087c56fa | 874 | (solar-sunrise-sunset-string date)) |
9f34a2a0 JB |
875 | |
876 | (defun diary-sabbath-candles () | |
877 | "Local time of candle lighting diary entry--applies if date is a Friday. | |
878 | No diary entry if there is no sunset on that date." | |
879 | (if (not (and calendar-latitude calendar-longitude calendar-time-zone)) | |
880 | (solar-setup)) | |
881 | (if (= (% (calendar-absolute-from-gregorian date) 7) 5);; Friday | |
087c56fa | 882 | (let* ((sunset (car (cdr (solar-sunrise-sunset date)))) |
c68488d2 | 883 | (light (if sunset |
75af4a4a | 884 | (dst-adjust-time |
c68488d2 | 885 | date |
087c56fa | 886 | (- (car sunset) (/ 18.0 60.0)))))) |
c68488d2 ER |
887 | (if (and light (calendar-date-equal date (car light))) |
888 | (format "%s Sabbath candle lighting" | |
889 | (apply 'solar-time-string (cdr light))))))) | |
9f34a2a0 | 890 | |
087c56fa ER |
891 | (defun solar-equinoxes/solstices (k year) |
892 | "Date of equinox/solstice K for YEAR. | |
893 | K=0, spring equinox; K=1, summer solstice; K=2, fall equinox; | |
894 | K=3, winter solstice. | |
895 | RESULT is a gregorian local date. | |
896 | ||
897 | Accurate to less than a minute between 1951 and 2050." | |
898 | (let* ((JDE0 (solar-mean-equinoxes/solstices k year)) | |
899 | (T (/ (- JDE0 2451545.0) 36525)) | |
900 | (W (- (* 35999.373 T) 2.47)) | |
901 | (Delta-lambda (+ 1 (* 0.0334 (solar-cosine-degrees W)) | |
902 | (* 0.0007 (solar-cosine-degrees (* 2 W))))) | |
903 | (S (apply '+ (mapcar '(lambda(x) | |
904 | (* (car x) (solar-cosine-degrees | |
905 | (+ (* (car (cdr (cdr x))) T) | |
906 | (car (cdr x)))))) | |
907 | solar-seasons-data))) | |
908 | (JDE (+ JDE0 (/ (* 0.00001 S) Delta-lambda))) | |
909 | (correction (+ 102.3 (* 123.5 T) (* 32.5 T T))) | |
910 | ; ephemeris time correction | |
911 | (JD (- JDE (/ correction 86400))) | |
912 | (date (calendar-gregorian-from-absolute (floor (- JD 1721424.5)))) | |
913 | (time (- (- JD 0.5) (floor (- JD 0.5)))) | |
914 | ) | |
915 | (list (car date) (+ (car (cdr date)) time | |
916 | (/ (/ calendar-time-zone 60.0) 24.0)) | |
917 | (car (cdr (cdr date)))))) | |
918 | ||
919 | ; from Meeus, 1991, page 166 | |
920 | (defun solar-mean-equinoxes/solstices (k year) | |
921 | "Julian day of mean equinox/solstice K for YEAR. | |
922 | K=0, spring equinox; K=1, summer solstice; K=2, fall equinox; K=3, winter | |
923 | solstice. These formulas are only to be used between 1000 BC and 3000 AD." | |
924 | (let ((y (/ year 1000.0)) | |
925 | (z (/ (- year 2000) 1000.0))) | |
926 | (if (< year 1000) ; actually between -1000 and 1000 | |
927 | (cond ((equal k 0) (+ 1721139.29189 | |
928 | (* 365242.13740 y) | |
929 | (* 0.06134 y y) | |
930 | (* 0.00111 y y y) | |
931 | (* -0.00071 y y y y))) | |
932 | ((equal k 1) (+ 1721233.25401 | |
933 | (* 365241.72562 y) | |
934 | (* -0.05323 y y) | |
935 | (* 0.00907 y y y) | |
936 | (* 0.00025 y y y y))) | |
937 | ((equal k 2) (+ 1721325.70455 | |
938 | (* 365242.49558 y) | |
939 | (* -0.11677 y y) | |
940 | (* -0.00297 y y y) | |
941 | (* 0.00074 y y y y))) | |
942 | ((equal k 3) (+ 1721414.39987 | |
943 | (* 365242.88257 y) | |
944 | (* -0.00769 y y) | |
945 | (* -0.00933 y y y) | |
946 | (* -0.00006 y y y y)))) | |
947 | ; actually between 1000 and 3000 | |
948 | (cond ((equal k 0) (+ 2451623.80984 | |
949 | (* 365242.37404 z) | |
950 | (* 0.05169 z z) | |
951 | (* -0.00411 z z z) | |
952 | (* -0.00057 z z z z))) | |
953 | ((equal k 1) (+ 2451716.56767 | |
954 | (* 365241.62603 z) | |
955 | (* 0.00325 z z) | |
956 | (* 0.00888 z z z) | |
957 | (* -0.00030 z z z z))) | |
958 | ((equal k 2) (+ 2451810.21715 | |
959 | (* 365242.01767 z) | |
960 | (* -0.11575 z z) | |
961 | (* 0.00337 z z z) | |
962 | (* 0.00078 z z z z))) | |
963 | ((equal k 3) (+ 2451900.05952 | |
964 | (* 365242.74049 z) | |
965 | (* -0.06223 z z) | |
966 | (* -0.00823 z z z) | |
967 | (* 0.00032 z z z z))))))) | |
968 | ||
969 | ; from Meeus, 1991, page 167 | |
970 | (defconst solar-seasons-data | |
971 | '((485 324.96 1934.136) | |
972 | (203 337.23 32964.467) | |
973 | (199 342.08 20.186) | |
974 | (182 27.85 445267.112) | |
975 | (156 73.14 45036.886) | |
976 | (136 171.52 22518.443) | |
977 | (77 222.54 65928.934) | |
978 | (74 296.72 3034.906) | |
979 | (70 243.58 9037.513) | |
980 | (58 119.81 33718.147) | |
981 | (52 297.17 150.678) | |
982 | (50 21.02 2281.226) | |
983 | (45 247.54 29929.562) | |
984 | (44 325.15 31555.956) | |
985 | (29 60.93 4443.417) | |
986 | (18 155.12 67555.328) | |
987 | (17 288.79 4562.452) | |
988 | (16 198.04 62894.029) | |
989 | (14 199.76 31436.921) | |
990 | (12 95.39 14577.848) | |
991 | (12 287.11 31931.756) | |
992 | (12 320.81 34777.259) | |
993 | (9 227.73 1222.114) | |
994 | (8 15.45 16859.074))) | |
995 | ||
ce4d3fff | 996 | ;;;###autoload |
a92ade89 | 997 | (defun solar-equinoxes-solstices () |
087c56fa | 998 | "*local* date and time of equinoxes and solstices, if visible in the calendar window. |
9f34a2a0 | 999 | Requires floating point." |
a92ade89 JB |
1000 | (let ((m displayed-month) |
1001 | (y displayed-year)) | |
9f34a2a0 JB |
1002 | (increment-calendar-month m y (cond ((= 1 (% m 3)) -1) |
1003 | ((= 2 (% m 3)) 1) | |
1004 | (t 0))) | |
1005 | (let* ((calendar-standard-time-zone-name | |
a92ade89 | 1006 | (if calendar-time-zone calendar-standard-time-zone-name "UTC")) |
9f34a2a0 JB |
1007 | (calendar-daylight-savings-starts |
1008 | (if calendar-time-zone calendar-daylight-savings-starts)) | |
1009 | (calendar-daylight-savings-ends | |
1010 | (if calendar-time-zone calendar-daylight-savings-ends)) | |
1011 | (calendar-time-zone (if calendar-time-zone calendar-time-zone 0)) | |
1012 | (k (1- (/ m 3))) | |
087c56fa ER |
1013 | (d0 (solar-equinoxes/solstices k y)) |
1014 | (d1 (list (car d0) (floor (car (cdr d0))) (car (cdr (cdr d0))))) | |
1015 | (h0 (* 24 (- (car (cdr d0)) (floor (car (cdr d0)))))) | |
1016 | (adj (dst-adjust-time d1 h0)) | |
1017 | (d (list (car d1) (+ (car (cdr d1)) | |
1018 | (/ (car (cdr adj)) 24.0)) | |
1019 | (car (cdr (cdr d1))))) | |
1020 | ; The following is nearly as accurate, but not quite: | |
1021 | ;(d0 (solar-date-next-longitude | |
1022 | ; (calendar-astro-from-absolute | |
1023 | ; (calendar-absolute-from-gregorian | |
1024 | ; (list (+ 3 (* k 3)) 15 y))) | |
1025 | ; 90)) | |
1026 | ;(abs-day (calendar-absolute-from-astro d))) | |
1027 | (abs-day (calendar-absolute-from-gregorian d))) | |
75af4a4a ER |
1028 | (list |
1029 | (list (calendar-gregorian-from-absolute (floor abs-day)) | |
1030 | (format "%s %s" | |
1031 | (nth k (if (and calendar-latitude | |
1032 | (< (calendar-latitude) 0)) | |
1033 | solar-s-hemi-seasons | |
1034 | solar-n-hemi-seasons)) | |
1035 | (solar-time-string | |
1036 | (* 24 (- abs-day (floor abs-day))) | |
1037 | (if (dst-in-effect abs-day) | |
1038 | calendar-daylight-time-zone-name | |
1039 | calendar-standard-time-zone-name)))))))) | |
1040 | ||
9f34a2a0 JB |
1041 | |
1042 | (provide 'solar) | |
1043 | ||
1044 | ;;; solar.el ends here |