MATLAB频域分析：频率响应曲线生成教程

MATLAB频域分析实战：手把手教你画出精准频率响应曲线

你有没有遇到过这样的情况？设计了一个滤波器，理论上通带很平、阻带衰减足够，但一上电测试却发现高频信号严重失真；或者调PID控制器时，系统总是莫名其妙地振荡——明明参数看起来“挺合理”。这时候，光看阶跃响应可能已经不够了，你需要的是从频率的角度重新审视这个系统。

这就是我们今天要讲的：如何用MATLAB真正“看懂”一个系统的频率性格。不是简单跑个bode(sys)就完事，而是搞清楚每一步背后的逻辑，让你不仅能画图，还能读懂图、改系统。

从传递函数开始：把微分方程变成可计算的对象

所有频域分析的第一步，都是建模。而对线性系统来说，最直观的模型就是传递函数。

比如你有一个典型的二阶系统，描述它的微分方程可能是：

$$
\ddot{y}(t) + 0.5\dot{y}(t) + y(t) = u(t)
$$

拉普拉斯变换后得到：

$$
G(s) = \frac{Y(s)}{U(s)} = \frac{1}{s^2 + 0.5s + 1}
$$

在MATLAB里，我们要做的就是把这个数学表达式“翻译”成计算机能处理的形式。核心函数是tf()：

num = [1]; % 分子系数：1（对应常数项） den = [1, 0.5, 1]; % 分母系数：s^2 + 0.5s + 1 sys = tf(num, den); disp(sys);

运行结果会显示：

1 ------------------- s^2 + 0.5 s + 1

别小看这三行代码，它已经把你抽象的系统封装成了一个“活”的对象。你可以对它求阶跃响应、零极点、状态空间表示，当然，最重要的——做频率分析。

经验提示：如果你的系统阶数很高（比如 > 6），直接用tf可能会因为数值精度问题导致不稳定。建议改用状态空间模型（ss）或零极点增益形式（zpk），它们在高阶情况下更稳健。

Bode图：你的系统“听力测试报告”

想象一下，给系统输入不同频率的正弦波，记录输出的幅度和相位变化。这就像给音响做频响测试——低音够不够沉？高音有没有刺耳？只不过我们现在测的是控制系统、滤波器或机械结构。

Bode图就是这份“听力报告”的标准格式：两张图，一张是幅频特性（dB vs log频率），一张是相频特性（度 vs log频率）。

一行代码出图？没错，但你要知道它干了啥

figure; bode(sys); grid on; title('二阶系统的Bode图');

就这么一句bode(sys)，MATLAB已经在后台完成了以下工作：

1. 自动选择合适的频率范围（通常覆盖系统动态主要区间）；
2. 在每个频率点 $ \omega $ 上计算复数响应 $ G(j\omega) $；
3. 提取幅值 $ |G(j\omega)| $ 并转换为分贝：$ 20\log_{10}|G(j\omega)| $；
4. 提取相位 $ \angle G(j\omega) $，单位转为角度；
5. 使用对数坐标绘图，突出关键频段细节。

生成的图像中你能看到：

• 低频段增益接近0 dB → 直流增益为1；
• 幅值在约1 rad/s附近略有凸起 → 存在轻微谐振；
• 相位从0°逐渐下降到-180° → 典型二阶系统的相移行为。

这些都不是随机的，而是由极点位置决定的。你可以用pole(sys)查看一下：

>> pole(sys) ans = -0.2500 + 0.9682i -0.2500 - 0.9682i

一对共轭复极点，实部负说明稳定，虚部决定了谐振频率。一切都能对应起来。

想自定义频率点？用freqresp精确控制采样

有时候你不满足于默认的频率划分。比如你想重点观察某个可疑频段是否有共振峰，或者需要将数据导出给嵌入式系统做查表补偿。

这时就得祭出底层函数：freqresp。

