/*-
 * SPDX-License-Identifier: BSD-3-Clause
 *
 * Copyright (c) 2019-2020 The DragonFly Project.  All rights reserved.
 *
 * This code is derived from software contributed to The DragonFly Project
 * by Aaron LI <aly@aaronly.me>
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 *    notice, this list of conditions and the following disclaimer.
 * 2. Redistributions in binary form must reproduce the above copyright
 *    notice, this list of conditions and the following disclaimer in
 *    the documentation and/or other materials provided with the
 *    distribution.
 * 3. Neither the name of The DragonFly Project nor the names of its
 *    contributors may be used to endorse or promote products derived
 *    from this software without specific, prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
 * FOR A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE
 * COPYRIGHT HOLDERS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
 * INCIDENTAL, SPECIAL, EXEMPLARY OR CONSEQUENTIAL DAMAGES (INCLUDING,
 * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED
 * AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
 * OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
 * OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
 * SUCH DAMAGE.
 *
 * Reference:
 * Calendrical Calculations, The Ultimate Edition (4th Edition)
 * by Edward M. Reingold and Nachum Dershowitz
 * 2018, Cambridge University Press
 */

/*
 * Implement the Chinese calendar of the 1645 version, established in the
 * second year of the Qīng dynasty.
 *
 * The winter solstice (dōngzhì; 冬至) always occurs during the eleventh
 * month of the year.  The winter-solstice-to-winter-solstice period is
 * called a suì (岁).
 *
 * The leap month of a 13-month winter-solstice-to-winter-solstice period
 * is the first month that does not contain a major solar term --- that is,
 * the first lunar month that is wholly within a solar month.
 */

#include <math.h>
#include <stdbool.h>
#include <stdio.h>

#include "calendar.h"
#include "basics.h"
#include "chinese.h"
#include "dates.h"
#include "gregorian.h"
#include "moon.h"
#include "sun.h"
#include "utils.h"

/*
 * Fixed date of the start of the Chinese calendar.
 * Ref: Sec.(19.3), Eq.(19.15)
 */
static const int epoch = -963099;  /* Gregorian: -2636, February, 15 */

/*
 * Timezone (in fraction of days) of Beijing adopted in Chinese calendar
 * calculations.
 * Ref: Sec.(19.1), Eq.(19.2)
 */
static double
chinese_zone(int rd)
{
	if (gregorian_year_from_fixed(rd) < 1929)
		return (1397.0 / 180.0 / 24.0);
	else
		return (8.0 / 24.0);
}

/*
 * The universal time (UT) of (clock) midnight at the start of the fixed date
 * $rd in China.
 * Ref: Sec.(19.1), Eq.(19.7)
 */
static double
midnight_in_china(int rd)
{
	return (double)rd - chinese_zone(rd);
}

/*
 * Calculate the last Chinese major solar term (zhōngqì) in range of [1, 12]
 * before the fixed date $rd.
 * Ref: Sec.(19.1), Eq.(19.1)
 */
static int
current_major_solar_term(int rd)
{
	double t_u = midnight_in_china(rd);
	double lon = solar_longitude(t_u);
	return mod1(2 + div_floor((int)lon, 30), 12);
}

/*
 * Calculate the fixed date (in China) of winter solstice on or before the
 * given fixed date $rd.
 * Ref: Sec.(19.1), Eq.(19.8)
 */
static int
chinese_winter_solstice_onbefore(int rd)
{
	/* longitude of Sun at winter solstice */
	const double winter = 270.0;

	/* approximate the time of winter solstice */
	double t = midnight_in_china(rd + 1);
	double approx = estimate_prior_solar_longitude(winter, t);

	/* search for the date of winter solstice */
	int day = (int)floor(approx) - 1;
	while (winter >= solar_longitude(midnight_in_china(day+1)))
		day++;

	return day;
}

/*
 * Calculate the fixed date (in China) of the first new moon on or after the
 * given fixed date $rd.
 * Ref: Sec.(19.2), Eq.(19.9)
 */
