控制系统辨识-扫频测试

嵌入式调参中有一个很常见的误区：能画出 Bode 图，就等于已经得到了被控对象模型。

实际上，扫频测试、频率响应分析和系统辨识是前后衔接、但并不相同的三件事。扫频负责激励系统，FRA 负责从输入输出中估计复频率响应，系统辨识则进一步把离散的频响点拟合成传递函数或状态空间模型。只有把这三层分清楚，后续的 Simulink 仿真、补偿器设计和嵌入式部署才不容易走偏。

这篇文章尝试把整条工程链路串起来：

定义测量通道 → 设计安全激励 → 板端同步采集 → 估计 FRF/Bode → 拟合并验证模型 → Simulink 设计 → 离散化部署 → 按正确通道复测

1. 扫频、FRA 与系统辨识分别在做什么？

名称	本质	主要产物
扫频测试（Sweep Test）	向系统注入频率随时间变化或逐点变化的已知激励	同步的输入、输出与状态数据
频率响应分析（FRA）	估计不同频率下输出相对输入的复数增益	FRF/FRD、Bode 图、相干性与置信信息
系统辨识（System Identification）	选择模型结构并利用数据估计参数，再用独立数据验证	传递函数、状态空间或其他参数模型

对于稳定、近似线性时不变的系统，单频正弦进入稳态后，输出仍以同一频率为主，变化体现在幅值和相位上：

$$ u(t)=A\sin(\omega t),\qquad y(t)=B\sin(\omega t+\varphi) $$

于是该频率下的响应可以写成：

$$ G(j\omega)=\frac{Y(j\omega)}{U(j\omega)} $$

扫过一组频率，便得到幅频和相频曲线。它们可以揭示有效带宽、相位滞后、共振与反共振、执行器或传感器滚降，以及不同工况下动态特性的变化。

但 FRF 仍然只是一组非参数频响点。它不天然对应唯一的传递函数，也不能直接完成普通的阶跃、饱和恢复等时域仿真。系统辨识的价值，就是用尽可能简单、又足以解释目标频段的模型，把频域证据转换成可供时域分析、观测器设计和控制器实现使用的数学对象。如果只做频域环路整形，也可以不急着拟合参数模型，直接使用 FRD 数据。

2. 什么样的系统适合扫频辨识？

经典的小信号 SISO FRA 比较适合以下情况：

条件	基本要求
工作点	测试期间基本固定，或变化速度远慢于扫频过程
线性度	小幅激励下近似线性，多次测试能够重复
稳定性	对象本身稳定，或能由已知闭环可靠稳定
可观测性	注入量和响应量都能测量，单位、符号与测量位置明确
信噪比	激励能够超过噪声、量化和背景扰动，又不触发非线性
安全性	具备独立的电流、速度、行程、温升、急停和异常退出限制
耦合性	单一通道占主导；强耦合对象应转向 MIMO FRA 或先做解耦

伺服轴、柔性传动、机器人关节、执行机构、电机试验台以及部分电力电子控制环路，在满足小信号、稳定性与可测量性条件时，常可采用 FRA 分析其局部动态。

相对地，以下对象不适合直接套用一次 SISO Chirp 小信号辨识：

饱和、死区、摩擦、迟滞、间隙或碰撞占主导；
工作点快速变化，参数在测试过程中明显漂移；
开环不稳定，又没有可靠的闭环保护与退出路径；
多通道耦合很强，却只激励和观察其中一个通道；
测试可能触发电气、机械或人员安全风险；
时间尺度极慢，逐频扫完的成本远高于阶跃、PRBS 或多正弦试验。

这里的“不适合”并不等于“完全不能辨识”。非线性或时变系统可以考虑多工作点局部线性化、增益调度、非线性模型或专门的闭环/MIMO 辨识方法，只是不能再把一条 Bode 曲线当成全工况下唯一的系统真相。

3. 扫频之前，先定义到底在测哪一段系统

注入点决定测量结果。参考输入到输出、执行器前扰动到输出、实际对象输入到输出，得到的并不是同一个传递关系。

扫频注入点、测量边界与适用条件

图 1：闭环中的扫频注入、实际输入与输出测量关系。图片为通用示意，不对应具体设备。

若控制器为 $C$、对象为 $P$、测量环节为 $H$，则有：

$$ L=CPH,\qquad S=\frac{1}{1+L},\qquad T=\frac{L}{1+L} $$

