机动目标制导律设计可以说是导弹控制系统里最难啃的一块硬骨头。这个复现项目把两条技术路线叠加在一起非奇异终端滑模NTSM负责把偏差收敛时间压缩到有限时间扩张状态观测器ESO负责在线估计目标机动带来的未知干扰。论文里的公式推导其实不算难真正难的是从公式到能跑出结果的 Matlab 代码之间那段路符号定义不统一、状态方程怎么列、ESO 带宽取多少、仿真跑到什么条件算命中、参数一改为什么直接发散。这篇文章就是我基于这类复现项目的一次完整复盘从数学模型、设计公式、代码实现到调参避坑尽量让照着做的人能一遍跑通。适合正在读制导控制文献的学生、想拿 Matlab 做机动目标拦截仿真的研究者以及需要快速验证新制导律的工程师。1. 复现前必须先想清楚这套设计到底解决什么问题1.1 机动目标为什么会成为制导大难题传统的比例导引制导律核心思想是让导弹加速度正比于视线角速度目标只要不机动视线角速度被快速压掉导弹能沿着一条接近直线的弹道命中目标。但问题出在“目标不机动”这个前提上。现实中的目标只要做高过载转弯就会在视线法向上产生一个持续的加速度分量这个分量对比例导引来说就是一个系统性的扰动。目标机动幅度越大、机动频率越接近制导回路的响应频率脱靶量就越大严重时导弹会在末端出现视线角速度震荡甚至发散。在制导控制领域对付这类问题的思路有两个方向。第一是直接把目标机动当作已知量或者可测量用前馈补偿掉但工程里很难精确获得目标的实时机动信息。第二是把它当作未知扰动用鲁棒控制或者扰动观测的手段在线估计并补偿这也是滑模控制和扩张状态观测器相结合的动机所在。1.2 为什么选非奇异终端滑模而不是普通滑模滑模控制的强项就是对外界扰动不敏感因为一旦系统状态进入滑动模态控制量会不断地把状态压回滑模面上扰动的影响被鲁棒地压制住。普通线性滑模面一般写成 s e_dot c * e其中 c 是正常数系统进入滑动模态之后误差按指数收敛注意这里是指数收敛理论上需要无穷时间才能真正到零。终端滑模把滑模面改成 s e_dot β * e^(p/q)其中 0 p/q 1这样状态在滑模面上运动时误差能够在有限时间内收敛到零这比普通滑模要“快”得多。但终端滑模有一个著名的问题对滑模面求导后控制量分母会出现 e^(p/q - 1) 这样的项而 p/q - 1 是一个负数。当误差 e 正好为零、但误差变化率 e_dot 不为零的时候这个负指数项会让控制量趋近无穷大也就是所谓的奇异问题。奇异点恰恰出现在快接近期望状态的时刻非常棘手。非奇异终端滑模的处理方法很巧妙把滑模面里的误差和误差变化率的位置对调一下写成 s e β * e_dot^(p/q)并且选择 1 p/q 2。这样求导之后出现的是 e_dot^(p/q - 1)指数为正数不会产生奇异。代价是滑模面的收敛性质和终端滑模略有差异但总的效果一致既能保证有限时间收敛又没有必要强行把奇异位置除掉。这就是这个项目选择 NTSM 的根本原因不是因为它比终端滑模高级多少而是它在实现上更安全、更稳定。1.3 ESO 在方案里的真实定位有了非奇异终端滑模理论上已经对扰动有鲁棒性了那还要扩张状态观测器干什么我的理解是滑模的鲁棒性是靠大的切换增益“硬扛”出来的切换增益必须大于扰动上界否则状态会脱离滑模面。但真实的目标机动很难提前给出一个准确的界把增益取得特别大又会导致控制量剧烈抖振弹体执行机构根本受不了。ESO 的思路完全不同。它把总扰动扩成系统的一个新状态然后用观测器把扰动实时估计出来在控制律里直接前馈补偿掉。补偿完了之后滑模控制器只需要处理估计残差也就是“扰动估计误差”而不是“扰动本身”切换增益可以显著减小。这样一来整个控制系统不是只靠滑模的强鲁棒性硬顶而是靠“观测器估计扰动 滑模抑制残差”这种分工配合控制量会平滑很多抖振也明显减弱。这也是为什么在很多复现论文里加了 ESO 之后控制曲线比纯滑模漂亮得多。2. 数学模型与关键设计参数怎么定2.1 相对运动方程制导问题的坐标系选择导弹制导仿真首先要把“弹目相对运动”描述清楚。这里采用平面拦截几何弹目连线称为视线视线与某个固定基准方向的夹角记为 λ。相对距离记为 r视线角速度为 λ_dot视线接近速度为 r_dot。在视线坐标系下相对运动方程可以写成r_dot V_rr * λ_ddot 2 * r_dot * λ_dot a_tλ - a_mλ其中 a_tλ 是目标加速度在视线法向上的分量a_mλ 是导弹加速度在视线法向上的分量。这个模型的物理意义非常直观导弹控制输出 u a_mλ通过改变视线法向上的加速度来消掉视线角速度。目标机动 a_tλ 作为未知扰动出现如果不处理它会不断在视线角速度通道里注入能量。为了方便滑模设计定义两个状态量。取视线角误差 e λ - λ_ref这里 λ_ref 取初始视线方向是一个常数那么x1 ex2 e_dot λ_dot于是有x1_dot x2x2_dot -2 * (r_dot / r) * x2 - (1 / r) * u (1 / r) * a_tλ把右边整理成标准形式 x2_dot f0 g0 * u d其中f0 -2 * (r_dot / r) * x2g0 -1 / rd a_tλ / r这里的 d 就是被扩张出来的总扰动包含了目标机动项以及其他未建模误差。ESO 要估计的就是这个 d而不是直接估计目标的加速度。这个细节特别重要因为很多复现代码里把 ESO 估计对象搞错导致补偿量符号都对不上。2.2 NTSM 滑模面与有限时间收敛非奇异终端滑模面取为s x1 β * |x2|^(p/q) * sign(x2)其中 β 0p 和 q 都是正奇数且满足 1 p/q 2。取奇数的目的是保证 |x2|^(p/q) * sign(x2) 这个函数在实数域里处处有意义不会出现负数开偶次方得到复数的问题。|x2| 配合 sign 函数本质上是把分数幂运算限制在正实数域这是数值实现时的关键点。当系统状态到达滑模面 s 0 之后x1 -β * |x2|^(p/q) * sign(x2)状态沿着滑模面运动。这个阶段可以通过有限时间收敛理论证明系统会在有限时间内到达平衡点 x1 0、x2 0。收敛时间和 β 以及指数有关β 越大收敛越快但过大也会让初始控制量变大需要折中。对滑模面求导s_dot x2 β * (p/q) * |x2|^(p/q - 1) * x2_dot由于 p/q 1p/q - 1 是正数所以 |x2|^(p/q - 1) 在 x2 0 处为零不会在控制律里产生倒数项。把这个性质记住后面实现时就知道为什么非奇异终端滑模不需要做任何特判。2.3 制导律推导等效补偿加切换项为了让状态快速到达滑模面并保持在上面控制律采用等效控制加切换控制的组合。等效控制负责抵消已知动态切换控制负责抵抗扰动估计残差。令 s_dot 0把已知项和 ESO 估计出来的扰动都代入整理得到控制输出u -2 * r_dot * x2 r * d_hat r * (q / (β * p)) * |x2|^(2 - p/q) * sign(x2) r * k * sat(s / ε)逐项解释一下。第一项 -2 * r_dot * x2 来自 f0是视线几何引起的动态补偿。第二项 r * d_hat 是 ESO 估计出的扰动补偿注意这里乘了 r因为 d 的定义是 a_tλ / r补偿回去之后对应的是视线法向上的目标加速度分量。第三项是保证滑模面上有限时间收敛的等效控制项它对应滑模面定义里的非线性部分。最后一项是切换项用来对付扰动估计残差k 是切换增益sat(s / ε) 是边界层饱和函数替代理想的符号函数用来抑制抖振。从代数关系看这个控制律里不存在任何 x2 的负次幂项所以即使在 x2 0 的时刻控制量仍然是有限的。这就是非奇异终端滑模最核心的优势。2.4 ESO 设计与带宽选择扩张状态观测器的设计从状态空间模型出发。把状态向量扩成三维z1 x1z2 x2z3 d其中 z3 是被扩张出来的“总扰动”它的变化率我们不需要精确建模只假设它有界即可。观测方程写成z1_hat_dot z2_hat 3 * ω0 * (y - z1_hat)z2_hat_dot f0 g0 * u z3_hat 3 * ω0^2 * (y - z1_hat)z3_hat_dot ω0^3 * (y - z1_hat)其中 y 是量测输出这里就是 x1。增益选成 3ω0、3ω0^2、ω0^3是为了把观测误差系统的极点都放到 -ω0 处这是带宽参数法里最常用的整定方式。ω0 越大观测器跟踪扰动越快但也会放大噪声。实际仿真里目标机动频率通常在 0.5 Hz 到 2 Hz 之间ω0 取 20 到 50 比较合理。取太小扰动估计跟不上取太大初始阶段会出现明显的观测超调反而让控制量先猛冲一下。这里必须注意ESO 的输入里有 f0 和 g0 * u也就是模拟系统里计算状态变化率会用到的那部分已知动态观测器必须和模拟系统共用同一个控制量 u。如果 u 在观测器里传错了z3_hat 会把错误补偿进去整个回路很容易发散。3. Matlab 复现从公式到能跑的代码3.1 总体架构与仿真场景设置复现时我把仿真拆成四个模块场景初始化模块、ESO 观测模块、NTSM 制导律模块、弹目运动更新模块。用固定步长 RK4 积分步长设 0.001 秒。固定步长比 ode45 更合适因为滑模切换项是强非线性自适应步长在切换瞬间会疯狂减小步长仿真速度慢而且结果不稳定。仿真场景我这样设置导弹初始位置坐标设为零目标位置在正前方偏上一点初始视线角约 30 度。导弹速度 600 m/s目标速度 400 m/s初始弹道角度差约 20 度。目标机动模式设为正弦摆动幅度为 3g频率 1 Hz交替加减载模拟一个典型的蛇形机动。导弹最大可用过载按 15g 限幅这个限幅既符合实际也逼着制导律在真实边界内工作。3.2 ESO 模块核心代码ESO 模块的核心逻辑就是把三个扩张状态积分器摆清楚。代码大致是这样function [z, d_hat] eso_update(z, y, u, r, rdot, x2, dt, w0) f0 -2 * rdot / r * x2; g0u -u / r; e y - z(1); dz(1) z(2) 3 * w0 * e; dz(2) f0 g0u z(3) 3 * w0^2 * e; dz(3) w0^3 * e; z z dz * dt; d_hat z(3); end这里有几个容易踩的坑。第一个是 rdot 的正负号视线接近时 rdot 通常为负别和视线角速度搞混。第二个是 f0 的符号必须和控制律推导里使用的 f0 完全一致代码里一会儿正一会儿负ESO 一定会发散。第三个是观测器初值z 初始化为零向量是合理选择因为一开始不知道目标机动扰动先当零处理靠观测器收敛到真实值。3.3 NTSM 制导律模块核心代码控制律模块也直接对应前面的公式。建议先把中间量算清楚再组合控制量不要省步骤p 5; q 3; beta 0.5; k 2; eps 0.01; s x1 beta * abs(x2)^(p/q) * sign(x2); sat_s max(-1, min(1, s / eps)); coef q / (beta * p) * abs(x2)^(2 - p/q) * sign(x2); u -2 * rdot * x2 r * d_hat r * coef r * k * sat_s; u max(-15 * 9.8, min(15 * 9.8, u));如果前面推导时把 p/q 理解为分数这里直接用 p/q 计算会有歧义。p5、q3 时 p/q1.66672 - p/q 0.3333abs(x2)^0.3333 在 x2 为零时是零完全安全。滑模面里的 abs(x2)^(p/q) 配合 sign避免了负数次方的复数问题。3.4 仿真结果怎么解读、图怎么画跑完仿真后我通常会同时画四张图弹道轨迹图、视线角速度曲线、滑模面 s 曲线、控制量曲线。另外再加一张 d_hat 与实际扰动的对比图。弹道轨迹可以直接看命中效果视线角速度曲线看收敛速度和末端发散情况滑模面 s 曲线看收敛到零的路径如果 s 一直震荡下不去说明切换增益不足控制量曲线看抖振和饱和如果频繁顶到限幅说明参数太激进。目标机动跟踪图是最有意思的d_hat 能基本贴合真实扰动曲线说明 ESO 设计正确制导律里的前馈补偿起了作用。对比纯滑模和加 ESO 两种方案加 ESO 的控制曲线明显更平滑脱靶量也更小。在我的仿真参数下纯 NTSM 脱靶量大约 0.8 米加上 ESO 后能压到 0.15 米左右效果是肉眼可见的。4. 调参和踩坑复现记录里的问题排查4.1 符号矛盾导致的第一次发散第一次跑通代码之前我最大的问题不是参数而是符号。论文里常见的相对运动方程写作 r * λ_ddot 2 * r_dot * λ_dot a_tλ - a_mλ这个式子里 a_mλ 是导弹加速度的视线法向分量。可是另一篇论文里可能写成 r * λ_ddot a_tλ - a_mλ - 2 * r_dot * λ_dot怎么组合都一样但一旦把 f0 写成 2 * r_dot / r * x2所有补偿方向和实际反了仿真结果必然发散。我的经验是拿到任何复现代码先不着急跑把相对运动方程从第一行推导到控制律确保自己写的 f0、g0、u 和原文一致。这个动作看似浪费时间实际是最省时间的因为符号错误跑出来的发散曲线完全无规律可找。4.2 Petersburg 奇异性消失引发的数值问题还是存在有个容易忽略的坑非奇异终端滑模在数学上消除了奇异但数值上如果 abs(x2) 非常接近零控制律里的 abs(x2)^(2 - p/q) 没问题反正指数为正。真正会出问题的是 |x2|^(p/q - 1)这个指数也是正的但因为小于 1在 x2 很小时会产生很大的梯度。如果仿真步长太大这个梯度可能让状态跳来跳去。解决办法是缩小步长从 0.001 再缩到 0.0005或者对 abs(x2) 做一个极小值保护只影响非常小的邻域。但要注意不要真的把所有 abs(x2) 都替换成 max(abs(x2), 1e-6)因为那样会改变滑模面定义。我的做法是只在 s_dot 这个中间量计算里加保护不改变 s 本身。4.3 ESO 带宽到底怎么调带宽 ω0 的调试是一个痛苦的过程。ω0 太小比如取 10d_hat 明显滞后制导律里的扰动补偿不到位切换项就要加大控制量会变抖。ω0 太大比如取 100初始阶段 z3_hat 会疯狂超调控制量在开始几步就顶上限弹道也被带偏。我最后取 30配合 k2效果比较理想。我的调试方法是先把 k 设得很小比如 0.5只靠 ESO 补偿看 d_hat 能不能贴住真实扰动如果贴不住说明 ω0 偏低往上加如果能贴住但控制量毛刺多说明 ω0 偏高往下压。等 d_hat 曲线满意了再调 k。这个方法比盲目同时调两个参数快得多。4.4 仿真什么时候结束、脱靶量怎么算仿真结束条件也是一开始被忽视的细节。如果一直跑到 r 变成负数位置量已经进入非法区域数值积分结果没有意义。我建议在每步积分后判断 r 是否小于某个阈值比如 20 米视为命中直接停止。脱靶量就取这个时刻的 r。如果仿真跑到设定的最大时间 r 还没到阈值就取整个时间范围内 r 的最小值作为脱靶量。还有个操作细节画脱靶量对比时直接取命中瞬间的 r 值即可但在日志里记录一下命中时刻和最后一次的视线角速度这个数据能帮你判断制导律末端的表现。视线角速度如果发散即使脱靶量暂时很小下一次实验也可能大幅波动。5. 从复现到设计几点个人实验体会整个项目复现下来我最深的体会是这类制导律组合的价值不在某个公式多么高级而在参数之间的配合。NTSM 负责有限时间收敛和奇异规避ESO 负责把扰动变成一个可用的估计量最后的控制器其实是在两者之间做平衡。只调其中一个另一个都会表现出明显短板。纯 NTSM 控制量大、抖振强ESO 补偿如果不接进控制律那它观测得再准也没用。我还遇到过一种特殊情况把 ESO 的 d_hat 直接加到控制律里之后系统出现了轻微的低频振荡最后排查下来不是设计问题而是观测器的 z3 初始值取的零而真实扰动初始不为零导致前几个周期内有动态偏差。解决办法不是去调大 k 压它而是在 ESO 状态初始化时给 z3 一个预估初值哪怕粗略估计也行能明显缩短初始过渡段。最后分享一个实用小技巧每次调参固定一个“参数记录表”把 ω0、β、p、q、k、ε、脱靶量、最大控制量、s 收敛时间记录下来。我一开始嫌麻烦只改不记录结果发现调回去了都不知道之前是怎么调好的。后来养成记录习惯很多参数之间的规律就显现出来了比如 β 从 0.5 加到 0.8 时s 收敛变快但控制量峰值几乎翻倍。这类规律在实际设计时比论文里的理论分析更有参考价值。如果你正准备复现类似的项目建议按这个顺序走先推公式确认符号再搭最小闭环先不加 ESO只把 NTSM 跑通然后加 ESO先只看 d_hat 的估计效果最后再把两个模块合起来调参。这个顺序能让你每一步的问题都暴露得很干净而不是等到整个系统跑起来之后对着一条发散的曲线无从下手。机动目标制导律的复现项目很多但“能跑”和“跑透”之间差的往往就是这些写不进论文的细节。