Published on

计算时间公式

Authors
  • Name
    lzray
    Twitter

基本定义与核心公式

**均时差(EoT)**定义为真太阳时平太阳时之差在钟表上的体现。按 Meeus / NREL-SPA,记:

  • L0L_0:太阳几何平均黄经(deg)
  • α\alpha:太阳视赤经(deg)
  • Δψ\Delta\psi黄经章动(deg)
  • ε\varepsilon真黄赤交角(deg)

则均时差(以分钟计)为

EoTmin=4(L00.0057183α+Δψcosε).\boxed{\mathrm{EoT}_{\text{min}}=4\big(L_0-0.0057183^\circ-\alpha+\Delta\psi\cos\varepsilon\big)}.

推导要点 太阳时由太阳时角决定:时角 == 恒星时 - 太阳赤经。以格林尼治为例,

(GASTα)真太阳时角(GMSTαm)平太阳时角=(GASTGMST)+(αmα).\underbrace{(\mathrm{GAST}-\alpha)}_{\text{真太阳时角}} -\underbrace{(\mathrm{GMST}-\alpha_{\rm m})}_{\text{平太阳时角}} =(\mathrm{GAST}-\mathrm{GMST})+(\alpha_{\rm m}-\alpha).

其中 GASTGMST=Δψcosε\mathrm{GAST}-\mathrm{GMST}=\Delta\psi\cos\varepsilon(章动沿赤道投影)。“平太阳赤经”近似取 αmL00.0057183\alpha_{\rm m}\approx L_0-0.0057183^\circ;这里的 0.00571830.0057183^\circ 来自太阳视黄经像差(约 20.59″)的常数处理。将角度差换成时间用 1=41^\circ=4 分钟(地球自转 36024360^\circ\to24 h),即得上式。号约定采用“真太阳超前为正”。


时间尺度与时间参数

所有高精度天文量以 TT(Terrestrial Time) 为自变量。若观测 UTC 为 UTC\mathrm{UTC}

TT=UTC+ΔT,ΔT=TTUT.\mathrm{TT}=\mathrm{UTC}+\Delta T,\qquad \Delta T=\mathrm{TT}-\mathrm{UT}.

ΔT\Delta T 可取 IERS 公布值或在 2005–2050 年用经验式

ΔT (s)62.92+0.32217(y2000)+0.005589(y2000)2,\Delta T~(\mathrm{s})\approx 62.92+0.32217\,(y-2000)+0.005589\,(y-2000)^2,

其中 yy 为历年(可含小数月份)。由 TT 得儒略日与世纪数:

JDTT,T=JDTT2451545.036525.\mathrm{JD}_{\mathrm{TT}},\qquad T=\frac{\mathrm{JD}_{\mathrm{TT}}-2451545.0}{36525}.

为什么要用 TT? 章动、岁差与行星理论的历书多以 TT 为自变量;用 TT 能与标准章动/岁差模型自洽,避免把不规则的地球自转(UT)抖动混入自变量。


太阳视位置的计算路径

1) 几何平均黄经、平近点角、轨道偏心率

L0=280.4664567+36000.76983T+0.0003032T2,M=357.52911+35999.05029T0.0001537T20.00000048T3,e=0.0167086340.000042037T0.0000001267T2.\begin{aligned} L_0 &= 280.4664567^\circ + 36000.76983^\circ\,T + 0.0003032^\circ\,T^2, \\ M &= 357.52911^\circ + 35999.05029^\circ\,T - 0.0001537^\circ\,T^2 - 0.00000048^\circ\,T^3, \\ e &= 0.016708634 - 0.000042037\,T - 0.0000001267\,T^2. \end{aligned}

(最终把所有角度归一化03600\sim360^\circ。)

这些多项式的来源 它们是把 VSOP/IAU 等理论中的长期项做截断近似并以 J2000 为基准后的标准表达:L0L_0 描述地球绕日的平均几何位置;MM 是平近点角;ee 是轨道偏心率随世纪缓变的拟合。

2) 中心差与真黄经

MM 的弧度为 MrM_r。太阳中心差(方程):

C=(1.9146020.004817T0.000014T2)sinMr+(0.0199930.000101T)sin(2Mr)+0.000289sin(3Mr)(deg),\begin{aligned} C &=(1.914602-0.004817\,T-0.000014\,T^2)\sin M_r\\ &\quad+(0.019993-0.000101\,T)\sin(2M_r)+0.000289\sin(3M_r)\quad(\text{deg}), \end{aligned}

真黄经λtrue=L0+C\lambda_{\text{true}}=L_0+C

中心差如何计算 解开普勒方程 EesinE=ME-e\sin E=M,由 tanν2=1+e1etanE2\tan\frac{\nu}{2}=\sqrt{\frac{1+e}{1-e}}\tan\frac{E}{2} 得真近点角 ν\nu,按 ee 展开:

νM=(2e14e3)sinM+(54e21124e4)sin2M+1312e3sin3M+\nu-M=(2e-\tfrac14e^3)\sin M+\big(\tfrac54 e^2-\tfrac{11}{24}e^4\big)\sin2M+\tfrac{13}{12}e^3\sin3M+\cdots