它不像bode那样直接画图，而是返回你在指定频率下的精确响应值。

% 定义频率向量：从0.01到100 rad/s，对数均匀分布 w = logspace(-2, 2, 1000); % 1000个点 % 计算频率响应（返回复数） H = freqresp(sys, 1i*w); % 注意：必须传 jω，即 1i*w % 拆解为幅值和相位 mag = abs(H); % 幅值 phase = angle(H) * 180/pi; % 弧度转角度 % 手动绘制 figure; subplot(2,1,1); semilogx(w, 20*log10(mag), 'b', 'LineWidth', 1.2); ylabel('Magnitude (dB)'); title('手动绘制的频率响应'); grid on; subplot(2,1,2); semilogx(w, phase, 'r', 'LineWidth', 1.2); ylabel('Phase (deg)'); xlabel('Frequency (rad/s)'); grid on;

这段代码的关键在于：

• freqresp(sys, 1i*w)中的1i*w是复频率 $ j\omega $，不能漏掉虚数单位；
• 输出H是一个三维数组（即使SISO系统也是[1,1,N]），实际使用时可以直接当作向量用；
• 你可以自由选择非均匀采样，比如在谐振区加密采样点，在平坦区稀疏些，节省计算资源。

实用技巧：如果你要做硬件部署前的数据预处理，可以把w,mag,phase导出为.csv或.h文件供C程序加载：

