4、坐标系与运动学：地面坐标系、机体坐标系、桨毂坐标系、倾转坐标系定义，欧拉角与四元数，坐标变换矩阵，运动学方程推导
倾转旋翼机的动力学建模之所以比常规固定翼或直升机复杂，核心在于其构型中存在可倾转的旋翼系统。当旋翼从垂直起降状态（直升机模式）向平飞状态（固定翼模式）过渡时，桨盘平面、旋翼拉力和力矩的方向都在连续变化。如果仅使用单一坐标系描述整机运动，将导致动力学方程中出现大量时变惯性耦合项，且物理意义不清晰。因此，建立一套层次分明、定义严谨的坐标系体系，并掌握其间的变换关系，是后续推导六自由度非线性动力学方程和进行飞控系统设计的前提。本章将系统阐述地面坐标系、机体坐标系、桨毂坐标系与倾转坐标系的定义，对比欧拉角与四元数两种姿态描述方法，推导坐标变换矩阵，并最终完成运动学方程的构建。
上图展示了倾转旋翼机建模中涉及的核心坐标系及其变换逻辑。地面坐标系（E系）作为惯性参考基准，机体坐标系（B系）用于描述飞机本体运动，桨毂坐标系（H系）固连于旋翼旋转轴，而倾转坐标系（T系）则反映了短舱倾转带来的旋翼方位变化。从E系到B系的变换由姿态角决定，从B系到T系由倾转角δ_n决定，从T系到H系则包含旋翼方位角。运动学方程的推导正是基于这些坐标系间的转换关系，将位置、速度、角速度等物理量在不同参考系中正确表达。
4.1 坐标系定义与物理内涵
4.1.1 地面坐标系（E系）
地面坐标系是惯性导航与飞行轨迹描述的基准。其原点O_E通常选取在地面某一点（如起飞点或机场参考点），X_E轴指向地理北向（或任务起始航向），Y_E轴指向东向（或按右手定则确定），Z_E轴垂直于地面向下（或向上，需在建模时统一约定）。在本课程中，为便于与航空惯例衔接，采用北-东-地（NED）顺序，即Z_E轴向下为正。该坐标系用于描述飞行器的位置（经度、纬度、高度）以及地速向量。在仿真中，E系被视为惯性系，忽略地球自转的影响，这对于中低速飞行的倾转旋翼机是合理的简化。
4.1.2 机体坐标系（B系）
机体坐标系与飞行器固连，原点O_B位于飞行器重心（CG）处。X_B轴沿机体结构纵轴指向机头前方，Y_B轴指向机体右侧（右翼方向），Z_B轴垂直于X_BOY_B平面指向机体下方。B系是定义飞行器转动惯量、角速度分量（p, q, r）以及气动力（升力、阻力、侧力）和力矩（滚转、俯仰、偏航）的基准。在倾转旋翼机中，B系还用于定义短舱倾转轴线以及机翼、尾翼等部件的安装位置。
4.1.3 桨毂坐标系（H系）
桨毂坐标系原点O_H位于旋翼桨毂中心（通常也是旋翼拉力作用点），其Z_H轴与旋翼旋转轴重合，指向桨盘上方（即拉力方向的反方向）。X_H轴和Y_H轴位于桨盘平面内，随旋翼旋转或固定于某一方位角（如初始方位角）。H系是计算旋翼气动力、力矩以及桨叶挥舞、摆振运动的自然坐标系。对于倾转旋翼机，左右两侧旋翼的H系定义需对称，且需明确旋转方向（通常左侧旋翼顺时针，右侧旋翼逆时针，从机头方向看）。
4.1.4 倾转坐标系（T系）
倾转坐标系是倾转旋翼机特有的坐标系。原点O_T同样位于桨毂中心，但该坐标系与短舱固连。当短舱倾转角δ_n = 0时（直升机模式），T系与H系重合（或仅差一个旋翼方位角）；当δ_n = 90°时（固定翼模式），T系的X_T轴旋转至与机体X_B轴平行。T系的引入将旋翼的拉力方向与机体姿态解耦，使得旋翼产生的力和力矩可以方便地通过倾转角δ_n投影到B系中。在动力学建模中，T系与B系之间的变换矩阵是倾转角δ_n的函数，这也是倾转旋翼机区别于常规飞行器的关键数学特征。
| 坐标系 | 原点 | X轴 | Y轴 | Z轴 | 主要用途 |
|---|---|---|---|---|---|
| 地面系 (E) | 地面固定点 | 北向 | 东向 | 向下 (NED) | 位置、速度、重力投影 |
| 机体系 (B) | 重心 | 机头方向 | 右翼方向 | 机体下方 | 惯量、角速度、气动力 |
| 桨毂系 (H) | 桨毂中心 | 桨盘平面内 | 桨盘平面内 | 旋翼轴向上 | 旋翼气动力、挥舞运动 |
| 倾转系 (T) | 桨毂中心 | 随短舱倾转 | 随短舱倾转 | 随短舱倾转 | 旋翼拉力方向解耦 |
4.2 姿态描述：欧拉角与四元数
4.2.1 欧拉角定义与旋转顺序
欧拉角是描述机体坐标系相对于地面坐标系姿态最直观的方法。本课程采用航空领域通用的Z-Y-X（偏航-俯仰-滚转）旋转顺序：首先绕Z_E轴旋转偏航角ψ，然后绕新的Y轴旋转俯仰角θ，最后绕新的X轴旋转滚转角φ。这一顺序对应的方向余弦矩阵（DCM）为：
L_EB = R_x(φ) · R_y(θ) · R_z(ψ)其中R_x、R_y、R_z分别为绕对应轴的基本旋转矩阵。欧拉角的物理意义明确，便于理解飞行姿态，但在θ = ±90°时存在万向节锁问题，导致方程奇异。在倾转旋翼机的过渡飞行阶段，俯仰角变化范围较大，可能接近90°，因此单纯依赖欧拉角进行全包线仿真存在风险。
4.2.2 四元数定义与优势
四元数是一种超复数，由一个实部和三个虚部组成：q = q0 + q1·i + q2·j + q3·k，满足约束条件q0² + q1² + q2² + q3² = 1。四元数描述姿态时无奇异点，且计算效率高（仅需乘法和加法，无三角函数）。从欧拉角到四元数的转换关系为：
q0 = cos(φ/2)cos(θ/2)cos(ψ/2) + sin(φ/2)sin(θ/2)sin(ψ/2)
q1 = sin(φ/2)cos(θ/2)cos(ψ/2) - cos(φ/2)sin(θ/2)sin(ψ/2)
q2 = cos(φ/2)sin(θ/2)cos(ψ/2) + sin(φ/2)cos(θ/2)sin(ψ/2)
q3 = cos(φ/2)cos(θ/2)sin(ψ/2) - sin(φ/2)sin(θ/2)cos(ψ/2)在仿真中，推荐使用四元数进行姿态更新，仅在需要输出可视化姿态或与飞控逻辑交互时，才将四元数转换为欧拉角。
q_dot = 0.5 · q ⊗ ω，其中ω为机体角速度四元数（实部为0）。
4.3 坐标变换矩阵构建
4.3.1 地面系到机体系（L_EB）
根据上述欧拉角旋转顺序，L_EB的完整表达式为：
L_EB =
[ 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θ ]该矩阵将E系中的向量变换到B系。其逆矩阵L_BE = L_EB^T，用于将B系中的向量（如重力在B系中的分量）变换回E系。
4.3.2 机体系到倾转系（L_BT）
倾转坐标系相对于机体坐标系仅绕Y_B轴旋转了倾转角δ_n（定义δ_n = 0为直升机模式，δ_n = 90°为固定翼模式）。因此变换矩阵为：
L_BT =
[ cosδ_n, 0, -sinδ_n ]
[ 0, 1, 0 ]
[ sinδ_n, 0, cosδ_n ]该矩阵将B系中的向量变换到T系。注意，当δ_n变化时，该矩阵是时变的，其导数在动力学方程中会产生附加耦合项。
4.3.3 倾转系到桨毂系（L_TH）
桨毂坐标系与倾转坐标系的差异在于旋翼的方位角（即旋翼旋转角）Ωt。若忽略旋翼旋转对平均气动力的影响（即采用叶素理论的平均值），在配平计算中可认为L_TH为单位阵。但在分析旋翼交变载荷或进行振动分析时，需引入方位角：
L_TH =
[ cosΩt, -sinΩt, 0 ]
[ sinΩt, cosΩt, 0 ]
[ 0, 0, 1 ]在整机动力学建模中，通常将旋翼力/力矩先在H系中计算，然后通过L_TH^T变换到T系，再通过L_BT^T变换到B系，最后通过L_EB^T变换到E系。
4.4 运动学方程推导
4.4.1 位置运动学方程
位置运动学描述了地面坐标系中飞行器位置向量P_E = [x_E, y_E, z_E]^T的变化率与机体速度V_B = [u, v, w]^T之间的关系：
P_dot_E = L_EB^T · V_B展开后得到：
x_dot = u·cosθ·cosψ + v·(sinφ·sinθ·cosψ - cosφ·sinψ) + w·(cosφ·sinθ·cosψ + sinφ·sinψ)
y_dot = u·cosθ·sinψ + v·(sinφ·sinθ·sinψ + cosφ·cosψ) + w·(cosφ·sinθ·sinψ - sinφ·cosψ)
z_dot = -u·sinθ + v·sinφ·cosθ + w·cosφ·cosθ注意，由于Z_E轴向下为正，z_dot为负值表示飞行器爬升。
4.4.2 姿态运动学方程（欧拉角形式）
姿态运动学描述了欧拉角变化率与机体角速度ω_B = [p, q, r]^T之间的关系：
[ φ_dot ] [ 1, sinφ·tanθ, cosφ·tanθ ] [ p ]
[ θ_dot ] = [ 0, cosφ, -sinφ ] [ q ]
[ ψ_dot ] [ 0, sinφ/cosθ, cosφ/cosθ ] [ r ]该方程在θ接近±90°时出现奇异（tanθ和1/cosθ趋于无穷），这正是欧拉角法的局限。
4.4.3 姿态运动学方程（四元数形式）
采用四元数时，姿态更新方程为：
q0_dot = -0.5·(q1·p + q2·q + q3·r)
q1_dot = 0.5·(q0·p + q2·r - q3·q)
q2_dot = 0.5·(q0·q + q3·p - q1·r)
q3_dot = 0.5·(q0·r + q1·q - q2·p)该方程无奇异，且为线性形式，便于数值积分。积分后需对四元数进行归一化处理，防止因数值误差导致模长偏离1。
4.5 案例拆解：倾转过渡段的姿态解算
背景与约束：某倾转旋翼机在从直升机模式向固定翼模式转换过程中，短舱倾转角δ_n从0°线性增加至90°，总耗时20秒。飞行器在此过程中保持高度不变，空速从10 m/s加速至60 m/s。飞控系统要求实时输出姿态角，且俯仰角不能超过±30°。
初始现象与目标：在过渡段初期（δ_n < 30°），旋翼拉力方向接近垂直，提供主要升力；随着δ_n增大，机翼升力逐渐承担更多重量。此时若仅使用欧拉角进行姿态解算，在俯仰角接近25°且伴随滚转扰动时，欧拉角速率方程中的tanθ项迅速增大，导致数值积分步长必须缩小，否则姿态发散。
关键输入数据：仿真步长Δt = 0.01s，初始姿态角（φ=0°, θ=5°, ψ=0°），初始角速度（p=0.5°/s, q=1°/s, r=0°/s），倾转角速率δ_dot = 4.5°/s。
分析判断过程：采用欧拉角法时，在t=10s（θ≈15°）时，方程中的tanθ≈0.268，尚可接受；但在t=15s（θ≈25°）时，tanθ≈0.466，且若同时存在滚转扰动φ=10°，则sinφ·tanθ≈0.081，耦合项开始显著影响俯仰角速率计算。若继续增大θ至30°，tanθ≈0.577，数值稳定性变差。而采用四元数法时，所有计算均为线性乘加，无奇异风险。
实际操作步骤：
- 在仿真初始化时，将初始欧拉角转换为四元数q0, q1, q2, q3。
- 在每个仿真步长内，读取机体角速度(p, q, r)（由动力学方程解算得出）。
- 使用四元数更新方程计算q_dot，并采用四阶龙格-库塔法积分得到新的四元数。
- 对四元数进行归一化：q = q / ||q||。
- 将归一化后的四元数转换为欧拉角用于飞控显示：φ = atan2(2(q0q1+q2q3), 1-2(q1²+q2²))，θ = asin(2(q0q2-q3q1))，ψ = atan2(2(q0q3+q1q2), 1-2(q2²+q3²))。
结果对比：在相同仿真条件下，欧拉角法在θ接近30°时出现高频振荡，需要将步长缩小至0.001s才能稳定；而四元数法在0.01s步长下全程稳定，姿态解算误差小于0.01°。此外，四元数法还避免了在θ=±90°时方程奇异导致的仿真中断风险。
风险与代价：四元数法的代价是物理直观性较差，飞控工程师难以直接通过四元数判断飞行姿态；且四元数积分后必须归一化，否则模长漂移会导致姿态误差累积。因此，在工程实践中，通常采用四元数进行内部解算，但输出端仍转换为欧拉角。
复盘经验：倾转旋翼机的过渡段飞行包线宽广，姿态变化剧烈，强烈建议采用四元数作为姿态解算的基准。同时，在动力学方程中，需将倾转角速率δ_dot作为状态变量引入，并考虑其对旋翼角动量的影响，否则在快速倾转时会产生未建模的陀螺力矩，导致仿真结果偏离实际。
4.6 常见误区与总结
- 误区一：混淆坐标系原点。桨毂坐标系原点在桨毂中心，而机体坐标系原点在重心。当旋翼拉力作用点与重心不重合时，会产生力矩。若将H系原点误设为重心，将导致力矩计算错误。
- 误区二：忽略倾转角速率耦合。在建立动力学方程时，若仅将L_BT视为常数矩阵，而忽略δ_dot的影响，则无法正确模拟短舱倾转过程中的惯性耦合效应。
- 误区三：欧拉角万向节锁。在倾转旋翼机进行大机动飞行或垂直起降阶段，俯仰角可能接近90°，此时欧拉角方程奇异。应使用四元数或旋转矢量法。
- 误区四：变换矩阵方向错误。L_EB表示将E系向量变换到B系，而L_BE是其转置。在推导速度方程时，应使用L_EB^T将B系速度变换到E系，方向反了会导致位置积分错误。