使用 Python 分析 phyphox 数据（真实飞机采样）

七月初的时候回了一趟内蒙古，在回来的路上，我打开了phyphox：

官网

官网

它是一个可以保存和读取可视化手机传感器的软件，我就保存了从滑行到爬升以及平稳的过程一些传感器数据，瞎几把分析一下，就当黑匣子分析了：

数据原貌

数据原貌

需要把文件转换成 xlsx 的

需要把文件转换成 xlsx 的

主要有：

加速度计 (Accelerometer)：三轴加速度（含重力分量）

陀螺仪 (Gyroscope)：三轴角速度

线性加速度 (Linear Acceleration)：已去除重力分量的三轴加速度

气压计 (Pressure)：环境气压随时间变化

先放个窗外的样子，俺们内蒙古还是好看的！

先放个窗外的样子，俺们内蒙古还是好看的！

提前声明，没分析出特别明显的好结果，因为手机是直接被我坐在屁股下面的；就略懂一点无人机知识，如果有问题，那我觉得也是没有问题滴（嘿嘿）

手机记录到的 飞机起飞时的线性加速度变化曲线（计算三轴合成）：

滑跑阶段：曲线在起始时接近零，然后逐渐出现明显上升，说明飞机在跑道上逐步加速。

持续推力：在滑跑后期，合加速度维持在一个较高的水平，表现出发动机全推力输出。

拉升阶段：起飞瞬间通常会看到合加速度的波动，随后进入较为稳定的水平（机体爬升、升力抵消部分重力）。

因为我忘了给数据打一个标，其实没太好能知道到底哪里是飞起来了。。。

0–120 s 放大图：验证早段主要是滑行/等待，不是起飞；

170–190 s 放大图：高亮“滑跑→抬头→爬升”的关键时刻与区间。

这次起飞发生在 t≈170–190 s 区间：

滑跑/全推力段（Detected Run）：169.9–189.0 s（约 19.2 s）——线性加速度（平滑后）持续高于基线 + 0.2 m/s²。

抬头（Rotation）：≈188.2 s——在滑跑末端附近出现明显的俯仰角速度峰值（由陀螺仪 y 轴自动检出）。

基线（0–120 s 的中位合加速度）：≈0.10 m/s²；滑跑判定阈值：≈0.30 m/s²。

在 170–189 s 区间，用线性加速度（已去重力）做了“地面滑跑速度”积分；用 气压→高度（标准大气公式）画了起飞前后高度曲线，并在“抬头”后 5–25 s 的窗口估算平均爬升率。

速度曲线：积分得到的速度量级只有 ≈ 1 m/s，并且符号为负（只是轴方向选取的符号问题）。不管符号，量级明显不合理—真实起飞速度通常几十 m/s。

高度曲线：在“抬头”附近，条形高度先下降后再上升，导致“平均爬升率”算出是 -0.69 m/s（负号）；这同样不符合物理直觉。

为什么会这样？

这是手机 IMU 做“惯导积分”的常见坑，主要有三点：

手机的“线性加速度”数据带有强滤波与慢漂移（厂商做了重力估计/滤波），基线误差只有 0.01–0.02 m/s²，20 s 积分就会累计成数 m/s 的速度误差。

手机坐标系与跑道前进方向未知。虽然我用 XY 平面做了 PCA 找“前向轴”，方向可用，但尺度/偏置仍决定性地影响积分结果。

气压高度在起飞阶段容易受空调引气/通风气流影响，短时间内会出现“压力上升→高度下降”的假象（尤其是舱门关闭、空调工况变化、地面滑行至加速段切换时），因此 抬头后 10–20 s 的条形高度并不一定单调上升。

直接用“线性加速度”在 170–190 s 内积分，得到的地面速度曲线不可信；抬头后的条形高度/爬升率在这 20 s 内也受机舱压力扰动，短窗统计会偏差很大。

更靠谱的算法

用“原始加速度 + 陀螺仪”做轻量级姿态解算（Mahony/Madgwick），把 重力分量去得更干净，得到机体坐标下的真实线加速度；重新做“前向轴”标定：

用滑跑段（t≈170–189 s）内的 水平面加速度做 PCA 得方向；在滑跑开始前 2–5 s 做“零速度约束”（ZUPT）校准加速度偏置；速度积分：从 t0（≈170 s）开始，v(t0)=0，梯形积分 + 轻度低通，并在抬头后短窗内限制纵向分量，避免把重力误差积进去；气压高度：改用更长窗的趋势（比如抬头后 30–120 s），并用 中值滤波 抑制空调/通风的短期扰动，再计算平均爬升率

