Matlab频谱与Bode图绘制：从FFT到频响验证的完整指南

简介：面向信号处理、通信工程与控制系统领域的初学者及工程师，提供一套实用的MATLAB频谱分析与Bode图绘制解决方案，专门解决两大高频需求：一是利用pwelch函数完成功率谱密度估计，涵盖数据读取、滤波去噪、窗函数选择与频率分辨率设定，最终输出清晰的频率-功率图谱；二是围绕bode及bodeplot函数展开系统频率响应分析，支持传递函数、状态空间或零极点增益模型的定义、频率范围调整以及幅频和相频响应的可视化。压缩包为zip格式，整体大小仅1KB，内含2个m脚本：一个演示从原始数据到频谱图的完整处理流程，另一个用于生成系统Bode图，并配有详细的参数设置说明与典型调用方式。资源轻量精炼，查阅方便，截至目前已有1648人学习下载，适用于课程实验、滤波器设计、噪声分析与控制系统稳定性验证等场景；通过运行和修改这两个脚本，读者可以快速掌握MATLAB频域分析的核心函数使用技巧，并直接迁移到自己的项目中。

1. 频谱绘制和 Bode 图绘制，Matlab 里其实是同一套功夫

拿到一段振动信号，先画频谱找峰值频率；再拿这个峰值去对比放大器或者滤波器的 Bode 图，看它是不是落在通带边缘。这两件事在 Matlab 里经常被分开写脚本，其实底层逻辑是同一个：把一个时域信号或系统模型从时间域换到频率域，然后按幅值和相位两个维度去读它。区别只在于，频谱图的横轴是物理频率，纵轴是幅值谱或功率谱；Bode 图的横轴是对数频率，纵轴是幅值（dB）和相位（度）。这篇内容覆盖从fft基础画法到bode函数调用，再到把实测频谱和理论 Bode 图叠在一起做验证的完整路径。适合正在做信号处理、控制系统分析、振动噪声测试的工程师，也适合刚把 Matlab 装好、想从第一行命令开始跑通全流程的初学者。读完可以把手头数据直接换成自己的采样率和信号长度。

2. 用 FFT 在 Matlab 里绘制频谱图：最小代码与频率轴映射

2.1 先理解fft输出的是什么

Matlab 里fft(x)返回的是离散傅里叶变换（DFT）的复数序列，长度和输入信号x相同。这个复数的模对应频率分量的幅值，幅角对应相位。但直接对fft结果画图毫无意义，因为频率轴需要自己构造，幅值也要经过换算。

Fs = 1000; % 采样率 1000 Hz T = 1/Fs; % 采样间隔 L = 2000; % 信号长度 2000 点 t = (0:L-1)*T; % 时间向量 x = 0.7*sin(2*pi*50*t) + 1.2*sin(2*pi*120*t); % 两个正弦叠加 Y = fft(x); P2 = abs(Y/L); % 双边频谱幅值 P1 = P2(1:L/2+1); % 取单边，只保留正频率部分 P1(2:end-1) = 2*P1(2:end-1); % 除了 0 频（直流）和 Nyquist 点，其他乘 2 f = Fs*(0:(L/2))/L; % 频率轴，单位 Hz plot(f, P1) xlabel('频率 (Hz)') ylabel('|X(f)|') grid on

这里Y = fft(x)得到复数谱，P2 = abs(Y/L)把复数的模除以信号长度，得到的是每个频率分量的真实幅值。因为是双边谱，正负频率各分一半能量，所以取单边后要把除直流和 Nyquist 频率外的幅值乘 2。f = Fs*(0:(L/2))/L构造了从 0 到 Fs/2 的频率向量，长度和P1完全一致，这样才能正确对应绘图。

实际工程里，Fs是采样设备直接给出的，L是你截取的数据长度，这两个参数决定了频率轴的刻度范围和频率分辨率。常见做法是先检查length(f)和length(P1)是否相等，不等会直接导致绘图报错。很多第一次做频谱分析的人直接把plot(abs(fft(x)))画出来，横轴是“点数”不是频率，这是最典型的错误。

2.2 频率分辨率、补零与窗函数

FFT 的频率分辨率由Δf = Fs / N决定，N 是参与 FFT 的点数。提高分辨率只有一条路：增加采样时长，也就是让L变大。补零到更大长度不会真正改善分辨率，只是对现有频谱做插值，曲线更平滑，但两个频率很接近的峰依然分不开。

如果信号不是整周期截断，频谱会泄漏。加窗是最常见的抑制手段。

