OPM：一种面向有界误差部署的紧凑主要天体星历表示及其 600 年密集验证

English title: OPM: a compact deployable ephemeris representation for major Solar-System bodies with 600-year dense validation

作者： Rz Liu

摘要

背景： JPL Development Ephemeris 系列为太阳系天体位置计算提供了高精度参考 [1]，但完整星历内核的数据体积和运行时依赖并不总是适合客户端、嵌入式或网络分发场景。将参考数值星历转换为紧凑的分段多项式表示，是降低星历分发和运行时成本的一条实用路线；Kammeyer (1988) [2] 的 Compressed Planetary and Lunar Ephemerides 是其中一个早期完整工程案例。现代应用仍然需要在文件大小、读取速度和误差控制之间取得平衡。

目标： 本文提出 OPM（Orbital Polynomial Model），一种面向运行时部署的紧凑主要天体星历表示。OPM 的目标不是替代动力学积分历表，而是在给定覆盖范围内提供小体积、快速读取、可校验且具有稳定尾部误差表现的位置重建模型。

方法： OPM 的设计借鉴了 Kammeyer (1988) [2] 在 Compressed Planetary and Lunar Ephemerides 中描述的若干方法：通过按天体配置的局部坐标和参考形状扣除移除主要轨道结构，存储残差切比雪夫系数，将其量化为整数，并按实际有效位宽进行按阶主序位级打包。OPM 使用分块容器、按天体区分的模型 ID、分节元数据和显式残差系数表。生成流程从 DE441 [1] 采样拟合初始分段模型，并进一步使用按天体类型划分的修正流程压低高分位和最坏情况误差。太阳使用原生千米度量修正；水星、金星和月球使用原生角度保护修正；以 SSB 存储的天体使用已经修正后的太阳模型作为锚点进行日心目标修正。验证时，本文以 DE441 为真值，将 OPM 和 Swiss Ephemeris [3] 在相同儒略日网格上转换到地心几何 J2000 ICRS 方向，并统计角误差分位数和最大值。

结果： 在以 J2000 为中心的 600 年范围内（约 JD 2378460—2597711，即约 1799 年 12 月至 2400 年 1 月；精确端点随天体覆盖交集略有不同），本文使用 OPM 段结构构造了每段 512 个切比雪夫节点加端点的密集确定性网格，10 个可由 Swiss Ephemeris 直接寻址的主要天体合计包含 6,810,404 个地心测试点。OPM 在这些天体上的第 99 百分位角误差最高约为 0.00119 角秒，其中 9 个天体低于 0.001 角秒。7 个天体的最坏情况误差低于 0.001 角秒；水星约为 0.00111 角秒，金星和火星分别约为 0.00309 和 0.00217 角秒。在同一确定性网格上，OPM 对全部 10 个天体的最坏情况地心角误差均低于 Swiss Ephemeris；最大的 OPM/Swiss 最坏值比约为 0.62。无外延消融显示，按天体设置的分段定义域外延可使太阳、地月质心、火星、木星、天王星、海王星和冥王星的解析速度最大残差降低约 41%—86%。该 600 年生产实例的 11 个 OPM 文件合计 1.634 MiB，折合约 278.93 KiB/世纪。

结论： OPM 的优势主要体现在误差尾部控制，而不是单纯追求中位误差最小。密集确定性验证显示，选定的 600 年生产实例可以在较小、可随机访问的表示中提供稳定的高分位误差表现，并且在该验证网格上全部天体的最坏情况地心角误差均低于 Swiss Ephemeris。剩余最大的尾部主要来自金星和火星的地心组合情形，这也提示后续可以继续发展组合目标修正和更长时间跨度生产验证。OPM 在本文验证范围内具备作为紧凑星历表示的实用潜力。

关键词： 星历 — 数值方法 — 天体力学 — 天体测量 — 软件：数据压缩

1. 引言

现代太阳系天体位置计算通常依赖高精度数值星历，例如 JPL Development Ephemeris 系列。DE440/DE441 等历表 [1] 在天体动力学建模、观测拟合和长期覆盖方面具有很高精度，但其完整数据文件体积较大，并且运行时通常需要读取和解释 SPK 内核。这对服务器端科学计算并非问题，但对客户端应用、网络分发、移动端、嵌入式环境以及需要快速冷启动的工程场景而言，完整数值星历并不总是最合适的表示。对于这类部署场景，应用往往更需要一种有限覆盖范围内的小体积、快速读取、跨平台一致且可独立验证的星历表示。

高精度星历的紧凑表示已有较长历史。Kammeyer (1988) [2] 在 Compressed Planetary and Lunar Ephemerides 中提出了一个早期而完整的工程方案：以 DE200 为参考，将 1801—2049 年的太阳、月球和行星位置压缩为约 830 KB 的数据文件，并通过 40 阶切比雪夫级数、局部坐标轴、内行星和月球的参考轨道扣除，以及整数量化系数的位打包，将位置误差控制在约 1 毫角秒（0.001 角秒）量级。该工作表明，在给定参考星历和有限时间范围内，从高精度数值星历派生紧凑切比雪夫表示是一条可行的工程路线。

OPM 的设计借鉴了 Kammeyer 论文 [2] 中描述的若干关键方法：将参考数值星历分段表示为切比雪夫级数，先通过局部坐标和参考轨道/参考形状扣除移除主要轨道结构，再将残差系数量化为整数并按有效位宽进行位级打包。由于 Kammeyer 原文主要描述系统结构和关键数值思想，许多工程细节需要在具体实现中明确指定。本文在同一思想框架下，将参考星历更新为 DE441 [1]，并进一步加入可随机访问的现代二进制封装、CRC 校验、按天体类型划分的误差修正目标，以及面向最坏情况误差的密集确定性验证。Swiss Ephemeris [3] 可视为这类紧凑星历思想在后续实际软件生态中的成熟工程代表之一。它通过压缩星历文件和统一运行时接口，为大量桌面、服务器和嵌入式应用提供了小体积星历能力。因此，本文将 Swiss Ephemeris 作为现实紧凑星历基线，而不是仅与原始 DE 内核比较。

本文提出 OPM（Orbital Polynomial Model），一种面向主要太阳系天体的紧凑二进制星历表示。OPM 不是新的动力学积分历表，而是从 DE441 派生出的分段多项式重建模型。它将主要太阳系天体的位置表示为分段切比雪夫多项式，并将多项式系数量化后写入可随机访问的二进制文件。与只关注平均误差或随机采样误差不同，OPM 的生成流程以高分位和最坏情况误差为主要修正目标。对于可部署星历，误差尾部是评价覆盖范围内稳定性的重要指标。用户查询时间不受验证采样分布约束；如果某些局部时间区域存在尾部尖峰，随机测试可能低估最坏情况误差。

为了评估该表示的实际精度，本文使用 DE441 作为真值来源，并将 OPM 与 Swiss Ephemeris 在同一地心几何 J2000 ICRS 口径下进行比较。验证范围为以 J2000 附近为中心的 600 年区间，约为 JD 2378460—2597711，即约从 1799 年 12 月到 2400 年 1 月；精确端点随天体覆盖交集略有不同。验证不采用单纯随机采样，而是基于 OPM 分段结构构造确定性密集网格：每个分段内使用 512 个切比雪夫节点加端点。OPM 和 Swiss Ephemeris 在完全相同的儒略日网格上计算地心向量，并与 DE441 派生向量比较方向角误差。这一验证方式的目标是尽可能暴露局部尾部误差，而不是估计某种用户查询分布下的平均表现。

在该 600 年范围内，10 个可由 Swiss Ephemeris 直接寻址的主要天体合计包含 6,810,404 个地心测试点。OPM 对多数主要天体具有稳定的亚毫角秒级尾部精度：10 个天体中有 9 个第 99 百分位低于 0.001 角秒，最高第 99 百分位约为 0.00119 角秒。7 个天体的最坏情况误差低于 0.001 角秒；水星最坏值约为 0.00111 角秒，金星和火星分别约为 0.00309 和 0.00217 角秒。在同一验证网格上，OPM 对全部 10 个天体的最坏情况地心角误差均低于 Swiss Ephemeris，最大的 OPM/Swiss 最坏值比约为 0.62。

本文其余部分结构如下。第 2 节介绍 OPM 的分段切比雪夫表示、坐标约定和二进制封装方式。第 3 节描述从 DE441 生成 OPM 的初始拟合与按天体类型划分的误差修正流程。第 4 节给出验证口径，包括地心重建、Swiss Ephemeris 标志位和密集网格构造。第 5 节报告 600 年密集比较结果。第 6 节讨论若干尝试过但不适合作为本文主线的候选压缩路线。第 7 节讨论 OPM 的适用场景、与 Swiss Ephemeris 及 Kammeyer (1988) 的关系和当前限制。第 8 节总结本文结论并给出后续工作方向。

2. OPM 表示

OPM 将一个天体在给定时间范围内的位置表示为一组连续时间段上的切比雪夫多项式。OPM 文件采用分节二进制容器：小型固定头部用于识别文件和节表，覆盖范围、网格与模型参数、参考形状数据、量化残差系数、元数据和完整性校验分别放在不同 section 中。运行时输入儒略日后，可通过段索引定位对应分段，再对该分段的残差切比雪夫系数反量化、求值，并与按天体定义的参考模型合成目标位置。

2.1 分段切比雪夫模型

对每个分段，OPM 将时间从儒略日区间 [t0, t1] 线性映射到 切比雪夫标准区间 [-1, 1]：

u=2t−t0t1−t0−1.u = 2\frac{t-t_0}{t_1-t_0}-1.

若直接表示完整位置向量，可写为：

r(u)=[x(u),y(u),z(u)],x(u)=∑iaiTi(u),y(u)=∑ibiTi(u),z(u)=∑iciTi(u).\begin{aligned} \mathbf r(u) &= [x(u), y(u), z(u)],\\ x(u) &= \sum_i a_i T_i(u),\\ y(u) &= \sum_i b_i T_i(u),\\ z(u) &= \sum_i c_i T_i(u). \end{aligned}