从姿态修正后的世界坐标线加速度看，滑跑持续到 t≈233.8 s 才达到“抬头（Rotation）”峰值；因此把积分窗口 扩展为 170–234 s 更合理。

地面速度曲线（170–234 s，已姿态修正）：速度单调上升，在 t≈233.8 s 标注了“抬头”。但积分得到的起飞瞬时速度 ≈ 4.64 m/s（图上已标），量级明显偏低（真实应为几十 m/s）。 这说明仅靠手机 IMU 做“零约束 + 前向投影”的绝对速度不可相信（偏置/姿态误差在 60 s 内已足够致命）。曲线形状仍可参考（相对变化），绝对值不可用。

气压→高度（平滑后）：抬头后 30–90 s 的窗口内，得到的平均“爬升率”≈ +0.39 m/s。这是客舱气压对应的“客舱高度”变化，受增压系统调度影响，并非飞机真实外界爬升率，因此仅能作为“舱内压力趋势”的佐证。

为什么绝对速度会偏离这么多？

加速度极小偏置（~0.01 m/s²）也会在几十秒积分成数 m/s 的速度误差；加速度计在强加速时参与姿态估计（我们用互补滤波），会把部分前向加速度“吃掉”；即便把滤波时常调长，仍不可避免；无磁力计/无 GPS 约束时，航向不可观测，前向轴用 PCA 估，尺度不受影响，但偏置/倾斜仍会漏进来；手机并未刚性固定，座椅/手持的柔顺性会吸收/滤去高频加速度。

第一张：俯仰角速度 gy（rad/s） 的时间曲线，蓝色虚线是检测到的峰值时刻。

第二张：用互补滤波（加计+陀螺）估出的俯仰角（°），同样标出峰值时刻，能看到抬头后的机头上仰趋势。

抬头峰值时刻：t ≈ 233.80 s

俯仰角速度峰值：≈ 0.079 rad/s（≈ 4.53 °/s）

半峰宽（|gy| ≥ 50% 峰值的持续时长）：≈ 0.094 s（比较尖锐的脉冲）

俯仰角变化（±5 s 窗口中位数差）：≈ +1.08 °

这些量化指标符合“起飞抬头”的典型特征：在滑跑末端出现一记短促而明显的俯仰角速度峰值，随后俯仰角进入更高的平台。

综合对齐图

综合对齐图

蓝色：前向加速度（推力）——在滑跑后期显著抬升。

橙色：俯仰角速度 gy（陀螺仪 y 轴）——在 t≈233.8 s 出现一个短促而明显的峰值（抬头）。

绿色：客舱气压换算的高度——在抬头后的几十秒内逐渐上升（受客舱压力控制影响，爬升率滞后于真实外界高度变化）。

黑色虚线：标记抬头峰值时刻。

从时序上可以清楚看到：

推力加速（蓝）在起飞前 20–30 s 已经持续增强；

俯仰角速度（橙）在 233.8 s 左右出现一个脉冲，明确对应“机头上仰”；

客舱高度（绿）在抬头后才逐步上升，且有一定延迟。

标准

标准

把三条曲线标准化（z-score）并叠加到同一坐标的对齐图（时间窗 200–280 s，已在 t≈233.8 s 标注抬头）：

前向加速度（蓝）：在 210 s 左右达到一波高峰，随后逐步回落到接近 0 的巡航加速度附近——对应推力主要用于地面滑跑加速，抬头前后加速度不再大幅攀升。

俯仰角速度 gy（橙）：在 233.8 s 附近出现尖峰，清晰对应“拉杆抬头”。该峰与蓝线高峰在时间上有明显错位（先加速、后抬头）。

客舱高度（绿）：为单调的上升趋势（经平滑），且滞后于抬头事件，符合“抬头后、离地后、客舱压力逐步调整”的时序特征。

{'fs_Hz': 10.58923580248911,
 'take_window_s': [200.0, 260.0],
 'cruise_window_s': [300.0, 360.0],
 'peak_take_Hz_dB': (0.6618272376555694, -37.38242013966849),
 'peak_cruise_Hz_dB': (0.6618272376555694, -45.340300291511326),
 'band_metrics': [{'band_Hz': '0.2-1.0',
   'RMS_take_rad_s': 0.010546424695589227,
   'RMS_cru_rad_s': 0.005790965600787031,
   'take_vs_cru_dB': np.float64(5.207085417040597)},
  {'band_Hz': '1.0-3.0',
   'RMS_take_rad_s': 0.008130079464038133,
   'RMS_cru_rad_s': 0.0034268772414093505,
   'take_vs_cru_dB': np.float64(7.50392486312456)},
  {'band_Hz': '3.0-8.0',
   'RMS_take_rad_s': 0.004788388290284797,
   'RMS_cru_rad_s': 0.003314200734378562,
   'take_vs_cru_dB': np.float64(3.1962110175779532)},
  {'band_Hz': '8.0-20.0',
   'RMS_take_rad_s': nan,
   'RMS_cru_rad_s': nan,
   'take_vs_cru_dB': np.float64(nan)}]
   }