matlab data = [w', 20*log10(mag)', phase']; writematrix(data, 'fr_data.csv');

多变量系统怎么分析？MIMO不是难题

现实中的系统往往不止一个输入一个输出。比如飞行器的姿态控制：三个角速率传感器、三个舵面执行机构，彼此之间还有强耦合。

这种多输入多输出（MIMO）系统的频率响应是一个矩阵：

$$
\mathbf{G}(j\omega) =
\begin{bmatrix}
G_{11}(j\omega) & G_{12}(j\omega) \
G_{21}(j\omega) & G_{22}(j\omega)
\end{bmatrix}
$$

每个元素代表某个输入对某个输出的影响。

在MATLAB中构建非常自然：

s = tf('s'); % 定义连续时间变量s % 构造四个通道的传递函数 G11 = 1 / (s + 1); G12 = 0.5 / (s + 2); G21 = 0.3 / (s + 1.5); G22 = 1 / (s^2 + s + 1); % 组合成2x2 MIMO系统 sys_mimo = [G11, G12; G21, G22];

然后直接调用bode(sys_mimo)，MATLAB会自动为你画出所有输入-输出组合的Bode图，一共 $ m \times n $ 组。

你会发现有些通道响应快、增益大，有些则慢且弱。如果交叉项（如G12）也很强，说明存在显著耦合，可能需要解耦控制策略。

工程洞察：当你看到某两个通道的相位曲线在关键频段相差接近180°，而幅值又相近时，就要警惕闭环下可能发生抵消或共振。这种情况在机械臂或多电机同步系统中尤为常见。

稳定吗？让 margin 告诉你离崩溃还有多远

画完Bode图，下一个问题是：这个开环系统闭环后稳不稳定？

答案藏在两个指标里：增益裕度（Gain Margin, GM）和相位裕度（Phase Margin, PM）。

它们的本质来自奈奎斯特判据，但在Bode图上更容易读取：

• 增益裕度：当相位降到 -180° 时，增益还差多少dB才到0dB；
• 相位裕度：当增益穿过0dB时，相位还差多少度才到 -180°。

经验值告诉我们：

• PM > 45°：响应平稳，超调小；
• PM < 30°：很可能出现大幅振荡；
• GM > 6 dB：允许增益波动而不失稳。

MATLAB提供了一个神器函数：margin()，它不仅能计算数值，还能在图上标出关键点：

figure; margin(sys); % 自动生成带裕度标注的Bode图 title('含稳定性裕度的Bode图');

同时获取具体数值用于分析或报告：

[Gm, Pm, Wcg, Wcp] = margin(sys); fprintf('增益裕度: %.2f dB (@ %.3f rad/s)\n', 20*log10(Gm), Wcg); fprintf('相位裕度: %.2f° (@ %.3f rad/s)\n', Pm, Wcp);

输出示例：

增益裕度: Inf dB (@ NaN rad/s) 相位裕度: 65.53° (@ 0.937 rad/s)

这里GM为无穷大，是因为该系统相位从未达到-180°，属于“超强稳定型”。但对于某些带延迟或高阶滤波的系统，这两个值往往会变得紧张，这时候你就得回头调整控制器了。

实际应用中的那些“坑”和应对策略

理论再完美，落地总有意外。以下是我在项目中踩过的几个典型坑，以及MATLAB里的解决办法：

❌ 坑1：频率范围没选好，漏掉了谐振峰

新手常犯的错误是依赖默认频率范围，结果在10kHz处有个尖锐共振，却被画图算法跳过了。

✅对策：手动指定频率向量，确保覆盖关心的频段：

w = logspace(0, 4, 2000); % 覆盖1 rad/s 到 10krad/s bode(sys, w);

❌ 坑2：单位混淆，Hz 和 rad/s 搞反了

很多传感器手册给的是Hz，但MATLAB默认用rad/s。忘了转换会导致整整 $ 2\pi $ 倍的偏差！

✅对策：明确标注单位，必要时转换：

f_Hz = logspace(1, 4, 1000); % 频率单位：Hz w_rad = 2*pi*f_Hz; % 转换为 rad/s H = freqresp(sys, 1i*w_rad);

❌ 坑3：离散系统当成连续系统分析

数字控制器本质是离散的，直接拿s域模型分析会误导。

✅对策：先用c2d转成离散模型再分析：

Ts = 0.01; % 采样周期10ms sys_d = c2d(sys, Ts, 'zoh'); % 零阶保持法离散化 bode(sys_d); % 分析离散系统频率响应

注意：离散系统的有效频率只到 $ \pi/T_s $（奈奎斯特频率），超过无意义。

❌ 坑4：非线性系统强行套用线性分析

饱和、死区、滞环……这些非线性环节会让频率响应随输入幅值变化。

✅对策：结合扫频仿真（chirp signal）+ FFT估算FRF（频率响应函数）：

% 生成扫频信号（chirp） t = 0:0.01:10; u = chirp(t, 0.1, 10, 10); % 从0.1Hz扫到10Hz % 仿真非线性系统响应（需建Simulink模型或自定义函数） % y = sim('nonlinear_system', t, u); % 对输入输出做FFT，估算FRF % H_est = fft(y)./fft(u);

这种方法更适合真实系统辨识，也适用于验证线性化模型的有效性。

写在最后：掌握频域思维，才能真正驾驭系统

这篇文章没有停留在“命令+截图”的层面，而是试图带你理解每一次bode()背后发生了什么，每一个裕度值意味着什么。

当你下次面对一个不稳定的控制系统时，不要急着调PID参数。先问自己几个问题：

• 它的带宽够吗？
• 在穿越频率附近，相位下降得太快了吗？
• 是否存在隐藏的谐振模式被激发？

这些问题的答案，都在频率响应图里。

随着AI与系统辨识的发展，MATLAB也在不断融合新的工具，比如基于实测数据的FRF估计、与Simulink联合的实时频率扫描等。未来甚至可能实现“自动诊断+推荐补偿器”的智能辅助设计。

但无论技术怎么变，理解频率响应的本质，始终是你作为工程师的核心竞争力。

如果你正在做滤波器设计、控制器调试或振动分析，不妨现在就打开MATLAB，把你手头的系统跑一遍bode，看看它的“声音”是什么样的。也许你会发现一些之前从未注意到的秘密。

欢迎在评论区分享你的频响图和发现！