static int
chinese_new_moon_onafter(int rd)
{
	double t = new_moon_atafter(midnight_in_china(rd));
	double t_s = t + chinese_zone(rd);  /* in standard time */

	return (int)floor(t_s);
}

/*
 * Calculate the fixed date (in China) of the first new moon before the
 * given fixed date $rd.
 * Ref: Sec.(19.2), Eq.(19.10)
 */
static int
chinese_new_moon_before(int rd)
{
	double t = new_moon_before(midnight_in_china(rd));
	double t_s = t + chinese_zone(rd);  /* in standard time */

	return (int)floor(t_s);
}

/*
 * Determine whether the Chinese lunar month starting on the given fixed
 * date $rd has no major solar term.
 * Ref: Sec.(19.2), Eq.(19.11)
 */
static bool
chinese_no_major_solar_term(int rd)
{
	int rd2 = chinese_new_moon_onafter(rd + 1);
	return (current_major_solar_term(rd) ==
		current_major_solar_term(rd2));
}

/*
 * Recursively determine whether there is a Chinese leap month on or
 * after the lunar month starting on fixed date $m1 and at or before
 * the lunar month starting on fixed date $m2.
 * Ref: Sec.(19.2), Eq.(19.12)
 */
static bool
chinese_prior_leap_month(int m1, int m2)
{
	int m2_prev = chinese_new_moon_before(m2);
	return ((m2 >= m1) &&
		(chinese_no_major_solar_term(m2) ||
		 chinese_prior_leap_month(m1, m2_prev)));
}

/*
 * Calculate the fixed date of Chinese New Year in suì containing the
 * given date $rd.
 * Ref: Sec.(19.2), Eq.(19.13)
 */
static int
chinese_new_year_in_sui(int rd)
{
	/* prior and following winter solstice */
	int s1 = chinese_winter_solstice_onbefore(rd);
	int s2 = chinese_winter_solstice_onbefore(s1 + 370);

	/* month after the 11th month: either 12 or leap 11 */
	int m12 = chinese_new_moon_onafter(s1 + 1);
	/* month after m12: either 12 (or leap 12) or 1 */
	int m13 = chinese_new_moon_onafter(m12 + 1);
	/* the next 11th month */
	int m11_next = chinese_new_moon_before(s2 + 1);

	bool leap_year = (lround((m11_next - m12) / mean_synodic_month) == 12);
	if (leap_year && (chinese_no_major_solar_term(m12) ||
			  chinese_no_major_solar_term(m13))) {
		/* either m12 or m13 is a leap month */
		return chinese_new_moon_onafter(m13 + 1);
	} else {
		return m13;
	}
}

/*
 * Calculate the fixed date of Chinese New Year on or before the given
 * fixed date $rd.
 * Ref: Sec.(19.2), Eq.(19.14)
 */
static int
chinese_new_year_onbefore(int rd)
{
	int newyear = chinese_new_year_in_sui(rd);
	if (rd >= newyear)
		return newyear;
	else
		return chinese_new_year_in_sui(rd - 180);
}

/*
 * Calculate the fixed date of Chinese New Year in Gregorian year $year.
 * Ref: Sec.(19.6), Eq.(19.26)
 */
int
chinese_new_year(int year)
{
	struct date date = { year, 7, 1 };
	int july1 = fixed_from_gregorian(&date);
	return chinese_new_year_onbefore(july1);
}

/*
 * Calculate the Chinese date (cycle, year, month, leap, day) corresponding
 * to the fixed date $rd.
 * Ref: Sec.(19.3), Eq.(19.16)
 */