其中 T_i 为第 i 阶切比雪夫多项式。切比雪夫表示在有限区间内具有良好的近似稳定性，也方便进行节点采样、误差分析和尾部控制。OPM 在此基础上进一步引入局部坐标和参考形状扣除，使实际存储对象从完整位置系数变为幅度更小的残差系数。

上式描述的是名义运行分段。本文 600 年生产实例在部分天体上会把多项式拟合域扩展为 [t0-Δ_L, t1+Δ_R]，但运行时查询区间仍然是 [t0,t1]。这种分段定义域外延是生成阶段策略，不是文件格式规则；它的主要作用是压低分段边界附近的解析导数尖峰。

2.2 局部坐标与参考形状扣除

Kammeyer (1988) [2] 的做法是在存储系数前，为月球和行星引入按天体配置的坐标轴，使运动尽量接近局部 XY 平面，并对内行星和月球扣除参考轨道。OPM 采用相同原则：在拟合残差切比雪夫系数前，先将 DE441 [1] 采样位置表示到按天体配置的局部坐标中，并扣除参考形状。

在局部坐标中，可以将目标位置写成：

rlocal(u)=rref(u)+δr(u).\mathbf r_{\mathrm{local}}(u)=\mathbf r_{\mathrm{ref}}(u)+\delta\mathbf r(u).

其中 r_ref(u) 为参考形状，δr(u) 为残差向量。OPM 文件主要存储 δr(u) 的残差切比雪夫系数，而不是直接存储完整位置多项式。由于参考形状已经吸收了大部分平滑轨道运动，残差系数的幅度通常显著小于原始坐标系下的完整系数。这会降低量化整数的动态范围，也减少后续位级打包所需的有效位宽。

需要强调的是，OPM 的“参考项 + 残差项”首先是一种统一的文件表示和读取思想，而不是要求所有天体在生成阶段使用同一种参考项构造方法。运行时读取器根据文件中的模型标识取得对应的参考信息和量化残差，将二者合成为当前分段可求值的切比雪夫表示，再按该模型的坐标约定得到全局位置。因此，读取端看到的是统一的模型标识、模型参数和残差系数载荷；它不需要关心这些数值在生成阶段是由哪一种拟合或搜索过程得到的。

生成阶段则可以按天体选择不同的参考构造路线。第一类是直接 XYZ 路线：对于太阳这类本身容易压缩、且不需要额外局部轨道结构的对象，可以直接拟合太阳系质心到目标体的三维切比雪夫系数，再对系数量化和打包。第二类是局部框架共享参考形状路线：对于水星、金星以及多数以质心存储的天体，生成器先把样本转到局部轨道框架，在位置空间中扣除共享参考形状，再拟合和存储残差。第三类是系数参考路线：对于月球这类周期结构很强但普通位置空间参考形状不够稳定的对象，生成器先在固定分段网格上得到每段的切比雪夫系数，再在系数空间中构造参考系数，文件中存储每段系数相对该参考系数的量化差值。

这三类路线的区别主要在参考项如何产生，而不是运行时是否要执行完全不同的系统。无论参考项退化为空参考、来自位置空间共享形状，还是来自系数空间参考系数，最终都被写成 OPM 容器中的模型参数和残差载荷，并由读取器按模型标识解释和求值。

这种做法与 Kammeyer 对内行星和月球使用参考轨道扣除的思想相通：先用稳定的轨道形状吸收主要运动，再存储残差。在 OPM 中，局部坐标和参考形状是用于数值压缩的低频基准项，由生成器根据天体配置和参考星历数据确定。

2.3 天体相关坐标约定

为了兼顾精度和重建便利性，OPM 对不同天体采用不同的自然坐标口径：

  • 太阳：太阳系质心到太阳的向量；
  • 水星 / 金星：以太阳为参考的原生日心向量；
  • 月球：地球到月球的地心向量；
  • 地月质心与外行星：太阳系质心到目标体或质心的向量。

在与地心参考数据比较时，OPM 会按相同规则重建 地球和目标体位置。例如 地球由地月质心和月球重建：

r⊕=rEMB−rMoon1+EMRAT.\mathbf r_{\oplus}=\mathbf r_{\mathrm{EMB}}-\frac{\mathbf r_{\mathrm{Moon}}}{1+\mathrm{EMRAT}}.

其中本文使用 DE441 的

EMRAT=81.300568221497215.\mathrm{EMRAT}=81.300568221497215.

反过来，如果外部星历接口提供地球和月球、但不直接提供 EMB，也可以用同一关系合成外部 EMB 参考：

rEMB=r⊕+rMoon1+EMRAT.\mathbf r_{\mathrm{EMB}}=\mathbf r_\oplus+\frac{\mathbf r_{\mathrm{Moon}}}{1+\mathrm{EMRAT}}.

2.4 量化与二进制封装

OPM 文件不直接存储浮点切比雪夫系数，而是先存储相对于参考项的残差系数，再将残差系数按固定量化单位转换为整数。整数系数使用 ZigZag 编码表示符号，并按坐标轴和阶数统计实际有效位宽后，进行按阶主序位级打包。OPM 容器中这些数据由显式的节（section）承载：覆盖范围和网格、模型标识与参数、量化步长表、位宽表、残差系数载荷、元数据以及带 CRC 的完整性字段相互分离，而不是被视为一个单一的系数数组。

OPM 使用残差存储、模型相关元数据和按阶主序的系数位级打包。这些表示选择独立于本文后续评估的 600 年生产实例；其他时间范围或天体子集也可以使用同一容器语义，只需替换对应的生成参数。

3. 生成与误差修正

OPM 的生成可以分成两步。第一步是从 DE441 采样，并为每个分段建立一个可压缩的初始模型；第二步是在不改变文件结构和位宽的前提下，对少量量化后的整数系数做局部修正，以压低高分位和最坏情况误差。前者决定模型的主体结构，后者主要处理有限阶数和整数量化带来的局部尾部尖峰。

3.1 从 DE441 采样到分段模型

生成器首先给定目标天体、覆盖时间范围和分段策略，然后把时间范围划分为连续分段。每个分段内使用归一化时间变量

u∈[−1,1].u\in[-1,1].

在一个分段 [t0,t1] 内，生成器在切比雪夫节点上调用 DE441，取得目标天体在对应原生坐标口径下的位置样本。这里的“原生坐标”指文件内部实际压缩的向量口径：例如太阳使用太阳系质心到太阳的向量，水星和金星使用以太阳为参考的向量，月球使用地月向量，外行星使用太阳系质心到目标体的向量。

从初始模型构造方式看，本文生产实例可以分为三条路线。第一条是直接裸拟合：太阳的太阳系质心向量直接拟合为三维切比雪夫系数，再进行量化和误差修正。第二条是行星和冥王星的局部参考形状路线：先把原生向量转入局部轨道框架，在位置空间扣除共享参考形状，再拟合残差。第三条是月球路线：月球也先在固定分段网格上对地月向量做裸切比雪夫拟合，但随后会在系数空间中扣除一组参考系数；这可以理解为在系数空间中扣除月球重复轨道形状的等价参考项。

分段长度、残差阶数、参考形状参数和量化参数均按天体配置给出。本文使用的 600 年数据集只是一个生产实例；OPM 格式本身不绑定固定年代、固定分段长度或固定天体集合。相同格式也可以用于其他时间范围，或按天体、按年代拆分为多个可随机访问的文件。

3.2 行星与冥王星的局部框架共享参考形状

除月球外，本文生产实例中的行星、地月质心和冥王星主要使用同一类压缩路线：先把原生向量转到随轨道缓慢变化的局部框架中，再在该框架里扣除共享参考形状，最后拟合和存储残差。这里的“行星”包括水星、金星、火星、木星、土星、天王星、海王星以及冥王星；地月质心（EMB）也使用同类质心向量路线。太阳的太阳系质心向量本身较容易直接拟合，因此不需要在本节展开额外的局部参考形状构造。

局部框架由三维单位法向量和一个平面内角定义：

n=(nx,ny,nz),∥n∥=1,α.\mathbf n=(n_x,n_y,n_z),\qquad \lVert\mathbf n\rVert=1,\qquad \alpha.

其中 n 是局部轨道平面的法向量，α 是该平面内近心方向相对基准轴的旋转角。二者不是动力学轨道元素，而是数值压缩坐标：它们只用于定义一个更容易压缩的位置表达方式。对于本节的行星、EMB 和冥王星路线，分段间的法向量和 α 不作为完全独立的自由表值存储，而是由低阶时间模型表示。生成器在各分段估计局部平面和近心方向后，用一阶时间模型描述这些框架参数在覆盖范围内的缓慢漂移。这里的一阶模型是生成和压缩层面的参数化，不表示物理轨道严格匀速进动。

生成时，先在每个分段内对 DE441 采样点拟合一个最佳平均平面。实现上可以把采样位置矩阵做主平面拟合，取最小方差方向作为平面法向量。由于 n 和 -n 表示同一个几何平面，生成器会在相邻分段之间选择连续的符号。如果发现

ns⋅ns−1<0,\mathbf n_s\cdot\mathbf n_{s-1}<0,

则将当前法向量翻转，并同时把平面内角平移半周：

ns←−ns,αs←αs+π.\mathbf n_s\leftarrow-\mathbf n_s,\qquad \alpha_s\leftarrow\alpha_s+\pi.

这样几何平面不变，但局部平面内的参考方向保持连续。随后再对 α 做相位展开（unwrap），避免跨越 2π 时产生跳变。

运行时，读取模型参数后，先得到当前分段的 (n_x,n_y,n_z) 和 α，再将法向量重新归一化。以 n 作为局部 z 轴，并构造与它正交的局部 x/y 基向量后，位置样本先投影到该局部三轴上。然后按 α 在局部平面内旋转：

x′=cos⁡α x+sin⁡α y,y′=−sin⁡α x+cos⁡α y,z′=z.\begin{aligned} x' &= \cos\alpha\,x+\sin\alpha\,y,\\ y' &=-\sin\alpha\,x+\cos\alpha\,y,\\ z' &=z. \end{aligned}