N = 4096; w = hann(L, 'periodic'); % 生成汉宁窗 xw = x(:) .* w; % 信号加窗 Yw = fft(xw, N); % 补零到 4096 点再做 FFT P2w = abs(Yw / sum(w)); % 注意，幅值修正用 sum(w) 而不是 L P1w = P2w(1:N/2+1); P1w(2:end-1) = 2*P1w(2:end-1); fw = Fs*(0:N/2)/N; plot(fw, P1w)

加窗后幅值修正公式从abs(Y/L)改为abs(Y / sum(w))，因为窗函数会让信号能量衰减。sum(w)在 Matlab 里可以直接算，不需要手查表格。hann(L, 'periodic')和hann(L)的区别在于周期性窗更适合 FFT，谱泄漏更小，这是我一般会优先采用'periodic'形式的原因。

2.3 分段频谱图：长时间信号的平均化处理

现场采集的一段振动数据往往长达几十秒，全段做一次 FFT 看不到随时间变化的频率特征。分段频谱图（也叫声谱图）是解决这个问题的常规做法。

segmentLength = 2048; overlap = 0.5; % 50% 重叠 nfft = 2048; [s, f_spec, t_spec] = spectrogram(x, hann(segmentLength), ... round(overlap*segmentLength), nfft, Fs); imagesc(t_spec, f_spec, 20*log10(abs(s) + eps)); axis xy; xlabel('时间 (s)'); ylabel('频率 (Hz)'); title('分段频谱图 (Spectrogram)'); colorbar;

spectrogram返回三个量：s是复数谱矩阵，f_spec是频率向量，t_spec是每个时间段的中心时刻。imagesc把矩阵画成彩色图，纵轴是频率，横轴是时间，颜色的深浅代表该时刻该频率上的能量。20*log10 是把幅值转成 dB，+eps防止 log 里出现 0。axis xy让 y 轴方向从下往上递增，符合直觉。

重叠率是这里最值得调的参数。重叠越多，时间轴上的变化越平滑，但计算量也越大。50% 重叠是默认偏好，既能看到瞬态变化，又不会让计算慢到影响交互式分析。segmentLength和nfft取相同值时，频率分辨率等于Fs/segmentLength，在喜欢分辨率的场合适当增大 segmentLength；在追求时间分辨率的冲击响应分析中，适当减小它。

3. Bode 图的 Matlab 绘制：从bode函数到自定义频响曲线

3.1 最小命令与基本设置

Bode 图是控制系统分析和滤波器设计里的标准工具，横轴是对数频率，包含幅频特性（dB）和相频特性（度）两个子图。Matlab 里只要有一个系统模型sys，就能直接画出来。

% 构造一个二阶低通滤波器模型 wn = 2*pi*100; % 自然频率 100 Hz zeta = 0.707; % 阻尼比 sys = tf(wn^2, [1 2*zeta*wn wn^2]); figure; bode(sys); grid on;

tf的分子分母都是按 s 的降幂排列的系数向量。分子wn^2保证直流增益为 0 dB，分母是标准的二阶系统特征多项式。bode(sys)不指定频率范围时，Matlab 会根据系统极点零点自动选择一个有意义的区间，从低于最低转折频率到高于最高转折频率，大约各扩展 10 倍。

如果你有一个传递函数，但没有tf的形式，只有实测的幅值和相位数组，也可以用bode的返回参数手动绘图。

w = logspace(1, 5, 500); % 从 10 Hz 到 100 kHz，对数间隔 500 点 [mag, phase] = bode(sys, w); % 返回值 mag 和 phase 的形状是 1x1x500，需要 squeeze 降维 mag_db = 20*log10(squeeze(mag)); phase_deg = squeeze(phase); figure; subplot(2,1,1); semilogx(w/(2*pi), mag_db); % 注意，输入频率单位是 rad/s，要转成 Hz ylabel('幅值 (dB)'); grid on; subplot(2,1,2); semilogx(w/(2*pi), phase_deg); ylabel('相位 (deg)'); xlabel('频率 (Hz)'); grid on;

bode(sys, w)的第二个参数要求频率点按 rad/s 为单位，因为是角频率。logspace(1, 5, 500)生成的是指数间隔的角频率点，覆盖范围宽且等对数间隔，是绘制 Bode 图的标准密度。squeeze是因为bode在 MIMO 系统或多输入多输出时会返回多维数组，单输入单输出时也有多余的维度。这里如果不 squeeze 也可以画图，但semilogx对 1x1x500 的三维数组处理方式不直观，有时会得到错误结果。

3.2 多个系统对比与频率范围选择