void
chinese_from_fixed(int rd, struct chinese_date *date)
{
	/* prior and following winter solstice */
	int s1 = chinese_winter_solstice_onbefore(rd);
	int s2 = chinese_winter_solstice_onbefore(s1 + 370);

	/* start of month containing the given date */
	int m = chinese_new_moon_before(rd + 1);
	/* start of the previous month */
	int m_prev = chinese_new_moon_before(m);
	/* month after the last 11th month: either 12 or leap 11 */
	int m12 = chinese_new_moon_onafter(s1 + 1);
	/* the next 11th month */
	int m11_next = chinese_new_moon_before(s2 + 1);

	bool leap_year = (lround((m11_next - m12) / mean_synodic_month) == 12);
	int month = (int)lround((m - m12) / mean_synodic_month);
	if (leap_year && chinese_prior_leap_month(m12, m))
		month--;
	month = mod1(month, 12);
	bool leap_month = (leap_year &&
			   chinese_no_major_solar_term(m) &&
			   !chinese_prior_leap_month(m12, m_prev));
	int elapsed_years = (int)floor(1.5 - month/12.0 +
				       (rd - epoch) / mean_tropical_year);

	date->cycle = div_floor(elapsed_years - 1, 60) + 1;
	date->year = mod1(elapsed_years, 60);
	date->month = month;
	date->leap = leap_month;
	date->day = rd - m + 1;
}

/*
 * Calculate the fixed date corresponding to the given Chinese date $date
 * (cycle, year, month, leap, day).
 * Ref: Sec.(19.3), Eq.(19.17)
 */
int
fixed_from_chinese(const struct chinese_date *date)
{
	int midyear = (int)floor(epoch + mean_tropical_year *
				 ((date->cycle - 1) * 60 + date->year - 0.5));
	int newyear = chinese_new_year_onbefore(midyear);

	/* new moon before the given date */
	int newmoon = chinese_new_moon_onafter(
			newyear + (date->month - 1) * 29);
	struct chinese_date date2;
	chinese_from_fixed(newmoon, &date2);
	if (date->month != date2.month || date->leap != date2.leap) {
		/* there was a prior leap month, so get the next month */
		newmoon = chinese_new_moon_onafter(newmoon + 1);
	}

	return newmoon + date->day - 1;
}

/*
 * Calculate the fixed date of Qīngmíng (清明) in Gregorian year $g_year.
 * Ref: Sec.(19.6), Eq.(19.28)
 */
int
chinese_qingming(int g_year)
{
	double lambda = 15.0;  /* Solar longitude of Qīngmíng */
	struct date date = { g_year, 4, 1 }; /* Qīngmíng is around April 5 */
	int rd = fixed_from_gregorian(&date);
	double zone = chinese_zone(rd);
	double t = solar_longitude_atafter(lambda, rd) + zone;
	return (int)floor(t);
}

/**************************************************************************/

/* celestial stem, tiāngān (天干) */
static const struct stem {
	const char *name;
	const char *zhname;
} STEMS[] = {
	{ "Jiǎ",  "甲" },
	{ "Yǐ",   "乙" },
	{ "Bǐng", "丙" },
	{ "Dīng", "丁" },
	{ "Wù",   "戊" },
	{ "Jǐ",   "己" },
	{ "Gēng", "庚" },
	{ "Xīn",  "辛" },
	{ "Rén",  "壬" },
	{ "Guǐ",  "癸" },
};

/* terrestrial branch, dìzhī (地支) */
static const struct branch {
	const char *name;
	const char *zhname;
	const char *zodiac;
	const char *zhzodiac;
} BRANCHES[] = {
	{ "Zǐ",   "子", "Rat",     "鼠" },
	{ "Chǒu", "丑", "Ox",      "牛" },
	{ "Yín",  "寅", "Tiger",   "虎" },
	{ "Mǎo",  "卯", "Rabbit",  "兔" },
	{ "Chén", "辰", "Dragon",  "龙" },
	{ "Sì",   "巳", "Snake",   "蛇" },
	{ "Wǔ",   "午", "Horse",   "马" },
	{ "Wèi",  "未", "Goat",    "羊" },
	{ "Shēn", "申", "Monkey",  "猴" },
	{ "Yǒu",  "酉", "Rooster", "鸡" },
	{ "Xū",   "戌", "Dog",     "狗" },
	{ "Hài",  "亥", "Pig",     "猪" },
};