这个旋转使不同分段中的近心方向具有更一致的相位。随后，生成器在局部坐标中拟合一个全文件共享的参考形状。本文生产路线只对主要平面内的 x' 和 y' 分量建立参考形状：

rref(u)=(Sx(u),Sy(u),0).\mathbf r_{\mathrm{ref}}(u)=\bigl(S_x(u),S_y(u),0\bigr).

从样本中扣除参考形状后得到残差：

δr(ui)=rlocal(ui)−rref(ui)=(x′(ui)−Sx(ui),  y′(ui)−Sy(ui),  z′(ui)).\begin{aligned} \delta\mathbf r(u_i) &=\mathbf r_{\mathrm{local}}(u_i)-\mathbf r_{\mathrm{ref}}(u_i)\\ &=\bigl(x'(u_i)-S_x(u_i),\; y'(u_i)-S_y(u_i),\; z'(u_i)\bigr). \end{aligned}

OPM 后续拟合、量化和打包的是 δr(u)，而不是完整位置函数。由于局部框架和参考形状已经吸收了主要平滑运动，残差系数的动态范围通常显著小于原始坐标中的完整系数。这是 OPM 能够缩小文件体积的主要原因之一。

3.3 直接系数模型：太阳与月球

太阳和月球都不使用上一节的普通位置空间共享参考形状。二者共同点是：生成器首先直接得到当前分段的切比雪夫系数，然后在系数层面进入量化、打包和后续误差修正流程。区别在于太阳的系数本身已经足够容易压缩，而月球还需要在系数空间中进一步扣除一组参考系数。

太阳路线最简单。对于太阳系质心到太阳的原生向量，生成器在每个分段上直接做三维切比雪夫裸拟合，得到该段的 XYZ 系数 C_s。这些系数不再扣除局部共享参考形状，而是直接作为待量化对象进入后续整数化和误差修正步骤。换言之，太阳可以看作参考项退化为空的直接系数模型。

月球路线则是在直接系数模型上再加一层系数参考。月球仍然使用地心月球向量和月球局部框架，但参考项不是在位置空间中的共享 S_x(u), S_y(u) 曲线，而是在系数空间中构造的。这样做的原因是月球周期短、扰动强，逐段位置曲线之间的形状变化比行星更复杂；直接在位置空间中扣除一条共享参考形状，不能像行星路线那样稳定地压低所有段的残差。图 1 给出了本文 600 年生产实例中月球局部框架参数 u 和 v 的逐段最佳估计散点，以及 OPM 实际使用的仿射项加 18.6 年周期三角项拟合曲线。散点本身已经显示出明显周期结构，因此月球框架参数采用周期形式表示，而不是沿用行星路线中的简单平滑进动模型。

图 1. 月球固定网格系数参考模型中的局部框架参数 u 和 v。半透明散点为逐段独立最佳框架参数，实线为 OPM 使用的仿射项加 18.6 年周期三角项拟合；上下两个面板分别给出 u 和 v，横轴为年份。

月球局部框架参数

具体地，月球路线先在固定分段网格上对每段地心月球向量做切比雪夫裸拟合，得到该段在局部月球框架中的三维系数 C_s。然后，生成器在这些逐段系数之间寻找一组参考系数 C_ref，使各段系数相对该参考系数的差值更小、更容易量化和打包。文件中实际存储的是

ΔCs=Cs−Cref,\Delta C_s = C_s-C_{\mathrm{ref}},

其中 s 是分段编号。读取时，先恢复量化后的 ΔC_s，再与 C_ref 相加得到该段实际用于求值的系数：

C^s=Cref+ΔCs^.\widehat C_s = C_{\mathrm{ref}}+\widehat{\Delta C_s}.

因此，太阳和月球都可以理解为“先得到分段切比雪夫系数”的系数路线；太阳直接量化这些系数，月球则先在系数空间扣除参考系数后再量化差值。它们和上一节行星路线的主要区别在于参考项所在的空间：行星路线在位置空间中扣共享参考形状，月球路线在系数空间中扣参考系数，而太阳路线不需要额外参考项。

3.4 残差切比雪夫拟合与量化

得到残差样本后，生成器对每个坐标分量分别拟合有限阶切比雪夫展开：

δra(u)≃∑kca,kTk(u).\delta r_a(u)\simeq \sum_k c_{a,k}T_k(u).

其中 a 表示坐标轴，T_k 为第 k 阶切比雪夫多项式，c_{a,k} 为残差系数。为了写入紧凑二进制文件，浮点系数不会直接保存，而是按预先配置的量化步长转换为整数：

na,k=round ⁣(ca,kqk).n_{a,k}=\mathrm{round}\!\left(\frac{c_{a,k}}{q_k}\right).

读取时再恢复为量化后的浮点系数：

c^a,k=na,kqk.\hat c_{a,k}=n_{a,k}q_k.

量化步长由基础步长 q_base 和按阶数变化的因子 m_k 给出：

qk=qbasemk,xk=kN.q_k=q_{\mathrm{base}}m_k,\qquad x_k=\frac{k}{N}.

其中 N 为残差最高阶。本文生产实例使用三种步长模式：

flat:mk=1,growth:gmk=gxk,linear:amk=1+axk.\begin{array}{ll} \text{flat:} & m_k=1,\\ \text{growth:}g & m_k=g^{x_k},\\ \text{linear:}a & m_k=1+a x_k. \end{array}

常量步长（flat）适合各阶系数范围差异不大的天体；指数增长步长（growth）允许高阶系数使用略粗的量化步长，以换取更小的整数幅度；线性增长步长（linear）则提供介于二者之间的折中。本文生产实例使用的 q_base 和量化模式列于附录 A.2。

OPM 在文件中保存的是量化步长表的浮点值。为了保证写入端验证（writer-side validation）和独立读取器（independent reader）得到完全一致的结果，生成器在写入前会先按文件实际存储精度重建一次 q_k，再用该重建值进行验证。

3.5 位宽统计与按阶主序位级打包

量化后的整数系数包含正负号。若直接使用补码表示，小的负数会带有大量前导 1，不利于按实际位宽压缩。OPM 因此先使用 ZigZag 编码：

zigzag(n)={2n,n≥0,−2n−1,n<0.\mathrm{zigzag}(n)= \begin{cases} 2n, & n\ge 0,\\ -2n-1, & n<0. \end{cases}

这样 0, -1, 1, -2, 2, ... 会映射为 0, 1, 2, 3, 4, ...，绝对值较小的正负整数都变成较小的无符号整数。

随后，生成器按坐标轴和阶数统计实际有效位宽。设 s 为分段编号，则第 a 个坐标轴、第 k 阶系数的位宽为：

wa,k=max⁡sbit_length ⁣(zigzag(ns,a,k)).w_{a,k}=\max_s \mathrm{bit\_length}\!\left(\mathrm{zigzag}(n_{s,a,k})\right).

文件保存这个 axis × degree 位宽表。实际载荷按坐标轴分流，并在每个坐标轴内部采用按阶主序（degree-major）的顺序写入：

for axis a:
  for degree k:
    for segment s:
      write_bits( zigzag(n_{s,a,k}), width = w_{a,k} )

这种顺序和普通的按分段主序浮点数组不同：它先把所有分段的同一阶系数连续写入，再进入下一阶。这样每个阶数都能使用自己的最小安全位宽，高阶小系数不会被低阶大系数迫使使用同样宽的字段。读取时使用同一位宽表和相同顺序反向解包，并进行 ZigZag 解码。

3.6 按天体划分的误差修正路由

初始量化后，OPM 对少量整数残差系数进行局部修正。这个步骤只发生在文件生成阶段；运行时读取 OPM 文件时不需要执行优化。修正过程也不改变文件格式、分段结构或位宽表，而是在已有量化整数附近尝试小幅调整，以降低目标验证网格上的高分位和最坏情况误差。

不同天体对最终地心方向误差的贡献方式不同，因此本文不使用单一误差度量处理所有天体，而是采用按天体类型划分的修正目标：

  1. 太阳. 太阳误差会直接进入地心重建，并且还会作为部分太阳系质心天体组合度量的锚点。因此太阳使用原生千米度量进行修正。
  2. 水星 / 金星 / 月球. 这些天体使用原生角度保护/细化最坏值修正，目标是在其原生几何口径下直接控制角误差尾部。
  3. 地月质心 / 火星 / 外行星. 这些天体的地心误差与地球/太阳锚点相关。因此修正时使用已修正太阳模型的地心组合度量，使太阳系质心向量的优化目标更接近最终地心误差。

这种分流避免了用同一尺度处理所有天体时的目标不匹配，也使生成阶段的优化更接近最终用户关心的地心视方向误差。

3.7 主动/保护网格与细化最坏值目标

误差修正阶段主要控制尾部误差，而不是优化 p50 或 RMS。为避免候选修改只在某一组采样点上变好，生成器把评估点分成两个集合：主动网格和保护网格。主动网格用于提出和排序候选修改；保护网格不参与候选生成，只用于检查候选是否在另一组相位错开的采样点上引入新的尾部突刺。只有同时改善主动网格、且不破坏保护网格表现的候选，才会被接受。

主动网格由三类点组成。第一类是切比雪夫中心节点，用于覆盖多项式误差容易振荡的位置；第二类是分段内的均匀线性采样点，用于补充 切比雪夫节点在时间轴上的非均匀分布；第三类是分段端点和端点邻域点，用于捕捉段边界附近的误差尖峰。保护网格则使用与主动网格错开的相位平移采样，并额外加入端点邻域节点，使候选不会只在主动网格上变好。

在每轮局部搜索中，生成器尝试对少量量化整数系数做 ±1 调整。每个候选先在主动网格上计算误差分布，并按带上限的词典序目标排序：优先降低超过最坏值软上限的峰值，其次降低实际最坏情况误差，再考虑高分位误差。对于主动网格上最接近峰值的若干局部区域，生成器还会做三点峰值细化（three-point peak refinement）：在峰值点左右加入邻近时间点重新估计局部峰顶，以降低真实最大值落在离散采样点之间的风险。

本文生产路线使用以下设置：

  • 主动网格：切比雪夫中心节点、均匀线性采样点和端点组合；
  • 保护网格：相位平移采样，并额外加入端点邻域节点；
  • 峰值细化：对最接近峰值的 3 个局部区域做三点细化；
  • 目标函数：带上限的词典序保护目标；
  • 最坏值软上限：0.00070 角秒；
  • 接受策略：优先降低最坏情况误差。

其中最坏值软上限是生成器侧的软目标和刹车线，不是文件格式参数，也不是严格数学保证。它的作用是引导搜索优先降低最坏情况误差，同时避免为了极小改善而牺牲稳定性或运行时间。最终误差仍以第 4 节描述的独立密集验证结果为准。

4. 验证设计

4.1 参考真值与比较系统

本文使用 DE441 [1] 作为真值来源。OPM 和 Swiss Ephemeris [3] 均在相同儒略日（Julian Date, JD）上计算目标体地心向量，并与 DE441 派生地心向量比较方向角误差。这里需要区分 Swiss Ephemeris 的两个调用接口：本文使用 swe.calc()，即星历时接口，而不是使用世界时接口 swe.calc_ut()。因此比较中不会引入由 ΔT 转换造成的世界时/星历时差异。DE441、OPM 和 Swiss Ephemeris 均使用同一组 JD 数值作为星历自变量；本文不在调用层面额外引入 TT/UT/TDB 转换，因此比较结果反映的是同一数值输入下各系统输出相对于 DE441 的差异。

Swiss Ephemeris [3] 通过 pyswisseph 调用；该 Python 扩展是 Swiss Ephemeris C 库的语言绑定（binding），而非独立重实现。本文使用的 Swiss Ephemeris 库版本为 2.10.03，Python 绑定版本为 20230604；C 库来源为官方 GitHub 仓库 aloistr/swisseph，本地校验提交为 ff04db0（2026-04-28）。Swiss Ephemeris 星历文件使用同一官方仓库随附的星历文件。调用时使用下列标志位，以输出几何 J2000 ICRS 赤道直角坐标：

FLG_SWIEPH | FLG_XYZ | FLG_EQUATORIAL | FLG_J2000 | FLG_TRUEPOS | FLG_ICRS

Swiss Ephemeris 输出的直角坐标以天文单位（AU）为单位，本文使用

1 AU = 149597870.7 km

转换到 km 后再与 DE441/OPM 比较。

4.2 地心重建

对每个天体，DE441、OPM 和 Swiss Ephemeris 都转换到同一地心向量口径：

rgeo(body)=rbary(body)−rbary(⊕).\mathbf r_{\mathrm{geo}}(\mathrm{body})= \mathbf r_{\mathrm{bary}}(\mathrm{body})-\mathbf r_{\mathrm{bary}}(\oplus).

OPM 的内部存储向量不是全部直接地心向量，因此验证时需要按存储口径组合。地球由地月质心（EMB）和月球向量重建：

r⊕=rEMB−rMoon1+EMRAT.\mathbf r_\oplus=\mathbf r_{\mathrm{EMB}}- \frac{\mathbf r_{\mathrm{Moon}}}{1+\mathrm{EMRAT}}.

随后各天体的地心向量按下列规则得到：

rgeo(Sun)=rSun−r⊕,rgeo(Mercury)=rSun+rMercury,helio−r⊕,rgeo(Venus)=rSun+rVenus,helio−r⊕,rgeo(Moon)=rMoon,rgeo(body)=rbody,SSB−r⊕.\begin{aligned} \mathbf r_{\mathrm{geo}}(\mathrm{Sun}) &=\mathbf r_{\mathrm{Sun}}-\mathbf r_\oplus,\\ \mathbf r_{\mathrm{geo}}(\mathrm{Mercury}) &=\mathbf r_{\mathrm{Sun}}+\mathbf r_{\mathrm{Mercury,helio}}-\mathbf r_\oplus,\\ \mathbf r_{\mathrm{geo}}(\mathrm{Venus}) &=\mathbf r_{\mathrm{Sun}}+\mathbf r_{\mathrm{Venus,helio}}-\mathbf r_\oplus,\\ \mathbf r_{\mathrm{geo}}(\mathrm{Moon}) &=\mathbf r_{\mathrm{Moon}},\\ \mathbf r_{\mathrm{geo}}(\mathrm{body}) &=\mathbf r_{\mathrm{body,SSB}}-\mathbf r_\oplus. \end{aligned}

最后一行用于火星、木星、土星、天王星、海王星和冥王星等以太阳系质心向量存储的天体。Swiss Ephemeris 和 DE441 的输出也按同一地心几何口径组合后再比较方向角误差。这样可以避免不同内部存储口径对比较结果造成偏差。

4.3 密集确定性网格

随机采样容易漏掉短时尾部尖峰，尤其当误差峰值集中在分段边界附近或某些局部振荡区域时。因此本文使用 OPM 分段结构构造密集确定性网格：

512 个切比雪夫节点/分段+分段端点.512\ \text{个切比雪夫节点/分段}+\text{分段端点}.

每个天体使用其对应 OPM 文件的分段边界，并限制在 OPM、DE441 和 Swiss Ephemeris 均可计算的覆盖交集内。OPM 和 Swiss Ephemeris 在完全相同的 JD 网格上计算误差。该设计的重点不是模拟真实用户查询分布，而是系统性扫描模型误差尾部。

角误差定义为：

err=atan2 ⁣(∥rtruth×rcandidate∥,rtruth⋅rcandidate).\mathrm{err}=\mathrm{atan2}\!\left( \lVert\mathbf r_{\mathrm{truth}}\times\mathbf r_{\mathrm{candidate}}\rVert, \mathbf r_{\mathrm{truth}}\cdot\mathbf r_{\mathrm{candidate}} \right).

并转换为角秒。

4.4 原生向量位置与速度残差

上述地心角误差是运行时应用最直接关心的口径；它回答的是从地球看天体方向相差多少。为了和 Kammeyer (1988) 的压缩星历图表采用更接近的口径，本文另外计算原生存储向量上的位置和速度残差。对每个 OPM 文件，原生口径按其存储向量定义：太阳为太阳系质心到太阳，水星和金星为日心向量，月球为地球到月球向量，EMB、火星、木星、土星、天王星、海王星和冥王星为太阳系质心到目标体或质心的向量。

与 Swiss Ephemeris 做原生向量对比时，太阳、火星和外行星使用 barycentric 标志，水星和金星使用 heliocentric 标志，月球使用地心月球向量。Swiss Ephemeris 不直接提供 EMB 对象时，EMB 参考由地球 SSB 状态加地月向量按 EMRAT 合成。

原生位置残差定义为：

Δr=∥rOPM,native−rDE441,native∥,\Delta r=\lVert\mathbf r_{\mathrm{OPM,native}}-\mathbf r_{\mathrm{DE441,native}}\rVert,

单位为 km。速度残差使用 OPM 切比雪夫表示的解析导数，而不是有限差分。若分段展开变量为 τ\tau，实际拟合区间为 [a,b][a,b]，则

drdJD=drdτ2b−a.\frac{d\mathbf r}{d\mathrm{JD}}=\frac{d\mathbf r}{d\tau}\frac{2}{b-a}.

对于带参考项和局部框架的天体，速度残差按文件读取器的逐段固定框架语义计算：先对参考项和残差项的切比雪夫多项式求导，再应用与位置重建相同的分段旋转。这个诊断用于评估压缩表示在段内的解析速度连续性和边界尾部行为；它不是对跨分段框架参数时间导数的独立动力学建模。DE441 真值速度由 SPK 段的 compute_and_differentiate() 给出，单位为 km/day。报告时速度残差同时换算为 mm/s，并保留最大值对应的 AU/day 形式，以便和 Kammeyer 文中速度误差单位对应。

5. 结果

5.1 与 Swiss Ephemeris 的密集比较

表 1 给出以 600 年生产实例覆盖范围为基础的 512 节点密集比较（约 JD 2378460—2597711；精确端点随天体覆盖交集略有不同）。表内误差单位为毫角秒（mas；1 mas = 0.001 角秒）。Swiss Ephemeris 在这里作为成熟、广泛使用的紧凑星历基线；下列比较只对应本文定义的几何地心口径、时间范围和密集 JD 网格，不代表对两套系统全部功能或全部时间范围的综合评价。

表 1. OPM 和 Swiss Ephemeris 在确定性密集网格上相对于 DE441 的地心角误差。

天体样本数Swiss p50Swiss p95Swiss p99Swiss maxOPM p50OPM p95OPM p99OPM maxmax 比值
太阳625,3480.453031.154731.57343.084710.07255090.2686310.3945180.6122020.198463
月球4,082,0400.1830820.3815540.4728060.9972960.0882570.1611440.1945630.3540620.355022
水星1,279,3750.4811561.299111.81314.016050.0964650.305520.4790461.109240.276202
金星502,1780.4623572.082193.9473910.90450.1202910.6336071.188613.089290.283304
火星164,4800.3865271.668672.708575.590260.1170220.5441120.9159972.166010.387461
木星37,9630.2448760.5805330.7618041.327450.169210.3471070.4284320.628950.473803
土星32,8330.2463080.5700780.8045691.330820.1747830.3459930.4267830.5934560.445931
天王星28,7290.1815870.3998190.521340.7443380.1300980.247630.2901370.4528780.60843
海王星28,7290.1759970.3693080.4425450.564350.135760.2459820.2781740.3480890.616796
冥王星28,7290.179380.3986750.5793520.9021870.1031390.2012140.2552120.329770.365523

密集测试点总数为 6,810,404。

图 1 给出同一密集网格上的 OPM 与 Swiss Ephemeris 地心角误差曲线。为避免主文过长，这里只放三个代表对象：水星代表内行星短周期高偏心率情形，月球代表高分段数和强摄动情形，海王星代表外行星长周期情形。完整 10 个天体的 SVG 曲线文件列于附录 D，并由 out/opm2-600y-selected-polished-swiss-geocentric-512.txt 对应的同一次密集验证生成。

图 1. 代表天体在 600 年密集网格上的地心角误差曲线。绿色为 OPM，橙色为 Swiss Ephemeris；纵轴单位为角秒。

（a）水星

水星密集误差曲线

（b）月球

月球密集误差曲线

（c）海王星

海王星密集误差曲线

5.2 原生向量位置与速度残差

表 2 给出与图 1 同一 600 年生产实例对应的原生向量残差诊断。该表不与 Swiss Ephemeris 作比较，也不替代前述地心角误差结果；它用于从压缩模型本身的坐标口径观察 OPM 相对于 DE441 的位置和速度表示误差。速度列为解析导数误差，因此也检验了量化切比雪夫系数在导数层面的连续性和尾部表现。

表 2. OPM 原生存储向量相对于 DE441 的位置和速度残差。位置单位为 km，速度单位为 mm/s；最后一列给出速度最大值的 AU/day 换算。

天体样本数位置 p99位置 max速度 p99速度 max速度 max (AU/day)
太阳626,5660.03686920.05931640.3479940.5793463.35e-10
EMB309,4280.3422820.5745720.8619841.489428.60e-10
月球4,089,8980.0003997910.0006538960.02334320.05207533.01e-11
水星1,281,4020.1168960.1735731.390144.513472.61e-09
金星502,1780.209270.3542751.403042.574931.49e-09
火星164,4800.6176661.103810.7935721.367647.90e-10
木星38,0361.774592.513660.6545381.03685.99e-10
土星32,8963.345074.56531.27461.926811.11e-09
天王星28,7844.524326.539121.033221.657579.57e-10
海王星28,7846.9871510.65282.028443.60562.08e-09
冥王星28,7847.1829210.30891.691462.605611.50e-09

图 2 给出两条代表曲线：土星位置残差用于和 Kammeyer 原文中以 km 表示的 Saturn 残差图作直观类比；水星速度残差则展示内行星短周期、高曲率情况下的导数误差。完整原生位置和速度残差曲线列于附录 E，并由 out/opm2-600y-selected-polished-native-residuals-512.txt 对应的同一次诊断生成。

图 2. 原生向量残差的代表曲线。上：土星原生位置残差，单位 km；下：水星原生速度残差，单位 mm/s。

（a）土星原生位置残差

土星原生位置残差

（b）水星原生速度残差

水星原生速度残差

5.3 分段外延对解析速度残差的影响

表 3 给出 OPM 配置中对解析速度行为最相关天体的无外延消融结果。外延不是 OPM 文件格式的固定规则，而是生成阶段的按天体策略；其主要收益不是单调降低位置残差，而是把运行查询区间的段端点移入拟合区间内部，从而显著压低解析导数在段边界处的尖峰。下表比较的是同一套 selected 600 年 raw 配置在启用外延时的结果，以及将 segment_domain_expansion_fraction = 0.0 后重新生成的 raw 配置；两者均使用相同的 512 节点原生残差诊断。Swiss Ephemeris 速度最大值作为同一原生向量口径下不变的外部参考一并列出；EMB 由 Swiss 地球 SSB 状态加地月向量按 DE441 质量比合成。

表 3. 分段定义域外延对解析速度最大残差的影响。速度单位为 mm/s。

天体selected 速度 max无外延速度 max降低幅度Swiss 速度 max
太阳0.5266872.30981177.198%28.673538
EMB1.2685228.72824085.466%20.603949
火星1.2108308.52717585.800%19.486554
木星0.8474873.95638078.579%5.271369
天王星1.3032325.18409074.861%4.720532
海王星2.9167345.51783147.140%5.464040
冥王星2.2499763.79215940.668%7.193522

5.4 按统计指标的对比概览

表 4 汇总了在本文验证网格和角误差定义下，两个系统分别在哪些统计量上给出较低误差。这里的“较低”仅指本文定义的比较口径，不表示对系统整体功能的排序。

表 4. 在本文验证集上，各统计量中误差较低的天体数量。

统计量OPM 较低Swiss Ephemeris 较低说明
p50100OPM 在全部列出天体上较低
p95100OPM 在全部列出天体上较低
p99100OPM 在全部列出天体上较低
max100OPM 在全部列出天体上较低

5.5 最坏情况误差相对比值

表 5 给出 Swiss Ephemeris 与 OPM 的最大误差比值，用于量化两者在误差尾部上的差异。比值大于 1 表示在本文验证口径下，Swiss Ephemeris 的最大误差高于 OPM。

表 5. 相对于 OPM 的最坏情况角误差比值。

天体Swiss Ephemeris 最大值 / OPM 最大值
太阳5.04×
月球2.82×
水星3.62×
金星3.53×
火星2.58×
木星2.11×
土星2.24×
天王星1.64×
海王星1.62×
冥王星2.74×

从分位数可以看出，两者的差异主要出现在误差分布尾部。更新后的 600 年生产实例在全部列出天体的 p50、p95、p99 和最大误差上均低于 Swiss Ephemeris；其中海王星和天王星仍是最接近的外行星个例，但其最坏情况误差仍分别约为 Swiss Ephemeris 的 62% 和 61%。

5.6 文件体积

表 6 给出本文 600 年生产实例的实际 OPM 文件大小。该大小为未经过 gzip 等外部通用压缩的文件大小，包含文件头、分段索引、参考形状、模型表、量化残差位流和校验信息。每世纪折算值仅为将 600 年文件大小除以 6 后得到的平均量级；实际生成更长时间范围时，分段配置、边界元数据和误差修正策略可能随时间范围调整。

表 6. 600 年 OPM 生产文件大小及每世纪折算值。

天体600 年大小（KiB）每世纪折算（KiB/世纪）
太阳99.7316.62
水星182.3030.38
金星86.7214.45
EMB75.5112.59
月球1158.91193.15
火星39.826.64
木星9.921.65
土星8.161.36
天王星4.750.79
海王星3.670.61
冥王星4.100.68
合计1673.60 KiB = 1.634 MiB278.93

月球是当前体积的主要来源，约占 600 年文件总量的 69%；水星、太阳、金星和地月质心次之。外行星由于运动较慢、分段数少，单体文件在 600 年范围内只有数 KiB 到十余 KiB。作为同一时间分片的文件体积参照，本地 DE441 版本 Swiss Ephemeris 中相对应的核心文件为 sepl_18.se1（太阳至冥王星）和 semo_18.se1（月球），二者合计 1,788,832 字节，即约 1.71 MiB；这里不计入小行星文件 seas_18.se1。若按表 6 的总每世纪体积线性外推到完整 DE441 时间跨度 JD -3100015—8000016（约 303.9 个儒略世纪），全套主要天体数据约为 82.9 MiB。该外推只用于量级比较，未经完整 DE441 范围生产验证。

6. 候选压缩路线与方法取舍

在收敛到本文采用的 OPM 表示之前，本文尝试过多种看似自然的紧凑星历构造方式。这些路线并非完全无效；许多方案在短时间范围、特定天体或特定运行时目标下可以得到可用结果。但从长期稳定性、文件体积、运行时成本、实现复杂度和最坏情况误差共同考虑，它们都不如本文采用的“参考形状扣除、局部框架、残差切比雪夫和位级打包”路线稳定。本节记录这些取舍，以说明 OPM 的设计目标不是简单地对完整位置向量作切比雪夫压缩。

6.1 VSOP87 与截断解析公式路线

最直接的思路是使用已有解析星历作为基础模型，再压缩它与 DE441/JPL 数值星历之间的差值。Swiss Ephemeris 文档 [3] 第 2.1.4 节 “The Swiss Ephemeris compression mechanism” 在说明其压缩机制时写道：

Instead of the positions we store the differences between JPL and the mean orbits of the analytical theory VSOP87. These differences are a lot smaller than the position values, wherefore they require less storage. They are stored in Chebyshew polynomials covering a period of an anomalistic cycle each.

这段描述提供了重要线索：压缩对象不是完整位置，而是相对于某种平均轨道的残差；分段也与近心周期有关。但若将其中 “JPL 与 VSOP87 平均轨道的差值” 按可执行实现理解为“直接计算完整 VSOP87A/VSOP87B 位置，再存储 JPL 减 VSOP87 的残差”，会高估 VSOP87 完整解析级数在压缩结构中的作用，并遇到明显的工程取舍问题。

在早期试验中，VSOP87A 加 DE441/JPL 残差可以把误差压到约 0.001" 量级，VSOP87B 加残差甚至可到约 0.0001" 量级。但完整 VSOP87B 项数较多，在客户端或嵌入式运行时中的计算成本偏高，并不适合作为小型运行时星历的主体。若截断 VSOP87 以换取速度，近现代范围内可以通过修正项得到可用结果，但远离现代历元后残差会明显恶化。内部试验中，截断 VSOP87B 加修正在约 ±4000 年范围内可达到约 0.8" 峰值，但再向外会出现发散；寿星万年历等软件中的截断公式在历法应用中足够实用，但在 ±3000 年以外也会退化到数角秒或更高。

进一步的变体包括合并接近频率和接近相位的 VSOP 项、按世纪单独调整截断项，或构造“截断 VSOP + 世纪修正 + 世纪内小段修正”的三层模型。这些方法可以在局部时间范围内降低误差，有时可达到约 0.05" 量级，但长期相位漂移仍会在残差中形成局部尖峰；同时，模型结构也变得不够统一。这些试验说明，解析级数路线的问题并不是“完全不准”，而是精度、速度和残差形态难以同时满足。高精度解析级数太重，强截断又留下长期相位误差；在其上叠加经验修正会增加格式复杂度，却不一定得到低熵、低位宽、易打包的残差。

6.2 Steve Moshier 的 PLAN404 三角级数参考模型

在 VSOP87 截断路线暴露出长期发散问题之后，另一个尝试是使用 Steve Moshier 发布的 PLAN404 数据包 [4] 作为参考模型。该数据包包含九大行星运动的三角级数展开，并按 JPL 的 DE404 Long 星历（约公元前 3000 年至公元 3000 年）作最佳拟合，输出日心黄道坐标。Moshier 原说明中给出的精度约从地球的 0.1" 到冥王星的 1"。在本文的方法探索阶段，采用 PLAN404 的目的不是把它作为最终星历，而是利用其长期相位行为较平滑的特点，缓解 VSOP87 截断残差在远离现代历元后不够平滑的问题。

Moshier/PLAN404 本身并不是最终高精度模型。内部试验中，原始 Moshier 相对于 DE441 的某些经纬残差在 ±5000 年范围内可达到数十角秒，例如地球相关试验中原始峰值约为 75"。但其误差随时间变化较平滑，经低阶分段多项式修正后可以在 ±5000 年范围内压到约 2"；同一试验还显示，在 ±7000 年附近仍未出现 VSOP87 截断修正那样的发散。其主要残差瓶颈表现为约 40 年周期的准周期项，远离现代历元时振幅增大，但这类结构比截断 VSOP 的长期相位漂移更容易用低阶分段模型处理。

因此，Moshier/PLAN404 类模型的价值不在于单独达到亚毫角秒精度，而在于它提供了一个更平滑、更不易远期发散的半解析参考。它说明参考模型不必本身就是最终星历；只要能吸收足够多的长期相位结构，使残差变得平滑且可预测，就有压缩价值。

不过，这一路线没有解决文件体积问题。以水星为例，一组经过调参的 PLAN404 加球面残差方案约为 116 KB/世纪，虽然比当时约 173 KB/世纪 的方案小约三分之一，且运行时成本最低，但仍明显大于后续类 Kammeyer 局部框架和参考形状路线希望达到的体积。换言之，Moshier/PLAN404 解决了“远期残差更平滑”的问题，却没有给出足够紧凑、统一的最终格式。

6.3 直接数值拟合与慢变天体例外

另一类路线是完全绕开解析星历，直接对 DE441/JPL 位置表作数值拟合。最简单的方案是把每个天体在固定时间段内的三维位置直接表示为切比雪夫多项式。这一路线实现清楚、验证方便，但对内行星和月球的压缩率不理想：如果不先扣除主要轨道结构，切比雪夫系数必须同时表示完整椭圆运动、平面变化和摄动项，低阶系数幅度大，高阶系数也难以压到很低位宽。

也曾尝试在拟合中加入速度约束，或引入空间特征、轨道特征辅助时间拟合。这些方法可以让拟合问题看起来更充分，但并没有从根本上降低待压缩函数的信息量。位置和速度共同拟合会增加目标约束；额外特征会提高运行时和生成器复杂度；若主轨道运动仍在残差中，最终系数仍然不够低熵。因此，直接数值拟合路线的教训是：切比雪夫多项式本身不是压缩的关键，关键是切比雪夫要拟合已经扣除主结构后的残差。

直接拟合并非在所有对象上都失败。对于冥王星或部分慢变小天体，长周期、低曲率的运动可以被直接切比雪夫拟合或固定主成分框架表示得很紧凑。早期冥王星直接拟合格式在随机测试中，第 99 百分位误差约为 0.00047"，文件体积也很小。这说明压缩策略应当按天体运动结构选择：慢变对象可以直接拟合，内行星和月球则需要近心相位、局部框架和参考形状扣除。

6.4 Kepler 椭圆、轨道根数与尖峰补丁

还有一类方案试图从低维轨道模型出发，例如按年或按世纪修正 Kepler 椭圆，或拟合一组随时间变化的轨道根数，再在局部坐标中压缩残差。这一路线看起来物理含义更强，但并不天然适合作为压缩坐标。

水星是最典型的反例。试验中曾尝试在一个世纪范围内拟合半长轴、偏心率、倾角、升交点黄经、近日点方向和平均近点角等根数，再通过 Kepler 方程重建参考椭圆并存储径向-切向-法向残差。但水星每世纪约有 415 个近心周期，摄动导致的瞬时轨道根数含有明显短周期振荡；世纪级低阶拟合会产生很大的结构误差。记录中的半长轴拟合误差可达约 431,000 km，平均近点角误差约 0.17 rad（约 10°），偏心率误差约 0.02。同时，该路线运行时需要多组切比雪夫求值、Kepler 方程迭代、大量三角函数和坐标旋转，三角函数成本估计可达 PLAN404 路线的 5--10 倍。

针对残差中的局部毛刺，也可以设计特殊补丁。例如在少数时间段增加局部修正项，或对特定分段采用不同的拟合策略。这类方法短期有效，但容易把格式变成一组例外规则。更重要的是，毛刺往往不是孤立异常点，而是参考相位、局部框架或参考形状与真实运动之间存在系统性错配。若不先修正相位和框架，补丁只能治标，不能降低残差整体复杂度。

6.5 月球与水星的特殊困难

月球和水星是压缩路线中最能暴露问题的对象。水星的问题主要来自短近心周期、高偏心率和长期进动；月球的问题则来自更短周期的地月运动、强摄动和复杂的节点/近地点变化。早期月球格式曾使用近地点分段、节点/轨道框架、混合位宽量化等策略，并在随机测试中达到约 0.0011" 的第 99 百分位角误差，但对应 30,000 年级别数据体积仍由月球主导。

围绕月球还尝试过多种专门路线：物理轨道框架、主成分分析参考形状、近地点/远地点对齐、固定旋转、等距轨道根数参数、低位宽尾部系数、容差扫描、阶数扫描和采样预算扫描等。其中一类尝试使用 IERS 公式构造基础月球框架，再只存储小的低位宽角度残差，目标是替代每段显式浮点框架元数据。这些试验说明，月球的主要困难不只是残差阶数，而是相位、框架和元数据三者如何共同降低低阶系数和高阶尾部。

水星的早期试验也给出类似结论：真正有效的参考不应是直接拟合瞬时轨道根数，而应是围绕平均近心相位和局部平面构造的压缩参考。本文的水星路线使用全局近心时钟修正、局部框架、参考形状和保护/细化误差修正，正是对这些试验结论的整理。月球路线则采用固定分段网格、局部框架和系数空间参考项：它与行星路线共享“参考项 + 残差项”的读取思想，但参考项是在裸拟合后的系数空间中构造的。

6.6 类 Kammeyer 的早期 OPM2/OPV2 路线

在上述路线之后，早期 OPM2/OPV2 是第一条在文件体积上接近 Swiss Ephemeris 核心星历文件量级的统一路线。该路线已经采用类 Kammeyer 的基本结构：每段局部框架、固定量化单位、混合位宽打包、水星参考形状扣除，以及面向冥王星或慢变对象的固定主成分框架。在部分早期测试中，该路线的第 99 百分位角误差约为 0.0008--0.0015"。Kammeyer 原作中还包含一种段边界控制思想：对水星、金星、地月质心、火星、外行星和太阳，切比雪夫展开区间在分段起止时间两侧各外延段长的 5%；月球的展开区间则与分段本身相同。这一点对后续 OPM 试验有启发，但最终并未作为固定规则照搬。

同一阶段也尝试过变长分段和自适应分段。直观上，误差大的时间区间可以切得更短，误差平滑的区间可以切得更长；这种策略在单个天体、单个误差阈值下往往有效，也适合离线压缩器用二分或动态规划搜索局部最优段长。但它没有成为本文 OPM 表示的主线，原因主要有三点。第一，变长分段需要在文件中显式保存每段起止时刻或持续时间，段数较多的内行星和月球会付出可见的索引/头部开销；若再要求快速随机访问，还需要额外的段查找结构。第二，不等长分段会削弱全局近心时钟、共享参考形状和按阶主序位宽打包的规则性，使生成器和运行时代码都更难保持简单。第三，局部缩短高误差分段虽然能降低某些峰值，但也可能把问题转移到相邻边界，最终仍需要密集验证和保护网格修正。因此，本文采用的是按天体选择的规则分段或相位分段，再通过参考形状、量化策略和最坏值修正控制尾部；变长分段保留为未来针对小天体或特殊时间区间的可选优化，而不是默认格式假设。

该阶段也留下了若干重要工程教训。第一，局部框架角度不能简单用单精度浮点（float32）保存；早期试验中，单精度框架角曾导致约 0.017" 的重建误差。第二，逐段自适应比例尺虽然看起来可以降低单段溢出风险，却会削弱参考形状扣除带来的共享收益；水星试验中逐段比例尺约为 87 KB/世纪，而固定量化单位加参考扣除约为 55 KB/世纪。第三，Kammeyer 原始量级的 10^-8 轨道半径量化对水星来说过粗，约 0.698 km 的单位会产生约 0.008" 误差；要达到约 0.001"，需要更细的量化单位。

因此，早期 OPM2/OPV2 已经证明了局部框架、参考扣除和位宽打包的方向是正确的；本文采用的 OPM 表示则进一步把该方向整理为更明确的生产格式：使用单位法向量加平面内角作为默认框架，采用按阶主序位级打包，按天体路由误差修正目标，并通过主动/保护网格与三点峰值细化来控制尾部误差。

6.7 局部框架编码的后续可能

在确定了残差压缩的主体结构以后，局部框架本身仍有进一步编码空间。本文主线使用单位法向量 n 加平面内角 α，因为它直接表示分段最优平面方向，不依赖单一投影图，并且对高倾角或快速进动的轨道族更稳健。

一个更紧凑的候选是用两个参数 p,q 表示分段最优平面。例如当局部轨道平面不接近竖直方向时，可以写成

z=px+qy.z=px+qy.

或等价地用未归一化法向量

n~=(−p,−q,1),n=n~∥n~∥.\tilde{\mathbf n}=(-p,-q,1),\qquad \mathbf n=\frac{\tilde{\mathbf n}}{\lVert\tilde{\mathbf n}\rVert}.

这样平面部分只需要两个参数，和单位法向量的真实自由度一致；对于长期停留在参考平面附近的低倾角主要行星，p,q 也可能更省。但它本质上仍是单图参数化：当轨道平面接近该图的奇异方向时，p,q 会变大、不连续或发散。类似的球极投影、多图投影加图标志以及四元数框架，都可以作为后续格式变体继续探索。

因此，本文采用的单位法向量方案是更保守、更通用的默认选择；p,q 等方案则更适合作为特定低倾角数据集的专用压缩优化。本文的 600 年生产结果只使用单位法向量加平面内角，不依赖这些候选编码。

7. 讨论

7.1 为什么第 99 百分位和最大误差重要

如果只看随机采样或中位误差，一个模型可能看起来非常稳定，但仍在少数时间区域产生较大误差。对运行时星历而言，这类尾部误差很重要，因为用户查询时间并不受采样分布约束。一个可部署星历模型应当尽量避免在覆盖范围内出现局部尖峰。

OPM 的误差修正目标正是控制这种尾部误差表现。本文 600 年密集验证显示，该目标在实际基线比较中有效：在本文定义的几何地心验证口径下，OPM 在列出的 10 个天体的中位数、高分位数和最大误差上均低于 Swiss Ephemeris。这说明 OPM 的优势不是偶然的单点最坏值调参，而是整体误差分布和尾部误差的改善。

7.2 与 Swiss Ephemeris 的比较

作为成熟、紧凑且广泛使用的系统，Swiss Ephemeris 是一个现实的强基线。本文的比较不应理解为否定其工程价值；相反，它代表了实际软件中会采用的压缩星历，而不是原始大型 DE 内核。

在这一基线下，OPM 的结果可以概括为：

  • 中位误差：OPM 在 10/10 个天体上更小；
  • 高分位误差：OPM 在第 95 百分位和第 99 百分位上均为 10/10 占优；
  • 最坏情况误差：OPM 在 10/10 天体上占优。

这表明 OPM 在本文的 DE441 几何地心压缩表示验证口径下具有更低的误差分布和尾部误差；Swiss Ephemeris 的优势则体现在成熟生态、长期使用、接口覆盖和完整软件系统能力。

7.3 密集确定性验证的作用

随机 10k JD 测试可以快速给出大致精度，但它可能漏掉局部尖峰。本文使用 基于分段的密集网格 的原因是 OPM 本身为分段模型，误差结构很可能与分段内部位置、边界和 切比雪夫振荡有关。通过在每段内使用 512 个切比雪夫节点加端点，可以更系统地扫过潜在尾部区域。

这种验证方式对 Swiss Ephemeris 也是公平的，因为 Swiss Ephemeris 和 OPM 使用完全相同的 JD 网格。换言之，本文比较的是同一组时刻下两个系统相对于 DE441 的误差，而不是各自选择有利采样点。

7.4 与 Kammeyer (1988) 的关系

Kammeyer (1988) [2] 的 Compressed Planetary and Lunar Ephemerides 为本文提供了重要历史参考：从参考数值星历派生分段切比雪夫表示，通过移除共同轨道结构降低系数幅度，再将量化后的整数字系数紧凑打包。原作以 DE200 为参考，覆盖约 1801—2049 年，数据文件约 830 KB，报告的位置精度约为 0.001 角秒 量级；其中残差图采用的是以 km 表示的位置残差幅值，文中还以 AU/day 报告速度误差。若只按覆盖时间粗略折算，其体积约为 330 KB/世纪；本文 600 年生产实例为 1.634 MiB，约 278.93 KiB/世纪，在更长覆盖、更现代参考星历和显式校验封装下仍处于同一量级并略小。这个比较只是量级参照，因为两者覆盖时间、参考星历、天体集合、误差口径和文件封装并不完全相同。

OPM 采用同一基本路线，并将其适配到 DE441 和现代运行时部署：参考形状与局部坐标由本文数据和天体配置确定，残差系数以按阶主序位流存储，文件头部显式记录 位宽表、量化参数和 CRC 校验。

本文相对于 Kammeyer 论文中描述的方法，主要扩展不在于改变切比雪夫压缩表示的基本思想，而在于生成和验证目标。Kammeyer 原作对水星、金星、地月质心、火星、外行星和太阳采用两端各 5% 的展开区间：实际拟合的切比雪夫区间略早于分段起点、略晚于分段终点，以减轻在段边界附近最容易出现的振荡和截断尾部；月球的展开区间则与分段本身相同。OPM 保留了“分段区间外延可以改善边界尾部”的思想，但把外延比例作为按天体验证的生成策略，而不是格式规则，也没有继承统一的 5% 取值。

最终生产实例显示，外延对位置精度并非单调收益：固定阶数的切比雪夫多项式需要覆盖更宽时间区间，过大的外延可能降低分段内部的位置近似效率。但外延对解析速度残差的收益非常明显，因为它把实际查询区间的端点移到拟合区间内部，显著降低由段边界端点效应放大的导数尖峰。本文因此将原生位置残差、地心角误差和解析速度残差同时纳入密集验证，并在地月质心、火星、木星、天王星、海王星和冥王星上采用 0.5%--1.0% 的小比例外延。表 3 显示，这一生成策略对不同天体收益不同：太阳、地月质心、火星、木星和天王星的解析速度最大残差降低约 75%—86%，海王星和冥王星也分别降低约 47% 和 41%。

OPM 在初始残差拟合 后加入按天体类型划分的 保护/细化误差修正，使优化目标更贴近地心角误差尾部；验证则使用 600 年范围内的确定性密集网格，而不是仅依赖稀疏残差图或随机采样。这样可以更清楚地评估部署型星历在覆盖范围内的高分位和最坏情况误差。

7.5 局限性

本文仍有若干限制：

  1. 本文比较集中于 600 年范围，尚未覆盖 DE441 的完整时间跨度；
  2. 本文比较为几何地心 J2000 ICRS，不包含光行时、视差、章动、岁差、引力偏折或大气折射等视位置与站心计算流水线；
  3. OPM 文件已经使用残差存储与位级打包，但文件体积仍可通过更强的参考形状模型、自适应阶数、更细粒度的位宽分配和共享低维轨道参数继续优化；
  4. 本文生产实现采用三维单位法向量加平面内角作为默认局部框架，避免了旧版单图球极投影编码的奇点问题；但局部框架仍是数值压缩策略，不是动力学轨道元素。对于人造卫星、高倾角轨道或轨道平面快速进动的对象，未来仍可继续引入带图选择标志的多图投影、四元数框架或其他局部框架策略；
  5. 海王星是本文数据集中 OPM 与 Swiss Ephemeris 误差分布最接近的情形之一，其中中位数和第 95 百分位仍有进一步优化空间；
  6. 最坏值软上限是生成器策略参数，不是数学证明意义上的全域误差界。

8. 结论

本文提出 OPM，一种面向运行时部署的紧凑主要天体星历表示。OPM 使用分段切比雪夫多项式、按天体选择的局部参考结构、残差系数量化和分节二进制封装，并通过按天体类型划分的保护/细化误差修正控制高分位和最坏情况误差。

在以 DE441 为真值的 600 年密集确定性验证中，OPM 在 6,810,404 个地心测试点上表现出稳定的误差尾部控制。在本文定义的几何地心 J2000 ICRS 口径下，OPM 对 10 个可与 Swiss Ephemeris 直接比较的主要天体，在中位数、第 95 百分位、第 99 百分位和最坏情况角误差上均给出更低数值；最大的 OPM 最坏值约为 Swiss Ephemeris 对应最坏值的 62%。该生产实例的 11 个 OPM 文件合计 1.634 MiB，约 278.93 KiB/世纪。

这些结果表明，OPM 是一种以稳定尾部误差和紧凑随机访问为目标的星历表示，尤其适合需要小体积、快速读取、跨平台校验和稳定最坏情况误差的应用。本文结果并不替代 Swiss Ephemeris 等成熟系统的完整功能，而是说明在 DE441 几何位置压缩这一明确口径下，OPM 可以用接近 Swiss 核心星历文件量级的体积提供更强的密集误差尾部表现。后续工作将集中在进一步压缩文件体积、扩展时间覆盖、补充 C/C++ 读取器的冷启动、随机访问和批量重建性能测试，以及将 OPM 接入视位置与站心计算流水线。

附录 A. 生产配置表（Production configuration）

本附录记录本文 600 年生产实例使用的主要生成参数。需要强调的是，这些参数属于本文数据集的生成配置，而不是 OPM 文件格式本身的限制；OPM 可按其他时间范围、分段长度、阶数和量化参数重新生成。

A.1 天体配置概览（Body configuration overview）

表 A.1 列出本文 600 年生产实例使用的主要天体配置。主文只讨论生成原则；具体分段长度、残差阶数、参考项配置和误差修正目标集中放在这里，避免方法部分被实现参数打断。表中的分段长度为本文生产实例使用的固定长度或固定全局周期；这些是数据集配置，不是 OPM 格式限制。除太阳的原始 XYZ 路线外，本文生产实例主要使用由局部轨道几何确定的压缩坐标；行星、EMB 和冥王星使用单位法向量加平面内角，月球使用单独的周期局部框架参数模型。

天体原生坐标口径分段长度 / 时钟残差阶数参考项配置误差修正目标
太阳太阳系质心（SSB）到太阳固定 180 d25无参考形状；原始 XYZ 切比雪夫表示原生千米误差
水星以太阳为中心的原生向量固定全局近点周期；P = 87.969349505206 d；8 阶切比雪夫时钟修正24平均拱点局部坐标；参考形状 40 阶原生角度保护/细化最坏值修正
金星以太阳为中心的原生向量固定全局近点周期；P = 224.700615924424 d24平均拱点局部坐标；参考形状 40 阶原生角度保护/细化最坏值修正
月球地球到月球向量固定全局近地点周期；P = 27.554538221087 d；世纪 int16 时钟表24固定分段网格；周期局部框架；系数空间参考项原生角度保护/细化最坏值修正
地月质心（EMB）太阳系质心到地月质心固定全局惯性相位；P = 365.256362982910 d28固定局部坐标；参考形状 22 阶以已修正太阳为锚点的地心组合度量
火星太阳系质心到火星固定全局近点相位；P = 686.996026060798 d28固定局部坐标；参考形状 22 阶以已修正太阳为锚点的地心组合度量
木星太阳系质心到木星固定 3000 d24固定局部坐标；参考形状 16 阶以已修正太阳为锚点的地心组合度量
土星太阳系质心到土星固定 3500 d24固定局部坐标；参考形状 16 阶以已修正太阳为锚点的地心组合度量
天王星太阳系质心到天王星固定 8000 d30固定局部坐标；参考形状 12 阶以已修正太阳为锚点的地心组合度量
海王星太阳系质心到海王星固定 10000 d30固定局部坐标；参考形状 12 阶以已修正太阳为锚点的地心组合度量
冥王星太阳系质心到冥王星固定 10000 d30固定局部坐标；参考形状 12 阶以已修正太阳为锚点的地心组合度量

A.2 量化与打包配置（Quantization and packing configuration）

表 A.2 列出本文生产实例的残差系数量化配置。base_km 为基础量化步长，单位为 km；flat 表示各阶使用相同步长，growth:x 或 linear:x 表示量化步长随阶数按对应模式变化。所有天体均采用 ZigZag 编码和按阶主序位级打包；位宽表由量化后的整数系数统计得到，并写入文件头部或相关元数据区。

天体基础量化步长 base_km量化模式备注
太阳0.01各阶相同（flat）太阳系质心到太阳的原始 XYZ 残差系数
水星0.032线性变化（linear:0.65）修正后的全局近点时钟；分段定义域扩展 0.01
金星0.06各阶相同（flat）全局近点时钟；分段定义域扩展 0.01
月球0.00025各阶相同（flat）世纪 int16 时钟表；分段定义域扩展 0.01
地月质心（EMB）0.02随阶数增长（growth:1.25）固定局部坐标参考形状；分段定义域扩展 0.01
火星0.04各阶相同（flat）固定局部坐标参考形状；分段定义域扩展 0.01
木星0.5随阶数增长（growth:1.25）固定局部坐标参考形状；分段定义域扩展 0.0075
土星1.0随阶数增长（growth:1.25）分段定义域扩展 0.01
天王星1.6线性变化（linear:0.5）固定局部坐标参考形状；分段定义域扩展 0.01
海王星3.5各阶相同（flat）固定局部坐标参考形状；分段定义域扩展 0.005
冥王星3.5随阶数增长（growth:1.25）固定局部坐标参考形状；分段定义域扩展 0.0075

A.3 误差修正配置

本文 600 年生产实例的误差修正采用以下统一保护设置：

主动网格：切比雪夫中心节点、均匀采样点和端点
保护网格：相位平移采样，并加入端点邻域节点
峰值细化区域数：3
目标函数：带上限的词典序保护目标
最坏值软上限：0.00070 角秒
接受策略：优先降低最坏情况误差

其中最坏值软上限是生成器侧软目标，不是文件格式参数，也不是严格数学误差界。

附录 B. 坐标重建细节（Coordinate reconstruction details）

B.1 DE441 质心坐标重建

DE441 中部分天体以质心和相对向量形式存储。例如地球和月球可由地月质心加相对向量重建：

rbary(⊕)=rbary(EMB)+rEMB→⊕,rbary(Moon)=rbary(EMB)+rEMB→Moon.\begin{aligned} \mathbf r_{\mathrm{bary}}(\oplus)&=\mathbf r_{\mathrm{bary}}(\mathrm{EMB})+\mathbf r_{\mathrm{EMB}\to\oplus},\\ \mathbf r_{\mathrm{bary}}(\mathrm{Moon})&=\mathbf r_{\mathrm{bary}}(\mathrm{EMB})+\mathbf r_{\mathrm{EMB}\to\mathrm{Moon}}. \end{aligned}

在 OPM 中，地球使用地月质心和月球地心向量重建：

rbary(⊕)=rbary(EMB)−rgeo(Moon)1+EMRAT.\mathbf r_{\mathrm{bary}}(\oplus)=\mathbf r_{\mathrm{bary}}(\mathrm{EMB})- \frac{\mathbf r_{\mathrm{geo}}(\mathrm{Moon})}{1+\mathrm{EMRAT}}.

B.2 角误差定义

方向角误差使用 cross/dot 的稳定形式：

errrad=atan2(∥a×b∥,a⋅b),errarcsec=errrad180π×3600.\begin{aligned} \mathrm{err}_{\mathrm{rad}}&=\mathrm{atan2}(\lVert\mathbf a\times\mathbf b\rVert,\mathbf a\cdot\mathbf b),\\ \mathrm{err}_{\mathrm{arcsec}}&=\mathrm{err}_{\mathrm{rad}}\frac{180}{\pi}\times3600. \end{aligned}

其中 a 为 DE441 派生真值向量，b 为 OPM 或 Swiss Ephemeris 向量。

附录 C. 生产生成与验证路线

本文 600 年生产实例的生成与验证路线可以概括为：

raw OPM fit from DE441
  -> polish 太阳 with native km metric
  -> polish 水星/金星/月球 with native angular guarded objective
  -> polish EMB and SSB-stored bodies with polished-Sun-anchor composite metric
  -> dense geocentric validation against DE441 and Swiss Ephemeris
  -> native position/analytic-velocity residual diagnostics

该路线把生成阶段优化和验证阶段诊断分开：误差修正目标用于降低本文生产实例的尾部误差，密集验证则独立报告地心角误差、原生位置残差和解析速度残差。若后续继续优化，应优先区分两类目标：一类是继续压缩文件体积，尤其是月球和水星；另一类是在不显著增大体积的前提下，进一步压低天王星、海王星等最接近基线天体的误差尾部。

附录 D. 密集误差曲线缩略图

本文 600 年生产实例的密集地心角误差验证输出见 out/opm2-600y-selected-polished-swiss-geocentric-512.txt，完整 SVG 曲线位于 out/opm2-600y-selected-polished-swiss-plots/angular/。下列缩略图均由同一次密集验证生成，使用每段 512 个切比雪夫节点加端点的确定性网格。

太阳月球
太阳密集误差曲线月球密集误差曲线
水星金星
水星密集误差曲线金星密集误差曲线
火星木星
火星密集误差曲线木星密集误差曲线
土星天王星
土星密集误差曲线天王星密集误差曲线
海王星冥王星
海王星密集误差曲线冥王星密集误差曲线

附录 E. 原生残差诊断曲线缩略图

本文 600 年生产实例的原生残差诊断输出见 out/opm2-600y-selected-polished-native-residuals-512.txt，完整 SVG 曲线位于 out/opm2-600y-selected-polished-native-plots/native/。下列缩略图均由同一次诊断生成，位置残差单位为 km，速度残差单位为 mm/s。

太阳位置太阳速度
太阳原生位置残差太阳原生速度残差
水星位置水星速度
水星原生位置残差水星原生速度残差
金星位置金星速度
金星原生位置残差金星原生速度残差
EMB 位置EMB 速度
EMB 原生位置残差EMB 原生速度残差
月球位置月球速度
月球原生位置残差月球原生速度残差
火星位置火星速度
火星原生位置残差火星原生速度残差
木星位置木星速度
木星原生位置残差木星原生速度残差
土星位置土星速度
土星原生位置残差土星原生速度残差
天王星位置天王星速度
天王星原生位置残差天王星原生速度残差
海王星位置海王星速度
海王星原生位置残差海王星原生速度残差
冥王星位置冥王星速度
冥王星原生位置残差冥王星原生速度残差

附录 F. OPM 与 Swiss Ephemeris 原生位置/速度残差对比曲线

本文 600 年生产实例的 OPM 与 Swiss Ephemeris 原生向量残差对比数据来自 out/opm2-600y-selected-polished-native-residuals-512.txt，完整 SVG 曲线位于 out/opm2-600y-selected-polished-swiss-plots/native-opm-vs-swiss/。下列曲线均以 DE441 为真值；绿色为 OPM，橙色为 Swiss Ephemeris。位置残差单位为 km，速度残差单位为 mm/s。EMB 的 Swiss 原生参考在表 3 中按 Earth SSB 与地月向量合成，缩略图部分只列出直接绘制的 Swiss 可寻址天体。

太阳位置太阳速度
太阳 OPM Swiss 原生位置残差对比太阳 OPM Swiss 原生速度残差对比
水星位置水星速度
水星 OPM Swiss 原生位置残差对比水星 OPM Swiss 原生速度残差对比
金星位置金星速度
金星 OPM Swiss 原生位置残差对比金星 OPM Swiss 原生速度残差对比
月球位置月球速度
月球 OPM Swiss 原生位置残差对比月球 OPM Swiss 原生速度残差对比
火星位置火星速度
火星 OPM Swiss 原生位置残差对比火星 OPM Swiss 原生速度残差对比
木星位置木星速度
木星 OPM Swiss 原生位置残差对比木星 OPM Swiss 原生速度残差对比
土星位置土星速度
土星 OPM Swiss 原生位置残差对比土星 OPM Swiss 原生速度残差对比
天王星位置天王星速度
天王星 OPM Swiss 原生位置残差对比天王星 OPM Swiss 原生速度残差对比
海王星位置海王星速度
海王星 OPM Swiss 原生位置残差对比海王星 OPM Swiss 原生速度残差对比
冥王星位置冥王星速度
冥王星 OPM Swiss 原生位置残差对比冥王星 OPM Swiss 原生速度残差对比

参考文献

[1] Park, R. S., Folkner, W. M., Williams, J. G., & Boggs, D. H. (2021). The JPL planetary and lunar ephemerides DE440 and DE441. The Astronomical Journal, 161(3), 105. https://doi.org/10.3847/1538-3881/abd414

[2] Kammeyer, P. (1988). Compressed planetary and lunar ephemerides. Celestial Mechanics, 45(1—3), 311—316.

[3] Koch, D., & Treindl, A. (1997—2022). Swiss Ephemeris — computer ephemeris for developers of astrological software. Astrodienst AG. https://www.astro.com/swisseph/swisseph.htm

[4] Moshier, S. L. PLAN404: The planets according to DE404. Steve Moshier’s numerical astronomy software archive. https://moshier.net/ ; data package: http://www.moshier.net/plan404.zip