做了 gy（俯仰角速度）的频谱对比：把起飞阶段 200–260 s和后段 300–360 s（等长的“巡航”窗口）进行 Welch PSD 分析。上图已给出对比曲线（对数频轴）。

采样率估计 fs ≈ 10.59 Hz（Nyquist ≈ 5.3 Hz），所以本次只能观察到 ≤ 5 Hz 的姿态微振动。

起飞段在 0.2–3 Hz 的能量显著高于巡航段：

峰值频率（0.5–20 Hz 搜索范围内）：≈ 0.66 Hz；起飞段峰值 -37.38 dB，巡航段 -45.34 dB（起飞更强）。

频带 RMS（单位：rad/s）与相对增益：

0.2–1 Hz：起飞 0.0105，巡航 0.00579，起飞比巡航 +5.21 dB

1–3 Hz：起飞 0.00813，巡航 0.00343，起飞比巡航 +7.50 dB

3–8 Hz：起飞 0.00479，巡航 0.00331，起飞比巡航 +3.20 dB

8–20 Hz：受限于采样率，不可观测

起飞滑跑/抬头阶段包含更强的低频姿态激励（0.5–1 Hz 一带），可能与地面滚动/起落架—机体耦合、推力/气动加载变化与座椅/人体传递叠加有关。进入巡航后，这些低频成分显著减弱。10.6 Hz 的采样率过低，无法分析 5–20 Hz 这类更“颤动感”的座椅/机体模态；建议把陀螺采样提升到 100–200 Hz，可同时看清 5–30 Hz 的结构/内饰件共振。

加速度分析

加速度分析

三轴加速度：

ax（假设是沿飞机前进方向）在起飞时有明显的正峰值（推力加速）。

ay 基本在 0 附近波动（横向加速度小）。

az 以 -9.8 m/s² 为基准（重力），但起飞和爬升阶段有小幅扰动。

速度（由 ax 积分得到）：积分后得到速度随时间上升，在滑跑 → 拉升阶段达到最大值，基本符合“起飞速度”；由于噪声和漂移，后期速度曲线会慢慢漂移。

高度（由气压推算）：地面阶段高度接近常数；起飞点以后高度迅速上升，进入稳定爬升。

加速度峰值 → 速度快速增加 → 高度开始上升，符合飞机起飞物理过程。

加入模长

加入模长

这里给出了三轴加速度（ax、ay、az）的时域曲线，以及它们合成的加速度模长 |a|：三个方向的加速度曲线可以看到运动过程中的不同分量变化；合成向量 |a| 代表了总的加速度强度，能直观反映整体运动强度（包括平移和震动）。

幅度谱（FFT 单边）

Welch 功率谱密度（PSD，log 频轴）

对 合成加速度向量 |a| = √(ax²+ay²+az²) 做了频谱分析，并给出两张图：

估计采样率 ≈ 10.0 Hz，Nyquist ≈ 5 Hz，因此这次频谱只能覆盖到 ≤ 5 Hz 的颠簸/振动。

|a| 的低频能量（<1 Hz）较强，这是合理的：起飞/爬升过程中的缓慢姿态/推力变化会在合成模长里体现为低频成分；在 1–3 Hz 一带也能见到一些能量凸起（具体主峰见“Top Peaks”表），可能对应跑道/座椅传递的低频振动。但由于 10 Hz 采样率有限，5–20 Hz 更“嗡嗡”的机体/内饰件模态无法观测。

使用模长 |a| 会把 重力 1g “折叠”进来（我已经用中位数去掉了DC），但它仍然是“总强度”，方向信息丢失。如果要区分“前向推力 vs 垂向过载”，需要回到三轴或做姿态解算后的分解。频谱峰值受窗长与采样率限制：本次 fs≈10 Hz，频率分辨率约 Δf ≈ 1/T_total，想更清晰地区分窄峰，需要更长记录或更高采样率。

时频图（spectrogram）

时频图（spectrogram）

这是合成加速度向量的时频图（spectrogram）：

横轴：时间（s），覆盖了整个起飞和飞行过程