/* Chinese names of months and days */
static const char *months[] = {
	"正", "二", "三", "四", "五", "六",
	"七", "八", "九", "十", "冬", "腊",
};
static const char *mdays[] = {
	"初一", "初二", "初三", "初四", "初五", "初六",
	"初七", "初八", "初九", "初十", "十一", "十二",
	"十三", "十四", "十五", "十六", "十七", "十八",
	"十九", "二十", "廿一", "廿二", "廿三", "廿四",
	"廿五", "廿六", "廿七", "廿八", "廿九", "三十",
};

/* 24 major and minor solar terms (节气) */
static const struct chinese_jieqi jieqis[] = {
	{ "Lìchūn",      "立春", false, 315 },
	{ "Yǔshuǐ",      "雨水", true,  330 },
	{ "Jīngzhé",     "惊蛰", false, 345 },
	{ "Chūnfēn",     "春分", true,  0   },
	{ "Qīngmíng",    "清明", false, 15  },
	{ "Gǔyǔ",        "谷雨", true,  30  },
	{ "Lìxià",       "立夏", false, 45  },
	{ "Xiǎomǎn",     "小满", true,  60  },
	{ "Mángzhòng",   "芒种", false, 75  },
	{ "Xiàzhì",      "夏至", true,  90  },
	{ "Xiǎoshǔ",     "小暑", false, 105 },
	{ "Dàshǔ",       "大暑", true,  120 },
	{ "Lìqiū",       "立秋", false, 135 },
	{ "Chǔshǔ",      "处暑", true,  150 },
	{ "Báilù",       "白露", false, 165 },
	{ "Qiūfēn",      "秋分", true,  180 },
	{ "Hánlù",       "寒露", false, 195 },
	{ "Shuāngjiàng", "霜降", true,  210 },
	{ "Lìdōng",      "立冬", false, 225 },
	{ "Xiǎoxuě",     "小雪", true,  240 },
	{ "Dàxuě",       "大雪", false, 255 },
	{ "Dōngzhì",     "冬至", true,  270 },
	{ "Xiǎohán",     "小寒", false, 285 },
	{ "Dàhán",       "大寒", true,  300 },
};

/*
 * Calculate the fixed date (in China) of the following jiéqì on or after
 * the given fixed date $rd, with the calculated jiéqì information stored
 * in $jieqi.
 * If $type is C_JIEQI_ALL, find either major or minor jiéqì;
 * If $type is C_JIEQI_MAJOR, find only major jiéqì;
 * If $type is C_JIEQI_MINOR, find only minor jiéqì.
 */
int
chinese_jieqi_onafter(int rd, int type, const struct chinese_jieqi **jieqi)
{
	const struct chinese_jieqi *jq1, *jq2;
	const double zone = chinese_zone(rd);
	double t_u = midnight_in_china(rd);
	double lon = solar_longitude(t_u);
	double lon1, lon2;
	size_t n = nitems(jieqis);

	for (size_t i = 0; i < n; i++) {
		jq1 = &jieqis[i];
		jq2 = &jieqis[(i+1) % n];
		if ((type == C_JIEQI_MINOR && jq2->is_major) ||
		    (type == C_JIEQI_MAJOR && !jq2->is_major)) {
			jq2 = &jieqis[(i+2) % n];
		}

		lon1 = jq1->longitude;
		lon2 = jq2->longitude;
		if (lon2 < lon1)
			lon2 += 360;

		if ((lon >= lon1 && lon < lon2) ||
		    (lon + 360 >= lon1 && lon + 360 < lon2)) {
			break;
		}
	}

	t_u = solar_longitude_atafter(jq2->longitude, floor(t_u));
	*jieqi = jq2;
	return (int)floor(t_u + zone);
}


/*
 * Format the Chinese date of the given fixed date $rd in $buf.
 * Return the formatted string length.
 */
int
chinese_format_date(char *buf, size_t size, int rd)
{
	struct chinese_date cdate;
	chinese_from_fixed(rd, &cdate);
	return snprintf(buf, size, "%s%s月%s",
			cdate.leap ? "闰" : "", months[cdate.month - 1],
			mdays[cdate.day - 1]);
}


