实历冬至合朔表（BC722-BC104）

利用python的pyephem模块（需先安装）可以很容易获得合朔及冬至时间。但需注意，该模块在历元2000.0年±3000年的平均误差在10分钟以内，但此外的误差是难以接受的。

关于朔旦冬至的判断：冬至月判断可能存在的几个误差。

现代测算误差：现代的天文计算，是近代天文学发展起来以后利用精密仪器测量天体的实际运行状况，由此构建的理论模型。但天体的运动状况非常复杂，当今的测量结果能难适用于久远前或久远后的宇宙，距离时间越远，误差越大。对于先秦时期的历表，应注意定朔或冬至时间在0-1时或23-24时段，可能几分钟的误差在进位时表现为1天。
古代测算误差：战国前后是以冬至所在的月份作为正月，但月朔是无法观测的，推朔会有误差，而冬至的测量也存在误差。以《左传》两次“日南至”的记载为例，皆与今测算差2日。

历表的排布：每个冬至合朔点有3项数据（干支、日期、时间），其中时间部分主要用于参考可能的误差问题，分三行排列。每年从冬至月起排，如果到下一个冬至月出现了13次合朔，该年则排出13个月份，但并不表示该年存在闰月，月名仅供参考。古代的置闰会因为不同的历元和置闰方法有变化。如果连续两个月朔日干支的天干部分相同，说明前一个月为大月，否则为小月。如果连续三个月朔日干支的天干部分相同，说明前两个月为连大月。古历的朔策较今值为大，大小月的情况也与实历不符。

wangnian部分的代码参见：https://blog.csdn.net/weixin_42763614/article/details/82955539

import ephem
from wangnian import *

nlrq = ["初一","初二","初三","初四","初五","初六","初七","初八","初九","初十","十一","十二","十三","十四","十五","十六","十七","十八","十九","二十","二一","二二","二三","二四","二五","二六","二七","二八","二九","三十"]

def jd2ce(JD):
	JD = JD + 0.5  # 以-4812年1月1日0时为历元
	Z = int(JD)
	F = JD - Z
	if Z < 2299161:  # 儒略历
		A = Z
	else:  # 格里历
		a = int((Z - 1867216.25) / 36524.25)
		A = Z + 1 + a - int(a / 4)
	B = A + 1524
	C = int((B - 122.1) / 365.25)
	D = int(365.25 * C)
	E = int((B - D) / 30.6001)
	day = B - D - int(30.6001 * E) + F
	if E < 14:
		month = E - 1
	elif E < 16:
		month = E - 13
	if month > 2:
		year = C - 4716
	elif month in [1, 2]:
		year = C - 4715
	d = day - int(day) #取出一日的小数部分，转换为时分秒
	h = int(d * 24)
	m = int((d * 24 - h) * 60)
	s = ((d *24 - h) * 60 - int((d *24 - h) * 60)) * 60
	if year <= 0: # 将天文年表示为公历年
		year -= 1
	return year, month, int(day),h,m,round(s)

def ut8(t): # 由UT0h转为UT8h
	d = ephem.Date(t)
	date = d.tuple() # 获得年月日时分秒
	jd0 = ephem.julian_date(t)
	jd = jd0 + 8/24
	ce = jd2ce(jd)
	return ce

def jd2gz(t):
	JD = ephem.julian_date(t)
	tiangan = ["甲", "乙", "丙", "丁", "戊", "己", "庚", "辛", "壬", "癸"]
	dizhi = ["子", "丑", "寅", "卯", "辰", "巳", "午", "未", "申", "酉", "戌", "亥"]
	gz = [''] * 60  # 六十甲子表
	for i in range(60):		gz[i] = tiangan[i % 10] + dizhi[i % 12]
	return gz[int(JD+0.5+8/24+49)%60]

def time_format(t,type):
	if type == 0: # 日期
		return "BC{}-{:02d}-{:02d}".format(-ut8(t)[0],ut8(t)[1],ut8(t)[2])
	elif type == 1: # 时间
		return "{:02d}:{:02d}:{:02d}".format(ut8(t)[3], ut8(t)[4], ut8(t)[5])

header = ["参考年","王公纪年","冬至",'冬至日期',"冬至月朔","二月（丑）","三月（寅）","四月","五月","六月","七月","八月","九月","十月","十一月","十二月","十三月"]
dzhsb = open("冬至合朔表.xls",'w')
for i in range(len(header)):
	dzhsb.write(header[i]+'\t')

for i in range(620): # BC722-BC104
	year = -723 + i
	time = str(year)+'/12'
	dz = ephem.next_solstice(time) # 冬至时间
	dz_jd = ephem.julian_date(dz)  # 冬至儒略日
	dzs0 = jd2ce(dz_jd - 29) # 冬至前一月
	dzs = ephem.next_new_moon(dzs0) # 预判的冬至朔
	dzs_jd = ephem.julian_date(dzs) # 在UT8h下比较冬至朔是否为实际的冬至所在月朔日
	if int(dz_jd+0.5+8/24) < int(dzs_jd+0.5+8/24):
		dzs = ephem.next_new_moon(jd2ce(dz_jd - 30))
		dzs_jd = ephem.julian_date(dzs)
	# 写入上一年历表
	if i != 0:
		if str(heshuo) == str(dzs):
			hsb.pop() # 上年最后一月朔为是年冬至朔，则删除
		dzhsb.write('\n')
		for m in range(3):
			dzhsb.write(hsb[0][0]+'\t')  # 参考年
			dzhsb.write(hsb[0][1]+'\t')  # 王公纪年年
			for n in range(len(hsb)): # 月表
				if n == 0 : continue
				if n == 2 : # 在冬至后插入冬至的日期
					if m == 1: dzhsb.write(dzrq)
					dzhsb.write('\t')
				dzhsb.write(hsb[n][m+0]+'\t')
			dzhsb.write('\n')
		hsb.clear()
	# 生成历表
	hsb = [[str(year + 1), wgjn[841 + year + 1]]]
	hsb.extend([[jd2gz(dz),time_format(dz,0),time_format(dz,1)]]) # 干支，年月日，时分秒
	hsb.extend([[jd2gz(dzs),time_format(dzs,0),time_format(dzs,1)]])
	dzrq = nlrq[int(dz_jd+0.5+8/24)-int(dzs_jd+0.5+8/24)]
	for j in range(12): # 冬至以后的每月合朔表
		hs = jd2ce(dzs_jd + 29*(j+1)) # 估测的合朔日期
		heshuo = ephem.next_new_moon(hs) # 实际合朔日期
		hsb.extend([[jd2gz(heshuo),time_format(heshuo,0),time_format(heshuo,1)]])

实历冬至合朔表（BC722-BC104）

猜你喜欢