纵轴：频率（Hz）

颜色：能量强度（dB）

在 170–190 s 起飞阶段，可以看到 低频能量（<2 Hz）显著增强，对应滑跑与抬头时的机体整体加速度变化；在进入巡航阶段后，能量逐渐集中在更低频，说明机体进入稳定状态。

合成加速度向量在起飞（170–190s）与巡航（400–420s）阶段的功率谱对比：

起飞阶段：低频（<2 Hz）能量显著增强，代表飞机起飞时较大的振动和姿态变化。

巡航阶段：谱能量整体下降，曲线更平稳，说明飞行进入稳定状态，主要是背景噪声和轻微震动。

Spectrogram：

起飞段（约 200–240 s）内，**<1 Hz 与 1–3 Hz** 的能量明显上涨，抬头附近尤为突出——说明当时整体姿态/推力变化带来了低频加速度增强；进入后段后，谱能量回落到更低水平，呈现更稳定的背景。

Band power：

0.2–1 Hz：在抬头前后出现清晰的脉冲式增强。

1–3 Hz：起飞阶段显著高于其他时段，常见于跑道/结构耦合与抬头机动。

3–5 Hz：也能看到起飞时的抬升，但幅度较低（受限于本次采样率 ~10 Hz 的 Nyquist）。

上面用 STFT 做了 |a|（合成加速度）的时频图，并且把能量按频带做成了随时间变化的曲线，还用陀螺仪自动检测的抬头时刻（≈233.8 s）做了标注。

gy 在 233.8 s 有一个短促而明显的峰值，典型的“拉杆抬头”信号；gx/gz 在该时刻没有同级别尖峰，说明主要是俯仰机动。

抬头附近低频能量瞬时上冲，然后回落到较稳态的低水平，符合“起飞→稳定爬升/巡航”的过渡。

先用前 120 s 的中位数当基线，给 gx/gy/gz 做了偏置校正；然后绘制三轴角速度（已去偏置）随时间的变化，并标记 抬头峰值时刻；对 gy 在抬头附近做积分，得到 俯仰角变化 Δpitch 的近似轨迹（只用陀螺，演示用，易漂移）；把 抬头前后 60 s（203.8–263.8 s） 和 较后段等长 60 s（323.8–383.8 s） 的 gy 频谱做对比。

时频图（spectrogram）：查看 gy 能量随时间的迁移，并在图上标出 抬头时刻。

采样率估计：fs ≈ 10.59 Hz

零偏（前 120 s 中位数）：gx −0.00153、gy +0.00061、gz −0.00046 rad/s

抬头峰值：t ≈ 233.80 s，gy ≈ +0.0785 rad/s（≈ 4.50 °/s）

PSD 窗口：起飞/抬头段 203.8–263.8 s vs 巡航段 323.8–383.8 s

接下来就是陀螺仪：

三轴（去偏置）+ |ω|

|ω| 的功率谱密度

|ω| 的能量随时间迁移

采样率估计：fs ≈ 10.59 Hz（Nyquist ≈ 5.3 Hz，当前仅能观察到 ≤5 Hz 的角振动）

零偏（前 120 s 中位数）：gx −0.00153、gy +0.00061、gz −0.00046 rad/s

时间序列：|ω| 在起飞/抬头附近相对增强，随后趋稳。

PSD：低频（约 0.2–2 Hz）能量较高，反映姿态慢变化和大尺度机动；受采样率限制，5–30 Hz 的座椅/机体模态无法观测。

时频图：在起飞阶段 |ω| 的低频能量有明显“冲高”，之后回落，符合“起飞 → 稳定段”的过渡。

这个是没 g 的数据，也是加速度，和上面的分析方法是一样的：

采样率估计：fs ≈ 10.59 Hz（Nyquist ≈ 5.29 Hz）

残余零偏（前 120 s 中位数）：

lx ≈ −3.88×10⁻⁴ m/s²

ly ≈ +3.20×10⁻³ m/s²

lz ≈ +1.80×10⁻³ m/s²

时间序列：起飞前后 |a| 明显抬升，随后趋于较低水平，符合“滑跑加速 → 抬头 → 稳定段”的变化。

PSD：能量主要集中在 0.2–2 Hz 的低频（受采样率限制，>5 Hz 不可观测）。

时频图：在起飞窗附近低频能量“冲高”，之后回落，和陀螺姿态时频的结论一致。

滑跑区间速度与高度对比

滑跑区间速度与高度对比

把 Linear Acceleration 的 ax 在滑跑窗口（170–190 s）积分得到前向速度，并与 Pressure 换算的高度叠加绘制：

