[bpt/emacs.git] / lisp / calendar / lunar.el

;;; lunar.el --- calendar functions for phases of the moon

;; Copyright (C) 1992, 1993, 1995, 1997, 2001, 2002, 2003, 2004, 2005,
;;   2006, 2007, 2008  Free Software Foundation, Inc.

;; Author: Edward M. Reingold <reingold@cs.uiuc.edu>
;; Maintainer: Glenn Morris <rgm@gnu.org>
;; Keywords: calendar
;; Human-Keywords: moon, lunar phases, calendar, diary

;; This file is part of GNU Emacs.

;; GNU Emacs is free software; you can redistribute it and/or modify
;; it under the terms of the GNU General Public License as published by
;; the Free Software Foundation; either version 3, or (at your option)
;; any later version.

;; GNU Emacs is distributed in the hope that it will be useful,
;; but WITHOUT ANY WARRANTY; without even the implied warranty of
;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
;; GNU General Public License for more details.

;; You should have received a copy of the GNU General Public License
;; along with GNU Emacs; see the file COPYING.  If not, write to the
;; Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
;; Boston, MA 02110-1301, USA.

;;; Commentary:

;; See calendar.el.

;; Based on ``Astronomical Formulae for Calculators,'' 3rd ed., by Jean Meeus,
;; Willmann-Bell, Inc., 1985 and ``Astronomical Algorithms'' by Jean Meeus,
;; Willmann-Bell, Inc., 1991.
;;
;; WARNING: The calculations will be accurate only to within a few minutes.

;; The author would be delighted to have an astronomically more sophisticated
;; person rewrite the code for the lunar calculations in this file!

;;; Code:

(require 'calendar)
(require 'solar)
(require 'cal-dst)
;; calendar-astro-to-absolute and v versa are cal-autoloads.
;;;(require 'cal-julian)

(defun lunar-phase (index)
  "Local date and time of lunar phase INDEX.
Integer below INDEX/4 gives the lunation number, counting from Jan 1, 1900;
remainder mod 4 gives the phase: 0 new moon, 1 first quarter, 2 full moon,
3 last quarter."
  (let* ((phase (mod index 4))
         (index (/ index 4.0))
         (time (/ index 1236.85))
         (date (+ (calendar-absolute-from-gregorian '(1 0.5 1900))
                  0.75933
                  (* 29.53058868 index)
                  (* 0.0001178 time time)
                  (* -0.000000155 time time time)
                  (* 0.00033
                     (solar-sin-degrees (+ 166.56
                                           (* 132.87 time)
                                           (* -0.009173 time time))))))
         (sun-anomaly (mod
                       (+ 359.2242
                          (* 29.105356 index)
                          (* -0.0000333 time time)
                          (* -0.00000347 time time time))
                       360.0))
         (moon-anomaly (mod
                        (+ 306.0253
                           (* 385.81691806 index)
                           (* 0.0107306 time time)
                           (* 0.00001236 time time time))
                        360.0))
         (moon-lat (mod
                    (+ 21.2964
                       (* 390.67050646 index)
                       (* -0.0016528 time time)
                       (* -0.00000239 time time time))
                    360.0))
         (adjustment
          (if (memq phase '(0 2))
              (+ (* (- 0.1734 (* 0.000393 time))
                    (solar-sin-degrees sun-anomaly))
                 (* 0.0021 (solar-sin-degrees (* 2 sun-anomaly)))
                 (* -0.4068 (solar-sin-degrees moon-anomaly))
                 (* 0.0161 (solar-sin-degrees (* 2 moon-anomaly)))
                 (* -0.0004 (solar-sin-degrees (* 3 moon-anomaly)))
                 (* 0.0104 (solar-sin-degrees (* 2 moon-lat)))
                 (* -0.0051 (solar-sin-degrees (+ sun-anomaly moon-anomaly)))
                 (* -0.0074 (solar-sin-degrees (- sun-anomaly moon-anomaly)))
                 (* 0.0004 (solar-sin-degrees (+ (* 2 moon-lat) sun-anomaly)))
                 (* -0.0004 (solar-sin-degrees (- (* 2 moon-lat) sun-anomaly)))
                 (* -0.0006 (solar-sin-degrees
                             (+ (* 2 moon-lat) moon-anomaly)))
                 (* 0.0010 (solar-sin-degrees (- (* 2 moon-lat) moon-anomaly)))
                 (* 0.0005 (solar-sin-degrees
                            (+ (* 2 moon-anomaly) sun-anomaly))))
            (+ (* (- 0.1721 (* 0.0004 time))
                  (solar-sin-degrees sun-anomaly))
               (* 0.0021 (solar-sin-degrees (* 2 sun-anomaly)))
               (* -0.6280 (solar-sin-degrees moon-anomaly))
               (* 0.0089 (solar-sin-degrees (* 2 moon-anomaly)))
               (* -0.0004 (solar-sin-degrees (* 3 moon-anomaly)))
               (* 0.0079 (solar-sin-degrees (* 2 moon-lat)))
               (* -0.0119 (solar-sin-degrees (+ sun-anomaly moon-anomaly)))
               (* -0.0047 (solar-sin-degrees (- sun-anomaly moon-anomaly)))
               (* 0.0003 (solar-sin-degrees (+ (* 2 moon-lat) sun-anomaly)))
               (* -0.0004 (solar-sin-degrees (- (* 2 moon-lat) sun-anomaly)))
               (* -0.0006 (solar-sin-degrees (+ (* 2 moon-lat) moon-anomaly)))
               (* 0.0021 (solar-sin-degrees (- (* 2 moon-lat) moon-anomaly)))
               (* 0.0003 (solar-sin-degrees
                          (+ (* 2 moon-anomaly) sun-anomaly)))
               (* 0.0004 (solar-sin-degrees
                          (- sun-anomaly (* 2 moon-anomaly))))
               (* -0.0003 (solar-sin-degrees
                           (+ (* 2 sun-anomaly) moon-anomaly))))))
         (adj (+ 0.0028
                 (* -0.0004 (solar-cosine-degrees
                             sun-anomaly))
                 (* 0.0003 (solar-cosine-degrees
                            moon-anomaly))))
         (adjustment (cond ((= phase 1) (+ adjustment adj))
                           ((= phase 2) (- adjustment adj))
                           (t adjustment)))
         (date (+ date adjustment))
         (date (+ date (/ (- calendar-time-zone
                             (solar-ephemeris-correction
                              (calendar-extract-year
                               (calendar-gregorian-from-absolute
                                (truncate date)))))
                          60.0 24.0)))
         (time (* 24 (- date (truncate date))))
         (date (calendar-gregorian-from-absolute (truncate date)))
         (adj (dst-adjust-time date time)))
    (list (car adj) (apply 'solar-time-string (cdr adj)) phase)))

(defun lunar-phase-list (month year)
  "List of lunar phases for three months starting with Gregorian MONTH, YEAR."
  (let* ((end-month month)
         (end-year year)
         (start-month month)
         (start-year year)
         (end-date (progn
                     (calendar-increment-month end-month end-year 3)
                     (list (list end-month 1 end-year))))
         (start-date (progn
                       (calendar-increment-month start-month start-year -1)
                       (list (list start-month
                                   (calendar-last-day-of-month
                                    start-month start-year)
                                   start-year))))
         (index (* 4 (truncate
                      (* 12.3685
                         (+ year
                            ( / (calendar-day-number (list month 1 year))
                                366.0)
                            -1900)))))
         (new-moon (lunar-phase index))
         list)
    (while (calendar-date-compare new-moon end-date)
      (if (calendar-date-compare start-date new-moon)
          (setq list (append list (list new-moon))))
      (setq index (1+ index)
            new-moon (lunar-phase index)))
    list))

(defun lunar-phase-name (phase)
  "Name of lunar PHASE.
0 = new moon, 1 = first quarter, 2 = full moon, 3 = last quarter."
  (cond ((= 0 phase) "New Moon")
        ((= 1 phase) "First Quarter Moon")
        ((= 2 phase) "Full Moon")
        ((= 3 phase) "Last Quarter Moon")))

(defvar displayed-month)                ; from calendar-generate
(defvar displayed-year)

;;;###cal-autoload
(defun calendar-phases-of-moon ()
  "Create a buffer with the lunar phases for the current calendar window."
  (interactive)
  (message "Computing phases of the moon...")
  (let ((m1 displayed-month)
        (y1 displayed-year)
        (m2 displayed-month)
        (y2 displayed-year))
    (calendar-increment-month m1 y1 -1)
    (calendar-increment-month m2 y2 1)
    (calendar-in-read-only-buffer lunar-phases-buffer
      (calendar-set-mode-line
       (if (= y1 y2)
           (format "Phases of the Moon from %s to %s, %d%%-"
                   (calendar-month-name m1) (calendar-month-name m2) y2)
         (format "Phases of the Moon from %s, %d to %s, %d%%-"
                 (calendar-month-name m1) y1 (calendar-month-name m2) y2)))
      (insert
       (mapconcat
        (lambda (x)
          (let ((date (car x))
                (time (cadr x))
                (phase (nth 2 x)))
            (concat (calendar-date-string date)
                    ": "
                    (lunar-phase-name phase)
                    " "
                    time)))
        (lunar-phase-list m1 y1) "\n")))
    (message "Computing phases of the moon...done")))

;; FIXME ?
;;;###autoload
(defun phases-of-moon (&optional arg)
  "Display the quarters of the moon for last month, this month, and next month.
If called with an optional prefix argument ARG, prompts for month and year.
This function is suitable for execution in a .emacs file."
  (interactive "P")
  (save-excursion
    (let* ((date (if arg (calendar-read-date t)
                   (calendar-current-date)))
           (displayed-month (calendar-extract-month date))
           (displayed-year (calendar-extract-year date)))
      (calendar-phases-of-moon))))

(defvar date)

;; To be called from diary-list-sexp-entries, where DATE is bound.

;;;###diary-autoload
(defun diary-phases-of-moon (&optional mark)
  "Moon phases diary entry.
An optional parameter MARK specifies a face or single-character string to
use when highlighting the day in the calendar."
  (let* ((index (* 4
                   (truncate
                    (* 12.3685
                       (+ (calendar-extract-year date)
                          ( / (calendar-day-number date)
                              366.0)
                          -1900)))))
         (phase (lunar-phase index)))
    (while (calendar-date-compare phase (list date))
      (setq index (1+ index)
            phase (lunar-phase index)))
    (if (calendar-date-equal (car phase) date)
        (cons mark (concat (lunar-phase-name (nth 2 phase)) " "
                           (cadr phase))))))

;; For the Chinese calendar the calculations for the new moon need to be more
;; accurate than those above, so we use more terms in the approximation.
(defun lunar-new-moon-time (k)
  "Astronomical (Julian) day number of K th new moon."
  (let* ((T (/ k 1236.85))
         (T2 (* T T))
         (T3 (* T T T))
         (T4 (* T2 T2))
         (JDE (+ 2451550.09765
                 (* 29.530588853 k)
                 (* 0.0001337 T2)
                 (* -0.000000150 T3)
                 (* 0.00000000073 T4)))
         (E (- 1 (* 0.002516 T) (* 0.0000074 T2)))
         (sun-anomaly (+ 2.5534
                         (* 29.10535669 k)
                         (* -0.0000218 T2)
                         (* -0.00000011 T3)))
         (moon-anomaly (+ 201.5643
                          (* 385.81693528 k)
                          (* 0.0107438 T2)
                          (* 0.00001239 T3)
                          (* -0.000000058 T4)))
         (moon-argument (+ 160.7108
                           (* 390.67050274 k)
                           (* -0.0016341 T2)
                           (* -0.00000227 T3)
                           (* 0.000000011 T4)))
         (omega (+ 124.7746
                   (* -1.56375580 k)
                   (* 0.0020691 T2)
                   (* 0.00000215 T3)))
         (A1  (+ 299.77 (*  0.107408 k) (* -0.009173 T2)))
         (A2  (+ 251.88 (*  0.016321 k)))
         (A3  (+ 251.83 (* 26.641886 k)))
         (A4  (+ 349.42 (* 36.412478 k)))
         (A5  (+  84.66 (* 18.206239 k)))
         (A6  (+ 141.74 (* 53.303771 k)))
         (A7  (+ 207.14 (*  2.453732 k)))
         (A8  (+ 154.84 (*  7.306860 k)))
         (A9  (+  34.52 (* 27.261239 k)))
         (A10 (+ 207.19 (*  0.121824 k)))
         (A11 (+ 291.34 (*  1.844379 k)))
         (A12 (+ 161.72 (* 24.198154 k)))
         (A13 (+ 239.56 (* 25.513099 k)))
         (A14 (+ 331.55 (*  3.592518 k)))
         (correction
          (+ (* -0.40720   (solar-sin-degrees moon-anomaly))
             (*  0.17241 E (solar-sin-degrees sun-anomaly))
             (*  0.01608   (solar-sin-degrees (* 2 moon-anomaly)))
             (*  0.01039   (solar-sin-degrees (* 2 moon-argument)))
             (*  0.00739 E (solar-sin-degrees (- moon-anomaly sun-anomaly)))
             (* -0.00514 E (solar-sin-degrees (+ moon-anomaly sun-anomaly)))
             (*  0.00208 E E (solar-sin-degrees (* 2 sun-anomaly)))
             (* -0.00111   (solar-sin-degrees
                            (- moon-anomaly (* 2 moon-argument))))
             (* -0.00057   (solar-sin-degrees
                            (+ moon-anomaly (* 2 moon-argument))))
             (*  0.00056 E (solar-sin-degrees
                            (+ (* 2 moon-anomaly) sun-anomaly)))
             (* -0.00042   (solar-sin-degrees (* 3 moon-anomaly)))
             (*  0.00042 E (solar-sin-degrees
                            (+ sun-anomaly (* 2 moon-argument))))
             (*  0.00038 E (solar-sin-degrees
                            (- sun-anomaly (* 2 moon-argument))))
             (* -0.00024 E (solar-sin-degrees
                            (- (* 2 moon-anomaly) sun-anomaly)))
             (* -0.00017   (solar-sin-degrees omega))
             (* -0.00007   (solar-sin-degrees
                            (+ moon-anomaly (* 2 sun-anomaly))))
             (*  0.00004   (solar-sin-degrees
                            (- (* 2 moon-anomaly) (* 2 moon-argument))))
             (*  0.00004   (solar-sin-degrees (* 3 sun-anomaly)))
             (*  0.00003   (solar-sin-degrees (+ moon-anomaly sun-anomaly
                                                 (* -2 moon-argument))))
             (*  0.00003   (solar-sin-degrees
                            (+ (* 2 moon-anomaly) (* 2 moon-argument))))
             (* -0.00003   (solar-sin-degrees (+ moon-anomaly sun-anomaly
                                                 (* 2 moon-argument))))
             (*  0.00003   (solar-sin-degrees (- moon-anomaly sun-anomaly
                                                 (* -2 moon-argument))))
             (* -0.00002   (solar-sin-degrees (- moon-anomaly sun-anomaly
                                                 (* 2 moon-argument))))
             (* -0.00002   (solar-sin-degrees
                            (+ (* 3 moon-anomaly) sun-anomaly)))
             (*  0.00002   (solar-sin-degrees (* 4 moon-anomaly)))))
         (additional
          (+ (* 0.000325 (solar-sin-degrees A1))
             (* 0.000165 (solar-sin-degrees A2))
             (* 0.000164 (solar-sin-degrees A3))
             (* 0.000126 (solar-sin-degrees A4))
             (* 0.000110 (solar-sin-degrees A5))
             (* 0.000062 (solar-sin-degrees A6))
             (* 0.000060 (solar-sin-degrees A7))
             (* 0.000056 (solar-sin-degrees A8))
             (* 0.000047 (solar-sin-degrees A9))
             (* 0.000042 (solar-sin-degrees A10))
             (* 0.000040 (solar-sin-degrees A11))
             (* 0.000037 (solar-sin-degrees A12))
             (* 0.000035 (solar-sin-degrees A13))
             (* 0.000023 (solar-sin-degrees A14))))
         (newJDE (+ JDE correction additional)))
    (+ newJDE
       (- (solar-ephemeris-correction
           (calendar-extract-year
            (calendar-gregorian-from-absolute
             (floor (calendar-astro-to-absolute newJDE))))))
       (/ calendar-time-zone 60.0 24.0))))

(defun lunar-new-moon-on-or-after (d)
  "Julian day number of first new moon on or after Julian day number D.
The fractional part is the time of day.

The date and time are local time, including any daylight saving rules,
as governed by the values of `calendar-daylight-savings-starts',
`calendar-daylight-savings-starts-time', `calendar-daylight-savings-ends',
`calendar-daylight-savings-ends-time', `calendar-daylight-time-offset', and
`calendar-time-zone'."
  (let* ((date (calendar-gregorian-from-absolute
                (floor (calendar-astro-to-absolute d))))
         (year (+ (calendar-extract-year date)
                  (/ (calendar-day-number date) 365.25)))
         (k (floor (* (- year 2000.0) 12.3685)))
         (date (lunar-new-moon-time k))
         (a-date (progn
                   (while (< date d)
                     (setq k (1+ k)
                           date (lunar-new-moon-time k)))
                   (calendar-astro-to-absolute date)))
         (time (* 24 (- a-date (truncate a-date))))
         (date (calendar-gregorian-from-absolute (truncate a-date)))
         (adj (dst-adjust-time date time)))
    (calendar-astro-from-absolute
     (+ (calendar-absolute-from-gregorian (car adj))
        (/ (cadr adj) 24.0)))))

(provide 'lunar)

;; arch-tag: 72f0b8a4-7bcc-4a1b-b67a-ff53c4a1d222
;;; lunar.el ends here
Commit	Line	Data
3afbc435	1	;;; lunar.el --- calendar functions for phases of the moon
7e1dae73	2
a20b3848	3	;; Copyright (C) 1992, 1993, 1995, 1997, 2001, 2002, 2003, 2004, 2005,
8b72699e	4	;; 2006, 2007, 2008 Free Software Foundation, Inc.
7e1dae73 JB	5
7e1dae73 JB	6	;; Author: Edward M. Reingold <reingold@cs.uiuc.edu>
dbfca9c4	7	;; Maintainer: Glenn Morris <rgm@gnu.org>
e9571d2a ER	8	;; Keywords: calendar
e9571d2a ER	9	;; Human-Keywords: moon, lunar phases, calendar, diary
7e1dae73 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
59243403 RS	14	;; it under the terms of the GNU General Public License as published by
075969b4	15	;; the Free Software Foundation; either version 3, or (at your option)
59243403 RS	16	;; any later version.
59243403 RS	17
7e1dae73	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	24	;; along with GNU Emacs; see the file COPYING. If not, write to the
3a35cf56 LK	25	;; Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
3a35cf56 LK	26	;; Boston, MA 02110-1301, USA.
7e1dae73 JB	27
	28	;;; Commentary:
	29
4bd7ad5f	30	;; See calendar.el.
7e1dae73 JB	31
7e1dae73 JB	32	;; Based on ``Astronomical Formulae for Calculators,'' 3rd ed., by Jean Meeus,
8a45b040 ER	33	;; Willmann-Bell, Inc., 1985 and ``Astronomical Algorithms'' by Jean Meeus,
8a45b040 ER	34	;; Willmann-Bell, Inc., 1991.
7e1dae73 JB	35	;;
	36	;; WARNING: The calculations will be accurate only to within a few minutes.
	37
	38	;; The author would be delighted to have an astronomically more sophisticated
	39	;; person rewrite the code for the lunar calculations in this file!
	40
7e1dae73 JB	41	;;; Code:
7e1dae73 JB	42
164e17e1	43	(require 'calendar)
7e1dae73	44	(require 'solar)
164e17e1	45	(require 'cal-dst)
5c645a20	46	;; calendar-astro-to-absolute and v versa are cal-autoloads.
164e17e1	47	;;;(require 'cal-julian)
7e1dae73	48
7e1dae73 JB	49	(defun lunar-phase (index)
	50	"Local date and time of lunar phase INDEX.
	51	Integer below INDEX/4 gives the lunation number, counting from Jan 1, 1900;
	52	remainder mod 4 gives the phase: 0 new moon, 1 first quarter, 2 full moon,
	53	3 last quarter."
b6c0aabf	54	(let* ((phase (mod index 4))
7e1dae73 JB	55	(index (/ index 4.0))
	56	(time (/ index 1236.85))
	57	(date (+ (calendar-absolute-from-gregorian '(1 0.5 1900))
	58	0.75933
	59	(* 29.53058868 index)
	60	(* 0.0001178 time time)
	61	(* -0.000000155 time time time)
	62	(* 0.00033
	63	(solar-sin-degrees (+ 166.56
	64	(* 132.87 time)
	65	(* -0.009173 time time))))))
e7bb3aeb	66	(sun-anomaly (mod
7e1dae73 JB	67	(+ 359.2242
	68	(* 29.105356 index)
	69	(* -0.0000333 time time)
	70	(* -0.00000347 time time time))
	71	360.0))
e7bb3aeb	72	(moon-anomaly (mod
7e1dae73 JB	73	(+ 306.0253
	74	(* 385.81691806 index)
	75	(* 0.0107306 time time)
	76	(* 0.00001236 time time time))
	77	360.0))
e7bb3aeb	78	(moon-lat (mod
7e1dae73 JB	79	(+ 21.2964
	80	(* 390.67050646 index)
	81	(* -0.0016528 time time)
	82	(* -0.00000239 time time time))
	83	360.0))
	84	(adjustment
	85	(if (memq phase '(0 2))
	86	(+ (* (- 0.1734 (* 0.000393 time))
	87	(solar-sin-degrees sun-anomaly))
	88	(* 0.0021 (solar-sin-degrees (* 2 sun-anomaly)))
	89	(* -0.4068 (solar-sin-degrees moon-anomaly))
	90	(* 0.0161 (solar-sin-degrees (* 2 moon-anomaly)))
	91	(* -0.0004 (solar-sin-degrees (* 3 moon-anomaly)))
	92	(* 0.0104 (solar-sin-degrees (* 2 moon-lat)))
	93	(* -0.0051 (solar-sin-degrees (+ sun-anomaly moon-anomaly)))
	94	(* -0.0074 (solar-sin-degrees (- sun-anomaly moon-anomaly)))
	95	(* 0.0004 (solar-sin-degrees (+ (* 2 moon-lat) sun-anomaly)))
	96	(* -0.0004 (solar-sin-degrees (- (* 2 moon-lat) sun-anomaly)))
	97	(* -0.0006 (solar-sin-degrees
	98	(+ (* 2 moon-lat) moon-anomaly)))
	99	(* 0.0010 (solar-sin-degrees (- (* 2 moon-lat) moon-anomaly)))
	100	(* 0.0005 (solar-sin-degrees
	101	(+ (* 2 moon-anomaly) sun-anomaly))))
	102	(+ (* (- 0.1721 (* 0.0004 time))
	103	(solar-sin-degrees sun-anomaly))
	104	(* 0.0021 (solar-sin-degrees (* 2 sun-anomaly)))
	105	(* -0.6280 (solar-sin-degrees moon-anomaly))
	106	(* 0.0089 (solar-sin-degrees (* 2 moon-anomaly)))
	107	(* -0.0004 (solar-sin-degrees (* 3 moon-anomaly)))
	108	(* 0.0079 (solar-sin-degrees (* 2 moon-lat)))
	109	(* -0.0119 (solar-sin-degrees (+ sun-anomaly moon-anomaly)))
	110	(* -0.0047 (solar-sin-degrees (- sun-anomaly moon-anomaly)))
	111	(* 0.0003 (solar-sin-degrees (+ (* 2 moon-lat) sun-anomaly)))
	112	(* -0.0004 (solar-sin-degrees (- (* 2 moon-lat) sun-anomaly)))
	113	(* -0.0006 (solar-sin-degrees (+ (* 2 moon-lat) moon-anomaly)))
	114	(* 0.0021 (solar-sin-degrees (- (* 2 moon-lat) moon-anomaly)))
	115	(* 0.0003 (solar-sin-degrees
	116	(+ (* 2 moon-anomaly) sun-anomaly)))
	117	(* 0.0004 (solar-sin-degrees
	118	(- sun-anomaly (* 2 moon-anomaly))))
	119	(* -0.0003 (solar-sin-degrees
9c0b9187	120	(+ (* 2 sun-anomaly) moon-anomaly))))))
7e1dae73 JB	121	(adj (+ 0.0028
	122	(* -0.0004 (solar-cosine-degrees
	123	sun-anomaly))
	124	(* 0.0003 (solar-cosine-degrees
	125	moon-anomaly))))
	126	(adjustment (cond ((= phase 1) (+ adjustment adj))
	127	((= phase 2) (- adjustment adj))
	128	(t adjustment)))
	129	(date (+ date adjustment))
71ea27ee GM	130	(date (+ date (/ (- calendar-time-zone
71ea27ee GM	131	(solar-ephemeris-correction
e803eab7	132	(calendar-extract-year
80e48f9f	133	(calendar-gregorian-from-absolute
238c5cc1	134	(truncate date)))))
71ea27ee	135	60.0 24.0)))
7e1dae73	136	(time (* 24 (- date (truncate date))))
71ea27ee	137	(date (calendar-gregorian-from-absolute (truncate date)))
8a45b040	138	(adj (dst-adjust-time date time)))
e483c805	139	(list (car adj) (apply 'solar-time-string (cdr adj)) phase)))
7e1dae73	140
9c0b9187 GM	141	(defun lunar-phase-list (month year)
9c0b9187 GM	142	"List of lunar phases for three months starting with Gregorian MONTH, YEAR."
4bd7ad5f GM	143	(let* ((end-month month)
	144	(end-year year)
	145	(start-month month)
	146	(start-year year)
	147	(end-date (progn
e803eab7	148	(calendar-increment-month end-month end-year 3)
4bd7ad5f GM	149	(list (list end-month 1 end-year))))
4bd7ad5f GM	150	(start-date (progn
e803eab7	151	(calendar-increment-month start-month start-year -1)
4bd7ad5f	152	(list (list start-month
9c0b9187 GM	153	(calendar-last-day-of-month
9c0b9187 GM	154	start-month start-year)
4bd7ad5f GM	155	start-year))))
4bd7ad5f GM	156	(index (* 4 (truncate
9c0b9187 GM	157	(* 12.3685
	158	(+ year
	159	( / (calendar-day-number (list month 1 year))
	160	366.0)
	161	-1900)))))
4bd7ad5f GM	162	(new-moon (lunar-phase index))
	163	list)
	164	(while (calendar-date-compare new-moon end-date)
	165	(if (calendar-date-compare start-date new-moon)
	166	(setq list (append list (list new-moon))))
	167	(setq index (1+ index)
	168	new-moon (lunar-phase index)))
	169	list))
9c0b9187	170
7e1dae73 JB	171	(defun lunar-phase-name (phase)
	172	"Name of lunar PHASE.
	173	0 = new moon, 1 = first quarter, 2 = full moon, 3 = last quarter."
	174	(cond ((= 0 phase) "New Moon")
	175	((= 1 phase) "First Quarter Moon")
	176	((= 2 phase) "Full Moon")
	177	((= 3 phase) "Last Quarter Moon")))
	178
e803eab7	179	(defvar displayed-month) ; from calendar-generate
9c0b9187 GM	180	(defvar displayed-year)
9c0b9187 GM	181
fb97eeac	182	;;;###cal-autoload
7e1dae73 JB	183	(defun calendar-phases-of-moon ()
	184	"Create a buffer with the lunar phases for the current calendar window."
	185	(interactive)
	186	(message "Computing phases of the moon...")
	187	(let ((m1 displayed-month)
	188	(y1 displayed-year)
	189	(m2 displayed-month)
cc4879d2	190	(y2 displayed-year))
e803eab7 GM	191	(calendar-increment-month m1 y1 -1)
e803eab7 GM	192	(calendar-increment-month m2 y2 1)
48844538 GM	193	(calendar-in-read-only-buffer lunar-phases-buffer
	194	(calendar-set-mode-line
	195	(if (= y1 y2)
	196	(format "Phases of the Moon from %s to %s, %d%%-"
	197	(calendar-month-name m1) (calendar-month-name m2) y2)
	198	(format "Phases of the Moon from %s, %d to %s, %d%%-"
	199	(calendar-month-name m1) y1 (calendar-month-name m2) y2)))
	200	(insert
	201	(mapconcat
	202	(lambda (x)
	203	(let ((date (car x))
	204	(time (cadr x))
	205	(phase (nth 2 x)))
	206	(concat (calendar-date-string date)
	207	": "
	208	(lunar-phase-name phase)
	209	" "
	210	time)))
	211	(lunar-phase-list m1 y1) "\n")))
7e1dae73 JB	212	(message "Computing phases of the moon...done")))
7e1dae73 JB	213
8c34d83e	214	;; FIXME ?
7e1dae73 JB	215	;;;###autoload
	216	(defun phases-of-moon (&optional arg)
	217	"Display the quarters of the moon for last month, this month, and next month.
f3924c53	218	If called with an optional prefix argument ARG, prompts for month and year.
7e1dae73 JB	219	This function is suitable for execution in a .emacs file."
	220	(interactive "P")
	221	(save-excursion
9c0b9187	222	(let* ((date (if arg (calendar-read-date t)
cc4879d2	223	(calendar-current-date)))
e803eab7 GM	224	(displayed-month (calendar-extract-month date))
e803eab7 GM	225	(displayed-year (calendar-extract-year date)))
7e1dae73 JB	226	(calendar-phases-of-moon))))
7e1dae73 JB	227
9c0b9187 GM	228	(defvar date)
9c0b9187 GM	229
8c34d83e	230	;; To be called from diary-list-sexp-entries, where DATE is bound.
9c0b9187	231
fb97eeac	232	;;;###diary-autoload
e0ab9b68	233	(defun diary-phases-of-moon (&optional mark)
9c0b9187	234	"Moon phases diary entry.
a1506d29	235	An optional parameter MARK specifies a face or single-character string to
e0ab9b68	236	use when highlighting the day in the calendar."
7e1dae73 JB	237	(let* ((index (* 4
	238	(truncate
	239	(* 12.3685
e803eab7	240	(+ (calendar-extract-year date)
7e1dae73 JB	241	( / (calendar-day-number date)
	242	366.0)
	243	-1900)))))
	244	(phase (lunar-phase index)))
	245	(while (calendar-date-compare phase (list date))
9c0b9187 GM	246	(setq index (1+ index)
9c0b9187 GM	247	phase (lunar-phase index)))
7e1dae73	248	(if (calendar-date-equal (car phase) date)
9c0b9187 GM	249	(cons mark (concat (lunar-phase-name (nth 2 phase)) " "
9c0b9187 GM	250	(cadr phase))))))
7e1dae73	251
9c0b9187 GM	252	;; For the Chinese calendar the calculations for the new moon need to be more
9c0b9187 GM	253	;; accurate than those above, so we use more terms in the approximation.
8a45b040 ER	254	(defun lunar-new-moon-time (k)
	255	"Astronomical (Julian) day number of K th new moon."
	256	(let* ((T (/ k 1236.85))
71ea27ee GM	257	(T2 (* T T))
	258	(T3 (* T T T))
	259	(T4 (* T2 T2))
	260	(JDE (+ 2451550.09765
	261	(* 29.530588853 k)
	262	(* 0.0001337 T2)
	263	(* -0.000000150 T3)
	264	(* 0.00000000073 T4)))
	265	(E (- 1 (* 0.002516 T) (* 0.0000074 T2)))
	266	(sun-anomaly (+ 2.5534
	267	(* 29.10535669 k)
	268	(* -0.0000218 T2)
	269	(* -0.00000011 T3)))
	270	(moon-anomaly (+ 201.5643
	271	(* 385.81693528 k)
	272	(* 0.0107438 T2)
	273	(* 0.00001239 T3)
	274	(* -0.000000058 T4)))
	275	(moon-argument (+ 160.7108
	276	(* 390.67050274 k)
	277	(* -0.0016341 T2)
	278	(* -0.00000227 T3)
	279	(* 0.000000011 T4)))
	280	(omega (+ 124.7746
	281	(* -1.56375580 k)
	282	(* 0.0020691 T2)
	283	(* 0.00000215 T3)))
	284	(A1 (+ 299.77 (* 0.107408 k) (* -0.009173 T2)))
	285	(A2 (+ 251.88 (* 0.016321 k)))
	286	(A3 (+ 251.83 (* 26.641886 k)))
	287	(A4 (+ 349.42 (* 36.412478 k)))
	288	(A5 (+ 84.66 (* 18.206239 k)))
	289	(A6 (+ 141.74 (* 53.303771 k)))
	290	(A7 (+ 207.14 (* 2.453732 k)))
	291	(A8 (+ 154.84 (* 7.306860 k)))
	292	(A9 (+ 34.52 (* 27.261239 k)))
	293	(A10 (+ 207.19 (* 0.121824 k)))
	294	(A11 (+ 291.34 (* 1.844379 k)))
	295	(A12 (+ 161.72 (* 24.198154 k)))
	296	(A13 (+ 239.56 (* 25.513099 k)))
	297	(A14 (+ 331.55 (* 3.592518 k)))
	298	(correction
9c0b9187 GM	299	(+ (* -0.40720 (solar-sin-degrees moon-anomaly))
	300	(* 0.17241 E (solar-sin-degrees sun-anomaly))
	301	(* 0.01608 (solar-sin-degrees (* 2 moon-anomaly)))
	302	(* 0.01039 (solar-sin-degrees (* 2 moon-argument)))
	303	(* 0.00739 E (solar-sin-degrees (- moon-anomaly sun-anomaly)))
	304	(* -0.00514 E (solar-sin-degrees (+ moon-anomaly sun-anomaly)))
	305	(* 0.00208 E E (solar-sin-degrees (* 2 sun-anomaly)))
	306	(* -0.00111 (solar-sin-degrees
	307	(- moon-anomaly (* 2 moon-argument))))
	308	(* -0.00057 (solar-sin-degrees
	309	(+ moon-anomaly (* 2 moon-argument))))
	310	(* 0.00056 E (solar-sin-degrees
	311	(+ (* 2 moon-anomaly) sun-anomaly)))
	312	(* -0.00042 (solar-sin-degrees (* 3 moon-anomaly)))
	313	(* 0.00042 E (solar-sin-degrees
	314	(+ sun-anomaly (* 2 moon-argument))))
	315	(* 0.00038 E (solar-sin-degrees
	316	(- sun-anomaly (* 2 moon-argument))))
	317	(* -0.00024 E (solar-sin-degrees
	318	(- (* 2 moon-anomaly) sun-anomaly)))
	319	(* -0.00017 (solar-sin-degrees omega))
	320	(* -0.00007 (solar-sin-degrees
	321	(+ moon-anomaly (* 2 sun-anomaly))))
	322	(* 0.00004 (solar-sin-degrees
	323	(- (* 2 moon-anomaly) (* 2 moon-argument))))
	324	(* 0.00004 (solar-sin-degrees (* 3 sun-anomaly)))
	325	(* 0.00003 (solar-sin-degrees (+ moon-anomaly sun-anomaly
	326	(* -2 moon-argument))))
	327	(* 0.00003 (solar-sin-degrees
	328	(+ (* 2 moon-anomaly) (* 2 moon-argument))))
	329	(* -0.00003 (solar-sin-degrees (+ moon-anomaly sun-anomaly
	330	(* 2 moon-argument))))
	331	(* 0.00003 (solar-sin-degrees (- moon-anomaly sun-anomaly
	332	(* -2 moon-argument))))
	333	(* -0.00002 (solar-sin-degrees (- moon-anomaly sun-anomaly
	334	(* 2 moon-argument))))
	335	(* -0.00002 (solar-sin-degrees
	336	(+ (* 3 moon-anomaly) sun-anomaly)))
	337	(* 0.00002 (solar-sin-degrees (* 4 moon-anomaly)))))
71ea27ee	338	(additional
9c0b9187 GM	339	(+ (* 0.000325 (solar-sin-degrees A1))
	340	(* 0.000165 (solar-sin-degrees A2))
	341	(* 0.000164 (solar-sin-degrees A3))
	342	(* 0.000126 (solar-sin-degrees A4))
	343	(* 0.000110 (solar-sin-degrees A5))
	344	(* 0.000062 (solar-sin-degrees A6))
	345	(* 0.000060 (solar-sin-degrees A7))
	346	(* 0.000056 (solar-sin-degrees A8))
	347	(* 0.000047 (solar-sin-degrees A9))
	348	(* 0.000042 (solar-sin-degrees A10))
	349	(* 0.000040 (solar-sin-degrees A11))
	350	(* 0.000037 (solar-sin-degrees A12))
	351	(* 0.000035 (solar-sin-degrees A13))
	352	(* 0.000023 (solar-sin-degrees A14))))
71ea27ee	353	(newJDE (+ JDE correction additional)))
8a45b040 ER	354	(+ newJDE
8a45b040 ER	355	(- (solar-ephemeris-correction
e803eab7	356	(calendar-extract-year
8a45b040	357	(calendar-gregorian-from-absolute
5c645a20	358	(floor (calendar-astro-to-absolute newJDE))))))
8a45b040 ER	359	(/ calendar-time-zone 60.0 24.0))))
	360
	361	(defun lunar-new-moon-on-or-after (d)
f3924c53 GM	362	"Julian day number of first new moon on or after Julian day number D.
f3924c53 GM	363	The fractional part is the time of day.
8a45b040	364
7b96378f	365	The date and time are local time, including any daylight saving rules,
9c0b9187 GM	366	as governed by the values of `calendar-daylight-savings-starts',
	367	`calendar-daylight-savings-starts-time', `calendar-daylight-savings-ends',
	368	`calendar-daylight-savings-ends-time', `calendar-daylight-time-offset', and
	369	`calendar-time-zone'."
8a45b040	370	(let* ((date (calendar-gregorian-from-absolute
5c645a20	371	(floor (calendar-astro-to-absolute d))))
e803eab7	372	(year (+ (calendar-extract-year date)
71ea27ee GM	373	(/ (calendar-day-number date) 365.25)))
71ea27ee GM	374	(k (floor (* (- year 2000.0) 12.3685)))
4bd7ad5f GM	375	(date (lunar-new-moon-time k))
	376	(a-date (progn
	377	(while (< date d)
	378	(setq k (1+ k)
	379	date (lunar-new-moon-time k)))
5c645a20	380	(calendar-astro-to-absolute date)))
4bd7ad5f GM	381	(time (* 24 (- a-date (truncate a-date))))
	382	(date (calendar-gregorian-from-absolute (truncate a-date)))
	383	(adj (dst-adjust-time date time)))
	384	(calendar-astro-from-absolute
	385	(+ (calendar-absolute-from-gregorian (car adj))
	386	(/ (cadr adj) 24.0)))))
8a45b040	387
7e1dae73 JB	388	(provide 'lunar)
7e1dae73 JB	389
b3a6c0ca	390	;; arch-tag: 72f0b8a4-7bcc-4a1b-b67a-ff53c4a1d222
7e1dae73	391	;;; lunar.el ends here