若把输出定义为传感器反馈量 $Hy$，则参考到反馈量为 $T=L/(1+L)$；若输出指对象的物理输出 $y$，则 $r\rightarrow y$ 为 $CP/(1+CPH)$；
在对象前注入扰动 $d$，从 $d$ 到输出得到的往往与 $P/(1+L)$ 有关，并不是裸对象 $P$；
若同步记录对象真正收到的输入 $u$ 与对象输出 $y$，才更接近对象边界上的 $Y/U$。

因此，嵌入式端最好同时保存 excitation_cmd、经过叠加和限幅后的 plant_input，以及对应的 plant_output。只记录“计划注入的正弦量”，一旦中途发生限幅或内环没有完全跟随，后续计算就会把错误的输入当成真实输入。

闭环数据还有一个额外问题：反馈可能使对象输入与噪声、扰动相关，简单的谱相除会产生偏差。专业的闭环辨识必须结合正确的注入点、已知环路结构和合适的估计方法，不能只看公式表面上的 $Y/U$。

4. 扫频信号如何选择？

激励类型	优点	局限	常见用途
分段定频正弦（Stepped Sine）	每个频点都能等待稳态，幅相精度高	测试时间较长	精细 FRA、共振附近复核
对数 Chirp	跨越多个数量级时频率分配合理	需要处理瞬态、泄漏和扫速	宽频段试探性扫描
线性 Chirp	生成逻辑简单	跨倍频程时低频和高频分配不均	频率跨度较窄的测试
多正弦 / PRBS	可同时激励多个频率，效率高	峰均比、谱设计与非线性判断更复杂	成熟的自动化测试系统

频率范围不应照搬一个固定数字：

$f_{min}$ 要覆盖最慢的关注动态，并给低频留出足够周期；
$f_{max}$ 要越过目标交越频率和关注模态，同时受执行器、传感器、抗混叠滤波器、采样率与数据链路限制；
Nyquist 频率只是数学上的绝对边界，不是推荐测试上限，实际可信频段通常要低得多；
弱阻尼模态需要建立响应的时间，扫得过快可能把峰值和相位“抹平”。

振幅选择可以概括成一句话：

大到能够超过噪声，小到不触发饱和、限位和显著非线性。

工程上应先做小幅、窄频、短时的 pilot sweep，再逐步扩展。若不同幅值下归一化后的 FRF 明显变化，说明局部线性假设已经开始失效。低相干性也不能简单归因于“振幅太小”：外部扰动、频谱泄漏、响应零点、时间变化、非线性和错误同步都可能造成相干性下降。

5. 嵌入式 C 代码怎样生成扫频？

连续扫频的核心不是直接计算 sin(2π × 当前频率 × 总时间)，而是对瞬时频率做相位积分：

$$ u_{exc}[k]=A(f_k)w[k]\sin\phi[k] $$

$$ \phi[k+1]=\phi[k]+2\pi f_kT_s $$

对数扫频可以在初始化阶段计算频率倍率：

$$ r=\left(\frac{f_1}{f_0}\right)^{1/(N-1)},\qquad f_{k+1}=rf_k $$

下面是一段只表达核心逻辑的 C 代码：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19


typedef struct {
    float phase;
    float frequency;
    float frequency_ratio;
    uint32_t index;
    uint32_t total;
} FraSweep;

float fra_step(FraSweep *s, float amplitude, float ts)
{
    const float envelope = fra_envelope(s->index, s->total);
    const float excitation = amplitude * envelope * sinf(s->phase);

    s->phase = fmodf(s->phase + 2.0f * PI * s->frequency * ts,
                     2.0f * PI);
    s->frequency *= s->frequency_ratio;
    s->index++;
    return excitation;
}

真正接入控制板时，还应补上这些工程细节：

用淡入、淡出包络避免开始和结束处的突变形成宽带冲击；
用 IDLE → RAMP_IN → SWEEP → RAMP_OUT → DONE/ABORT 状态机管理过程；
powf、expf 等尽量在初始化阶段完成，算力紧张时可使用 DDS 或查表；
采样、控制计算和数据入队保持确定周期；
不要在高优先级控制中断里逐点 printf，应使用环形缓冲区，并在低优先级任务或 DMA 中发送；
每帧至少保留 seq、板端时间戳、频率、请求激励、实际输入、输出和状态标志；
原型阶段可以使用 CSV 文本帧，高数据率时应换成带版本、序号和校验的二进制协议。