蓝线：前向速度（由 ax 积分得到），清晰展现了滑跑加速过程。

红线：高度（由压力推算），在滑跑阶段几乎保持水平，起飞后才快速上升。

绘制出 巡航阶段 (t > 300 s) 的前向速度曲线（由 ax 积分得到）与气压高度曲线 对比图。蓝色表示速度，红色表示高度，可以用来分析在巡航阶段速度的稳定性与高度保持情况。

集合

集合

把整段飞行（从 t=0 开始）的关系画在了一张图上：

左轴（线）：前向速度（由 Linear Acceleration 的 lx 在滑跑窗 170–234 s 内积分而得；窗外不积分，避免漂移）。

右轴（线）：高度（由 Pressure 换算的条形高度，整段时间）。

灰色区块：滑跑阶段（170–234 s）。

虚线：抬头时刻 ≈ 233.8 s。

可以看到：滑跑阶段速度快速上升；在抬头附近后，高度开始明显上升，随后进入爬升/巡航段。

扩展速度曲线（抑制漂移）

扩展速度曲线（抑制漂移）

以前窗的 5–20 s 做零偏，再对前向加速度做 ~20 s 滑动均值高通 抑制低频漂移，然后积分得到全程的相对前向速度（蓝线）；在滑跑开始时重新锚定 v(t0)=0，避免一开始的累积误差。

自动相位边界检测并标注

Ground roll（灰色区）：从加速度模长高于基线的最长片段自动找出 → ≈ 168.83–191.78 s。

Rotation（黑虚线）：在滑跑末段用 gy 峰值自动检测 → ≈ 188.19 s。

Liftoff（绿虚线）：自 Rotation 后，寻找高度 5 s 连续上升的首时刻 → ≈ 196.97 s。

Cruise start（橙虚线）：在后段，寻找 |爬升率| < 0.2 m/s 持续 ≥ 60 s 的最早时刻 → ≈ 405.67 s。

右轴为条形高度（由气压换算，整段时间）；可以清楚看到速度（蓝）与高度（浅蓝）的先后关系与转折点。

这里把 滑跑 → 起飞 → 爬升 → 巡航 四个阶段的关键参数（速度、高度、加速度、俯仰角速度）统一绘制在一张图里，可以完整看到：

滑跑阶段 (0–50s)：速度线性上升，高度保持 0，加速度有小抖动，角速度基本为 0。

起飞抬头 (50–80s)：速度继续增加，高度开始上升，加速度逐渐减弱，俯仰角速度出现峰值。

爬升 (80–150s)：速度平缓增长，高度快速上升，加速度和角速度有小扰动。

巡航 (≥150s)：速度和高度基本稳定，只有微小扰动。

IMU 与气压计在时间上的一致性与相位关系

IMU 与气压计在时间上的一致性与相位关系

在同一张图上叠加了前向速度、条形高度、以及“由高度反算得到的等效气压（并与原始气压对比）”，并标注了 抬头时刻 ≈ 233.8 s。这样你能同时看到。

绿色（高度）在抬头后单调上升；

红/橙（原始气压 vs 根据高度反算气压）高度一致，说明换算与原始数据互相印证；

蓝色（速度，来自线性加速度经高通抑制后积分）在滑跑段前后出现明显相对变化，但由于手机放置与采样率限制，绝对量纲只能作为相对指标。

滑跑 / 离地+爬升 / 巡航三段阴影区块与标注

滑跑 / 离地+爬升 / 巡航三段阴影区块与标注

灰色：Ground roll（自动检测的滑跑段）

金色：Takeoff + Climb（自“Rotation”后到进入平稳爬升结束）

青色：Cruise（持续低爬升率区段）

同图叠加了：前向速度（蓝，左轴）、条形高度（绿，右轴）、原始气压与由高度反算气压（红/橙，右次轴），以及三条关键竖线：Rotation / Liftoff / Cruise start。

Forward Velocity、Altitude、加速度、角速度、线性加速度

Forward Velocity、Altitude、加速度、角速度、线性加速度

标注了 滑跑、起飞爬升、巡航 三个阶段的阴影分区。

把 加速度模量、角速度模量、速度、气压高度 全部叠加显示在 4 个子图中，并且在 X 轴上加入了三个飞行阶段的阴影分区：

滑跑（0–20s）

起飞/爬升（20–60s）

巡航（60–120s）；每个阶段都用灰色阴影标注，并在图上写了文字。

这东西早就分析完了，我也不知道分析的怎么样，如果感兴趣的人多了，相关分析代码也可以开源出来。