工程里最常见的需求是把实测的频率响应和理论模型对比，或者把不同设计参数的滤波器放在同一张图里。

sys1 = tf(wn^2, [1 2*0.707*wn wn^2]); % 阻尼比 0.707 sys2 = tf(wn^2, [1 2*0.3*wn wn^2]); % 阻尼比 0.3 figure; bode(sys1, 'b-', sys2, 'r--', {2*pi*10, 2*pi*10000}); legend('zeta=0.707', 'zeta=0.3'); grid on;

bode支持直接传入多个系统，后面用线型字符串区分，频率范围用元胞数组{wmin, wmax}限定。这里的单位同样是 rad/s。阻尼比从 0.707 降低到 0.3，幅频特性在自然频率附近会出现明显的谐振峰，并且相位变化更陡峭。这正好提醒你，Bode 图上看到尖峰不一定是系统不稳定，要看阻尼比和相位裕度结合判断。

需要注意bode对连续系统和离散系统的处理方式不同。离散系统用c2d转换后，Bode 图的频率范围上限是 Nyquist 频率。这个容易出错的地方是在tf后直接 bode，没有意识到还要关心采样率。

3.3 用bodeplot控制显示细节

bode画图后要精调坐标轴或者隐藏某个子图，直接操作句柄比较麻烦。bodeplot返回句柄，配合setoptions可以更精细地控制。

sys = tf(wn^2, [1 2*zeta*wn wn^2]); h = bodeplot(sys); setoptions(h, 'FreqUnits', 'Hz', 'MagUnits', 'abs', 'PhaseWrapping', 'on');

FreqUnits设为'Hz'可以避免手动把 rad/s 转成 Hz。MagUnits设为'abs'时纵轴是线性幅值不是 dB。PhaseWrapping在需要看相位卷绕时会用到，一般默认关掉。bode函数本身不带这些选项，多数情况下bode够用，但如果你的报告需要特定单位、需要导出图例、需要把网格样式改掉，bodeplot是更合适的入口。

4. 实测频谱与 Bode 图叠加分析：验证系统幅频特性的完整流程

4.1 用 FFT 计算系统的实测频率响应

Bode 图是理论模型给出的频率响应，而实测频率响应可以通过输入输出信号的 FFT 比值获得。原理很简单：对输入信号 x(t) 做 FFT 得 X(f)，对输出信号 y(t) 做 FFT 得 Y(f)，系统的频率响应 H(f) = Y(f) / X(f)。幅值就是abs(H)，相位是angle(H)。

% 输入是扫频信号或白噪声，输出是同长度的响应信号 Fs = 5000; t = (0:length(x)-1)/Fs; X = fft(x); Y = fft(y); % 计算频率响应 H = Y ./ X; % 只取单边 half = floor(length(H)/2); H = H(1:half); f = Fs*(0:half-1)/length(H*2); % 正确频率轴构造 % 幅值与相位 H_db = 20*log10(abs(H) + eps); H_phase = unwrap(angle(H)) * 180/pi;

这里最容易出问题的是X在某些频率点上接近 0，导致H出现极大尖峰。常见做法是用输入信号的幅值阈值做掩码，只保留信噪比高的频段。eps是防止对数计算溢出。unwrap是相位展开函数，避免相位角在 ±180 度之间跳变。如果信号里有噪声，实测 H 的曲线会很毛糙，通常配合medfilt1做平滑。

H_smooth = medfilt1(H_db, 21); % 21 点中值滤波

medfilt1的窗口宽度要试，越大曲线越平滑，但会抹掉真实的窄带特征。21 点对于 2048 点 FFT 来说是适中的选择。

4.2 与理论 Bode 图叠在同一张图里

实测频响算出来后，自然的目标是和理论模型画在一起。

w = 2*pi*f; % 实测频率向量转成 rad/s 与 bode 对齐 [mag_th, phase_th] = bode(sys, w); mag_th_db = 20*log10(squeeze(mag_th) + eps); phase_th_deg = squeeze(phase_th); figure; subplot(2,1,1); semilogx(f, H_db, 'b.', 'MarkerSize', 3); hold on; semilogx(w/(2*pi), mag_th_db, 'r-', 'LineWidth', 1.2); xlabel('频率 (Hz)'); ylabel('幅值 (dB)'); legend('实测', '理论'); grid on; subplot(2,1,2); semilogx(f, H_phase, 'b.', 'MarkerSize', 3); hold on; semilogx(w/(2*pi), phase_th_deg, 'r-', 'LineWidth', 1.2); xlabel('频率 (Hz)'); ylabel('相位 (deg)'); grid on;