电脑收到串口数据的时间不能替代板端采样时间。网络、USB 和串口缓存都会引入不确定延迟，真正用于幅相计算的时间基准必须来自控制板。

6. 扫频数据怎样形成 Bode 图？

从扫频数据到Bode图的处理流程

图 2：从激励设计、同步采集到 FRF/Bode 与质量信息的处理链。图片为通用示意。

对于适合宽带谱估计的数据，常用 H1 估计器：

$$ \hat G_{H1}(f)=\frac{S_{yu}(f)}{S_{uu}(f)},\qquad S_{yu}(f)=\mathbb{E}!\left[Y(f)U^*(f)\right] $$

这里采用 $S_{yu}=\mathbb{E}[YU^*]$ 的记号。SciPy 的 csd(x, y) 按 $\mathbb{E}[X^*Y]$ 定义，因此代码中仍调用 csd(u, y)。

然后计算：

$$ M(f)=20\log_{10}|\hat G(f)|,\qquad \Phi(f)=\operatorname{unwrap}(\angle\hat G(f)) $$

幅值平方相干性可作为数据质量指标之一：

$$ \gamma^2(f)=\frac{|S_{yu}(f)|^2}{S_{uu}(f)S_{yy}(f)} $$

MATLAB 中可以使用 tfestimate 与 mscohere：

1
2
3
4
5


[H, f] = tfestimate(u, y, win, overlap, nfft, Fs);
coh = mscohere(u, y, win, overlap, nfft, Fs);

mag_db = mag2db(abs(H));
phase_deg = rad2deg(unwrap(angle(H)));

Python/SciPy 中可以使用互功率谱和 Welch 功率谱：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13


from scipy import signal
import numpy as np

f, Syu = signal.csd(u, y, fs=fs, window="hann",
                    nperseg=nseg, noverlap=noverlap)
_, Suu = signal.welch(u, fs=fs, window="hann",
                      nperseg=nseg, noverlap=noverlap)

H = Syu / Suu
_, coherence = signal.coherence(u, y, fs=fs,
                                nperseg=nseg, noverlap=noverlap)
mag_db = 20 * np.log10(np.maximum(np.abs(H), 1e-12))
phase_deg = np.rad2deg(np.unwrap(np.angle(H)))

连续 Chirp 是非平稳信号，不能机械地认为“做一次 FFT 再相除”就一定可靠。精细测试可以采用分频处理，或者使用分段定频正弦，在每个频点舍弃过渡段后用同步解调、正余弦最小二乘估计输入与输出的复相量。

一份可信的分析结果至少应同时保留：原始数据、采样率与单位、频率响应复数值、幅值、相位、相干性、重复性、异常频点和处理配置。曲线能画出来，与曲线值得相信，仍然是两回事。

7. 从 FRF 到 Simulink 被控对象模型

FRD/IDFRD 可以直接用于 Bode 分析、环路整形和参数模型拟合，但要做阶跃、扰动、饱和、状态观测等普通时域分析，仍然需要拟合传递函数或状态空间模型。

一个比较稳妥的过程是：

用一组数据拟合尽可能低阶的模型；
保留另一组扫频或时域数据作为验证集；
同时比较幅值、相位、时域输出、残差和延迟；
只要求模型在控制设计关注的频段内可信，不追求用高阶模型拟合全频段噪声；
若测量结果已经包含传感器滤波和延迟，不要在 Simulink 中重复叠加。

MATLAB 中可以把复频响整理成 idfrd，再用 tfest 或 ssest 估计模型：

1
2


frData = idfrd(reshape(H, 1, 1, []), 2*pi*f, Ts);
plant = tfest(frData, pole_count, zero_count);

Simulink 中建议至少保留这些环节：参考输入与预滤波、控制器、补偿环节、饱和与抗积分饱和、被控对象、传感器/采样/计算延迟、输入扰动、测量噪声，以及控制量和状态记录。

时域分析也不应只看一个“好看的阶跃”。还要检查轨迹跟踪、低频扰动、窄带扰动、控制量峰值、饱和恢复、参数变化和噪声放大。

8. 根据频域证据选择补偿结构