将其换算为“度”并代入地球的 e(T)e(T)、截到常用阶,数值系数即化为上面的 CC。因此有“真黄经 = 平黄经 + 中心差”。

3) 章动与黄赤交角

Ω=125.041934.136T,ε0=232621.44846.8150T0.00059T2+0.001813T3.\Omega=125.04^\circ-1934.136^\circ\,T,\qquad \varepsilon_0=23^\circ26'21.448''-46.8150''\,T-0.00059''\,T^2+0.001813''\,T^3.

真黄赤交角(含主章动):

ε=ε0+0.00256cosΩ.\varepsilon=\varepsilon_0+0.00256^\circ\cos\Omega.

黄经章动(SPA 精简项):

Δψ0.00478sinΩ.\Delta\psi\approx-0.00478^\circ\sin\Omega.

Ω 是什么?为何是主项? Ω\Omega月球升交点黄经的简化表达,主导 18.6 年章动周期;其 sinΩ/cosΩ\sin\Omega/\cos\Omega 项给出黄经与黄赤交角的最大振幅主项,足以把 EoT 误差压到秒量级。

4) 视黄经与视赤经

λ=λtrue0.005690.00478sinΩ+Δψ.\lambda=\lambda_{\text{true}}-0.00569^\circ-0.00478^\circ\sin\Omega+\Delta\psi.

λ,ε\lambda,\varepsilon 变换得

α=atan2(cosεsinλ, cosλ)(deg,归一到 0360),δ=arcsin(sinεsinλ).\alpha=\operatorname{atan2}\big(\cos\varepsilon\,\sin\lambda,\ \cos\lambda\big)\quad(\text{deg,归一到 }0\sim360^\circ), \quad \delta=\arcsin(\sin\varepsilon\,\sin\lambda).

小贴士

  • 这里的 0.00569-0.00569^\circ 是像差的常数处理,0.00478sinΩ-0.00478^\circ\sin\Omega 是章动对视黄经的主修正,与前述 Δψ\Delta\psi 项配套。
  • 计算 atan2\operatorname{atan2} 后要做象限归一化;赤经单位用“度”时再参与 EoT 公式。

计算流程

  1. 输入: 公历日期与时刻(UTC),或本地时并知时区。
  2. ΔT\Delta T,由 UTC \to TT,得 JDTT\mathrm{JD}_{\mathrm{TT}}TT
  3. 按上式计算 L0,M,e,C,λtrueL_0,\,M,\,e,\,C,\,\lambda_{\text{true}}
  4. Ω,ε0,ε,Δψ\Omega,\,\varepsilon_0,\,\varepsilon,\,\Delta\psi(精简或长级数)。
  5. 视黄经 λ\lambda,再求视赤经 α\alpha
  6. 代入 核心式,输出 EoTmin\mathrm{EoT}_{\text{min}}(分钟)。

实现要点

  • 角度单位: 三角函数用弧度;多项式常数若以度给出,先转弧度再取 sin,cos\sin,\cos
  • 归一化: 所有经度/赤经建议归一到 03600\sim360^\circ(或 180 ⁣ ⁣180-180^\circ\!\sim\!180^\circ)以避免溢出。
  • 角度→时间: 1=41^\circ=4 分钟;弧度直接换分钟用系数 229.18=14402π229.18=\dfrac{1440}{2\pi}
  • 精度取舍: 仅用 Ω\Omega 主项的 Δψ,Δε\Delta\psi,\Delta\varepsilon 通常把 EoT 误差压到 \sim1 s;若需 <0.1\lt0.1 s,使用 IAU 2000/2006 章动长级数与更精确的 ΔT\Delta T

近似式(NOAA 谐波)

便于快速估算,常见 NOAA 近似把 EoT 写成年分角的谐波:

EoTmin229.18 ⁣(0.000075+0.001868cosγ0.032077sinγ0.014615cos2γ0.040849sin2γ),γ=2πN1365,\mathrm{EoT}_{\text{min}}\approx 229.18\!\left( 0.000075+0.001868\cos\gamma-0.032077\sin\gamma -0.014615\cos2\gamma-0.040849\sin2\gamma \right), \quad \gamma=2\pi\frac{N-1}{365},

其中 NN 为年内日序。

这些系数的简单物理含义 把 EoT 拆成两部分: 椭圆项(轨道偏心 ee 导致日速不匀)与倾斜项(地轴倾角 ε\varepsilon 导致黄道投影到赤道的非线性)。令 y=tan2(ε/2)0.043y=\tan^2(\varepsilon/2)\approx0.043,对小量 e0.0167e\simeq0.0167yy 进行低阶展开,可得

EoT229.18 ⁣(ysin2L02esinM+4eysinMcos2L012y2sin4L054e2sin2M).\mathrm{EoT}\approx 229.18\!\left( y\sin2L_0-2e\sin M+4ey\sin M\cos2L_0-\tfrac12 y^2\sin4L_0-\tfrac54 e^2\sin2M \right).

再把 L0,ML_0,M 用近似线性相位改写为统一的年分角 γ\gamma,合并相位与振幅,便得到 NOAA 的 sinγ,cosγ,sin2γ,cos2γ\sin\gamma,\cos\gamma,\sin2\gamma,\cos2\gamma 结构与上面的系数。常数 229.18229.18 正是“弧度→分钟”的换算因子。