3. 参考坐标系（下）：欧拉角、方向余弦矩阵、四元数基础
好，咱们接着聊坐标系。上一节我们把几个常用的坐标系定义清楚了，这一节要解决一个核心问题：怎么在两个坐标系之间来回切换？
说白了，飞行器在空中的姿态，就是机体坐标系相对于地面坐标系（或导航坐标系）的转动关系。描述这个转动关系，数学上主要有三种工具：欧拉角、方向余弦矩阵、四元数。我做了这么多年飞控，这三种方法都用过，各有各的脾气。
核心思想：任何两个坐标系之间的相对方位，都可以通过三次连续旋转来描述。这三次旋转的角度，就是欧拉角。
3.1 欧拉角：最直观，但有个坑
欧拉角的概念其实很简单。你想想看，一架飞机在空中，它先绕自己的Z轴转一个角度（偏航角 ψ），再绕新的Y轴转一个角度（俯仰角 θ），最后绕新的X轴转一个角度（滚转角 φ）。三次旋转，搞定。
我习惯用“Z-Y-X”这个旋转顺序，这也是航空航天领域最常用的。当然还有其他顺序，但咱们搞飞行器设计的，用这个就够了。
| 欧拉角 | 符号 | 旋转轴 | 范围 |
|---|---|---|---|
| 偏航角 | ψ | Z轴（初始） | [-π, π] |
| 俯仰角 | θ | Y轴（新） | [-π/2, π/2] |
| 滚转角 | φ | X轴（新） | [-π, π] |
优点很明显：物理意义清晰，飞行员一看就懂。“当前俯仰10度，滚转5度”，多直观。
但坑也在这里。我曾经在做一个全姿态控制算法时，飞机做大机动动作，俯仰角接近90度，结果姿态解算直接炸了——这就是著名的“万向节死锁”问题。当俯仰角为±90°时，滚转和偏航的旋转轴重合，丢失了一个自由度。
⚠ 避坑指南：我曾经在仿真中遇到过，当俯仰角接近90度时，欧拉角微分方程会出现奇异。如果你做的是全姿态飞行（比如特技飞行、导弹），千万别只用欧拉角。我当时花了两天排查这个问题，最后换成了四元数。
3.2 方向余弦矩阵：数学上完美，但计算量大
方向余弦矩阵（DCM），说白了就是一个3×3的矩阵，它能把一个向量从一个坐标系变换到另一个坐标系。
假设我们有一个向量在机体坐标系下的坐标 vb，想求它在地面坐标系下的坐标 vn，那么：
v_n = C_n^b · v_b其中 C_n^b 就是方向余弦矩阵。它的每个元素，其实就是两个坐标系对应轴之间的夹角余弦值。
用欧拉角表示的话，这个矩阵长这样（Z-Y-X顺序）：
C_n^b =
[ cosθ·cosψ, cosθ·sinψ, -sinθ ]
[ sinφ·sinθ·cosψ - cosφ·sinψ, sinφ·sinθ·sinψ + cosφ·cosψ, sinφ·cosθ ]
[ cosφ·sinθ·cosψ + sinφ·sinψ, cosφ·sinθ·sinψ - sinφ·cosψ, cosφ·cosθ ]嗯，看着确实有点吓人。但别怕，实际编程时我们一般不会手写这个矩阵，而是用现成的库函数。
DCM的优点：没有奇异点，可以描述任意姿态。而且矩阵运算在数学上非常优雅，组合旋转就是矩阵乘法。
缺点也很明显：9个元素，计算量大。而且随着时间推移，数值误差会导致矩阵不再正交，需要定期做正交化修正。我在做捷联惯导系统时，每迭代几步就要做一次Gram-Schmidt正交化，挺麻烦的。
💡 我的经验：如果你做的是低动态、小角度的应用（比如固定翼巡航），用DCM完全没问题。但如果是高动态、大机动，我建议你直接上四元数。
3.3 四元数：计算高效，但抽象
四元数这东西，我第一次接触时也觉得挺玄乎的。它其实就是一个四维的超复数：
q = q0 + q1·i + q2·j + q3·k其中 i² = j² = k² = i·j·k = -1。
用四元数表示旋转，核心思想是：任何三维空间中的旋转，都可以表示为绕某个单位轴旋转某个角度。
假设旋转轴的单位向量为 u = (ux, uy, uz)，旋转角度为 θ，那么对应的四元数为：
q = [cos(θ/2), u_x·sin(θ/2), u_y·sin(θ/2), u_z·sin(θ/2)]你看，只有4个参数，比DCM的9个少了一半多。而且四元数乘法可以组合旋转，没有奇异点，数值稳定性也好。
我个人的习惯：在飞控算法内部，全部用四元数做姿态解算和更新。只在需要输出给地面站显示时，才把四元数转换成欧拉角。
四元数 vs DCM vs 欧拉角：
- 欧拉角：最直观，但有死锁，适合人机交互
- DCM：数学完美，但计算量大，需要正交化
- 四元数：计算高效，无奇异，但抽象难理解
3.4 三种方法的转换
实际工程中，这三种方法经常需要互相转换。比如从传感器读出来的是四元数，但你要显示给飞行员看欧拉角。
从四元数到欧拉角的转换公式（Z-Y-X顺序）：
φ = atan2(2(q0·q1 + q2·q3), 1 - 2(q1² + q2²))
θ = asin(2(q0·q2 - q3·q1))
ψ = atan2(2(q0·q3 + q1·q2), 1 - 2(q2² + q3²))注意，这里俯仰角用了 asin，当 θ 接近 ±90° 时，依然会有数值问题。所以实际工程中，我一般会用 atan2 的变体来处理边界情况。
⚠ 注意：不同旋转顺序的转换公式不一样。如果你用的是“Z-X-Y”或其他顺序，公式要重新推导。我见过有人直接把别人的代码拷过来用，结果姿态全反了——这种错误很隐蔽，排查起来特别痛苦。
3.5 知识体系总览
下面这张图，是我自己梳理的本章知识结构。你可以看到，三种描述方法各有优劣，但最终都服务于同一个目标：描述飞行器的空间姿态。
嗯，这张图基本把本章的核心逻辑串起来了。三种方法不是互斥的，而是互补的。你在实际项目中，很可能三种方法都会用到。
我个人建议刚入门的同学，先把欧拉角搞透，因为它最直观。等你对姿态有了感性认识，再深入四元数。至于DCM，可以作为一个过渡理解工具。
💡 一个小技巧：在代码里，我习惯把姿态相关的变量都封装成一个结构体，里面同时保留四元数和欧拉角，用标志位标记哪个是最新的。这样既方便计算，又方便调试。
公众号:蓝海资料掘金营，微信deep3321