 | This module is rated as pre-alpha. It is incomplete and may or may not be in active development. Do not use it in article namespace pages. A module remains in pre-alpha until its developer, or another editor who adopts it if it is abandoned for some time, considers the basic structure complete. |
This module is part of an ongoing project to translate the code in pages 361–439 of Calendrical Calculations by Nachum Dershowitz and Edward Reingold from Lisp to Lua. When finished, it will allow any date to be converted between any two of the roughly 30 calendars in the book.
Usage
editI will write a complete manual on how to use this module once I finish making it
Source
edit- Dershowitz, Nachum; Reingold, Edward M. (2008). Calendrical Calculations (Third ed.). Cambridge: Cambridge University Press. ISBN 9780521885409. OCLC 144768713.
-- This code is from pages 361-439 of Calendrical Calculations by Dershowitz and Reingold, translated from Lisp to Lua. Most comments in the code are verbatim from the book
local getArgs = require('Module:Arguments').getArgs
local p = {}
local M = {} -- This namespace is to get around Lua's local variable limit of 200
----- B.1 Lisp Preliminaries -----
M.SUNDAY = 0
M.MONDAY = 1
M.TUESDAY = 2
M.WEDNESDAY = 3
M.THURSDAY = 4
M.FRIDAY = 5
M.SATURDAY = 6
M.BOGUS = "bogus" -- Used to denote nonexistent dates
local function make_public(func, do_unpack)
-- This is not in the book, but is necessary to make the functions work with #invoke on Wikipedia
if do_unpack == nil then do_unpack = true end
return function(frame)
local args = getArgs(frame)
-- Convert arguments from strings to numbers, because it breaks otherwise
local numArgs = {}
for i, v in ipairs(args) do
numArgs[i] = tonumber(v) or v
end
local result
if do_unpack then
-- Pass arguments individually: func(arg1, arg2, ...)
result = func(unpack(numArgs))
else
-- Pass argumemts as a single table: func({arg1, arg2, ...})
result = func(numArgs)
end
if type(result) == "boolean" then
-- Return 1 for true and 0 for false
return result and "1" or "0"
elseif type(result) == "table" then
-- Handle functions that return tables
local idx = args.index
if idx then
-- If the user provided an 'index' parameter in the template, like {{#invoke:Calendrical Calculations|clock_from_moment|739737.12|index=1}}, then use it
local val = result[tonumber(idx) or idx]
-- Apply boolean check to table elements as well in case they contain booleans
if type(val) == "boolean" then return val and "1" or "0" end
return tostring(val)
else
-- Otherwise return the whole table joined by a separator. So {2026, 4, 30} would become "2026-4-30"
return table.concat(result, "-")
end
else
-- Submit normal results as a string
return tostring(result)
end
end
end
-- Identity function for fixed dates/moments. If internal timekeeping is shifted, change epoch to be RD date of origin of internal count. Epoch should be an integer.
function M.rd(tee)
local epoch = 0
return tee - epoch
end
function M.day_of_week_from_fixed(date)
-- The residue class of the day of the week of date
return (date - M.rd(0) - M.SUNDAY) % 7
end
p.day_of_week_from_fixed = make_public(M.day_of_week_from_fixed)
M.JD_EPOCH = M.rd(-1721424.5) -- Fixed time of start of the julian day number
function M.moment_from_jd(jd)
-- Moment of julian day number jd
return jd + M.JD_EPOCH
end
p.moment_from_jd = make_public(M.moment_from_jd)
function M.jd_from_moment(tee)
-- Julian day number of moment tee
return tee - M.JD_EPOCH
end
p.jd_from_moment = make_public(M.jd_from_moment)
function M.fixed_from_jd(jd)
-- Fixed date of julian day number jd
return math.floor(M.moment_from_jd(jd))
end
p.fixed_from_jd = make_public(M.fixed_from_jd)
function M.jd_from_fixed(date)
-- Fixed date of julian day number jd
return M.jd_from_moment(date)
end
p.jd_from_fixed = make_public(M.jd_from_fixed)
M.MJD_EPOCH = M.rd(678576) -- Fixed time of start of the modified julian day number
function M.fixed_from_mjd(mjd)
-- Fixed date of modified julian day number mjd
return mjd + M.MJD_EPOCH
end
p.fixed_from_mjd = make_public(M.fixed_from_mjd)
function M.mjd_from_fixed(date)
-- Modified julian day number of fixed date
return date - M.MJD_EPOCH
end
p.mjd_from_fixed = make_public(M.mjd_from_fixed)
function M.quotient(m, n)
-- Whole part of m/n
return math.floor(m / n)
end
function M.amod(x, y)
-- The value of (x mod y) with y instead of 0
return y + (x % (-y))
end
p.amod = make_public(M.amod)
function M.mod3(x, a, b)
-- The value of x shifted into the range [a..b). Returns x if a=b
if a == b then
return x
else
return a + ((x - a) % (b - a))
end
end
function M.sign(y)
if y < 0 then return -1
elseif y > 0 then return 1
else return 0 end
end
function M.sum(f, i, p)
-- Sum expression for i, i+1, i+2, and so on, as long as condition holds
local s = 0
while p(i) do
s = s + f(i)
i = i + 1
end
return s
end
function M.sigma(list, f)
-- Sum of body for indices il..in running simultaneously thru lists ll..ln.
local s = 0
for _, item in ipairs(list) do
s = s + f(item)
end
return s
end
function M.poly(x, a)
-- Sum powers of x with coefficients (from order 0 up) in list a
local sum = 0
for i = #a, 1, -1 do
sum = (sum * x) + a[i]
end
return sum
end
function M.next(i, p)
-- First integer greater or equal to initial such that condition holds
while not p(i) do
i = i + 1
end
return i
end
function M.final(i, p)
-- Last integer greater or equal to initial such that condition holds
while p(i) do
i = i + 1
end
return i - 1
end
function M.binary_search(lo, hi, test)
-- Bisection search for x in [lo, hi] such that end holds. test determines when to go left
local l = lo
local h = hi
local x
while true do
x = (l + h) / 2
local is_left, is_end = test(x, l, h)
if is_end then
return x
elseif is_left then
h = x
else
l = x
end
end
end
function M.invert_angular(f, y, a, b)
-- Use bisection to find inverse of angular function f at y within interval [a,b]
local varespsilon = 0.00001 -- Desired accuracy
return M.binary_search(a, b, function(x, low, high)
local is_left = ( (f(x) - y) % 360 ) < 180
local is_end = (high - low) < varespsilon
return is_left, is_end
end)
end
----- B.2 Basic Code -----
function M.fixed_from_moment(tee)
-- Fixed-date from moment tee
return math.floor(tee)
end
p.fixed_from_moment = make_public(M.fixed_from_moment)
function M.time_from_moment(tee)
-- Time from moment tee
return tee % 1
end
p.time_from_moment = make_public(M.time_from_moment)
function M.in_range(tee, range)
-- True if tee is in range
return tee >= range[1] and tee <= range[2]
end
function M.list_range(ell, range)
-- Those moments in list ell that occur in range
local results = {}
for _, tee in ipairs(ell) do
if M.in_range(tee, range) then
table.insert(results, tee)
end
end
return results
end
function M.time_from_clock(hms)
-- Time of day from hms = (hour minute second)
return (1/24) * (hms[1] + (hms[2] + (hms[3] / 60)) / 60)
end
p.time_from_clock = make_public(M.time_from_clock, false)
function M.clock_from_moment(tee)
-- Clock time hour:minute: second from moment tee
local time = M.time_from_moment(tee)
local hour = math.floor(time * 24)
local minute = math.floor((time * 24 * 60) % 60)
local second = (time * 24 * 60 * 60) % 60
return {hour, minute, second}
end
p.clock_from_moment = make_public(M.clock_from_moment)
function M.angle_from_degrees(alpha)
-- List of degrees-arcminutes-arcseconds from angle alpha in degrees
local d = math.floor(alpha)
local m = math.floor(60 * (alpha % 1))
local s = (alpha * 60 * 60) % 60
return {d, m, s}
end
p.angle_from_degrees = make_public(M.angle_from_degrees)
----- B.3 The Egyptian/Armenian Calendars -----
M.EGYPTIAN_EPOCH = M.fixed_from_jd(1448638) -- Fixed date of start of the Egyptian (Nabonasser) calendar. JD 1448638 = February 26, 747 BCE (Julian)
function M.fixed_from_egyptian(e_date)
-- Fixed date of Egyptian date e-date. Not to be confused with the kind of e-date that happens on Discord
local year, month, day = e_date[1], e_date[2], e_date[3]
return M.EGYPTIAN_EPOCH -- Days before start of calendar
+ (365 * (year - 1)) -- Days in prior years
+ (30 * (month - 1)) -- Days in prior months this year
+ (day - 1) -- Days so far this month
end
p.fixed_from_egyptian = make_public(M.fixed_from_egyptian, false)
function M.egyptian_from_fixed(date)
-- Egyptian equivalent of fixed date
local days = date - M.EGYPTIAN_EPOCH -- Elapsed days since epoch
local year = M.quotient(days, 365) + 1 -- Years since epoch
local month = 1 + M.quotient((days % 365), 30) -- Calculate the month by division
local day = days - (365 * (year - 1)) - (30 * (month - 1)) + 1 -- Calculate the day by subtraction
return {year, month, day}
end
p.egyptian_from_fixed = make_public(M.egyptian_from_fixed)
M.ARMENIAN_EPOCH = M.rd(201443) -- Fixed date of start of the Armenian calendar. = July 11, 552 CE (Julian)
function M.fixed_from_armenian(a_date)
-- Fixed date of Armenian date a-date
return M.ARMENIAN_EPOCH + M.fixed_from_egyptian(a_date) - M.EGYPTIAN_EPOCH
end
p.fixed_from_armenian = make_public(M.fixed_from_armenian, false)
function M.armenian_from_fixed(date)
-- Armenian equivalent of fixed date
return M.egyptian_from_fixed(date + (M.EGYPTIAN_EPOCH - M.ARMENIAN_EPOCH))
end
p.armenian_from_fixed = make_public(M.armenian_from_fixed)
----- B.4 Cycles of Days -----
function M.kday_on_or_before(k, date)
-- Fixed date of the k-day on or before fixed date
-- k=0 means M.Sunday, k=1 means M.Monday, and so on
return date - M.day_of_week_from_fixed(date - k)
end
p.kday_on_or_before = make_public(M.kday_on_or_before)
function M.kday_on_or_after(k, date)
-- Fixed date of the k-day on or after fixed date
-- k=0 means M.Sunday, k=1 means M.Monday, and so on
return M.kday_on_or_before(k, date + 6)
end
p.kday_on_or_after = make_public(M.kday_on_or_after)
function M.kday_nearest(k, date)
-- Fixed date of the k-day nearest fixed date
-- k=0 means M.Sunday, k=1 means M.Monday, and so on
return M.kday_on_or_before(k, date + 3)
end
p.kday_nearest = make_public(M.kday_nearest)
function M.kday_after(k, date)
-- Fixed date of the k-day after fixed date
-- k=0 means M.Sunday, k=1 means M.Monday, and so on
return M.kday_on_or_before(k, date + 7)
end
p.kday_after = make_public(M.kday_after)
function M.kday_before(k, date)
-- Fixed date of the k-day before fixed date
-- k=0 means M.Sunday, k=1 means M.Monday, and so on
return M.kday_on_or_before(k, date - 1)
end
p.kday_before = make_public(M.kday_before)
----- B.5 The Gregorian Calendar -----
M.GREGORIAN_EPOCH = M.rd(1) -- Fixed date of start of the (proleptic) Gregorian calendar
M.JANUARY = 1 -- January on Julian/Gregorian calendar
M.FEBRUARY = 2 -- February on Julian/Gregorian calendar
M.MARCH = 3 -- March on Julian/Gregorian calendar
M.APRIL = 4 -- April on Julian/Gregorian calendar
M.MAY = 5 -- May on Julian/Gregorian calendar
M.JUNE = 6 -- June on Julian/Gregorian calendar
M.JULY = 7 -- July on Julian/Gregorian calendar
M.AUGUST = 8 -- August on Julian/Gregorian calendar
M.SEPTEMBER = 9 -- September on Julian/Gregorian calendar
M.OCTOBER = 10 -- October on Julian/Gregorian calendar
M.NOVEMBER = 11 -- November on Julian/Gregorian calendar
M.DECEMBER = 12 -- December on Julian/Gregorian calendar
function M.gregorian_leap_year_p(g_year)
-- True if g-year is a leap year on the Gregorian calendar
local mod400 = g_year % 400
return (g_year % 4 == 0) and not (mod400 == 100 or mod400 == 200 or mod400 == 300)
end
p.gregorian_leap_year_p = make_public(M.gregorian_leap_year_p)
function M.fixed_from_gregorian(g_date)
-- Fixed date equivalent to the Gregorian date g-date
local year, month, day = g_date[1], g_date[2], g_date[3]
local correction
if month <= 2 then
correction = 0
elseif M.gregorian_leap_year_p(year) then
correction = -1
else
correction = -2
end
return M.GREGORIAN_EPOCH - 1 -- Days before start of calendar
+ (365 * (year - 1)) -- Ordinary days since epoch
+ M.quotient(year - 1, 4) -- Julian leap days since epoch...
- M.quotient(year - 1, 100) -- ...minus century years since epoch..
+ M.quotient(year - 1, 400) -- ..plus years since epoch divisible by 400
+ M.quotient(367 * month - 362, 12) -- Days in prior months this year assuming 30-day Feb
+ correction -- Correct for 28- or 29-day Feb
+ day -- Days so far this month
end
p.fixed_from_gregorian = make_public(M.fixed_from_gregorian, false)
function M.gregorian_new_year(g_year)
-- Fixed date of January 1 in g-year
return M.fixed_from_gregorian({g_year, M.JANUARY, 1})
end
p.gregorian_new_year = make_public(M.gregorian_new_year)
function M.gregorian_year_end(g_year)
-- Fixed date of December 31 in g-year
return M.fixed_from_gregorian({g_year, M.DECEMBER, 31})
end
p.gregorian_year_end = make_public(M.gregorian_year_end)
function M.gregorian_year_from_fixed(date)
-- Gregorian year corresponding to the fixed date
local d0 = date - M.GREGORIAN_EPOCH -- Prior days
local n400 = M.quotient(d0, 146097) -- Completed 400-year cycles
local d1 = d0 % 146097 -- Prior days not in n400
local n100 = M.quotient(d1, 36524) -- 100-year cycles not in n400
local d2 = d1 % 36524 -- Prior days not in n400 or n100
local n4 = M.quotient(d2, 1461) -- 4-year cycles not in n400 or n100
local d3 = d2 % 1461 -- Prior days not in n400, n100, or n4
local n1 = M.quotient(d3, 365) -- Years not in n400, n100, or n4
local year = (400 * n400) + (100 * n100) + (4 * n4) + n1
if (n100 == 4 or n1 == 4) then
return year -- Date is day 366 in a leap year
else
return year + 1 -- Date is ordinal day ((d3 % 365) + 1) in (year + 1)
end
end
p.gregorian_year_from_fixed = make_public(M.gregorian_year_from_fixed)
function M.gregorian_year_range(g_year)
-- The range of moments in Gregorian year g-year
return {M.gregorian_new_year(g_year), M.gregorian_year_end(g_year)}
end
function M.gregorian_from_fixed(date)
-- Gregorian {year, month, day} corresponding to the fixed date
local year = M.gregorian_year_from_fixed(date)
local prior_days = date - M.gregorian_new_year(year) -- This year
local correction -- To simulate a 30-day Feb
if date < M.fixed_from_gregorian({year, M.MARCH, 1}) then
correction = 0
elseif M.gregorian_leap_year_p(year) then
correction = 1
else
correction = 2
end
local month = M.quotient(12 * (prior_days + correction) + 373, 367) -- Assuming a 30-day Feb
local day = date - M.fixed_from_gregorian({year, month, 1}) + 1 -- Calculate the day by subtraction
return {year, month, day}
end
p.gregorian_from_fixed = make_public(M.gregorian_from_fixed)
function M.gregorian_date_difference(g_date1, g_date2)
-- Number of days from Gregorian date g-date1 until g-date2
return M.fixed_from_gregorian(g_date2) - M.fixed_from_gregorian(g_date1)
end
function M.gregorian_date_difference_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.gregorian_date_difference({table[1], table[2], table[3]}, {table[4], table[5], table[6]})
end
p.gregorian_date_difference = make_public(M.gregorian_date_difference_for_wikipedia, false)
function M.day_number(g_date)
-- Day number in year of Gregorian date g-date
return M.gregorian_date_difference({g_date[1] - 1, M.DECEMBER, 31}, g_date)
end
p.day_number = make_public(M.day_number, false)
function M.days_remaining(g_date)
-- Days remaining in year after Gregorian date g-date
return M.gregorian_date_difference(g_date, {g_date[1], M.DECEMBER, 31})
end
p.days_remaining = make_public(M.days_remaining, false)
function M.alt_fixed_from_gregorian(g_date)
-- Alternative calculation of fixed date equivalent to the Gregorian date g-date
local year, month, day = g_date[1], g_date[2], g_date[3]
local m = M.amod(month - 2, 12)
local y = year + M.quotient(month + 9, 12)
return M.GREGORIAN_EPOCH - 1 - 306 -- Days in March...December
+ (365 * (y - 1)) -- Ordinary days since epoch
+ M.uotient(y - 1, 4) -- Julian leap days since epoch...
- M.quotient(y - 1, 100) -- ...minus century years since epoch...
+ M.quotient(y - 1, 400) -- ...plus years since epoch divisible by 400
+ M.quotient(3 * m - 1, 5) -- Days in prior months this year
+ (30 * (m - 1))
+ day -- Days so far this month
end
p.alt_fixed_from_gregorian = make_public(M.alt_fixed_from_gregorian, false)
function M.alt_gregorian_from_fixed(date)
-- Alternative calculation of Gregorian {year, month, day} corresponding to fixed date
local y = M.gregorian_year_from_fixed((M.GREGORIAN_EPOCH - 1) + date + 306)
local prior_days = date - M.fixed_from_gregorian({y - 1, M.MARCH, 1})
local month = M.amod(M.quotient((5 * prior_days) + 2, 153) + 3, 12)
local year = y - M.quotient(month + 9, 12)
local day = (date - M.fixed_from_gregorian({year, month, 1})) + 1
return {year, month, day}
end
p.alt_gregorian_from_fixed = make_public(M.alt_gregorian_from_fixed)
function M.alt_gregorian_year_from_fixed(date)
-- Gregorian year corresponding to the fixed date
local approx = M.quotient(date - M.GREGORIAN_EPOCH + 2, 146097/400) -- approximate year
local start = M.GREGORIAN_EPOCH + (365 * approx) + M.quotient(approx, 4) - M.quotient(approx, 100) + M.quotient(approx, 400) -- start of next year
if date < start then
return approx
else
return approx + 1
end
end
p.alt_gregorian_year_from_fixed = make_public(M.alt_gregorian_year_from_fixed)
function M.independence_day(g_year)
-- Fixed date of United States Independence Day in Gregorian year g-yaer [sic]
return M.fixed_from_gregorian({g_year, M.JULY, 4})
end
p.independence_day = make_public(M.independence_day)
function M.nth_kday(n, k, g_date)
-- Fixed date of n-th k-day after/before Gregorian date g-date. If n>0, return the n-th k-day on or after g-date. If n<0, return the n-th k-day on or before g-date
-- A k-day of 0 means M.Sunday, 1 means M.Monday, and so on.
if n > 0 then
return (7 * n) + M.kday_before(k, M.fixed_from_gregorian(g_date))
else
return (7 * n) + M.kday_after(k, M.fixed_from_gregorian(g_date))
end
end
function M.nth_kday_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.nth_kday(table[1], table[2], {table[3], table[4], table[5]})
end
p.nth_kday = make_public(M.nth_kday_for_wikipedia, false)
function M.first_kday(k, g_date)
-- Fixed date of first k-day on or after Gregorian date g-date
-- A k-day of 0 means M.Sunday, 1 means M.Monday, and so on
return M.nth_kday(1, k, g_date)
end
function M.first_kday_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.first_kday(table[1], {table[2], table[3], table[4]})
end
p.first_kday = make_public(M.first_kday_for_wikipedia, false)
function M.last_kday(k, g_date)
-- Fixed date of last k-day on or before Gregorian date g-date
-- A k-day of 0 means M.Sunday, 1 means M.Monday, and so on
return M.nth_kday(-1, k, g_date)
end
function M.last_kday_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.last_kday(table[1], {table[2], table[3], table[4]})
end
p.last_kday = make_public(M.last_kday_for_wikipedia, false)
function M.labor_day(g_year)
-- Fixed date of United States Labor Day in Gregorian year g-year (the first Monday in September)
return M.first_kday(M.MONDAY, {g_year, M.SEPTEMBER, 1})
end
p.labor_day = make_public(M.labor_day)
function M.memorial_day(g_year)
-- Fixed date of United States Memorial Day in Gregorian year g-year (the last Monday in May)
return M.last_kday(M.MONDAY, {g_year, M.MAY, 31})
end
p.memorial_day = make_public(M.memorial_day)
function M.election_day(g_year)
-- Fixed date of United States Election Day in Gregorian year g-year (the Tuesday after the first Monday in November)
return M.first_kday(M.TUESDAY, {g_year, M.NOVEMBER, 2})
end
p.election_day = make_public(M.election_day)
function M.daylight_saving_start(g_year)
-- Fixed date of the start of United States daylight saving time in Gregorian year g-year (the second Sunday in March)
return M.nth_kday(2, M.SUNDAY, {g_year, M.MARCH, 1})
end
p.daylight_saving_start = make_public(M.daylight_saving_start)
function M.daylight_saving_end(g_year)
-- Fixed date of of the end of United States daylight saving time in Gregorian year g-year (the first Sunday in November)
return M.first_kday(M.SUNDAY, {g_year, M.NOVEMBER, 1})
end
p.daylight_saving_end = make_public(M.daylight_saving_end)
function M.christmas(g_year)
-- Fixed date of Christmas in Gregorian year g-year
return M.fixed_from_gregorian({g_year, M.DECEMBER, 25})
end
p.christmas = make_public(M.christmas)
function M.advent(g_year)
-- Fixed date of Advent in Gregorian year g-year (the Sunday closest to November 30)
return M.kday_nearest(M.SUNDAY, M.fixed_from_gregorian({g_year, M.NOVEMBER, 30}))
end
p.advent = make_public(M.advent)
function M.epiphany(g_year)
-- Fixed date of Epiphany in U.S. in Gregorian year g-year (the first Sunday after January 1)
return M.first_kday(M.SUNDAY, {g_year, M.JANUARY, 2})
end
p.epiphany = make_public(M.epiphany)
----- B.6 The Julian Calendar -----
function M.julian_leap_year_p(j_year)
-- True if j-year is a leap year on the Julian calendar.
if j_year > 0 then
return (j_year % 4) == 0
else
return (j_year % 4) == 3
end
end
p.julian_leap_year_p = make_public(M.julian_leap_year_p)
M.JULIAN_EPOCH = M.fixed_from_gregorian({0, M.DECEMBER, 30}) -- Fixed date of start of the Julian calendar
function M.fixed_from_julian(j_date)
-- Fixed date equivalent to the Julian date j-date
local year, month, day = j_date[1], j_date[2], j_date[3]
local y
if year < 0 then
y = year + 1 -- No year zero
else
y = year
end
local correction
if month <= 2 then
correction = 0
elseif M.julian_leap_year_p(year) then
correction = -1
else
correction = -2
end
return (M.JULIAN_EPOCH - 1) -- Days before start of calendar
+ (365 * (y - 1)) -- Ordinary days since epoch
+ M.quotient(y - 1, 4) -- Leap days since epoch...
+ M.quotient(367 * month - 362, 12) -- Days in prior months this year assuming 30-day Feb
+ correction -- Correct for 28- or 29-day Feb
+ day -- Days so far this month
end
p.fixed_from_julian = make_public(M.fixed_from_julian, false)
function M.julian_from_fixed(date)
-- Julian {year, month, day} corresponding to fixed date.
local approx = M.quotient((4 * (date - M.JULIAN_EPOCH)) + 1464, 1461) -- Nominal year
local year
if approx <= 0 then
year = approx - 1 -- No year 0
else
year = approx
end
local prior_days = date - M.fixed_from_julian({year, M.JANUARY, 1}) -- This year
local correction -- To simulate a 30-day Feb
if date < M.fixed_from_julian({year, M.MARCH, 1}) then
correction = 0
elseif M.julian_leap_year_p(year) then
correction = 1
else
correction = 2
end
local month = M.quotient((12 * (prior_days + correction)) + 373, 367) -- Assuming a 30-day Feb
local day = 1 + (date - M.fixed_from_julian({year, month, 1})) -- Calculate the day by subtraction
return {year, month, day}
end
p.julian_from_fixed = make_public(M.julian_from_fixed)
M.YEAR_ROME_FOUNDED = -753 -- Year on the Julian calendar of the founding of Rome.
function M.julian_year_from_auc_year(year)
-- Julian year equivalent to AUC year
if 1 <= year and year <= -M.YEAR_ROME_FOUNDED then
return year + M.YEAR_ROME_FOUNDED - 1
else
return year + M.YEAR_ROME_FOUNDED
end
end
p.julian_year_from_auc_year = make_public(M.julian_year_from_auc_year)
function M.auc_year_from_julian_year(year)
-- AUC equivalent to Julian year
if M.YEAR_ROME_FOUNDED <= year and year <= - 1 then
return year - M.YEAR_ROME_FOUNDED + 1
else
return year - M.YEAR_ROME_FOUNDED
end
end
p.auc_year_from_julian_year = make_public(M.auc_year_from_julian_year)
function M.julian_in_gregorian(j_month, j_day, g_year)
-- List of the fixed dates of Julian month j-month, day j-day that occur in Gregorian year g-year
local jan1 = M.gregorian_new_year(g_year)
local y = M.julian_from_fixed(jan1)[1]
local y_prime
if y == -1 then
y_prime = 1
else
y_prime = y + 1
end
-- The possible occurrences in one year are
local date1 = M.fixed_from_julian({y, j_month, j_day})
local date2 = M.fixed_from_julian({y_prime, j_month, j_day})
return M.list_range({date1, date2}, M.gregorian_year_range(g_year))
end
p.julian_in_gregorian = make_public(M.julian_in_gregorian)
-- Since we can only receive numerical arguments from Wikipedia, I made Kalends = 1, Nones = 2 and Ides = 3
M.KALENDS = 1 -- Class of Kalends
M.NONES = 2 -- Class of Nones
M.IDES = 3 -- Class of Ides
function M.ides_of_month(month)
-- Date of Ides in Roman month
local long_months = { [M.MARCH]=true, [M.MAY]=true, [M.JULY]=true, [M.OCTOBER]=true }
if long_months[month] then
return 15
else
return 13
end
end
p.ides_of_month = make_public(M.ides_of_month)
function M.nones_of_month(month)
-- Date of Nones in Roman month
return M.ides_of_month(month) - 8
end
p.nones_of_month = make_public(M.nones_of_month)
function M.fixed_from_roman(r_date)
-- Fixed date for Roman name r-date
local year = r_date[1]
local month = r_date[2]
local event = r_date[3]
local count = r_date[4]
local leap = r_date[5]
if event == M.KALENDS then
local correction1
local correction2
if (month == M.MARCH) and (M.julian_leap_year_p(year)) and 16 >= count and count >= 6 then
correction1 = 0 -- After Ides until leap day
else
correction1 = 1 -- Otherwise
end
if leap then
correction2 = 1 -- Leap day
else
correction2 = 0 -- Non-leap day
end
return (M.fixed_from_julian({year, month, 1}) - count) + correction1 + correction2
elseif event == M.NONES then
return (M.fixed_from_julian({year, month, M.nones_of_month(month)}) - count) + 1
else -- event == M.IDES
return (M.fixed_from_julian({year, month, M.ides_of_month(month)}) - count) + 1
end
end
function M.fixed_from_roman_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.fixed_from_roman({table[1], table[2], table[3], table[4], table[5] == 1})
end
p.fixed_from_roman = make_public(M.fixed_from_roman_for_wikipedia, false)
function M.roman_from_fixed(date)
-- Roman name for fixed date
local j_date = M.julian_from_fixed(date)
local year, month, day = j_date[1], j_date[2], j_date[3]
local month_prime = M.amod(month + 1, 12)
local year_prime
if month_prime ~= 1 then
year_prime = year
elseif year ~= -1 then
year_prime = year + 1
else
year_prime = 1
end
local kalends1 = M.fixed_from_roman({year_prime, month_prime, M.KALENDS, 1, false})
if day == 1 then
return {year, month, M.KALENDS, 1, false}
elseif day <= M.nones_of_month(month) then
return {year, month, M.NONES, (M.nones_of_month(month) - day) + 1, false}
elseif day <= M.ides_of_month(month) then
return {year, month, M.IDES, (M.ides_of_month(month) - day) + 1, false}
elseif (month ~= M.FEBRUARY) or (not M.julian_leap_year_p(year)) then
-- After the Ides, in a month that is not February of a leap year
return {year_prime, month_prime, M.KALENDS, (kalends1 - date) + 1, false}
elseif day < 25 then -- February of a leap year, before leap day
return {year_prime, month_prime, M.KALENDS, 30 - day, false}
else -- February of a leap year, on or after leap day
local count = 31 - day
return {year_prime, month_prime, M.KALENDS, count, (day == 25)}
end
end
p.roman_from_fixed = make_public(M.roman_from_fixed)
----- B.7 The Coptic and Ethiopic Calendars -----
M.COPTIC_EPOCH = M.fixed_from_julian({284, M.AUGUST, 29}) -- Fixed date of start of the Coptic calendar
function M.coptic_leap_year_p(c_year)
-- True if c-year is a leap year on the Coptic calendar
return c_year % 4 == 3
end
p.coptic_leap_year_p = make_public(M.coptic_leap_year_p)
function M.fixed_from_coptic(c_date)
-- Fixed date of Coptic date c-date
local year, month, day = c_date[1], c_date[2], c_date[3]
return M.COPTIC_EPOCH - 1 -- Days before start of calendar
+ 365 * (year - 1) -- Ordinary days in prior years
+ M.quotient(year, 4) -- Leap days in prior years
+ 30 * (month - 1) -- Days in prior months this year
+ day -- Days so far this month
end
p.fixed_from_coptic = make_public(M.fixed_from_coptic, false)
function M.coptic_from_fixed(date)
-- Coptic equivalent of fixed date
local year = M.quotient(4 * (date - M.COPTIC_EPOCH) + 1463, 1461) -- Calculate the year by cycle-of-years formula
local month = M.quotient(date - M.fixed_from_coptic({year, 1, 1}), 30) + 1 -- Calculate the month by division
local day = date + 1 - M.fixed_from_coptic({year, month, 1}) -- Calculate the day by subtraction
return {year, month, day}
end
p.coptic_from_fixed = make_public(M.coptic_from_fixed)
M.ETHIOPIC_EPOCH = M.fixed_from_julian({8, M.AUGUST, 29}) -- Fixed date of start of the Ethiopic calendar
function M.fixed_from_ethiopic(e_date)
-- Fixed date of Ethiopic date e-date
return M.ETHIOPIC_EPOCH + M.fixed_from_coptic(e_date) - M.COPTIC_EPOCH
end
p.fixed_from_ethiopic = make_public(M.fixed_from_ethiopic, false)
function M.ethiopic_from_fixed(date)
-- Ethiopic equivalent of fixed date
return M.coptic_from_fixed(date + (M.COPTIC_EPOCH - M.ETHIOPIC_EPOCH))
end
p.ethiopic_from_fixed = make_public(M.ethiopic_from_fixed)
function M.coptic_in_gregorian(c_month, c_day, g_year)
-- List of the fixed dates of Coptic month c-month, day c-day that occur in Gregorian year g-year
local jan1 = M.gregorian_new_year(g_year)
local y = M.coptic_from_fixed(jan1)[1]
-- The possible occurrences in one year are
local date0 = M.fixed_from_coptic({y, c_month, c_day})
local date1 = M.fixed_from_coptic({y + 1, c_month, c_day})
return M.list_range({date0, date1}, M.gregorian_year_range(g_year))
end
p.coptic_in_gregorian = make_public(M.coptic_in_gregorian)
function M.coptic_christmas(g_year)
-- List of zero or one fixed dates of Coptic Christmas in Gregorian year g-year
return M.coptic_in_gregorian(4, 29, g_year)
end
p.coptic_christmas = make_public(M.coptic_christmas)
----- B.8 The ISO Calendar -----
function M.fixed_from_iso(i_date)
-- Fixed date equivalent to ISO i-date
local year, week, day = i_date[1], i_date[2], i_date[3]
-- Add fixed date of Sunday preceding date plus day in week
return M.nth_kday(week, M.SUNDAY, {year - 1, M.DECEMBER, 28}) + day
end
p.fixed_from_iso = make_public(M.fixed_from_iso, false)
function M.iso_from_fixed(date)
-- ISO (year week day) corresponding to the fixed date
local approx = M.gregorian_year_from_fixed(date - 3) -- Year may be one too small
local year
if date >= M.fixed_from_iso({approx + 1, 1, 1}) then
year = approx + 1
else
year = approx
end
local week = M.quotient(date - M.fixed_from_iso({year, 1, 1}), 7) + 1
local day = M.amod(date - M.rd(0), 7)
return {year, week, day}
end
p.iso_from_fixed = make_public(M.iso_from_fixed)
function M.iso_long_year_p(i_year)
-- True if i-year is a long (53-week) year
local jan1 = M.day_of_week_from_fixed(M.gregorian_new_year(i_year))
local dec31 = M.day_of_week_from_fixed(M.gregorian_year_end(i_year))
return jan1 == M.THURSDAY or dec31 == M.THURSDAY
end
p.iso_long_year_p = make_public(M.iso_long_year_p)
----- B.9 The Islamic Calendar -----
----- B.10 The Hebrew Calendar -----
----- B.11 The Ecclesiastical Calendars -----
function M.eastern_orthodox_christmas(g_year)
-- List of zero or one fixed dates of Eastern Orthodox Christmas in Gregorian year g-year
return M.julian_in_gregorian(M.DECEMBER, 25, g_year)
end
p.eastern_orthodox_christmas = make_public(M.eastern_orthodox_christmas)
function M.orthodox_easter(g_year)
-- Fixed date of Orthodox Easter in Gregorian year g-year.
local shifted_epact = (14 + 11 * (g_year % 19)) % 30 -- Age of moon for April 5
local j_year -- Julian year number
if g_year > 0 then j_year = g_year else j_year = g_year - 1 end
local paschal_moon = M.fixed_from_julian({j_year, M.APRIL, 19}) - shifted_epact -- Day after full moon on or after March 21
-- Return the Sunday following the Paschal moon
return M.kday_after(M.SUNDAY, paschal_moon)
end
p.orthodox_easter = make_public(M.orthodox_easter)
function M.alt_orthodox_easter(g_year)
-- Alternative calculation of fixed date of Orthodox Easter in Gregorian year g-year
local paschal_moon = 354 * g_year -- Day after full moon on or after March 21
+ 30 * M.quotient(7 * g_year + 8, 19)
+ M.quotient(g_year, 4)
- M.quotient(g_year, 19)
- 273
+ M.GREGORIAN_EPOCH
-- Return the Sunday following the Paschal moon
return M.kday_after(M.SUNDAY, paschal_moon)
end
p.alt_orthodox_easter = make_public(M.alt_orthodox_easter)
function M.easter(g_year)
-- Fixed date of Easter in Gregorian year g-year
local century = M.quotient(g_year, 100) + 1
local shifted_epact = (14 + 11 * (g_year % 19) -- Age of moon for April 5 by Nicaean rule
- M.quotient(3 * century, 4) -- ...corrected for the Gregorian century rule
+ M.quotient(5 + 8 * century, 25)) % 30 -- ...corrected for Metonic cycle inaccuracy
local adjusted_epact -- Adjust for 29.5 day month
if shifted_epact == 0 or (shifted_epact == 1 and (g_year % 19) > 10) then
adjusted_epact = shifted_epact + 1
else
adjusted_epact = shifted_epact
end
local paschal_moon = M.fixed_from_gregorian({g_year, M.APRIL, 19}) - adjusted_epact -- Day after full moon on or after March 21
-- Return the Sunday following the Paschal moon
return M.kday_after(M.SUNDAY, paschal_moon)
end
p.easter = make_public(M.easter)
function M.pentecost(g_year)
-- Fixed date of Pentecost in Gregorian year g-year
return M.easter(g_year) + 49
end
p.pentecost = make_public(M.pentecost)
----- B.12 The Old Hindu Calendars -----
M.HINDU_EPOCH = M.fixed_from_julian({-3102, M.FEBRUARY, 18}) -- Fixed date of start of the Hindu calendar (Kali Yuga)
function M.hindu_day_count(date)
return date - M.HINDU_EPOCH -- Elapsed days (Ahargana) to date since Hindu epoch (KY)
end
M.ARYA_JOVIAN_PERIOD = 1577917500/364224 -- Number of days in one revolution of Jupiter around the Sun
function M.jovian_year(date)
-- Year of Jupiter cycle at fixed date
return M.amod(27 + M.quotient(M.hindu_day_count(date), M.ARYA_JOVIAN_PERIOD / 12), 60)
end
M.ARYA_SOLAR_YEAR = 1577917500/4320000 -- Length of Old Hindu solar year
M.ARYA_SOLAR_MONTH = M.ARYA_SOLAR_YEAR / 12 -- Length of Old Hindu solar month
function M.old_hindu_solar_from_fixed(date)
-- Old Hindu solar date equivalent to fixed date
local sun = M.hindu_day_count(date) + M.hr(6) -- Sunrise on Hindu date
local year = M.quotient(sun, M.ARYA_SOLAR_YEAR) -- Elapsed years
local month = M.quotient(sun, M.ARYA_SOLAR_MONTH) % 12 + 1
local day = math.floor(sun % M.ARYA_SOLAR_MONTH) + 1
return {year, month, day}
end
p.old_hindu_solar_from_fixed = make_public(M.old_hindu_solar_from_fixed)
function M.fixed_from_old_hindu_solar(s_date)
-- Fixed date corresponding to Old Hindu solar date s-date
local year, month, day = s_date[1], s_date[2], s_date[3]
return math.ceil(M.HINDU_EPOCH -- Since start of era
+ year * M.ARYA_SOLAR_YEAR -- Days in elapsed years
+ (month - 1) * M.ARYA_SOLAR_MONTH -- ...in months
+ day + M.hr(-30)) -- Midnight of day
end
p.fixed_from_old_hindu_solar = make_public(M.fixed_from_old_hindu_solar, false)
M.ARYA_LUNAR_MONTH = 1577917500/53433336 -- Length of Old Hindu lunar month
M.ARYA_LUNAR_DAY = M.ARYA_LUNAR_MONTH / 30 -- Length of Old Hindu lunar day
function M.old_hindu_lunar_leap_year(l_year)
-- True if l-year is a leap year on the old Hindu calendar
return ((l_year * M.ARYA_SOLAR_YEAR - M.ARYA_SOLAR_MONTH) % M.ARYA_LUNAR_MONTH) >= 23902504679/1282400064
end
function M.old_hindu_lunar_from_fixed(date)
-- Old Hindu lunar date equivalent to fixed date
local sun = M.hindu_day_count(date) + M.hr(6) -- Sunrise on Hindu date
local new_moon = sun - (sun % M.ARYA_LUNAR_MONTH) -- Beginning of lunar month
local leap = (M.ARYA_SOLAR_MONTH - M.ARYA_LUNAR_MONTH >= (new_moon % M.ARYA_SOLAR_MONTH))
and ((new_moon % M.ARYA_SOLAR_MONTH) > 0) -- If lunar contained in solar
local month = math.ceil(new_moon / M.ARYA_SOLAR_MONTH) % 12 + 1 -- Next solar month's name
local day = M.quotient(sun, M.ARYA_LUNAR_DAY) % 30 + 1 -- Lunar days since beginning of lunar month
local year = math.ceil((new_moon + M.ARYA_SOLAR_MONTH) / M.ARYA_SOLAR_YEAR) - 1 -- Solar year at end of lunar month(s)
return {year, month, leap, day}
end
p.old_hindu_lunar_from_fixed = make_public(M.old_hindu_lunar_from_fixed)
function M.fixed_from_old_hindu_lunar(l_date)
-- Fixed date corresponding to Old Hindu lunar date l-date
local year = M.old_hindu_lunar_year(l_date)
local month = M.old_hindu_lunar_month(l_date)
local leap = M.old_hindu_lunar_leap(l_date)
local day = M.old_hindu_lunar_day(l_date)
local mina = (12 * year - 1) * M.ARYA_SOLAR_MONTH -- One solar month before solar new year
local lunar_new_year = M.ARYA_LUNAR_MONTH * (M.quotient(mina, M.ARYA_LUNAR_MONTH) + 1) -- New moon after mina
local m
if (not leap) and math.ceil((lunar_new_year - mina) / (M.ARYA_SOLAR_MONTH - M.ARYA_LUNAR_MONTH)) <= month then
m = month
else
m = month - 1
end
return math.ceil(M.HINDU_EPOCH
+ lunar_new_year
+ M.ARYA_LUNAR_MONTH * m
+ (day - 1) * M.ARYA_LUNAR_DAY -- Lunar days
+ M.hr(-6)) -- Subtract 1 if phase begins before sunrise
end
----- B.13 The Mayan Calendars -----
-- Fixed date of start of the Mayan calendar, according to the Goodman-Martinez-Thompson correlation.
-- That is, August 11, -3113
M.MAYAN_EPOCH = M.fixed_from_jd(584283)
function M.fixed_from_mayan_long_count(count)
-- Fixed date corresponding to the Mayan long count, which is a list (baktun katun tun uinal kin)
local baktun, katun, tun, uinal, kin = count[1], count[2], count[3], count[4], count[5]
return M.MAYAN_EPOCH + ((((baktun * 20) + katun) * 20 + tun) * 18 + uinal) * 20 + kin
end
p.fixed_from_mayan_long_count = make_public(M.fixed_from_mayan_long_count, false)
function M.mayan_long_count_from_fixed(date)
-- Mayan long count date of fixed date
local long_count = date - M.MAYAN_EPOCH
local baktun = M.quotient(long_count, 144000)
local day_of_baktun = long_count % 144000
local katun = M.quotient(day_of_baktun, 7200)
local day_of_katun = day_of_baktun % 7200
local tun = M.quotient(day_of_katun, 360)
local day_of_tun = day_of_katun % 360
local uinal = M.quotient(day_of_tun, 20)
local kin = day_of_tun % 20
return {baktun, katun, tun, uinal, kin}
end
p.mayan_long_count_from_fixed = make_public(M.mayan_long_count_from_fixed)
function M.mayan_haab_ordinal(h_date)
-- Number of days into haab cycle of Mayan haab date h-date
local day, month = h_date[2], h_date[1]
return (month - 1) * 20 + day
end
M.MAYAN_HAAB_EPOCH = M.MAYAN_EPOCH - M.mayan_haab_ordinal({18, 8}) -- Fixed date of start of haab cycle
function M.mayan_haab_from_fixed(date)
-- Mayan haab date of fixed date
local count = (date - M.MAYAN_HAAB_EPOCH) % 365
local day = count % 20
local month = M.quotient(count, 20) + 1
return {month, day}
end
p.mayan_haab_from_fixed = make_public(M.mayan_haab_from_fixed)
function M.mayan_haab_on_or_before(haab, date)
-- Fixed date of latest date on or before fixed date that is Mayan haab date haab
return M.mod3(M.mayan_haab_ordinal(haab) + M.MAYAN_HAAB_EPOCH, date - 365, date)
end
function M.mayan_tzolkin_ordinal(t_date)
-- Number of days into Mayan tzolkin cycle of t-date
local number, name = t_date[1], t_date[2]
return (number - 1 + 39 * (number - name)) % 260
end
M.MAYAN_TZOLKIN_EPOCH = M.MAYAN_EPOCH - M.mayan_tzolkin_ordinal({4, 20}) -- Start of tzolkin date cycle
function M.mayan_tzolkin_from_fixed(date)
-- Mayan tzolkin date of fixed date
local count = date - M.MAYAN_TZOLKIN_EPOCH + 1
local number = M.amod(count, 13)
local name = M.amod(count, 20)
return {number, name}
end
p.mayan_tzolkin_from_fixed = make_public(M.mayan_tzolkin_from_fixed)
function M.mayan_tzolkin_on_or_before(tzolkin, date)
-- Fixed date of latest date on or before fixed date that is Mayan tzolkin date tzolkin
return M.mod3(M.mayan_tzolkin_ordinal(tzolkin) + M.MAYAN_TZOLKIN_EPOCH, date - 260, date)
end
function M.mayan_year_bearer_from_fixed(date)
-- Year bearer of year containing fixed date. Returns BOGUS for uayeb
local x = M.mayan_haab_on_or_before({1, 0}, date)
if M.mayan_haab_month(M.mayan_haab_from_fixed(date)) == 19 then
return M.BOGUS
else
return M.mayan_tzolkin_from_fixed(x)[2]
end
end
p.mayan_year_bearer_from_fixed = make_public(M.mayan_year_bearer_from_fixed)
function M.mayan_calendar_round_on_or_before(haab, tzolkin, date)
-- Fixed date of latest date on or before date, that is Mayan haab date haab and an tzolkin date tzolkin.
-- Returns bogus for impossible combinations
local haab_count = M.mayan_haab_ordinal(haab) + M.MAYAN_HAAB_EPOCH
local tzolkin_count = M.mayan_tzolkin_ordinal(tzolkin) + M.MAYAN_TZOLKIN_EPOCH
local diff = tzolkin_count - haab_count
if diff % 5 == 0 then
return M.mod3(haab_count + 365 * diff, date - 18980, date)
else
return M.BOGUS -- haab-tzolkin combination is impossible
end
end
-- Aztec calendars
M.AZTEC_CORRELATION = M.fixed_from_julian({1521, M.AUGUST, 13}) -- Known date of Aztec cycles (Caso's correlation)
function M.aztec_xihuitl_ordinal(x_date) -- Number of elapsed days into cycle of Aztec xihuitl x-date
local day, month = x_date[2], x_date[1]
return (month - 1) * 20 + (day - 1)
end
M.AZTEC_XIHUITL_CORRELATION = M.AZTEC_CORRELATION - M.aztec_xihuitl_ordinal({11, 2}) -- Start of a xihuitl cycle
function M.aztec_xihuitl_from_fixed(date)
-- Aztec xihuitl date of fixed date
local count = (date - M.AZTEC_XIHUITL_CORRELATION) % 365
local day = (count % 20) + 1
local month = M.quotient(count, 20) + 1
return {month, day}
end
p.aztec_xihuitl_from_fixed = make_public(M.aztec_xihuitl_from_fixed)
function M.aztec_xihuitl_on_or_before(xihuitl, date)
-- Fixed date of latest date on or before date that is Aztec xihuitl xihuitl
return M.mod3(M.AZTEC_XIHUITL_CORRELATION + M.aztec_xihuitl_ordinal(xihuitl), date - 365, date)
end
function M.aztec_tonalpohualli_ordinal(t_date)
-- Number of days into Aztec tonalpohualli cycle of t-date
local number, name = t_date[1], t_date[2]
return (number - 1 + 39 * (number - name)) % 260
end
M.AZTEC_TONALPOHUALLI_CORRELATION = M.AZTEC_CORRELATION - M.aztec_tonalpohualli_ordinal({1, 5}) -- Start of a tonalpohualli date cycle
function M.aztec_tonalpohualli_from_fixed(date)
-- Aztec tonalpohualli date of fixed date
local count = date - M.AZTEC_TONALPOHUALLI_CORRELATION + 1
local number = M.amod(count, 13)
local name = M.amod(count, 20)
return {number, name}
end
p.aztec_tonalpohualli_from_fixed = make_public(M.aztec_tonalpohualli_from_fixed)
function M.aztec_tonalpohualli_on_or_before(tonalpohualli, date)
-- Fixed date of latest date on or before date that is Aztec tonalpohualli tonalpohualli
return M.mod3(M.AZTEC_TONALPOHUALLI_CORRELATION + M.aztec_tonalpohualli_ordinal(tonalpohualli),
date - 260, date)
end
function M.aztec_xiuhmolpilli_from_fixed(date)
-- Designation of year containing fixed date. Returns BOGUS for nemontemi
local x = M.aztec_xihuitl_on_or_before({18, 20}, date + 364)
local month = M.aztec_xihuitl_from_fixed(date)[1]
if month == 19 then return M.BOGUS
else return M.aztec_tonalpohualli_from_fixed(x) end
end
function M.aztec_xihuitl_tonalpohualli_on_or_before(xihuitl, tonalpohualli, date)
-- Fixed date of latest xihuitl-tonalpohualli combination on or before date. That is the date on or before date that is Aztec xihuitl date xihuitl and tonalpohualli date tonalpohualli
-- Returns bogus for impossible combinations.
local xihuitl_count = M.aztec_xihuitl_ordinal(xihuitl) + M.AZTEC_XIHUITL_CORRELATION
local tonalpohualli_count = M.aztec_tonalpohualli_ordinal(tonalpohualli) + M.AZTEC_TONALPOHUALLI_CORRELATION
local diff = tonalpohualli_count - xihuitl_count
if diff % 5 == 0 then
return M.mod3(xihuitl_count + 365 * diff, date - 18980, date)
else
return M.BOGUS -- xihuitl-tonalpohualli combination is impossible
end
end
----- B.14 The Balinese Pawukon Calendar -----
----- B.15 Time and Astronomy -----
-- I had to rearrange some of the functions here, because the code in the book often calls constants and functions before defining them, which is not possible in Lua
function M.degrees(theta)
-- Normalize angle theta to range [0,360) degrees
return theta % 360
end
p.degrees = make_public(M.degrees)
function M.radians_from_degrees(theta)
-- Convert angle theta from degrees to radians
return M.degrees(theta) * math.pi / 180
end
p.radians_from_degrees = make_public(M.radians_from_degrees)
function M.degrees_from_radians(theta)
-- Convert angle theta from degrees to radians
return M.degrees(theta * 180 / math.pi)
end
p.degrees_from_radians = make_public(M.degrees_from_radians)
function M.sin_degrees(theta) return math.sin(M.radians_from_degrees(theta)) end -- Sine of theta (given in degrees)
function M.cosine_degrees(theta) return math.cos(M.radians_from_degrees(theta)) end -- Cosine of theta (given in degrees)
function M.tangent_degrees(theta) return math.tan(M.radians_from_degrees(theta)) end -- Tangent of theta (given in degrees)
function M.arctan_degrees(y, x)
-- Arctangent of y/x in degrees
local result
if x == 0 and y ~= 0 then
result = M.sign(y) * 90
else
local alpha = M.degrees_from_radians(math.atan(y / x))
if x >= 0 then result = alpha
else result = alpha + 180 end
end
return result % 360
end
function M.arcsin_degrees(x) return M.degrees_from_radians(math.asin(x)) end -- Arcsine of x in degrees
function M.arccos_degrees(x) return M.degrees_from_radians(math.acos(x)) end -- Arccosine of x in degrees
function M.hr(x) return x / 24 end -- x hours
function M.mn(x) return x / (24 * 60) end -- x minutes
function M.sec(x) return x / (24 * 60 * 60) end -- x seconds
function M.mins(x) return x / 60 end -- x arcminutes -> degrees
function M.secs(x) return x / 3600 end -- x arcseconds -> degrees
function M.angle(d, m, s)
-- d degrees, m arcminutes, s arcseconds
return d + (m + s/60) / 60
end
p.angle = make_public(M.angle)
M.MECCA = {M.angle(21, 25, 24), M.angle(39, 49, 24), 298, M.hr(2)}
M.JERUSALEM = {31.8, 35.2, 298, M.hr(2)}
function M.direction(locale, focus)
-- Angle (clockwise from North) to face focus when standing in locale. Subject to errors near focus and its antipode
local phi = locale[1]
local phi_prime = focus[1]
local psi = locale[2]
local psi_prime = focus[2]
local y = M.sin_degrees(psi_prime - psi)
local x = M.cosine_degrees(phi) * M.tangent_degrees(phi_prime)
- M.sin_degrees(phi) * M.cosine_degrees(psi - psi_prime)
if (x == 0 and y == 0) or phi_prime == 90 then return 0
elseif phi_prime == -90 then return 180
else return M.arctan_degrees(y, x) end
end
function M.direction_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.direction({table[1], table[2], table[3], table[4], table[5]}, {table[6], table[7], table[8], table[9], table[10]})
end
p.direction = make_public(M.direction_for_wikipedia, false)
function M.zone_from_longitude(phi)
-- Difference between UT and local mean time at longitude phi (fraction of day)
return phi / 360
end
function M.universal_from_local(tee_ell, locale)
-- Universal time from local tee_ell at locale
return tee_ell - M.zone_from_longitude(locale[2])
end
function M.local_from_universal(tee_rom_u, locale)
-- Local time from universal tee_rom-u at locale
return tee_rom_u + M.zone_from_longitude(locale[2])
end
function M.standard_from_universal(tee_rom_u, locale)
-- Standard time from tee_rom-u in universal time at locale
return tee_rom_u + locale[4]
end
function M.universal_from_standard(tee_rom_s, locale)
-- Universal time from tee_rom-s in standard time at locale
return tee_rom_s - locale[4]
end
function M.standard_from_local(tee_ell, locale)
-- Standard time from local tee_ell at locale
return M.standard_from_universal(M.universal_from_local(tee_ell, locale), locale)
end
function M.local_from_standard(tee_rom_s, locale)
-- Local time from standard tee_rom-s at locale
return M.local_from_universal(M.universal_from_standard(tee_rom_s, locale), locale)
end
M.ephemeris_correction = function(tee)
-- Dynamical Time minus Universal Time (in days) for moment tee
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 1991
local year = M.gregorian_year_from_fixed(math.floor(tee))
local c = M.gregorian_date_difference({1900, M.JANUARY, 1}, {year, M.JULY, 1}) / 36525
local x = M.hr(12) + M.gregorian_date_difference({1810, M.JANUARY, 1}, {year, M.JANUARY, 1})
if 1988 <= year and year <= 2019 then
return (1/86400) * (year - 1933)
elseif 1900 <= year and year <= 1987 then
return M.poly(c, {-0.00002, 0.000297, 0.025184, -0.181133, 0.553040, -0.861938, 0.677066, -0.212591})
elseif 1800 <= year and year <= 1899 then
return M.poly(c, {-0.000009, 0.003844, 0.083563, 0.865736, 4.867575, 15.845535, 31.332267, 38.291999, 28.316289, 11.636204, 2.043794})
elseif 1700 <= year and year <= 1799 then
return (1/86400) * M.poly((year - 1700), {8.118780842, -0.005092142, 0.003336121, -0.0000266484})
elseif 1620 <= year and year <= 1699 then
return (1/86400) * M.poly((year - 1600), {196.58333, -4.0675, 0.0219167})
else
return (1/86400) * ((x * x) / 41048480 - 15)
end
end
function M.dynamical_from_universal(tee)
-- Dynamical time at Universal moment tee
return tee + M.ephemeris_correction(tee)
end
function M.universal_from_dynamical(tee)
-- Universal moment from Dynamical time tee
return tee - M.ephemeris_correction(tee)
end
M.J2000 = M.hr(12) + M.gregorian_new_year(2000) -- Noon at start of Gregorian year 2000
function M.julian_centuries(tee)
-- Julian centuries since 2000 at moment tee
return (M.dynamical_from_universal(tee) - M.J2000) / 36525
end
p.julian_centuries = make_public(M.julian_centuries)
function M.obliquity(tee)
-- Obliquity of ecliptic at moment tee
local c = M.julian_centuries(tee)
return M.angle(23, 26, 21.448)
+ M.poly(c, {0, M.angle(0, 0, -46.8150), M.angle(0, 0, -0.00059), M.angle(0, 0, 0.001813)})
end
function M.declination(tee, beta, lambda)
-- Declination at moment UT tee of object at latitude beta and longitude lambda
local varespsilon = M.obliquity(tee)
return M.arcsin_degrees(M.sin_degrees(beta) * M.cosine_degrees(varespsilon)
+ M.cosine_degrees(beta) * M.sin_degrees(varespsilon) * M.sin_degrees(lambda))
end
function M.right_ascension(tee, beta, lambda)
-- Right ascension at moment UT tee of object at latitude beta and longitude lambda
-- Note: The comment in the book confusingly says "latitude lambda and longitude beta" here. But this seems to be a typo, as it doesn't match with page 186
local varespsilon = M.obliquity(tee)
return M.arctan_degrees(
M.sin_degrees(lambda) * M.cosine_degrees(varespsilon) - M.tangent_degrees(beta) * M.sin_degrees(varespsilon),
M.cosine_degrees(lambda))
end
M.equation_of_time = function(tee)
-- Equation of time (as fraction of day) for moment tee.
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 1991
local c = M.julian_centuries(tee)
local lambda = M.poly(c, {280.46645, 36000.76983, 0.0003032})
local anomaly = M.poly(c, {357.52910, 35999.05030, -0.0001559, -0.00000048})
local eccentricity = M.poly(c, {0.016708617, -0.000042037, -0.0000001236})
local varespsilon = M.obliquity(tee)
local y = M.tangent_degrees(varespsilon / 2) ^ 2
local equation = (1 / (2 * math.pi)) * (
y * M.sin_degrees(2 * lambda)
- 2 * eccentricity * M.sin_degrees(anomaly)
+ 4 * eccentricity * y * M.sin_degrees(anomaly) * M.cosine_degrees(2 * lambda)
- 0.5 * y * y * M.sin_degrees(4 * lambda)
- 1.25 * eccentricity * eccentricity * M.sin_degrees(2 * anomaly))
return M.sign(equation) * math.min(math.abs(equation), M.hr(12))
end
p.equation_of_time = make_public(M.equation_of_time)
function M.apparent_from_local(tee, locale)
-- Sundial time at local time tee at locale
return tee + M.equation_of_time(M.universal_from_local(tee, locale))
end
function M.local_from_apparent(tee, locale)
-- Local time from sundial time tee at locale
return tee - M.equation_of_time(M.universal_from_local(tee, locale))
end
function M.midday(date, locale)
-- Standard time on fixed date of midday at locale
return M.standard_from_local(M.local_from_apparent(date + M.hr(12), locale), locale)
end
function M.midnight(date, locale)
-- Standard time on fixed date of true (apparent) midnight at locale
return M.standard_from_local(M.local_from_apparent(date, locale), locale)
end
function M.sidereal_from_moment(tee)
-- Mean sidereal time of day from moment tee expressed as hour angle
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 1991
local c = (tee - M.J2000) / 36525
return M.poly(c, {280.46061837, 36525 * 360.98564736629, 0.000387933, -1/38710000}) % 360
end
M.nutation = function(tee)
-- Longitudinal nutation at moment tee
local c = M.julian_centuries(tee) -- moment in Julian centuries
local cap_A = M.poly(c, {124.90, -1934.134, 0.002063})
local cap_B = M.poly(c, {201.11, 72001.5377, 0.00057})
return -0.004778 * M.sin_degrees(cap_A) - 0.0003667 * M.sin_degrees(cap_B)
end
M.aberration = function(tee)
-- Aberration at moment tee
local c = M.julian_centuries(tee) -- moment in Julian centuries
return 0.0000974 * M.cosine_degrees(177.63 + 35999.01848 * c) - 0.005575
end
M.MEAN_TROPICAL_YEAR_AST = 365.242189
M.MEAN_SIDEREAL_YEAR = 365.25636
M.solar_longitude = function(tee)
-- Longitude of sun at moment tee
-- Adapted from "Planetary Programs and Tables from -4000 to +2800" by Pierre Bretagnon and Jean-Louis Simon, Willmann-Bell, Inc., 1986
local c = M.julian_centuries(tee) -- moment in Julian centuries
local coefficients = {
403406, 195207, 119433, 112392, 3891, 2819, 1721,
660, 350, 334, 314, 268, 242, 234, 158, 132, 129, 114,
99, 93, 86, 78, 72, 68, 64, 46, 38, 37, 32, 29, 28, 27, 27,
25, 24, 21, 21, 20, 18, 17, 14, 13, 13, 13, 12, 10, 10, 10, 10
}
local multipliers = {
0.9287892, 35999.1376958, 35999.4089666,
35998.7287385, 71998.20261, 71998.4403,
36000.35726, 71997.4812, 32964.4678,
-19.4410, 445267.1117, 45036.8840, 3.1008,
22518.4434, -19.9739, 65928.9345,
9038.0293, 3034.7684, 33718.148, 3034.448,
-2280.773, 29929.992, 31556.493, 149.588,
9037.750, 107997.405, -4444.176, 151.771,
67555.316, 31556.080, -4561.540,
107996.706, 1221.655, 62894.167,
31437.369, 14578.298, -31931.757,
34777.243, 1221.999, 62894.511,
-4442.039, 107997.909, 119.066, 16859.071,
-4.578, 26895.292, -39.127, 12297.536, 90073.778
}
local addends = {
270.54861, 340.19128, 63.91854, 331.26220,
317.843, 86.631, 240.052, 310.26, 247.23,
260.87, 297.82, 343.14, 166.79, 81.53,
3.50, 132.75, 182.95, 162.03, 29.8,
266.4, 249.2, 157.6, 257.8, 185.1, 69.9,
8.0, 197.1, 250.4, 65.3, 162.7, 341.5,
291.6, 98.5, 146.7, 110.0, 5.2, 342.6,
230.9, 256.1, 45.3, 242.9, 115.2, 151.8,
285.3, 53.3, 126.6, 205.7, 85.9, 146.1
}
local s = 0
for i = 1, #coefficients do
s = s + coefficients[i] * M.sin_degrees(addends[i] + multipliers[i] * c)
end
local lambda = 282.7771834 + 36000.76953744 * c + 0.000005729577951308232 * s
return (lambda + M.aberration(tee) + M.nutation(tee)) % 360
end
M.SPRING = 0 -- Longitude of sun at vernal equinox
M.SUMMER = 90 -- Longitude of sun at summer solstice
M.AUTUMN = 180 -- Longitude of sun at autumnal equinox
M.WINTER = 270 -- Longitude of sun at winter solstice
function M.solar_longitude_after(lambda, tee)
-- Moment UT of first time at or after tee when the solar longitude will be lambda degrees
local rate = M.MEAN_TROPICAL_YEAR_AST / 360 -- Mean days for 1 degree change
local tau = tee + rate * ((lambda - M.solar_longitude(tee)) % 360) -- Estimate within 5 days
local a = math.max(tee, tau - 5) -- At or after tee
local b = tau + 5
return M.invert_angular(M.solar_longitude, lambda, a, b)
end
function M.precession(tee)
-- Precession at moment tee using 0,0 as J2000 coordinates
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 1991
local c = M.julian_centuries(tee)
local eta = M.poly(c, {0, M.secs(47.0029), M.secs(-0.03302), M.secs(0.000060)}) % 360
local cap_P = M.poly(c, {174.876384, M.secs(-869.8089), M.secs(0.03536)}) % 360
local p = M.poly(c, {0, M.secs(5029.0966), M.secs(1.11113), M.secs(0.000006)}) % 360
local cap_A = M.cosine_degrees(eta) * M.sin_degrees(cap_P)
local cap_B = M.cosine_degrees(cap_P)
local arg = M.arctan_degrees(cap_A, cap_B)
return ((p + cap_P) - arg) % 360
end
-- Skip sidereal_solar_longitude for now, as I need to do the Hindu calendars for it to work
function M.estimate_prior_solar_longitude(lambda, tee)
-- Approximate moment at or before tee when solar longitude just exceeded lambda degrees
local rate = M.MEAN_TROPICAL_YEAR_AST / 360 -- Mean change per degree
local tau = tee - rate * ((M.solar_longitude(tee) - lambda) % 360) -- First approximation
local cap_Delta = M.mod3(M.solar_longitude(tau) - lambda, -180, 180) -- Difference in longitude
return math.min(tee, tau - rate * cap_Delta)
end
M.MEAN_SYNODIC_MONTH = 29.530588853
function M.nth_new_moon(n)
-- Moment of n-th new moon after (or before) the new moon of January 11, 1
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
local n0 = 24724 -- Months from RD 0 until j2000
local k = n - n0 -- Months since j2000
local c = k / 1236.85 -- Julian centuries
local approx = M.J2000 + M.poly(c, {5.09766, M.MEAN_SYNODIC_MONTH * 1236.85,
0.0001337, -0.000000150, 0.00000000073})
local cap_E = M.poly(c, {1, -0.002516, -0.0000074})
local solar_anomaly = M.poly(c, {2.5534, 1236.85 * 29.10535669, -0.0000218, -0.00000011})
local lunar_anomaly = M.poly(c, {201.5643, 385.81693528 * 1236.85, 0.0107438, 0.00001239, -0.000000058})
local moon_argument = M.poly(c, {160.7108, 390.67050274 * 1236.85, -0.0016341, -0.00000227, 0.000000011}) -- Moon's argument of latitude
local cap_omega = M.poly(c, {124.7746, -1.56375588 * 1236.85, 0.0020672, 0.00000215}) -- Longitude of ascending node
local E_factor = {0, 1, 0, 0, 1, 1, 2, 0, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0}
local solar_coeff = {0, 1, 0, 0, -1, 1, 2, 0, 0, 1, 0, 1, 1, -1, 2, 0, 3, 1, 0, 1, -1, -1, 1, 0}
local lunar_coeff = {1, 0, 2, 0, 1, 1, 0, 1, 1, 2, 3, 0, 0, 2, 1, 2, 0, 1, 2, 1, 1, 1, 3, 4}
local moon_coeff = {0, 0, 0, 2, 0, 0, 0, -2, 2, 0, 0, 2, -2, 0, 0, -2, 0, -2, 2, 2, 2, -2, 0, 0}
local sine_coeff = {
-0.40720, 0.17241, 0.01608, 0.01039,
0.00739, -0.00514, 0.00208,
-0.00111, -0.00057, 0.00056,
-0.00042, 0.00042, 0.00038,
-0.00024, -0.00007, 0.00004,
0.00004, 0.00003, 0.00003,
-0.00003, 0.00003, -0.00002,
-0.00002, 0.00002}
local correction = -0.00017 * M.sin_degrees(cap_omega)
for i = 1, #sine_coeff do
correction = correction + sine_coeff[i] * (cap_E ^ E_factor[i])
* M.sin_degrees(solar_coeff[i] * solar_anomaly
+ lunar_coeff[i] * lunar_anomaly
+ moon_coeff[i] * moon_argument)
end
local add_const = {251.88, 251.83, 349.42, 84.66, 141.74, 207.14, 154.84,
34.52, 207.19, 291.34, 161.72, 239.56, 331.55}
local add_coeff = {0.016321, 26.651886, 36.412478, 18.206239, 53.303771,
2.453732, 7.306860, 27.261239, 0.121824, 1.844379,
24.198154, 25.513099, 3.592518}
local add_factor = {0.000165, 0.000164, 0.000126, 0.000110, 0.000062,
0.000060, 0.000056, 0.000047, 0.000042, 0.000040,
0.000037, 0.000035, 0.000023}
local extra = 0.000325 * M.sin_degrees(M.poly(c, {299.77, 132.8475848, -0.009173}))
local additional = 0
for i = 1, #add_const do
additional = additional + add_factor[i] * M.sin_degrees(add_const[i] + add_coeff[i] * k)
end
return M.universal_from_dynamical(approx + correction + extra + additional)
end
function M.mean_lunar_longitude(c)
-- Mean longitude of moon (in degrees) at moment given in Julian centuries c
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
return M.poly(c, {218.3164477, 481267.88123421, -0.0015786, 1/538841, -1/65194000}) % 360
end
function M.lunar_elongation(c)
-- Elongation of moon (in degrees) at moment given in Julian centuries c
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
return M.poly(c, {297.8501921, 445267.1114034, -0.0018819, 1/545868, -1/113065000}) % 360
end
function M.solar_anomaly(c)
-- Mean anomaly of sun (in degrees) at moment given in Julian centuries c
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
return M.poly(c, {357.5291092, 35999.0502909, -0.0001536, 1/24490000}) % 360
end
function M.lunar_anomaly(c)
-- Mean anomaly of moon (in degrees) at moment given in Julian centuries c
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
return M.poly(c, {134.9633964, 477198.8675055, 0.0087414, 1/69699, -1/14712000}) % 360
end
function M.moon_node(c)
-- Moon's argument of latitude (in degrees) at moment given in Julian centuries c
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
return M.poly(c, {93.2720950, 483202.0175233, -0.0036539, -1/3526000, 1/863310000}) % 360
end
M.lunar_longitude_func = function(tee)
-- Longitude of moon (degrees) at moment tee
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
local c = M.julian_centuries(tee)
local cap_L_prime = M.mean_lunar_longitude(c)
local cap_D = M.lunar_elongation(c)
local cap_M = M.solar_anomaly(c)
local cap_M_prime = M.lunar_anomaly(c)
local cap_F = M.moon_node(c)
local cap_E = M.poly(c, {1, -0.002516, -0.0000074})
local args_lunar_elongation = {
0, 2, 2, 0, 0, 0, 2, 2, 2, 2, 0, 1, 0, 2, 0, 0, 4, 0, 4, 2, 2, 1,
1, 2, 2, 4, 2, 0, 2, 2, 1, 2, 0, 0, 2, 2, 2, 4, 0, 3, 2, 4, 0, 2,
2, 2, 4, 0, 4, 1, 2, 0, 1, 3, 4, 2, 0, 1, 2}
local args_solar_anomaly = {
0, 0, 0, 0, 1, 0, 0, -1, 0, -1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1,
0, 1, -1, 0, 0, 0, 1, 0, -1, 0, -2, 1, 2, -2, 0, 0, -1, 0, 0, 1,
-1, 2, 2, 1, -1, 0, 0, -1, 0, 1, 0, 1, 0, 0, -1, 2, 1, 0}
local args_lunar_anomaly = {
1, -1, 0, 2, 0, 0, -2, -1, 1, 0, -1, 0, 1, 0, 1, 1, -1, 3, -2,
-1, 0, -1, 0, 1, 2, 0, -3, -2, -1, -2, 1, 0, 2, 0, -1, 1, 0,
-1, 2, -1, 1, -2, -1, -1, -2, 0, 1, 4, 0, -2, 0, 2, 1, -2, -3,
2, 1, -1, 3}
local args_moon_node = {
0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, -2, 2, -2, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, -2, 2, 0, 2, 0, 0, 0, 0,
0, 0, -2, 0, 0, 0, 0, -2, -2, 0, 0, 0, 0, 0, 0, 0}
local sine_coeff = {
6288774, 1274027, 658314, 213618, -185116, -114332,
58793, 57066, 53322, 45758, -40923, -34720, -30383,
15327, -12528, 10980, 10675, 10034, 8548, -7888,
-6766, -5163, 4987, 4036, 3994, 3861, 3665, -2689,
-2602, 2390, -2348, 2236, -2120, -2069, 2048, -1773,
-1595, 1215, -1110, -892, -810, 759, -713, -700, 691,
596, 549, 537, 520, -487, -399, -381, 351, -340, 330,
327, -323, 299, 294}
local correction = 0
for i = 1, #sine_coeff do
correction = correction + sine_coeff[i]
* (cap_E ^ math.abs(args_solar_anomaly[i]))
* M.sin_degrees(args_lunar_elongation[i] * cap_D
+ args_solar_anomaly[i] * cap_M
+ args_lunar_anomaly[i] * cap_M_prime
+ args_moon_node[i] * cap_F)
end
correction = correction / 1000000
local venus = (3958/1000000) * M.sin_degrees(119.75 + c * 131.849)
local jupiter = (318/1000000) * M.sin_degrees(53.09 + c * 479264.29)
local flat_earth = (1962/1000000) * M.sin_degrees(cap_L_prime - cap_F)
return (cap_L_prime + correction + venus + jupiter + flat_earth + M.nutation(tee)) % 360
end
function M.lunar_longitude(tee) return M.lunar_longitude_func(tee) end
M.lunar_phase_fn = function(tee)
-- Lunar phase, as an angle in degrees, at moment tee.
-- An angle of 0 means a new moon, 90 degrees means the first quarter, 180 means a full moon, and 270 degrees means the last quarter
local phi = (M.lunar_longitude(tee) - M.solar_longitude(tee)) % 360
local t0 = M.nth_new_moon(0)
local n = math.floor(0.5 + (tee - t0) / M.MEAN_SYNODIC_MONTH)
local phi_prime = 360 * (((tee - M.nth_new_moon(n)) / M.MEAN_SYNODIC_MONTH) % 1)
if math.abs(phi - phi_prime) > 180 then return phi_prime else return phi end
end
function M.lunar_phase(tee) return M.lunar_phase_fn(tee) end
function M.new_moon_before(tee)
-- Moment UT of last new moon before tee
local t0 = M.nth_new_moon(0)
local phi = M.lunar_phase(tee)
local n = math.floor(0.5 + (tee - t0) / M.MEAN_SYNODIC_MONTH - phi / 360)
local k = n - 1
while M.nth_new_moon(k) < tee do k = k + 1 end
return M.nth_new_moon(k - 1)
end
function M.new_moon_at_or_after(tee)
-- Moment UT of first new moon at or after tee
local t0 = M.nth_new_moon(0)
local phi = M.lunar_phase(tee)
local n = math.floor(0.5 + (tee - t0) / M.MEAN_SYNODIC_MONTH - phi / 360)
local k = n
while M.nth_new_moon(k) < tee do k = k + 1 end
return M.nth_new_moon(k)
end
function M.lunar_phase_at_or_before(phi, tee)
-- Moment UT of the last time at or before tee when the lunar-phase was phi degrees
local tau = tee - M.MEAN_SYNODIC_MONTH * (1/360) * ((M.lunar_phase(tee) - phi) % 360) -- Estimate
local a = tau - 2
local b = math.min(tee, tau + 2) -- At or before tee
return M.invert_angular(M.lunar_phase, phi, a, b)
end
function M.lunar_phase_at_or_after(phi, tee)
-- Moment UT of the next time at or after tee when the lunar-phase is phi degrees
local tau = tee + M.MEAN_SYNODIC_MONTH * (1/360) * ((phi - M.lunar_phase(tee)) % 360) -- Estimate
local a = math.max(tee, tau - 2) -- At or after tee
local b = tau + 2
return M.invert_angular(M.lunar_phase, phi, a, b)
end
M.NEW = 0
M.FULL = 180
M.FIRST_QUARTER = 90
M.LAST_QUARTER_PHASE = 270
function M.lunar_latitude(tee)
-- Latitude of moon (in degrees) at moment tee
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
local c = M.julian_centuries(tee)
local cap_L_prime = M.mean_lunar_longitude(c)
local cap_D = M.lunar_elongation(c)
local cap_M = M.solar_anomaly(c)
local cap_M_prime = M.lunar_anomaly(c)
local cap_F = M.moon_node(c)
local cap_E = M.poly(c, {1, -0.002516, -0.0000074})
local args_lunar_elongation = {
0, 0, 0, 2, 2, 2, 2, 0, 2, 0, 2, 2, 2, 2, 2, 2, 2, 0, 4, 0, 0, 0,
1, 0, 0, 0, 1, 0, 4, 4, 0, 4, 2, 2, 2, 2, 0, 2, 2, 2, 2, 4, 2, 2,
0, 2, 1, 1, 0, 2, 1, 2, 0, 4, 4, 1, 4, 1, 4, 2}
local args_solar_anomaly = {
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 1, -1, -1, -1, 1, 0, 1,
0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 1, 1,
0, -1, -2, 0, 1, 1, 1, 1, 1, 0, -1, 1, 0, -1, 0, 0, 0, -1, -2}
local args_lunar_anomaly = {
0, 1, 1, 0, -1, -1, 0, 2, 1, 2, 0, -2, 1, 0, -1, 0, -1, -1, -1,
0, 0, -1, 0, 1, 1, 0, 0, 3, 0, -1, 1, -2, 0, 2, 1, -2, 3, 2, -3,
-1, 0, 0, 1, 0, 1, 1, 0, 0, -2, -1, 1, -2, 2, -2, -1, 1, 1, -1,
0, 0}
local args_moon_node = {
1, 1, -1, -1, 1, -1, 1, 1, -1, -1, -1, -1, 1, -1, 1, 1, -1, -1,
-1, 1, 3, 1, 1, 1, -1, -1, -1, 1, -1, 1, -3, 1, -3, -1, -1, 1,
-1, 1, -1, 1, 1, 1, 1, -1, 3, -1, -1, 1, -1, -1, 1, -1, 1, -1,
-1, -1, -1, -1, -1, 1}
local sine_coeff = {
5128122, 280602, 277693, 173237, 55413, 46271, 32573,
17198, 9266, 8822, 8216, 4324, 4200, -3359, 2463, 2211,
2065, -1870, 1828, -1794, -1749, -1565, -1491, -1475,
-1410, -1344, -1335, 1107, 1021, 833, 777, 671, 607,
596, 491, -451, 439, 422, 421, -366, -351, 331, 315,
302, -283, -229, 223, 223, -220, -220, -185, 181,
-177, 176, 166, -164, 132, -119, 115, 107}
local beta = 0
for i = 1, #sine_coeff do
beta = beta + sine_coeff[i]
* (cap_E ^ math.abs(args_solar_anomaly[i]))
* M.sin_degrees(args_lunar_elongation[i] * cap_D
+ args_solar_anomaly[i] * cap_M
+ args_lunar_anomaly[i] * cap_M_prime
+ args_moon_node[i] * cap_F)
end
beta = beta / 1000000
local venus = (175/1000000) * (
M.sin_degrees(119.75 + c * 131.849 + cap_F)
+ M.sin_degrees(119.75 + c * 131.849 - cap_F))
local flat_earth = (-2235/1000000) * M.sin_degrees(cap_L_prime)
+ (127/1000000) * M.sin_degrees(cap_L_prime - cap_M_prime)
+ (-115/1000000) * M.sin_degrees(cap_L_prime + cap_M_prime)
local extra = (382/1000000) * M.sin_degrees(313.45 + c * 481266.484)
return beta + venus + flat_earth + extra
end
function M.lunar_altitude(tee, locale)
-- Geocentric altitude of moon at tee at locale, as a small positive/negative angle in degrees, ignoring parallax and refraction
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed., 1998
local phi = locale[1] -- Local latitude
local psi = locale[2] -- Local longitude
local lambda = M.lunar_longitude(tee) -- Lunar longitude
local beta = M.lunar_latitude(tee) -- Lunar latitude
local alpha = M.right_ascension(tee, beta, lambda) -- Lunar right ascension
local delta = M.declination(tee, beta, lambda) -- Lunar declination
local theta0 = M.sidereal_from_moment(tee) -- Sidereal time
local cap_H = (theta0 + psi - alpha) % 360 -- Local hour angle
local altitude = M.arcsin_degrees(M.sin_degrees(phi) * M.sin_degrees(delta)
+ M.cosine_degrees(phi) * M.cosine_degrees(delta) * M.cosine_degrees(cap_H))
return M.mod3(altitude, -180, 180)
end
function M.lunar_distance(tee)
-- Distance to moon (in meters) at moment tee
-- Adapted Adapted from from "Astronomical Astronomica Algorithms" by Jean Meeus, Willmann-Bell, Inc., 2nd ed.
local c = M.julian_centuries(tee)
local cap_D = M.lunar_elongation(c)
local cap_M = M.solar_anomaly(c)
local cap_M_prime = M.lunar_anomaly(c)
local cap_F = M.moon_node(c)
local cap_E = M.poly(c, {1, -0.002516, -0.0000074})
local args_lunar_elongation = {
0, 2, 2, 0, 0, 0, 2, 2, 2, 2, 0, 1, 0, 2, 0, 0, 4, 0, 4, 2, 2, 1,
1, 2, 2, 4, 2, 0, 2, 2, 1, 2, 0, 0, 2, 2, 2, 4, 0, 3, 2, 4, 0, 2,
2, 2, 4, 0, 4, 1, 2, 0, 1, 3, 4, 2, 0, 1, 2, 2}
local args_solar_anomaly = {
0, 0, 0, 0, 1, 0, 0, -1, 0, -1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1,
0, 1, -1, 0, 0, 0, 1, 0, -1, 0, -2, 1, 2, -2, 0, 0, -1, 0, 0, 1,
-1, 2, 2, 1, -1, 0, 0, -1, 0, 1, 0, 1, 0, 0, -1, 2, 1, 0, 0}
local args_lunar_anomaly = {
1, -1, 0, 2, 0, 0, -2, -1, 1, 0, -1, 0, 1, 0, 1, 1, -1, 3, -2,
-1, 0, -1, 0, 1, 2, 0, -3, -2, -1, -2, 1, 0, 2, 0, -1, 1, 0,
-1, 2, -1, 1, -2, -1, -1, -2, 0, 1, 4, 0, -2, 0, 2, 1, -2, -3,
2, 1, -1, 3, -1}
local args_moon_node = {
0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, -2, 2, -2, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, -2, 2, 0, 2, 0, 0, 0, 0,
0, 0, -2, 0, 0, 0, 0, -2, -2, 0, 0, 0, 0, 0, 0, 0, -2}
local cosine_coeff = {
-20905355, -3699111, -2955968, -569925, 48888, -3149,
246158, -152138, -170733, -204586, -129620, 108743,
104755, 10321, 0, 79661, -34782, -23210, -21636, 24208,
30824, -8379, -16675, -12831, -10445, -11650, 14403,
-7003, 0, 10056, 6322, -9884, 5751, 0, -4950, 4130, 0,
-3958, 0, 3258, 2616, -1897, -2117, 2354, 0, 0, -1423,
-1117, -1571, -1739, 0, -4421, 0, 0, 0, 0, 1165, 0, 0,
8752}
local correction = 0
for i = 1, #cosine_coeff do
correction = correction + cosine_coeff[i]
* (cap_E ^ math.abs(args_solar_anomaly[i]))
* M.cosine_degrees(args_lunar_elongation[i] * cap_D
+ args_solar_anomaly[i] * cap_M
+ args_lunar_anomaly[i] * cap_M_prime
+ args_moon_node[i] * cap_F)
end
return 385000560 + correction
end
function M.lunar_parallax(tee, locale)
-- Parallax of moon at tee at locale
-- Adapted from "Astronomical Algorithms" by Jean Meeus, Willmann-Bell, Inc., 1998
local geo = M.lunar_altitude(tee, locale)
local cap_Delta = M.lunar_distance(tee)
local alt = 6378140 / cap_Delta
return M.arcsin_degrees(alt * M.cosine_degrees(geo))
end
function M.topocentric_lunar_altitude(tee, locale)
-- Topocentric altitude of moon at tee at locale, as a small positive/negative angle in degrees, ignoring refraction
return M.lunar_altitude(tee, locale) - M.lunar_parallax(tee, locale)
end
function M.sine_offset(tee, locale, alpha)
-- Sine of angle between position of sun at local time tee and when its depression is alpha at locale. Out of range when it does not occur
local phi = locale[1]
local tee_prime = M.universal_from_local(tee, locale)
local delta = M.declination(tee_prime, 0, M.solar_longitude(tee_prime)) -- Declination of sun
return M.tangent_degrees(phi) * M.tangent_degrees(delta)
+ M.sin_degrees(alpha) / (M.cosine_degrees(delta) * M.cosine_degrees(phi))
end
function M.approx_moment_of_depression(tee, locale, alpha, early)
-- Moment in local time near tee when depression angle of sun is alpha (negative if above horizon) at locale; early? is true when morning event is sought and false for evening. Returns bogus if depression angle is not reached
local try = M.sine_offset(tee, locale, alpha)
local date = M.fixed_from_moment(tee)
local alt
if alpha >= 0 then
if early then alt = date else alt = date + 1 end
else
alt = date + M.hr(12)
end
local value
if math.abs(try) > 1 then value = M.sine_offset(alt, locale, alpha) else value = try end
if math.abs(value) <= 1 then -- Event occurs
local offset = M.mod3(M.arcsin_degrees(value) / 360, M.hr(-12), M.hr(12))
local t
if early then t = date + (M.hr(6) - offset) else t = date + (M.hr(18) + offset) end
return M.local_from_apparent(t, locale)
else
return M.BOGUS
end
end
function M.moment_of_depression(approx, locale, alpha, early)
-- Moment in local time near approx when depression angle of sun is alpha (negative if above horizon) locale; early? is true when morning event is sought, and false for evening
-- Returns bogus if depression angle is not reached
local tee = M.approx_moment_of_depression(approx, locale, alpha, early)
if tee == M.BOGUS then return M.BOGUS end
if math.abs(approx - tee) < M.sec(30) then return tee
else return M.moment_of_depression(tee, locale, alpha, early) end
end
M.MORNING = true -- Signifies morning
M.EVENING = false -- Signifies evening
function M.dawn(date, locale, alpha)
-- Standard time in morning on fixed date at locale when depression angle of sun is alpha. Returns bogus if there is no dawn on date
local result = M.moment_of_depression(date + M.hr(6), locale, alpha, M.MORNING)
if result == M.BOGUS then return M.BOGUS end
return M.standard_from_local(result, locale)
end
function M.dusk(date, locale, alpha)
-- Standard time in evening on fixed date at locale when depression angle of sun is alpha. Returns bogus if there is no dusk on date
local result = M.moment_of_depression(date + M.hr(18), locale, alpha, M.EVENING)
if result == M.BOGUS then return M.BOGUS end
return M.standard_from_local(result, locale)
end
function M.refraction(locale)
-- Refraction angle at locale
local h = math.max(0, locale[3])
local cap_R = 6.372e6 -- Radius of Earth
local dip = M.arccos_degrees(cap_R / (cap_R + h)) -- Depression of visible horizon
return M.mins(34) + dip + M.secs(19) * math.sqrt(h)
end
function M.sunrise(date, locale)
-- Standard time of sunrise on fixed date at locale
local alpha = M.refraction(locale) + M.mins(16)
return M.dawn(date, locale, alpha)
end
function M.sunset(date, locale)
-- Standard time of sunset on fixed date at locale
local alpha = M.refraction(locale) + M.mins(16)
return M.dusk(date, locale, alpha)
end
function M.jewish_dusk(date, locale)
-- Standard time of Jewish dusk on fixed date at locale (as per Vilna Gaon)
return M.dusk(date, locale, M.angle(4, 40, 0))
end
function M.jewish_sabbath_ends(date, locale)
-- Standard time of end of Jewish sabbath on fixed date at locale (as per Berthold Cohn)
return M.dusk(date, locale, M.angle(7, 5, 0))
end
function M.daytime_temporal_hour(date, locale)
-- Length of daytime temporal hour on fixed date at locale
-- Returns bogus if there no sunrise or sunset on date
local sunrise, sunset = M.sunrise(date, locale), M.sunset(date, locale)
if sunrise == M.BOGUS or sunset == M.BOGUS then return M.BOGUS end
return (sunset - sunrise) / 12
end
function M.nighttime_temporal_hour(date, locale)
-- Length of nighttime temporal hour on fixed date at locale
-- Returns bogus if there no sunrise or sunset on date
local sunrise_next = M.sunrise(date + 1, locale)
local sunset = M.sunset(date, locale)
if sunrise_next == M.BOGUS or sunset == M.BOGUS then return M.BOGUS end
return (sunrise_next - sunset) / 12
end
function M.standard_from_sundial(tee, locale)
-- Standard time of temporal moment tee at locale
-- Returns bogus if temporal hour is undefined that day
local date = M.fixed_from_moment(tee)
local hour = 24 * (tee % 1)
local h
if 6 <= hour and hour <= 18 then h = M.daytime_temporal_hour(date, locale) -- daytime today
elseif hour < 6 then h = M.nighttime_temporal_hour(date - 1, locale) -- early this morning
else h = M.nighttime_temporal_hour(date, locale) end -- this evening
if h == M.BOGUS then return M.BOGUS end
if 6 <= hour and hour <= 18 then return M.sunrise(date, locale) + (hour - 6) * h -- daytime today
elseif hour < 6 then return M.sunset(date - 1, locale) + (hour + 6) * h -- early this morning
else return M.sunset(date, locale) + (hour - 18) * h end -- this evening
end
function M.jewish_morning_end(date, locale)
-- Standard time on fixed date at locale of end of morning according to Jewish ritual
return M.standard_from_sundial(date + M.hr(10), locale)
end
function M.asr(date, locale)
-- Standard time of asr on fixed date at locale
local noon = M.midday(date, locale) -- Time when sun nearest zenith
local phi = locale[1]
local delta = M.declination(noon, 0, M.solar_longitude(noon)) -- Solar declination at noon
local altitude = M.arcsin_degrees(M.cos_degrees(delta) * M.cos_degrees(phi) -- Solar altitude at noon
+ M.sin_degrees(delta) * M.sin_degrees(phi))
local h = M.mod3(M.arctan_degrees(M.tan_degrees(altitude), -- Sun's altitude when shadow increases by double its length
1 + 2 * M.tan_degrees(altitude)), -90, 90)
-- For Shafii use instead:
-- 1 + M.tan_degrees(altitude)), -90, 90)
if altitude <= 0 then return M.BOGUS end
return M.dusk(date, locale, -h)
end
----- B.16 The Persian Calendar -----
----- B.18 The French Revolutionary Calendar -----
M.FRENCH_EPOCH = M.fixed_from_gregorian({1792, M.SEPTEMBER, 22}) -- Fixed date of start of the French Revolutionary calendar
M.PARIS = {M.angle(48, 50, 11), M.angle(2, 20, 15), 27, M.hr(1)} -- Location of Paris Observatory. Longitude corresponds to difference of 9m 21s between Paris time zone and Universal Time
function M.midnight_in_paris(date)
-- Universal time of true midnight at end of fixed date in Paris
return M.midnight(date + 1, M.PARIS)
end
function M.french_new_year_on_or_before(date)
-- Fixed date of French Revolutionary New Year on or before fixed date
local approx = M.estimate_prior_solar_longitude(M.AUTUMN, M.midnight_in_paris(date)) -- Approximate time of solstice
local day = math.floor(approx) - 1
while not (M.AUTUMN <= M.solar_longitude(M.midnight_in_paris(day))) do
day = day + 1
end
return day
end
function M.fixed_from_french(f_date)
-- Fixed date of French Revolutionary date
local year, month, day = f_date[1], f_date[2], f_date[3]
local new_year = M.french_new_year_on_or_before(
math.floor(M.FRENCH_EPOCH + 180 + M.MEAN_TROPICAL_YEAR_AST * (year - 1)))
return new_year - 1 + 30 * (month - 1) + day
end
p.fixed_from_french = make_public(M.fixed_from_french, false)
function M.french_from_fixed(date)
-- French Revolutionary date of fixed date
local new_year = M.french_new_year_on_or_before(date)
local year = math.floor(0.5 + (new_year - M.FRENCH_EPOCH) / M.MEAN_TROPICAL_YEAR_AST) + 1
local month = M.quotient(date - new_year, 30) + 1
local day = (date - new_year) % 30 + 1
local decade = M.quotient(day - 1, 10) + 1 -- I added decade and weekday, which are not in the book, for use on Wikipedia
local weekday = M.amod(day, 10)
return {year, month, decade, weekday}
end
p.french_from_fixed = make_public(M.french_from_fixed)
function M.french_leap_year(f_year)
-- True if f-year is a leap year on the French Revolutionary calendar
return (M.fixed_from_french({f_year + 1, 1, 1}) - M.fixed_from_french({f_year, 1, 1})) > 365
end
function M.arithmetic_french_leap_year(f_year)
-- True if f-year is a leap year on the Arithmetic French Revolutionary calendar
local m400 = f_year % 400
return (f_year % 4 == 0)
and not (m400 == 100 or m400 == 200 or m400 == 300)
and not (f_year % 4000 == 0)
end
function M.fixed_from_arithmetic_french(f_date)
-- Fixed date of Arithmetic French Revolutionary date f-date
local year, month, day = f_date[1], f_date[2], f_date[3]
return M.FRENCH_EPOCH - 1 -- Days before start of calendar
+ 365 * (year - 1) -- Ordinary days in prior years
+ M.quotient(year - 1, 4) -- Leap days in prior years
- M.quotient(year - 1, 100)
+ M.quotient(year - 1, 400)
- M.quotient(year - 1, 4000)
+ 30 * (month - 1) -- Days in prior months this year
+ day -- Days this month
end
p.fixed_from_arithmetic_french = make_public(M.fixed_from_arithmetic_french, false)
function M.arithmetic_french_from_fixed(date)
-- Arithmetic French Revolutionary (year month day) of fixed date
local approx = M.quotient(date - M.FRENCH_EPOCH + 2, 1460969/4000) + 1 -- Approximate year (may be off by 1)
local year
if date < M.fixed_from_arithmetic_french({approx, 1, 1}) then year = approx - 1
else year = approx end
local month = M.quotient(date - M.fixed_from_arithmetic_french({year, 1, 1}), 30) + 1 -- Calculate the month by division
local day = date - M.fixed_from_arithmetic_french({year, month, 1}) + 1 -- Calculate the day by subtraction
return {year, month, day}
end
p.arithmetic_french_from_fixed = make_public(M.arithmetic_french_from_fixed)
----- B.19 The Chinese Calendar -----
function M.chinese_location(tee)
-- Location of Beijing; time zone varies with tee
local year = M.gregorian_year_from_fixed(math.floor(tee))
if year < 1929 then
return {M.angle(39, 55, 0), M.angle(116, 25, 0), 43.5, M.hr(1397/180)}
else
return {M.angle(39, 55, 0), M.angle(116, 25, 0), 43.5, M.hr(8)}
end
end
function M.midnight_in_china(date)
-- Universal time of (clock) midnight at start of fixed date in China
return M.universal_from_standard(date, M.chinese_location(date))
end
function M.chinese_solar_longitude_on_or_after(lambda, date)
-- Moment (Beijing time) of the first date on or after fixed date (Beijing time) when the solar longitude will be lambda degrees
local tee = M.solar_longitude_after(lambda, M.universal_from_standard(date, M.chinese_location(date)))
return M.standard_from_universal(tee, M.chinese_location(tee))
end
function M.chinese_current_major_solar_term(date)
-- Last Chinese major solar term (zhongqi) before fixed date
local s = M.solar_longitude(M.universal_from_standard(date, M.chinese_location(date)))
return M.amod(2 + M.quotient(s, 30), 12)
end
function M.chinese_major_solar_term_on_or_after(date)
-- Moment (in Beijing) of first Chinese major solar term (zhongqi) on or after fixed date. The major terms begin when the sun's longitude is a multiple of 30 degrees
local s = M.solar_longitude(M.midnight_in_china(date))
local l = (30 * math.ceil(s / 30)) % 360
return M.chinese_solar_longitude_on_or_after(l, date)
end
function M.chinese_current_minor_solar_term(date)
-- Last Chinese minor solar term (jieqi) before date
local s = M.solar_longitude(M.universal_from_standard(date, M.chinese_location(date)))
return M.amod(3 + M.quotient(s - 15, 30), 12)
end
function M.chinese_minor_solar_term_on_or_after(date)
-- Moment (in Beijing) of the first Chinese minor solar term (jieqi) on or after fixed date. The minor terms begin when the sun's longitude is an odd multiple of 15 degrees
local s = M.solar_longitude(M.midnight_in_china(date))
local l = (30 * math.ceil((s - 15) / 30) + 15) % 360
return M.chinese_solar_longitude_on_or_after(l, date)
end
function M.chinese_solar_term_from_fixed(fixed)
-- This one is not in the book, but I added it so that Template:Currect Sekki can use it
local s = M.solar_longitude(M.universal_from_standard(fixed, M.chinese_location(fixed)))
local l = 15 * (M.quotient(s, 15) % 24) -- Last solar term
local next_l = (l + 15) % 360 -- Next solar term
local next_moment = math.floor(M.chinese_solar_longitude_on_or_after(next_l, fixed)) -- Calculate day when the next solar term begins
local days_until = next_moment - fixed
return {l, next_l, days_until}
end
p.chinese_solar_term_from_fixed = make_public(M.chinese_solar_term_from_fixed)
function M.chinese_new_moon_before(date)
-- Fixed date (Beijing) of first new moon before fixed date
local tee = M.new_moon_before(M.midnight_in_china(date))
return math.floor(M.standard_from_universal(tee, M.chinese_location(tee)))
end
function M.chinese_new_moon_on_or_after(date)
-- Fixed date (Beijing) of first new moon on or after fixed date
local tee = M.new_moon_at_or_after(M.midnight_in_china(date))
return math.floor(M.standard_from_universal(tee, M.chinese_location(tee)))
end
M.CHINESE_EPOCH = M.fixed_from_gregorian({-2636, M.FEBRUARY, 15}) -- Fixed date of start of the Chinese calendar
function M.chinese_no_major_solar_term(date)
-- True if Chinese lunar month starting on date has no major solar term
return M.chinese_current_major_solar_term(date)
== M.chinese_current_major_solar_term(M.chinese_new_moon_on_or_after(date + 1))
end
function M.chinese_winter_solstice_on_or_before(date)
-- Fixed date, in the Chinese zone, of winter solstice on or before fixed date
local approx = M.estimate_prior_solar_longitude(M.WINTER, M.midnight_in_china(date + 1)) -- Approximate time of solstice
local day = math.floor(approx) - 1
while not (M.WINTER < M.solar_longitude(M.midnight_in_china(day + 1))) do
day = day + 1
end
return day
end
function M.chinese_new_year_in_sui(date)
-- Fixed date of Chinese New Year in sui (period from solstice to solstice) containing date
local s1 = M.chinese_winter_solstice_on_or_before(date) -- prior solstice
local s2 = M.chinese_winter_solstice_on_or_before(s1 + 370) -- following solstice
local m12 = M.chinese_new_moon_on_or_after(s1 + 1) -- month after 11th month; either 12 or leaр 11
local m13 = M.chinese_new_moon_on_or_after(m12 + 1) -- month after m12; either 12 (or leap 12) or 1
local next_m11 = M.chinese_new_moon_before(s2 + 1) -- next 11th month
-- Either m12 or m13 is a leap month if there are 13 new moons (12 full lunar months) and either m12 or m13 has no major solar term
if math.floor(0.5 + (next_m11 - m12) / M.MEAN_SYNODIC_MONTH) == 12
and (M.chinese_no_major_solar_term(m12) or M.chinese_no_major_solar_term(m13)) then
return M.chinese_new_moon_on_or_after(m13 + 1)
else
return m13
end
end
function M.chinese_new_year_on_or_before(date)
-- Fixed date of Chinese New Year on or before fixed date
local new_year = M.chinese_new_year_in_sui(date)
if date >= new_year then return new_year
-- Got the New Year after; this happens if date is after the solstice but before the new year. So, go back half a year
else return M.chinese_new_year_in_sui(date - 180) end
end
function M.chinese_new_year(g_year)
-- Fixed date of Chinese New Year in Gregorian year g-year
return M.chinese_new_year_on_or_before(M.fixed_from_gregorian({g_year, M.JULY, 1}))
end
p.chinese_new_year = make_public(M.chinese_new_year)
M.chinese_prior_leap_month = function(m_prime, m)
-- True if there is a Chinese leap month on or after lunar month starting m-prime and at or before lunar month starting m
return m >= m_prime
and (M.chinese_no_major_solar_term(m)
or M.chinese_prior_leap_month(m_prime, M.chinese_new_moon_before(m)))
end
function M.chinese_from_fixed(date)
-- Chinese date (cycle year month leap day) of fixed date
local s1 = M.chinese_winter_solstice_on_or_before(date) -- Prior solstice
local s2 = M.chinese_winter_solstice_on_or_before(s1 + 370) -- Following solstice
local m12 = M.chinese_new_moon_on_or_after(s1 + 1) -- month after last 11th month
local next_m11 = M.chinese_new_moon_before(s2 + 1) -- next 11th month
local m = M.chinese_new_moon_before(date + 1) -- start of month containing date
local leap_year = math.floor(0.5 + (next_m11 - m12) / M.MEAN_SYNODIC_MONTH) == 12 -- if there are 13 new moons (12 full lunar months)
local prior = (leap_year and M.chinese_prior_leap_month(m12, m)) and 1 or 0 -- subtract 1 from ordinal position of month in year during or after a leap month
local month = M.amod(math.floor(0.5 + (m - m12) / M.MEAN_SYNODIC_MONTH) - prior, 12) -- month number
local leap_month = leap_year -- it's a leap month if there are 13 months
and M.chinese_no_major_solar_term(m) -- no major solar term
and (not M.chinese_prior_leap_month(m12, M.chinese_new_moon_before(m))) -- and no prior leap month
local elapsed_years = math.floor(1.5 - month/12 + (date - M.CHINESE_EPOCH) / M.MEAN_TROPICAL_YEAR_AST) -- Approximate since the epoch
local cycle = M.quotient(elapsed_years - 1, 60) + 1
local year = M.amod(elapsed_years, 60)
local day = date - m + 1
return {cycle, year, month, leap_month, day}
end
p.chinese_from_fixed = make_public(M.chinese_from_fixed)
function M.fixed_from_chinese(c_date)
-- Fixed date of Chinese date c-date.
local cycle, year, month, leap, day = c_date[1], c_date[2], c_date[3], c_date[4], c_date[5]
local mid_year = math.floor(M.CHINESE_EPOCH
+ ((cycle - 1) * 60 + (year - 1) + 1/2) * M.MEAN_TROPICAL_YEAR_AST)
local new_year = M.chinese_new_year_on_or_before(mid_year)
local p = M.chinese_new_moon_on_or_after(new_year + (month - 1) * 29)
local d = M.chinese_from_fixed(p)
local prior_new_moon
if month == d[3] and leap == d[4] then
prior_new_moon = p
else
prior_new_moon = M.chinese_new_moon_on_or_after(p + 1)
end
return prior_new_moon + day - 1
end
function M.fixed_from_chinese_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.fixed_from_chinese({table[1], table[2], table[3], table[4] == 1, table[5]})
end
p.fixed_from_chinese = make_public(M.fixed_from_chinese_for_wikipedia, false)
function M.chinese_sexagesimal_name(n)
-- The n-th name of the Chinese sexagesimal cycle
return {M.amod(n, 10), M.amod(n, 12)}
end
function M.chinese_name_difference(c_name1, c_name2)
-- Number of names from Chinese name c-name1 to next occurrence of c-name2
local stem1, stem2 = c_name1[1], c_name2[2]
local branch1, branch2 = c_name1[1], c_name2[2]
local stem_difference = stem2 - stem1
local branch_difference = branch2 - branch1
return ((stem_difference - 1 + 25 * (branch_difference - stem_difference)) % 60) + 1
end
function M.chinese_name_of_year(year)
-- Sexagesimal name for Chinese year of any cycle
return M.chinese_sexagesimal_name(year)
end
M.CHINESE_MONTH_NAME_EPOCH = 57 -- Elapsed months at start of Chinese sexagesimal month cycle
function M.chinese_name_of_month(month, year)
-- Sexagesimal name for month of Chinese year
local elapsed_months = 12 * (year - 1) + (month - 1)
return M.chinese_sexagesimal_name(elapsed_months - M.CHINESE_MONTH_NAME_EPOCH)
end
M.CHINESE_DAY_NAME_EPOCH = M.rd(45) -- RD date of a start of Chinese sexagesimal day cycle
function M.chinese_name_of_day(date)
-- Chinese sexagesimal name for date
return M.chinese_sexagesimal_name(date - M.CHINESE_DAY_NAME_EPOCH)
end
function M.chinese_day_name_on_or_before(name, date)
-- Fixed date of latest date on or before fixed date that has Chinese name
return M.mod3(M.chinese_name_difference(M.chinese_day_name(0), name), date - 60, date)
end
function M.dragon_festival(g_year)
-- Fixed date of the Dragon Festival in Gregorian year g-year
local elapsed_years = g_year - M.gregorian_year_from_fixed(M.CHINESE_EPOCH) + 1
local cycle = M.quotient(elapsed_years - 1, 60) + 1
local year = M.amod(elapsed_years, 60)
return M.fixed_from_chinese({cycle, year, 5, false, 5})
end
function M.qing_ming(g_year)
-- Fixed date of Qingming in Gregorian year g-year
return math.floor(M.chinese_minor_solar_term_on_or_after(M.fixed_from_gregorian({g_year, M.MARCH, 30})))
end
function M.chinese_age(birthdate, date)
-- Age at fixed date, given Chinese birthdate, according to the Chinese custom. Returns bogus if date is before birthdate
local today = M.chinese_from_fixed(date)
if date >= M.fixed_from_chinese(birthdate) then
return 60 * (M.chinese_cycle(today) - M.chinese_cycle(birthdate))
+ (M.chinese_year(today) - M.chinese_year(birthdate)) + 1
else
return M.BOGUS
end
end
M.WIDOW_AUGURY = 0 -- Lichun does not occur (double-blind year)
M.BLIND_AUGURY = 1 -- Lichun occurs once at the end
M.BRIGHT_AUGURY = 2 -- Lichun occurs once at the start
M.DOUBLE_BRIGHT_AUGURY = 3 -- Lichun occurs twice (double-happiness)
function M.chinese_year_marriage_augury(cycle, year)
-- The marriage augury type of Chinese year in cycle
local new_year = M.fixed_from_chinese({cycle, year, 1, false, 1})
local c, y
if year == 60 then c = cycle + 1; y = 1
else c = cycle; y = year + 1 end
local next_new_year = M.fixed_from_chinese({c, y, 1, false, 1})
local first_minor_term = M.chinese_current_minor_solar_term(new_year)
local next_first_minor_term = M.chinese_current_minor_solar_term(next_new_year)
if first_minor_term == 1 and next_first_minor_term == 12 then
return M.WIDOW_AUGURY
elseif first_minor_term == 1 and next_first_minor_term ~= 12 then
return M.BLIND_AUGURY
elseif first_minor_term ~= 1 and next_first_minor_term == 12 then
return M.BRIGHT_AUGURY
else
return M.DOUBLE_BRIGHT_AUGURY
end
end
function M.korean_location(tee)
-- Location for Korean calendar; varies with tee
-- Seoul city hall at a varying time zone
local z
if tee < M.fixed_from_gregorian({1908, M.APRIL, 1}) then z = 3809/450 -- mean local mean time for longitude 126 deg 58 min
elseif tee < M.fixed_from_gregorian({1912, M.JANUARY, 1}) then z = 8.5
elseif tee < M.fixed_from_gregorian({1954, M.MARCH, 21}) then z = 9
elseif tee < M.fixed_from_gregorian({1961, M.AUGUST, 10}) then z = 8.5
else z = 9 end
return M.location(M.angle(37, 34, 0), M.angle(126, 58, 0), 0, M.hr(z))
end
function M.korean_year(cycle, year)
-- Equivalent Korean year to Chinese cycle and year
return 60 * cycle + year - 364
end
function M.vietnamese_location(tee)
-- Location for Vietnamese calendar is Hanoi; varies with tee. Time zone has changed over the years
local z
if tee < M.gregorian_new_year(1968) then z = 8 else z = 7 end
return M.location(M.angle(21, 2, 0), M.angle(105, 51, 0), 0, M.hr(z))
end
----- Japanese calendar -----
function M.japanese_location(tee)
-- Location for Japanese calendar; varies with tee
local year = M.gregorian_year_from_fixed(math.floor(tee))
if year < 1888 then -- Tokyo (139 deg 46 min east) local time
return {35.7, M.angle(139, 46, 0), 24, M.hr(9 + 143/450)}
else -- Longitude 135 time zone
return {35, 135, 0, M.hr(9)}
end
end
function M.midnight_in_japan(date)
-- Universal time of (clock) midnight at start of fixed date in Japan
return M.universal_from_standard(date, M.japanese_location(date))
end
function M.japanese_solar_longitude_on_or_after(lambda, date)
-- Moment (Tokyo time) of the first date on or after fixed date (Tokyo time) when the solar longitude will be lambda degrees
local tee = M.solar_longitude_after(lambda, M.universal_from_standard(date, M.japanese_location(date)))
return M.standard_from_universal(tee, M.japanese_location(tee))
end
function M.japanese_current_major_solar_term(date)
-- Last Japanese major solar term (chuki) before fixed date
local s = M.solar_longitude(M.universal_from_standard(date, M.japanese_location(date)))
return M.amod(2 + M.quotient(s, 30), 12)
end
function M.japanese_major_solar_term_on_or_after(date)
-- Moment (in Tokyo) of first Japanese major solar term (chuki) on or after fixed date. The major terms begin when the sun's longitude is a multiple of 30 degrees
local s = M.solar_longitude(M.midnight_in_japan(date))
local l = (30 * math.ceil(s / 30)) % 360
return M.japanese_solar_longitude_on_or_after(l, date)
end
function M.japanese_current_minor_solar_term(date)
-- Last Japanese minor solar term (sekki) before date
local s = M.solar_longitude(M.universal_from_standard(date, M.japanese_location(date)))
return M.amod(3 + M.quotient(s - 15, 30), 12)
end
function M.japanese_minor_solar_term_on_or_after(date)
-- Moment (in Tokyo) of the first Japanese minor solar term (sekki) on or after fixed date. The minor terms begin when the sun's longitude is an odd multiple of 15 degrees
local s = M.solar_longitude(M.midnight_in_japan(date))
local l = (30 * math.ceil((s - 15) / 30) + 15) % 360
return M.japanese_solar_longitude_on_or_after(l, date)
end
function M.japanese_solar_term_from_fixed(fixed)
-- This one is not in the book, but I added it so that Template:Currect Sekki can use it
local s = M.solar_longitude(M.universal_from_standard(fixed - M.hr(12), M.japanese_location(fixed)))
local l = 15 * (M.quotient(s, 15) % 24) -- Last solar term
local next_l = (l + 15) % 360 -- Next solar term
local next_moment = math.floor(M.japanese_solar_longitude_on_or_after(next_l, fixed)) -- Calculate day when the next solar term begins
local days_until = next_moment - fixed
return {l, next_l, days_until}
end
p.japanese_solar_term_from_fixed = make_public(M.japanese_solar_term_from_fixed)
function M.japanese_new_moon_before(date)
-- Fixed date (Tokyo) of first new moon before fixed date
local tee = M.new_moon_before(M.midnight_in_japan(date))
return math.floor(M.standard_from_universal(tee, M.japanese_location(tee)))
end
function M.japanese_new_moon_on_or_after(date)
-- Fixed date (Tokyo) of first new moon on or after fixed date
local tee = M.new_moon_at_or_after(M.midnight_in_japan(date))
return math.floor(M.standard_from_universal(tee, M.japanese_location(tee)))
end
function M.japanese_no_major_solar_term(date)
-- True if Japanese lunar month starting on date has no major solar term
return M.japanese_current_major_solar_term(date)
== M.japanese_current_major_solar_term(M.japanese_new_moon_on_or_after(date + 1))
end
function M.japanese_winter_solstice_on_or_before(date)
-- Fixed date, in the Japanese zone, of winter solstice on or before fixed date
local approx = M.estimate_prior_solar_longitude(M.WINTER, M.midnight_in_japan(date + 1)) -- Approximate time of solstice
local day = math.floor(approx) - 1
while not (M.WINTER < M.solar_longitude(M.midnight_in_japan(day + 1))) do
day = day + 1
end
return day
end
function M.japanese_new_year_in_sui(date)
-- Fixed date of Japanese New Year in sui (period from solstice to solstice) containing date
local s1 = M.japanese_winter_solstice_on_or_before(date) -- prior solstice
local s2 = M.japanese_winter_solstice_on_or_before(s1 + 370) -- following solstice
local m12 = M.japanese_new_moon_on_or_after(s1 + 1) -- month after 11th month; either 12 or leaр 11
local m13 = M.japanese_new_moon_on_or_after(m12 + 1) -- month after m12; either 12 (or leap 12) or 1
local next_m11 = M.japanese_new_moon_before(s2 + 1) -- next 11th month
-- Either m12 or m13 is a leap month if there are 13 new moons (12 full lunar months) and either m12 or m13 has no major solar term
if math.floor(0.5 + (next_m11 - m12) / M.MEAN_SYNODIC_MONTH) == 12
and (M.japanese_no_major_solar_term(m12) or M.japanese_no_major_solar_term(m13)) then
return M.japanese_new_moon_on_or_after(m13 + 1)
else
return m13
end
end
function M.japanese_new_year_on_or_before(date)
-- Fixed date of Japanese New Year on or before fixed date
local new_year = M.japanese_new_year_in_sui(date)
if date >= new_year then return new_year
-- Got the New Year after; this happens if date is after the solstice but before the new year. So, go back half a year
else return M.japanese_new_year_in_sui(date - 180) end
end
M.japanese_prior_leap_month = function(m_prime, m)
-- True if there is a Japanese leap month on or after lunar month starting m-prime and at or before lunar month starting m
return m >= m_prime
and (M.japanese_no_major_solar_term(m)
or M.japanese_prior_leap_month(m_prime, M.japanese_new_moon_before(m)))
end
function M.japanese_from_fixed(date)
-- Japanese date (cycle year month leap day) of fixed date
local s1 = M.japanese_winter_solstice_on_or_before(date) -- Prior solstice
local s2 = M.japanese_winter_solstice_on_or_before(s1 + 370) -- Following solstice
local m12 = M.japanese_new_moon_on_or_after(s1 + 1) -- month after last 11th month
local next_m11 = M.japanese_new_moon_before(s2 + 1) -- next 11th month
local m = M.japanese_new_moon_before(date + 1) -- start of month containing date
local leap_year = math.floor(0.5 + (next_m11 - m12) / M.MEAN_SYNODIC_MONTH) == 12 -- if there are 13 new moons (12 full lunar months)
local prior = (leap_year and M.japanese_prior_leap_month(m12, m)) and 1 or 0 -- subtract 1 from ordinal position of month in year during or after a leap month
local month = M.amod(math.floor(0.5 + (m - m12) / M.MEAN_SYNODIC_MONTH) - prior, 12) -- month number
local leap_month = leap_year -- it's a leap month if there are 13 months
and M.japanese_no_major_solar_term(m) -- no major solar term
and (not M.japanese_prior_leap_month(m12, M.japanese_new_moon_before(m))) -- and no prior leap month
local elapsed_years = math.floor(1.5 - month/12 + (date - M.CHINESE_EPOCH) / M.MEAN_TROPICAL_YEAR_AST) -- Approximate since the epoch
local cycle = M.quotient(elapsed_years - 1, 60) + 1
local year = M.amod(elapsed_years, 60)
local day = date - m + 1
return {cycle, year, month, leap_month, day}
end
p.japanese_from_fixed = make_public(M.japanese_from_fixed)
function M.fixed_from_japanese(j_date)
-- Fixed date of Japanese date j-date.
local cycle, year, month, leap, day = j_date[1], j_date[2], j_date[3], j_date[4], j_date[5]
local mid_year = math.floor(M.CHINESE_EPOCH
+ ((cycle - 1) * 60 + (year - 1) + 1/2) * M.MEAN_TROPICAL_YEAR_AST)
local new_year = M.japanese_new_year_on_or_before(mid_year)
local p = M.japanese_new_moon_on_or_after(new_year + (month - 1) * 29)
local d = M.japanese_from_fixed(p)
local prior_new_moon
if month == d[3] and leap == d[4] then
prior_new_moon = p
else
prior_new_moon = M.japanese_new_moon_on_or_after(p + 1)
end
return prior_new_moon + day - 1
end
function M.fixed_from_japanese_for_wikipedia(table)
-- This isn't in the book, but I have to add it so the function would work from Wikipedia
return M.fixed_from_japanese({table[1], table[2], table[3], table[4] == 1, table[5]})
end
p.fixed_from_japanese = make_public(M.fixed_from_japanese_for_wikipedia, false)
----- B.20 The Modern Hindu Calendars -----
----- B.21 The Tibetan Calendar -----
----- B.22 Astronomical Lunar Calendars -----
return p