/*
 * Find days of the specified Chinese month ($month) and day ($day).
 * If month $month < 0, then month is ignored (i.e., all months).
 * (NOTE: The year $year is ignored.)
 */
int
chinese_find_days_ymd(int year __unused, int month, int day,
		      struct cal_day **dayp, char **edp __unused)
{
	struct cal_day *dp;
	struct chinese_date cdate;
	int rd, rd_begin, rd_end;
	int count = 0;

	rd_begin = chinese_new_moon_before(Options.day_begin);
	rd_end = chinese_new_moon_onafter(Options.day_end);
	rd = rd_begin;
	while (rd <= rd_end) {
		chinese_from_fixed(rd, &cdate);
		if (cdate.month == month || month < 0) {
			cdate.day = day;
			rd = fixed_from_chinese(&cdate);
			if ((dp = find_rd(rd, 0)) != NULL) {
				if (count >= CAL_MAX_REPEAT) {
					warnx("%s: too many repeats",
					      __func__);
					return count;
				}
				dayp[count++] = dp;
			}
		}

		/* go to next month */
		rd = chinese_new_moon_onafter(rd + (day == 0 ? 2 : 1));
	}

	return count;
}

/*
 * Find days of the specified Chinese day of month ($dom) of all months.
 */
int
chinese_find_days_dom(int dom, struct cal_day **dayp, char **edp __unused)
{
	return chinese_find_days_ymd(-1, -1, dom, dayp, NULL);
}


/*
 * Print the Chinese calendar of the given date $rd and events of the year.
 */
void
show_chinese_calendar(int rd)
{
	struct chinese_date date;
	chinese_from_fixed(rd, &date);
	/* Ref: Sec.(19.3), Eq.(19.18) */
	struct stem stem = STEMS[mod1(date.year, 10) - 1];
	struct branch branch = BRANCHES[mod1(date.year, 12) - 1];

	struct date gdate;
	gregorian_from_fixed(rd, &gdate);
	int g_year = gdate.year;
	printf("公历 (Gregorian): %d-%02d-%02d\n",
	       gdate.year, gdate.month, gdate.day);

	printf("农历: %s%s年 [%s], %s%s月%s\n",
	       stem.zhname, branch.zhname, branch.zhzodiac,
	       date.leap ? "闰" : "", months[date.month - 1],
	       mdays[date.day - 1]);
	printf("Chinese calendar: year %s%s [%s], %smonth %d, day %d\n",
	       stem.name, branch.name, branch.zodiac,
	       date.leap ? "leap " : "", date.month, date.day);

	/* the following Chinese New Year */
	int newyear = chinese_new_year_onbefore(
			chinese_new_year_onbefore(rd) + 30*13);
	gregorian_from_fixed(newyear, &gdate);
	printf("新年 (Chinese New Year): %d-%02d-%02d\n",
	       gdate.year, gdate.month, gdate.day);

	const double zone = chinese_zone(rd);
	const struct chinese_jieqi *jq;
	double t_jq;
	char buf_time[32], buf_zone[32];
	int lambda;

	/* 1st solar term (Lìchūn) is generally around February 4 */
	date_set(&gdate, g_year, 2, 1);
	t_jq = (double)fixed_from_gregorian(&gdate);

	format_zone(buf_zone, sizeof(buf_zone), zone);

	printf("\n二十四节气 (solar terms):\n");
	for (size_t i = 0; i < nitems(jieqis); i++) {
		jq = &jieqis[i];
		lambda = jq->longitude;
		t_jq = solar_longitude_atafter(lambda, floor(t_jq)) + zone;
		gregorian_from_fixed(floor(t_jq), &gdate);
		format_time(buf_time, sizeof(buf_time), t_jq);

		printf("%s (%-13s): %3d°, %d-%02d-%02d %s %s\n",
		       jq->zhname, jq->name, lambda,
		       gdate.year, gdate.month, gdate.day,
		       buf_time, buf_zone);
	}
}