在信号处理这块，Matlab加FFT的组合我用了差不多十年，从本科课程设计做到实际工程项目，几乎每个阶段都在跟它打交道。频谱分析说到底就是把时域信号拆解成不同频率的成分，FFT是目前最高效的工具；滤波则是频谱分析的直接延伸——先搞清信号里哪些频率是想要的、哪些是干扰，再设计滤波器把它们分开。这篇内容就是我做FFT频谱分析与滤波的一次完整梳理，除了算法原理和可复现的代码，还会把我踩过的坑、调参的经验一并整理出来。适合信号处理入门到进阶的读者，也适合正在做课程设计或需要处理实测信号的工程师参考。
1. 从时域到频域：FFT频谱分析到底在分析什么
1.1 先理解“频率成分”这件事
很多人第一次接触FFT时，容易被那一堆复数运算吓到，其实核心思想特别朴素。你手头有一段波形，比如一段录音、一段振动信号，或者一段电流采样，直接看时域波形往往只能看到振幅在上下跳，看不出规律。但这段波形背后可能藏着多个不同频率的周期性成分，它们叠加在一起，看起来就是一锅粥。
FFT做的事情就是把这一锅粥重新分开。它把信号分解成一系列不同频率、不同幅度、不同相位的正弦波之和，输出结果里横轴是频率，纵轴是幅度（或功率），这就是频谱。举个例子，一段50Hz和120Hz叠加的波形，时域上看是扭曲的波浪线，频谱上一眼就能看到两个峰，一个在50Hz，一个在120Hz。做设备故障诊断、语音降噪、振动分析，本质上都是先靠这一步把信号里的“成分”摸清楚。
Matlab里做FFT不需要自己写算法，直接调fft()函数就行。函数底层用的是FFTW库，效率和稳定性都没话说。但这里有个关键点：fft()算出来的是双边频谱（包含正负频率），工程上通常只看单边频谱，而且幅度要根据点数做归一化。这些细节决定了你最终读出来的幅值准不准，后面会详细讲。
1.2 采样率、点数与频率分辨率：三个参数先定下来
做频谱分析之前，有个三角关系必须心里有数：采样率fs、FFT点数N、频率分辨率Δf。它们的关系是：
Δf = fs / N
这个式子意思是，频率轴上相邻两个点之间的间隔。如果你用1000Hz采样率采了1000个点，那Δf就是1Hz，也就是说你能区分相差1Hz的两个频率成分。如果两个峰的频率差小于Δf，频谱上就会糊成一个峰，看不出来是两个独立分量。
采样率还直接决定了你能分析的最高频率。根据奈奎斯特采样定理，采样率至少要是信号最高频率的两倍，频谱分析中能显示的最高频率就是fs/2，也就是奈奎斯特频率。超过这个频率的信号叫混叠，会伪装成低频信号出现在频谱里，非常坑人。所以做频谱分析前的第一步不是写代码，而是确认采样率合规。
| 采样率 fs (Hz) | FFT点数 N | 频率分辨率 Δf (Hz) | 可分析最高频率 (Hz) |
|---|---|---|---|
| 1000 | 1000 | 1 | 500 |
| 1000 | 2000 | 0.5 | 500 |
| 5000 | 4096 | 约1.22 | 2500 |
| 20000 | 16384 | 约1.22 | 10000 |
表格里的规律很直观，栅栏效应说的是信号真实频率可能落在两个频点之间，导致谱线幅度偏低。另外一个高频参数是补零（zero padding），很多人误以为补零能提高分辨率。实际上补零只是让频谱曲线看起来更平滑，并没有增加真实信息，两个靠得近的峰该分不开还是分不开，想真正提高分辨率只能增加采样时长或增大N。这个坑我见太多人踩过，提前说一声。
2. 滤波方案：先把“滤什么”想明白再动手
2.1 频谱分析是诊断，滤波是手术
频谱分析好比体检报告，告诉你哪里有病灶；滤波就是根据报告做手术，把不要的成分去掉，把保留的成分留下来。两者是天然搭配，我在实操中最常见的流程是：先对原始信号做FFT，找出干扰频段和有用信号的频段范围，再据此设计滤波器。比如采样一段电机振动信号，频谱上发现50Hz工频干扰很突出，目标信号在200Hz，那就做一个带阻或高通滤波器把工频附近压掉。
滤波器的类型选择要结合信号的频带特征。分三大类：低通（保留低频、滤掉高频）、高通（保留高频、滤掉低频）、带通/带阻（保留或剔除指定频段）。Matlab里设计滤波器非常成熟，一行fir1()或butter()就能拿到系数，但关键在于参数怎么定——截止频率、阶数、过渡带宽度，这些都得结合频谱分析的结果来选，拍脑袋随便填往往效果很差。
滤波器还有个重要的硬件概念是阶数。阶数越高，过渡带越窄、阻带衰减越深，但计算量也随之增大，还更容易带来数值不稳定和时间延迟。对实时系统和离线分析的选择策略是完全不同的。离线分析用高阶数没问题，反正Matlab一秒算完；但如果你要把滤波器部署到单片机上，就得分外小心，一次算几十个乘加运算对低成本MCU来说压力不小。这些取舍要在动手之前想清楚。
2.2 FIR与IIR的选型对比
设计滤波器时，Matlab提供了两条主要技术路线：FIR（有限冲激响应）和IIR（无限冲激响应）。两者各有适用场景，没有绝对好坏，我直接在项目里对比过很多次，总结如下：
| 对比维度 | FIR | IIR |
|---|---|---|
| 相位特性 | 可做到严格线性相位，波形不失真 | 非线性相位，可能引起波形畸变 |
| 计算效率 | 需要较高阶数才能达到同样衰减，计算量大 | 阶数低，计算量小，效率高 |
| 稳定性 | 始终稳定 | 设计不当时可能不稳定 |
| 适用场景 | 对相位敏感的信号（如生物电、数传） | 对相位不敏感的场合（如音频EQ、功率滤波） |
实际工作中我的选型习惯是：如果是离线数据分析，优先FIR滤波器，因为线性相位能保住波形形态，处理振动、心电、语音这类信号时不会因为相位扭曲导致后续特征提取出错。用Matlab的fir1()设计时，指定阶数和截止频率即可，阶数通常取N=2*fs/fc的1.5~3倍，过渡带能控制得比较理想。
如果是开发实时处理系统，或者噪声比较大但对相位不敏感的工程场景，IIR滤波器是更经济的选择。butter()设计巴特沃斯滤波器是最常用的，通带平坦，阻带单调衰减，没有纹波，非常适合做数据预处理。需要注意butter()默认的截止频率是归一化频率（用fs/2做分母），很多初学者直接填实际Hz值，结果滤波器完全不对，这是我看到的高频错误之一。
3. 完整实操：从一段含噪信号开始
3.1 构造测试信号并做基础频谱分析
纸上谈兵永远不如直接动手。这里我给出一个可复现的完整例子，信号由50Hz和120Hz两个正弦波叠加，再加上随机噪声，模拟一个很典型的“有效信号+干扰”场景：
fs = 1000; % 采样率 1000 Hz N = 2000; % 采样点数 2000，共2秒数据 t = (0:N-1)/fs; % 时间向量 f1 = 50; % 目标信号频率1 f2 = 120; % 目标信号频率2 x = 1.5*sin(2*pi*f1*t) + 0.8*sin(2*pi*f2*t) + 0.3*randn(1,N);这段代码里randn生成的是高斯白噪声，幅度系数0.3决定噪声大小。这样构造的好处是：我们知道真实的频率分量和幅值，后面频谱分析和滤波的效果好坏可以直接对照评估。用plot(t, x)画出来的时域波形已经比较混乱，如果不做频谱分析，根本看不出里面有两个明显频率成分。
接下来做单边频谱分析：
X = fft(x); P2 = abs(X/N); % 双边频谱归一化 P1 = P2(1:N/2+1); % 取单边 P1(2:end-1) = 2*P1(2:end-1); % 除直流外，单边频谱翻倍 f = fs*(0:(N/2))/N; % 频率轴 plot(f, P1); xlabel('频率 (Hz)'); ylabel('幅度');这里最关键的是幅度归一化。FFT结果的绝对大小和N成正比，直接abs(fft(x))会得到很大的数，不是真实幅值，所以要先除以N。又因为双边频谱的能量分成了正负两边，单边谱要把负频率那半的幅度加回来，因此正频率幅度乘以2（直流分量除外）。
运行这段代码后，频谱中应该能看到50Hz和120Hz两个清晰的谱峰，幅度分别接近1.5和0.8，和原始设定吻合。噪声则表现为一条不平整的“基底”，在0~500Hz范围内均匀铺开。这一步跑通了，说明频谱分析的框架已经没问题了。
3.2 加窗处理：解决频谱泄露问题
我前面构造的信号是整周期截断的，FFT点数刚好包含整数个信号周期，算出来谱线干干净净。但实测信号基本不可能这么凑巧，截断时间多半不是周期的整数倍，此时频谱上会出现“泄洪”一般的能量扩散，原本尖锐的谱峰变得粗胖，旁边还拖出一堆假的旁瓣，这个现象叫频谱泄露，本质是矩形窗的频域旁瓣在作怪。
解决办法是加窗。Matlab里内置了hann、hamming、blackman、kaiser等多种窗函数，选择标准要看你的需求：
| 窗函数 | 主瓣宽度 | 旁瓣衰减 | 适用场景 |
|---|---|---|---|
| 矩形窗 | 最窄 | 最差（约-13dB） | 整周期截断或纯频率测量 |
| 汉宁窗 | 稍宽 | 较好（约-31dB） | 一般工程信号，最常用 |
| 海明窗 | 稍宽 | 好（约-41dB） | 语音、窄带信号检测 |
| 布莱克曼窗 | 较宽 | 很好（约-58dB） | 强调幅值精度近旁瓣抑制 |
实际处理中我最常用的是汉宁窗，推荐度最高。引入窗函数的代码很简单：
w = hann(N)'; xw = x .* w; Xw = fft(xw); P2w = abs(Xw/N); P1w = P2w(1:N/2+1); P1w(2:end-1) = 2*P1w(2:end-1); plot(f, [P1', P1w']); legend('未加窗', '汉宁窗');不过有个细节要注意：加窗虽然抑制了旁瓣，但主瓣会变宽，如果两个频率离得很近（比如50Hz和52Hz），未加窗时勉强能分辨，加窗后可能糊成一个峰。所以这里是一个取舍，要根据信号特征选择是否加窗。此外，加窗会导致信号能量略有变化，幅值读数和真实值之间有偏差，工程上可以引入幅度恢复系数做补偿，汉宁窗约为1.63左右，具体看应用精度要求。
3.3 滤波器的设计与对比分析
频谱分析确认了噪声分布之后，开始设计滤波。假设我们的目的是保留50Hz和120Hz成分，抑制高频段噪声，那么低通滤波器的截止频率大概定在150Hz比较合适，太高则噪声压制不够，太低会把120Hz的有用信号削掉。
先用IIR方案实现：
fc = 150; % 截止频率 [b, a] = butter(4, fc/(fs/2), 'low'); y_iir = filtfilt(b, a, x); % 零相位滤波这里窝了个大坑：我用的不是filter()，而是filtfilt()。原因在于IIR滤波器存在非线性相位，直接用filter()会让输出波形发生相位偏移，这个偏移随频率变化，可能会导致波形严重变形。filtfilt()对信号做正向和反向两次滤波，相位抵消，输出零相位失真，代价是计算量加倍，且信号首尾的处理逻辑不同。离线分析和Matlab原型验证阶段用filtfilt()几乎是无脑正确选择。
再看FIR方案：
M = 128; % 阶数，需自行测试调整 b_fir = fir1(M, fc/(fs/2), 'low'); y_fir = filtfilt(b_fir, 1, x);FIR的filtfilt不是必须的，因为fir1()设计的是线性相位滤波器，直接用filter()就不会引起相位畸变。我习惯统一用filtfilt()纯粹是为了变量对比时排除相位因素的干扰。FIR滤波器的优势在信号形态要求高的场景非常明显。
设计完成后，我想验证滤波效果，最直观的方式是再做一次FFT，把滤波前后的频谱叠在一起看：
Y_iir = fft(y_iir); P1y = abs(Y_iir/N); P1y = P1y(1:N/2+1); P1y(2:end-1) = 2*P1y(2:end-1); plot(f, [P1', P1y']);正常情况下你应该看到高频基底被压平，50Hz和120Hz谱峰保持完好。
3.4 滤波效果的定量评估：别只靠眼睛看
把滤波前后的波形和频谱画在一起对比，是初学者最常用的方式，但眼睛看出来的“好”不一定靠得住，特别是噪声比较顽固时，主观判断容易骗人。我在实际项目中通常用两个定量指标：信噪比（SNR）和均方根误差（RMSE）。
因为仿真信号里我们“知道”真实信号是什么（不含噪声的纯净版），所以可以准确地计算滤波算法带来的改善：
% 构造纯净参考信号 x_clean = 1.5*sin(2*pi*f1*t) + 0.8*sin(2*pi*f2*t); % 滤波前信噪比 snr_before = 10*log10(sum(x_clean.^2) / sum((x - x_clean).^2)); % 滤波后信噪比 snr_after = 10*log10(sum(x_clean.^2) / sum((y_iir - x_clean).^2)); % 均方根误差 rmse_before = sqrt(mean((x - x_clean).^2)); rmse_after = sqrt(mean((y_iir - x_clean).^2));这个思路在真实数据上无法完全套用，因为你不可能得到纯净参考信号。但可以退而求其次，把滤波前后的功率谱积分做对比，或者统计频带外剩余幅度的均值。我习惯导出滤波前后的功率谱数据，做一个频带外总能量对比，简单又好记录到报告里。
从实际项目经验来说，滤波后的SNR提升在10~20dB属于正常工作范围，如果低于这个幅度，要怀疑截止频率或阶数选得不合适。这套评估方法也为后续调参提供了客观依据，而不是凭感觉叠参数。
4. 实战中躲不开的坑与排查实录
4.1 谱峰幅度对不上：归一化与窗函数恢复
很多人用Matlab做完FFT，兴致勃勃看谱峰幅度，发现和信号真实幅值差了一截。这不是函数的问题，是归一化没做对。我之前讲过要除以N和单边乘2，但还有一个被忽略的因素是窗函数的衰减。加窗后信号总能量被削弱，和未加窗直接算出来的幅值自然不同。
如果你确实关心谱峰幅值的绝对值，需要按窗函数能量恢复。汉宁窗的幅度恢复系数近似为2（实际约1.97），hamming窗接近1.85，矩形窗就是1，不需要修正。具体公式为：
A_corrected = A_measured / mean(w)
其中mean(w)是窗函数的平均值。如果只关心相对频率分布，不在乎绝对幅值，这一项可以不考虑。做故障诊断系统时，我们更关注的是特定频段能量相对变化趋势，所以很多算法在特征提取阶段不做绝对值修正，反而更稳定。
4.2 频率测不准：栅栏效应与频率细化
FFT输出的频率是离散的，栅栏效应导致真实信号频率不一定恰好落在谱线上。比如信号是50.3Hz，fs=1000、N=2000时频率分辨率是0.5Hz，你的谱峰大概率不出现在50.5Hz，而是在49.5Hz或50Hz上，读出来的幅值也偏低。
对于这种情况，最直接的办法是增加FFT点数N，即增加采样时长，分辨率自然提高。但如果数据长度受限，做不了更长采样，可以尝试谱峰重心法：在峰附近取几条谱线做加权平均，估算真实频率。Matlab里可以很方便地提取局部峰值区间，用sum(f.*P)/sum(P)算加权中心，精度通常能比谱线分辨率高一个数量级。
更进阶的是Chirp-Z变换或Zoom FFT，在指定的窄频段内加密计算。Matlab的信号处理工具箱仍然有czt()函数，可以直接调用。不过对大部分课程设计和中等精度工程需求，谱峰重心法已经足够用，不必上太复杂的工具。
4.3 输出波形首尾异常：边界效应处理
用filter()或filtfilt()滤波时，很多人会发现输出波形的前几十毫秒或最后几十毫秒有明显抖动甚至突变，中间段反而是正常的。这是数字滤波器的初始状态产生的瞬态响应（或边界效应），所有有限长度信号的滤波都会遇到，只是严重程度不同。
最实用的缓解策略是拓延数据：滤波前在信号两端各补一段数据（比如镜像翻转、零延拓或线性预测），滤波后再裁掉补的部分。Matlab的filtfilt()内部已经做了类似的边界延拓优化，所以首尾畸变比filter()小得多，这也是我在离线分析时更偏好filtfilt()的原因之一。
如果你的场景必须用实时滤波器（如在线监测），就只能接受初始的几个采样点波形不理想，做法是给系统留一个预热时间，抛弃前M个点的输出再做后续分析。M的值可以视滤波器阶数而定，阶数越高预热时间越长。
4.4 常见问题速查表
我在多个项目里汇总了一份排查思路，几乎所有入门者和部分进阶者遇到的问题都可以在里面找对应解法：
| 症状 | 可能原因 | 解决方法 |
|---|---|---|
| 频谱出现意料之外的峰 | 混叠；直流偏置；工频干扰 | 检查采样率；信号去直流后再做FFT；级联带阻滤波 |
| 幅值偏小/偏大 | 未归一化；未乘2；窗函数衰减 | 除以N；单边谱乘2；用窗函数能量修正 |
| 两个频峰分不开 | 频率分辨率不足 | 增加采样点数N；提高采样时长 |
| 滤波后波形抖动剧烈 | IIR相位非线性；边界效应；阶数过高 | 改用filtfilt；拓延数据；重置阶数 |
| 滤波器截止频率与预期不符 | 归一化频率错误 | 确保使用 fc/(fs/2) 作为截止频率参数 |
| 频谱图底部噪声不平坦 | 意外混入趋势项 | 先做去趋势处理 detrend() |
这张表里有一个经常被忽略的：直流偏置。实测信号往往带有直流分量，FFT结果中0Hz处会有一个非常高的峰，导致其他频率分量被压缩看不清楚。处理办法极简单，先x = x - mean(x)，或者用detrend(x)，0Hz的大峰就消失了，其余频谱细节立刻清晰起来。
最后说一点我自己的体会。频谱分析和滤波这两件事，工具层面的操作其实不难，Matlab里的函数都是一行行调用，难的是对信号的判断——比如你到底该保留哪个频段、该选什么窗函数、滤波器阶数定多少。这些问题没有标准答案，只能在理解原理的基础上多做实验。我个人习惯是每次先把原始信号的特征彻底摸清楚，再谈滤波，直接套滤波器往往会把有意义的信息一起滤掉。这篇内容里的代码和应用思路，后续你可以很自然地沿用到语音去噪、振动特征提取、生物电信号处理甚至水声通信领域，核心方法都是同一套：先看清频谱，再决定滤波方案。
再分享一个小技巧：在实际工作中，我会把FFT和滤波的代码封装成两个独立的函数，输入参数只有信号、采样率和目标频段，输出频谱数据和处理后的波形。这样无论是处理电机振动数据还是分析客户发来的噪声录音，都能直接复用，不用每次从零写一遍。做信号处理，积累自己的工具箱，比收集一百份现成脚本更管用。