频域现象或目标	候选措施	需要付出的代价
可重复的窄带柔性模态	优先结构排查、降低带宽或加入陷波器	频率漂移、附加相位与中心偏差
恒值或慢变低频扰动	积分环节、扰动观测器	相位裕度、积分饱和、模型误差
交越附近相位不足	超前环节或 PID 重新整形	高频噪声和控制量增加
希望提高低频增益并压低交越	滞后环节	带宽下降、瞬态响应变慢
高频测量噪声	低通滤波	新增延迟、降低可达带宽

陷波器不是修复机械松动、装配和结构刚度问题的“橡皮擦”；积分环节提高低频增益的同时会增加相位滞后，并需要抗积分饱和。扰动观测器也不是简单串联一个滤波器，它通常由名义模型逆、$Q$ 滤波器和扰动估计反馈组成，其有效带宽不能超过模型可信的范围。

因此，不应把 Notch、DOB、积分和低通一股脑堆进系统，而应从跟踪、扰动、噪声、控制量和鲁棒性目标反推所需结构。

9. 离散化、C 代码部署与复测

从辨识模型到嵌入式固件的闭环流程

图 3：从 FRF、参数辨识到控制设计、离散化、C 实现与正确通道复测。图片为通用示意。

连续域并不是必经步骤。连续模型便于解释质量、阻尼、固有频率等物理含义；直接辨识离散模型则更贴近数字控制器实际看到的采样系统。若从连续域切回离散域，应明确采样周期和离散化方法：物理对象常结合 ZOH，控制器和滤波器可考虑 Tustin 或关键频率预畸变，但都不能当成固定规则。

以二阶数字环节为例：

$$ H(z)=\frac{b_0+b_1z^{-1}+b_2z^{-2}} {1+a_1z^{-1}+a_2z^{-2}} $$

对应差分方程为：

$$ y[k]=b_0x[k]+b_1x[k-1]+b_2x[k-2]-a_1y[k-1]-a_2y[k-2] $$

部署前需要逐项核对：

设计时的 $T_s$ 与真实控制任务周期一致；
计算、PWM、传感器与通信延迟已经纳入分析；
$a_0$ 归一化、分母符号和系数顺序与 C 实现一致；
float32 或定点量化后，极点与频响仍然满足要求；
状态初始化、复位、限幅、溢出和异常退出行为正确；
使用已知测试向量逐样本对比 MATLAB/Python 与 C 输出。

最后还要注意“复测的对象”。如果前后都用实际对象输入 $u$ 到输出 $y$ 计算 FRF，那么只更换控制器并不会改变物理对象 $P$。要验证补偿是否有效，应比较正确的环路传递、闭环 $T/S$、扰动响应、时域指标、噪声和控制量，而不是看到“对象 Bode 变平了”就宣布成功。

10. Agent_SFRA_Skills：把重复工作串成流水线

本文配套提供一个初版开源项目：

GitHub：a233a2/Agent_SFRA_Skills

它目前以 Codex Skill 的形式组织，当前初版提供：

无 HAL 依赖的嵌入式 C 扫频与遥测模板；
串口采集、文本回放、丢帧与状态标志检查；
Python 的 stepped-sine 复相量估计与 Welch H1 分析；
fra.csv、质量报告和 Bode 图输出；
基础陷波器、泄漏积分器、超前/滞后环节的离散系数与 C 头文件导出；
MATLAB/Simulink 被控对象拟合与闭环模型脚手架；
基准与候选结果的复扫对比。

1

git clone https://github.com/a233a2/Agent_SFRA_Skills.git

这个初版的目标，是让 AI 逐步承担“生成代码—采集数据—估计 FRF—拟合模型—提出补偿—离散化—回灌验证”中的重复工作，同时让每次运行都能留下配置、代码差异、原始数据和分析报告。

但“自动化”不等于让 AI 在没有边界的情况下直接驱动实物。激励范围、烧录、上电、使能和系数回写仍应由操作者确认；硬件限幅、急停、看门狗和异常退出也必须独立于 Agent。当前版本已经用通用合成数据跑通基础链路，真实硬件、不同控制架构和不同工具链还需要后续逐步验证。

写在最后

扫频测试真正有价值的地方，不只是找到一两个共振峰，而是建立一条可复现的证据链：知道测了什么、数据是否可信、模型在哪个频段有效、补偿为什么这样设计，以及回到 MCU 后该如何验证。

参数并不是脱离采样周期、工作点和测量边界的几个孤立数字。把这些上下文一起记录下来，频域调参才能从一次性的经验操作，逐渐变成可以复用、比较和自动化的工程流程。