状态：草稿 (2026-10-18)，实测。 数字是通过笔记本电脑（Apple M1 Max，macOS 26.5）上的 Python 3.9.6 和 numpy 2.0.2 直接计算出的模拟值。 原始数据与代码：工作目录
measure/robot-arm-basics/的common.py,m06_flight.py,m06b_constants.py,out/m06_flight.json,out/m06b_constants.json. 「真实轨迹」也是通过相同模型生成的，因此模型与实际球体之间的差异并未包含在本文的数字中（与实际球体轨迹的对比未实测）。 球的初始状态采用 DeepMind 数据集（CC-BY 4.0）。参考资料查阅于 2026-10-18。 编写规则：以实测数据为基础，无法测量的标注为「未实测」，无法确认的标注为「未确认」。
目录
1. 直觉：轻盈且高速旋转的球
乒乓球重 2.7 克，直径 40 毫米（ITTF 规章 2026，2.3.1·2.3.2）。它比同尺寸的铁球轻得多，因此空气推力相对于重量而言非常大。 以 5 m/s 飞行时，阻力约为球重的三分之一；若为 10 m/s，则大于球重。此外，运动员击出的球每秒旋转数十圈。旋转的球会把一侧的空气带走更多，从而受到侧向力。这就是马格尼斯效应，也是上旋球急剧下坠、下旋球向上飘的原因。
因此，乒乓球机器人的预测器不能仅依靠真空中的抛物线。然而，模型越复杂，需要了解的常数（阻力系数、马格尼斯系数、反弹与摩擦）和需要观测的状态（旋转）就越多。本文将用数字来权衡这种取舍。
2. 运动方程
m · dv/dt = m g − ½ ρ C_D A |v| v + ½ ρ A r C_M (ω × v)
dp/dt = v
m = 2.7 g, r = 0.02 m, A = π r², ρ = 1.225 kg/m³
- 第二项是阻力（drag）。它与速度的平方成正比，方向与速度相反。
- 第三项是马格尼斯力。它是旋转轴 ω 与速度 v 的外积方向，使球向侧面偏转。对于向机器人方向（−y）飞来的球，ω_x > 0 表示上旋，ω × v 指向下方，将球向下拉拽。
- 假设飞行过程中的旋转保持不变（第 1 部分的共同假设）。
球台反弹被建模为瞬时冲量模型。法向速度通过反弹系数 e 反转，接触点的滑动速度通过库仑摩擦 μ 减小，当达到滚动状态（接触点速度为 0）时停止。薄壳球（I = ⅔ m r²）在消除滑动速度的五分之二时即达到滚动状态。
def bounce(v, w, e=0.876, mu=0.25):
ux, uy = v[0] - R*w[1], v[1] + R*w[0] # 접점(공 아래쪽)의 미끄럼 속도
us = hypot(ux, uy)
jn = (1 + e) * abs(v[2]) # 질량당 법선 충격량
dv = min(mu * jn, 0.4 * us) # 마찰은 굴림에 닿으면 멈춤
jx, jy = -dv*ux/us, -dv*uy/us
k = M / I_BALL # = 1.5 / R²
return (v[0]+jx, v[1]+jy, -e*v[2]), (w[0] + k*R*jy, w[1] - k*R*jx, w[2])
2.1 常数从何而来
| 常数 | 本系列 | 文献中确认的值和条件 |
|---|---|---|
| 阻力系数 C_D | 0.47 | Sony AI 2026 表 I：2.5 m/s 时为 0.55，17.5 m/s 时根据旋转比为 0.37~0.51。Google 2023 的 PyBullet 设置为 0.47。将 Koç 等人 2018 表 A.3 中的组合常数 0.141 1/m 按空气密度 1.2 换算约为 0.50（本计算） |
| 马格尼斯系数（每单位质量） | 3.4×10⁻³ | Koç 等人 2018 表 A.3 中的升力常数 1.0×10⁻³ 1/rad。在同一项目的 MuJoCo 模拟器（DeepMind 系数）中测得的有效值为 7.6×10⁻³。三个值相差超过 7 倍 |
| 球台反弹 e | 0.876 | 在无空气情况下对 ITTF 2.1.3（30 cm → 约 23 cm）进行换算。Nakashima 等人 2010 实测值为 0.93，Koç 等人 2018 估计值为 0.883 |
| 球台摩擦 μ | 0.25 | Nakashima 等人 2010 测得的滑动摩擦值为 0.25（Sony AI 2026 采用相同值）。Koç 等人 2018 估计值为 0.102 |
文献中马格尼斯系数差异巨大，是因为定义（将什么组合并称为常数）和测量条件（速度、旋转比）不同。Sony AI 2026 拟合了随旋转比呈线性及二次变化的 C_M。本文将在第 5 节中直接测量这种不确定性对预测的影响程度。
3. 欧拉法与 RK4
要在计算机上求解微分方程，需要将时间细分为小的时间步长，一步步前进。最简单的欧拉法是利用当前的斜率沿 dt 直行。RK4（四阶龙格-库塔法）在单步内测量四次斜率并进行加权平均。
오일러: v ← v + a(v)·dt, p ← p + v·dt (걸음당 기울기 1번, 오차 ∝ dt)
RK4: k1 = f(x), k2 = f(x + k1·dt/2), k3 = f(x + k2·dt/2), k4 = f(x + k3·dt)
x ← x + (k1 + 2k2 + 2k3 + k4)·dt/6 (걸음당 4번, 오차 ∝ dt⁴)
3.1 实测：0.4秒后的位置误差
初始位置 (0, 1.2, 0.25) m，速度 (0.3, −7.0, 1.2) m/s（速率 7.1 m/s），上旋 80 rad/s 的球在无反弹情况下飞行了 0.4 秒。基准为 RK4 dt = 0.01 ms。
| 方法 | dt | 0.4秒后的误差 | 单条轨迹计算时间 |
|---|---|---|---|
| 欧拉法 | 0.5 ms | 0.96 mm | 808 µs |
| 欧拉法 | 1 ms | 1.92 mm | 406 µs |
| 欧拉法 | 10 ms | 19.3 mm | 42 µs |
| 半隐式欧拉法 | 1 ms | 2.81 mm | 430 µs |
| RK4 | 1 ms | 9.4×10⁻¹¹ mm | 3,296 µs |
| RK4 | 5 ms | 5.7×10⁻⁸ mm | 658 µs |
| RK4 | 10 ms | 9.1×10⁻⁷ mm | 339 µs |
- 对于欧拉法，dt 减半则误差也减半（一阶）。RK4 则减小为十六分之一（四阶）。
- RK4 在 dt = 10 ms 时比欧拉法 dt = 0.5 ms 更快，且精度高出 100 万倍。 乒乓球飞行是一个平滑的微分方程，因此高阶方法具有压倒性优势。
- 半隐式欧拉法（先修正速度，再用新速度更新位置）在此问题中表现略差。由于它是一种用于保守系统的振动问题的方法，对受阻力影响的球没有优势。
- 但是，诸如反弹和球拍碰撞等 存在不连续性的事件如果在较大的 dt 下容易错过时间点。本系列将包含碰撞判定的飞行过程以 1 ms 积分，而无碰撞的区间则缩短至 5 ms（第 7 部分）。
4. 旋转与反弹改变了什么
4.1 旋转改变落点
比较了从对方半场 (0, 1.2, 0.25) m 以速度 (0, −6, 1) m/s 出发的球首次接触机器人一侧球台的位置，仅改变旋转进行对比。
| 旋转 | 50 rad/s | 100 rad/s | 150 rad/s |
|---|---|---|---|
| 上旋：落点向球网方向偏移 | +10.3 cm | +19.4 cm | +27.3 cm |
| 下旋：落点向底线方向偏移 | −11.8 cm | −25.7 cm | −41.7 cm |
若无旋转，球会在 0.344 秒时落在 y = −0.62 m 处。由于第 1 部分数据集显示击球后旋转的中位数为 49 rad/s，95 百分位数为 100 rad/s，因此即使是一般业余水平的球，旋转也会使落点偏移 10~25 厘米。
4.2 反弹检查：ITTF 跌落测试
将球底部从 30 cm 高度落下并测量回弹高度。若 e = 0.876，则在无空气阻力时应为 23 cm，但在考虑阻力的模型中为 21.5 cm 。这是因为在下落和上升过程中阻力消耗了能量。若要将阻力考虑在内并使结果达到 23 cm，则需要 e = 0.907（二分法）。该值介于 Nakashima 等人 2010 实测的 0.93 与 Koç 等人 2018 的 0.883 之间。仅反弹系数一项，也会因「忽略了什么并如何拟合」而产生 3% 的差异。
4.3 摩擦系数与旋转
假设在底线后 23 cm 的平面（y = −1.6 m）接球，当改变球台摩擦 μ 时，穿过该平面的高度（球台上方的距离）变化如下。
| 旋转 | μ = 0.10 | μ = 0.25 | μ = 0.40 |
|---|---|---|---|
| 下旋 100 rad/s | 0.19 m | 0.18 m | 0.15 m |
| 无旋转 | 0.17 m | 0.12 m | 0.02 m |
| 上旋 100 rad/s | 0.10 m | 0.005 m | 0.005 m |
上旋球在 μ 分别为 0.25 和 0.40 时结果相同。这是因为在反弹过程中，接触点达到了滚动状态，摩擦力不再发挥作用。无旋转的球对摩擦最敏感。 旋转和摩擦不能单独拟合，必须同时拟合。
5. 常数的不确定性对预测的影响有多大
测量了用数据集中的 1,000 个球以默认常数飞行的结果，与仅更改一个常数后飞行的结果在通过击球平面（y = −1.6 m）时的位置和时间差（m06b_constants.py).
| 更改的常数 | 位置差中位数 / 90%分位数 | 时间差中位数 | 未能到达击球平面的球 |
|---|---|---|---|
| C_D 0.47 → 0.40 | 3.1 / 12.4 cm | 30 ms | 21 |
| C_D 0.47 → 0.55 | 3.9 / 14.8 cm | 35 ms | 35 |
| 马格尼斯 → Koç 值 (1.0×10⁻³) | 5.4 / 11.5 cm | 6.7 ms | 43 |
| 马格尼斯 → MuJoCo 有效值 (7.6×10⁻³) | 10.4 / 22.3 cm | 11.3 ms | 53 |
| 无马格尼斯效应 | 7.4 / 15.5 cm | 9.6 ms | 72 |
阻力系数对 时间的影响大于位置。0.07 的差异导致了 30 ms 的偏差，这相当于第 1 部分时间预算中的四帧相机画面。马格尼斯系数对位置的影响很大。仅在文献值之间切换，击球平面的位置就会产生 5~10 cm 的变化。在制作好机器人后，这些常数必须 针对我们的球和我们的球台 重新拟合 。Sony AI 2026 报告称，他们利用 9 台相机收集的大约 5.8 万条轨迹重新拟合了系数，将着陆点中心误差从 0.37 m 降低到了 0.15 m。
6. 卡尔曼滤波：从带噪声的观测中提取状态
相机仅提供位置，且位置本身包含噪声。速度需要对位置求导，而求导会放大噪声。卡尔曼滤波（Kalman filter）根据「模型预测的状态」各自的不确定性与「新的观测值」进行融合，共同估计位置和速度。如果模型是非线性的（如阻力），则会变成在每一步将模型线性化的扩展卡尔曼滤波（EKF）。
sp = f_state(s, dt) # 예측: RK4 로 한 프레임 앞으로
F = finite_difference_jacobian(f_state, s, dt) # 선형화
P = F @ P @ F.T + Q # 불확실성도 함께 전파
K = P @ H.T @ inv(H @ P @ H.T + R) # 칼만 이득: 모델과 관측의 신뢰 비율
s = sp + K @ (z - H @ sp) # 갱신: 관측 z(위치 3개)로 고침
P = (I - K @ H) @ P
6.1 实测：什么能改善预测
我们从数据集中挑选了 200 个飞向击球面的实际球轨迹，生成了真实轨迹（空气阻力 + 实际旋转 + 弹跳）。模拟了 60、120、240 fps 的相机观测，各轴高斯噪声 $\sigma$ 为 3、10 mm。在击球后观测 0.1、0.2、0.3 秒后，使用滤波器状态预测了穿过击球面 $y = -1.6\text{ m}$ 的位置。模型共有三种。
- A 仅重力: 线性卡尔曼滤波。预测结果也是真空抛物线 + 弹跳。
- B 重力 + 空气阻力: EKF。假设旋转为 0。
- C 重力 + 空气阻力 + 真实旋转: 假设准确知道旋转的 EKF。这是旋转估计完美时的下限。
这是 120 fps，$\sigma = 3\text{ mm}$ 的结果（200 个球的中位数）。
| 观测时间 | A 位置误差 | B 位置误差 | C 位置误差 | B 时间误差 | C 时间误差 |
|---|---|---|---|---|---|
| 0.1 s | 15.3 cm | 7.0 cm | 1.3 cm | 14 ms | 2.9 ms |
| 0.2 s | 13.3 cm | 5.6 cm | 0.5 cm | 13 ms | 1.3 ms |
| 0.3 s | 10.5 cm | 4.8 cm | 0.5 cm | 10 ms | 1.2 ms |
这是更改帧率和噪声后的结果（观测 0.2 秒，位置误差中位数）。
| 条件 | B（未知旋转） | C（已知旋转） |
|---|---|---|
| 60 fps, $\sigma$ 3 mm | 5.9 cm | 0.77 cm |
| 120 fps, $\sigma$ 3 mm | 5.6 cm | 0.51 cm |
| 240 fps, $\sigma$ 3 mm | 5.9 cm | 0.39 cm |
| 60 fps, $\sigma$ 10 mm | 6.2 cm | 2.0 cm |
| 240 fps, $\sigma$ 10 mm | 6.0 cm | 1.2 cm |
- 对于不考虑旋转的模型（B），无论观察多久、频率多高，误差都无法降到 5 cm 以下。 因为误差的大部分不是来自噪声，而是来自模型误差。
- 如果已知旋转（C），误差会减少到十分之一，从那时起帧率和噪声才具有意义。当噪声较大时，240 fps 比 60 fps 将误差减少了 40%。
- 真空模型（A）在 200 个球中有 28 到 62 个球甚至无法预测击球面的通过情况（弹跳位置定位错误）。
- 处理时间方面，Python EKF 每帧耗时 126 µs，预测到达击球面一次耗时 11 ms。
Mülling 等人（2013）使用了 4 台 60 Hz 相机和 EKF，由于旋转观测困难，在 EKF 模型中去掉了旋转（3.1 节）。Koç 等人（2018）也使用了 4 台 60 Hz 相机和 EKF，并在飞行模型中加入了升力项，但未能确认滤波器如何处理旋转（未确认）。这些实测表明，这种选择的代价是几厘米。由于球拍面半径在 7~8 cm 左右，这是一个可以承受的大小，但正如第 7 部分所示，它会直接蔓延到着陆点精度上。
7. 局限性
- 真实轨迹和滤波器使用相同的方程（仅常数不同）。实际的球具有旋转衰减、每个球不同的重心、风等模型之外的影响，因此误差可能会更大（未实测）。
- 观测噪声作为各轴独立的高斯噪声引入。实际的三角测量噪声在相机视线方向呈细长形，并且还存在球检测失败和延迟（未实测）。
- "已知"旋转的情况（C）仅仅是一个下限。不涉及从轨迹弯曲中估计旋转或通过球表面图案测量旋转的方法。
- 乒乓球专用风洞实验原文未能打开（未确认）。$C_D$ 范围是 Sony AI 2026 飞行轨迹拟合的结果。
8. 总结
- 乒乓球飞行可以通过重力、空气阻力、Magnus 三个力和反弹、摩擦弹跳来很好地解释。如果积分使用 RK4，5~10 ms 的间隔就足够了。
- 100 rad/s 的旋转会将弹跳位置移动 20 cm 以上，而 0.07 的空气阻力系数差异会改变到达时间 30 ms。
- 为了减少预测误差， 模型（尤其是旋转）比帧率更重要。如果不知道旋转，误差为 5 cm；如果知道，误差为 0.5 cm。
- 因此，在选择硬件时，继"这是几 fps 的相机"之后的下一个问题应该是"如何确定旋转"。
TT 点评。 即使将相机换成快四倍的相机，如果不知道旋转，球还是会飞到其他地方去。首要任务是弄清楚球是如何旋转的！
参考
所有资料均于 2026-10-18 查阅。
- ITTF, Statutes 2026, 2.1.3·2.3.1·2.3.2.
- Conti 等（Sony AI）, Physics Models for Sim-to-Real Transfer in Professional-Level Robot Table Tennis, arXiv:2606.28805, III-B·III-C 节，表 I。
- Koç, Maeda, Peters, Online optimal trajectory generation for robot table tennis, RAS 105 (2018), 作者公开版, 公式 9–16，表 A.3，附录 B。
- Nakashima 等, Modeling of Rebound Phenomenon of a Rigid Ball with Friction and Elastic Effects, ACC 2010, 机构存储库, IV 节。
- Mülling, Kober, Kroemer, Peters, Learning to select and generalize striking movements in robot table tennis, IJRR 32(3) (2013), 作者公开版, 3.1 节。
- D'Ambrosio 等, arXiv:2309.03315, 附录 D（模拟器设置）。
- Google DeepMind, competitive_robot_table_tennis, CC-BY 4.0。