这里频率轴的单位要保持一致。实测的f是 Hz，理论数据中w/(2*pi)是 Hz，两者匹配。对比时会发现实测幅频曲线在中高频段与理论偏离明显，这通常是测量噪声、传感器自身频响、或系统非线性导致的。实际工作中遇到这种情况，先别怀疑模型不对，检查一下输入信号的频谱在关注频段是否足够平坦。

4.3 随机噪声下频响估计的 3 个参数选择

实测频响估计有三种常见策略，各自对应不同的参数选择。

平均互功率谱法实际上是 Welch 估计的变体。把输入输出各分段加窗做 FFT，然后计算互功率谱Pxy = conj(X).*Y和自功率谱Pxx = abs(X).^2，再对多段做平均，最后相除。

segmentLen = 1024; noverlap = 512; [Pxx, fw] = pwelch(x, hann(segmentLen), noverlap, segmentLen, Fs); [Pxy, ~] = cpsd(x, y, hann(segmentLen), noverlap, segmentLen, Fs); H_welch = Pxy ./ Pxx; H_welch_db = 20*log10(abs(H_welch) + eps);

pwelch输出负频率到正频率的双边功率谱，cpsd输出互功率谱。实际绘图中只要取前一半即可，横轴对应频率上限是 Fs/2。这里的segmentLen决定了频率分辨率和平均次数之间的权衡——段越短，平均次数越多，方差越小，但频率分辨率越差。

5. 让 Bode 图分析更省力的三个细节技巧

5.1 用findpeaks自动定位谐振峰并换算 dB

手动画几条垂直参考线标注谐振峰比较费事。findpeaks可以直接返回峰值位置和幅值。

[peaks, locs] = findpeaks(H_welch_db, 'MinPeakHeight', -20, 'MinPeakDistance', 10); peak_freqs = fw(locs); % 峰值对应的频率 for k = 1:length(locs) text(peak_freqs(k), peaks(k)+2, sprintf('%.2f Hz', peak_freqs(k))); end

MinPeakDistance需要按点数设置，而不是频率值。如果频率分辨率是 1 Hz，设置 10 就表示两个峰至少隔 10 个点。MinPeakHeight用来滤掉低于 -20 dB 的微小波动。这里得到的是峰值点数和索引，要换算成频率必须用频率向量fw去索引定位。

5.2 频率轴单位换算：Hz、kHz 与 rad/s 的表驱动转换

Bode 图的横轴单位总是造成混乱。tf中默认的算子变量是角频率 rad/s，bode(sys)的横轴如果没设置FreqUnits，默认也是 rad/s。实测 FFT 得到的是物理频率 Hz。做理论实测对比前，先统一单位。

factor = 1; % 1 表示用 Hz % factor = 1000; % 改用 kHz 时，频率轴除以 1000 f_display = freq_vals / factor;

用bodeplot的setoptions(h, 'FreqUnits', 'Hz')是最干净的做法。手动semilogx绘图时，注意freq_vals需要提前确认自己到底用的什么单位，在代码开头定义fact是很好的习惯，避免后面复制粘贴时忘记换算。

5.3 不知道采样率时，只能从相对关系读图谱

网上很多人问“不知道采样率怎么求频率频谱”，实际上，FFT 谱图的频率轴完全依赖采样率。不知道 Fs，你只能读出一个比例关系：某个峰的位置除以总点数，得到的归一化频率。用这个归一化频率去对照系统的特征频率，比如已知电机转频，就能反推其他峰相对它的倍数关系。这是一个仅适用于“找不到采样率但信号本身有已知参考频率”场景的技巧，用于估计频谱结构，无法做到绝对标定。

% 如果有已知参考频率 f_ref 且它在频谱的 bin 位置 k_ref f_unknown = f_ref * (k ./ k_ref); % 将其他峰 k 换算成绝对频率

这段代码实现了一个简单的线性映射：k_ref是参考峰在 FFT 中的 bin 序号，k是任意峰的位置。只要有一个准确已知频率，整个频率轴就被标定出来了。实际项目里，我通常会用 50 Hz 工频干扰或者设备本身的转频作为参考。最稳的办法仍是回到采集端找到采样率配置，这个方法只是排查过程中拿来兜底的。

最后的落脚点：一个含有findpeaks自动标注和频率轴归一化的模板，能大幅拉低频响曲线分析的时间成本。把你自己的数据替换进去，两分钟就能得到一张带标注的对比图。

本文还有配套的精品资源，点击获取