1. 从一个采样率翻车的夜晚说起三年前调试一套电机振动监测系统加速度传感器输出接进采集卡采样率设定为 10 kHz理论上能覆盖到 5 kHz 以下的信号。可抓出来的频谱图在 2.8 kHz 附近总有一个莫名其妙的尖峰把轴承内圈故障频率的边带完全淹没了。我们换了传感器、换了采集卡、甚至怀疑是机械共振折腾到凌晨两点才发现真正的问题出在没做抗混叠滤波4.8 kHz 附近的一个高频谐波被折叠到了 2.8 kHz 处。这件事让我意识到一个很朴素的道理——离散傅里叶变换和频谱分析这两个词在教科书里通常是连着出现的但真正做过实际信号处理的人都知道DFT 本身并不可怕可怕的是你在把连续信号喂给 DFT 之前忽略的那些环节。采样定理、频谱泄露、窗函数选择、频率分辨率、幅度谱归一化任何一个环节没想清楚你看到的频谱都可能是一张谎言地图。这篇博文想写给那些已经在用 NumPy、MATLAB 或者嵌入式 DSP 库做频谱分析但总觉得结果看着对、用着虚的同学。我会从 DFT 的数学骨架讲起一步步过渡到工程实现再把踩过的坑整理成可直接查的排查表。如果你正在做振动分析、音频处理、电力谐波检测、雷达信号处理或者任何需要把时间波形变成频率柱状图的活儿下面的内容应该能帮你省掉几个通宵。2. 为什么连续傅里叶变换不够用2.1 计算机不认识连续函数连续傅里叶变换的定义是一个积分$$X(f) \int_{-\infty}^{\infty} x(t) e^{-j2\pi ft} dt$$这个公式在数学上很漂亮但计算机做不了两件事第一它无法存储无限长的时间信号第二它无法计算连续积分。所以工程上必须做两步离散化把时间轴切成 $N$ 个等间隔的点把频率轴也切成 $N$ 个等间隔的点。这一刀切下去就有了离散傅里叶变换。理解 DFT 最直观的方式是把它当成对连续频谱的一次等间隔采样。假设原始信号 $x(t)$ 的频谱是 $X(f)$采样周期 $T_s$采样点数 $N$那么 DFT 输出的第 $k$ 个频率点对应的物理频率是 $k \cdot f_s / N$其中 $f_s 1/T_s$。这个 $f_s / N$ 就是频率分辨率它决定了你能否把两个靠得很近的谱峰分开。很多人第一次看到 DFT 公式会觉得别扭因为它不是对称的$$X[k] \sum_{n0}^{N-1} x[n] e^{-j2\pi kn/N}$$其实可以这样理解$e^{-j2\pi kn/N}$ 是一组正交基DFT 就是在问信号 $x[n]$ 里含有多少这个频率的成分。正交性保证了每个频率点的投影互不干扰这是 DFT 可逆、能无损重构的前提。2.2 DFT 和 FFT 的关系别把两者混为一谈这是我带新人的时候必问的一个问题DFT 和 FFT 有什么区别正确答案是——DFT 是变换的定义FFT 是计算 DFT 的一种快速算法。DFT 的暴力计算复杂度是 $O(N^2)$FFT 通过分治把复杂度降到 $O(N \log N)$。当 $N 1024$ 时FFT 大约比暴力 DFT 快 100 倍当 $N 1$ 048 576 时差距超过 5 万倍。这意味着什么意味着你在实际工程中调用np.fft.fft()的时候底层走的是 FFT 算法但输出的结果在数学上与 DFT 定义完全一致。你不需要自己写蝶形运算但你必须知道 FFT 对点数有偏好——绝大多数实现要求点数是 2 的幂或者可以分解为小素数的乘积。如果你的数据长度是质数FFT 库可能会退化成慢速算法甚至直接报错。实操提示做 FFT 之前先把数据补齐到最接近的 2 的幂长度。补零不会增加真实频率分辨率但能让 FFT 跑得快、让频谱曲线更平滑。2.3 频率分辨率的真相物理分辨率 vs 视觉分辨率这是最容易混淆的一个点。假设采样率 1000 Hz采集 1 秒数据得到 1000 个点做 1000 点 FFT频率分辨率是 1 Hz。如果你把这 1000 个点补零到 8192 点再做 FFT频率轴上的点距变成 0.122 Hz看起来分辨率提高了但实际上两个原本间隔 1 Hz 以内的信号依然无法被分辨。原因很简单物理频率分辨率只取决于实际观测时长 T即 $\Delta f 1/T$。补零只是在原有频谱上做插值让曲线更光滑不会凭空变出新的频率信息。我把这个规律总结成一句话给团队新人采样率决定你能看到多高的频率观测时长决定你能分清多近的频率。参数决定因素典型影响最高分析频率采样率 $f_s$只能看到 $f_s/2$ 以下的信号物理频率分辨率观测时长 $T$间隔小于 $1/T$ 的谱峰无法分开频率轴点距FFT 点数 $N_{\text{fft}}$补零可减小点距但不提升物理分辨率幅度精度窗函数与归一化方式直接影响谱峰高度的可信度3. 从采样到频谱工程实现的关键环节3.1 采样定理不是大于两倍就万事大吉奈奎斯特采样定理说采样率要大于信号最高频率的两倍但实际工程中我建议至少取 2.5 倍甚至 5 倍。为什么因为真实信号几乎不可能是严格带限的总会有高频噪声、谐波或者瞬态成分。如果你的采样率刚好卡在 2 倍任何高于 $f_s/2$ 的成分都会折叠回低频造成混叠。抗混叠滤波器是必须的。模拟前端要放一个截止频率略低于 $f_s/2$ 的低通滤波器把高频成分在采样之前就砍掉。我见过太多项目为了省一个运放和几个电容结果在频谱上花了十倍时间做后处理得不偿失。判断是否发生混叠有一个简单方法改变采样率看频谱峰的位置是否跟着变。如果某个峰的位置固定不变那很可能是真实信号如果峰的位置随采样率变化而移动那基本可以确定是混叠产物。3.2 窗函数不是可选项是必选项对有限长数据做 FFT本质上是对无限长信号乘了一个矩形窗。矩形窗的频谱有较大的旁瓣会导致强信号的旁瓣淹没附近的弱信号这就是频谱泄露。不同的窗函数在主瓣宽度和旁瓣衰减之间做取舍窗类型主瓣宽度bin旁瓣衰减dB适用场景矩形窗0.89-13瞬态信号、整周期采样汉宁窗1.44-31连续信号通用分析汉明窗1.30-43音频、语音分析布莱克曼窗1.68-58强弱信号共存场景平顶窗2.94-70幅度精确测量选窗的核心逻辑是如果你关心频率定位精度选主瓣窄的如果你关心幅度精度选旁瓣低的。汉宁窗是最通用的默认选择我大概 70% 的场合都用它。注意加窗之后信号的幅度会被衰减必须在频域做补偿。常见的补偿系数是窗函数时域序列的均值汉宁窗的幅度补偿系数约为 2.0汉明窗约为 1.85。3.3 单边谱与双边谱为什么你看到的幅度是两倍FFT 输出的是双边谱频率范围从 $-f_s/2$ 到 $f_s/2$。对于实信号正负频率的幅度是对称的所以我们通常只看正半轴把负半轴的能量合并过来这就是单边谱。单边谱的幅度换算规则很简单直流分量$k0$和奈奎斯特分量$kN/2$保持不变其余频率点的幅度乘以 2。如果你用np.abs(fft_result)直接画图看到的幅度会比真实物理幅度小一半对非直流分量而言。归一化还要除以 FFT 点数 $N$。完整的单边幅度谱计算公式$$A[k] \frac{2}{N \cdot C} \left| X[k] \right|, \quad k 1, 2, \ldots, N/2-1$$其中 $C$ 是窗函数的幅度补偿系数。这个公式我建议每个做频谱分析的人都亲手推导一遍否则你永远不知道自己画的谱对不对。3.4 频率轴的构建别再手算 bin 号了频率轴第 $k$ 个点对应的物理频率$$f_k k \cdot \frac{f_s}{N_{\text{fft}}}$$其中 $N_{\text{fft}}$ 是实际做 FFT 的点数补零后的点数不是原始数据长度。这一点很容易搞混。如果你用 1000 个原始点补零到 8192 点做 FFT频率轴要按 8192 来算但物理分辨率仍然由 1000 个点的观测时长决定。代码层面我通常这样构建频率轴import numpy as np fs 10000 # 采样率 N_original 5000 # 原始点数 N_fft 8192 # 补零后 FFT 点数 freq_axis np.fft.rfftfreq(N_fft, d1/fs) # 或者手动构建 freq_axis_manual np.arange(N_fft // 2 1) * fs / N_fftnp.fft.rfftfreq返回的是单边谱的频率轴长度是 $N/21$直接对应np.fft.rfft的输出。用这个组合可以避免拼接负频率的麻烦。4. 频谱分析实操全流程4.1 完整代码框架从原始数据到可读频谱下面这套流程是我在实际项目中反复使用并逐步固化的。它包含去均值、加窗、FFT、幅度归一化、单边谱提取五个步骤适用于大多数振动、音频、电力信号分析场景。import numpy as np import matplotlib.pyplot as plt def spectrum_analyze(x, fs, windowhann, n_fftNone): 输入: x: 一维时域信号 fs: 采样率 window: 窗函数类型 n_fft: FFT 点数默认为最接近数据长度的 2 的幂 输出: freq: 频率轴 (Hz) amp: 单边幅度谱 N len(x) if n_fft is None: n_fft 2 ** int(np.ceil(np.log2(N))) # 1. 去均值消除直流偏置 x x - np.mean(x) # 2. 加窗 if window hann: win np.hanning(N) cg 0.5 # 汉宁窗幅度补偿系数 elif window hamming: win np.hamming(N) cg 0.54 else: win np.ones(N) cg 1.0 x_win x * win # 3. FFT X np.fft.rfft(x_win, nn_fft) # 4. 幅度归一化 amp np.abs(X) / (N * cg) amp[1:-1] * 2 # 单边谱非直流分量乘 2 # 5. 频率轴 freq np.fft.rfftfreq(n_fft, d1/fs) return freq, amp这段代码里有几个细节值得单独拎出来说。x x - np.mean(x)这一步很多人会忽略但直流分量在频谱上就是一个巨大的 0 Hz 峰它的旁瓣可能污染低频段。归一化除以的是 $N$ 而不是 $n_{\text{fft}}$因为窗函数只加在原始数据上补零部分没有信号能量。amp[1:-1] * 2跳过了直流和奈奎斯特点这两个频率点不应该乘 2。4.2 参数选择的计算过程假设你要分析一个电机振动信号关注的故障频率在 100 Hz 到 2000 Hz 之间其中有两个特征频率间隔约 8 Hz需要把它们分开。第一步确定采样率。关注最高频率 2000 Hz考虑 2.5 倍以上余量选 $f_s 10000$ Hz。这样最高分析频率 5000 Hz有足够余量。第二步确定观测时长。要分辨 8 Hz 间隔物理分辨率 $\Delta f \leq 8$ Hz所以观测时长 $T \geq 1/8 0.125$ 秒。实际工程中取 1 秒更稳妥于是原始点数 $N f_s \times T 10000$ 点。第三步确定 FFT 点数。10000 点不是 2 的幂最接近的是 16384。补零到 16384 点后频率轴点距变成 $10000 / 16384 \approx 0.61$ Hz曲线更光滑。但物理分辨率仍然是 $1/1 1$ Hz。第四步确定窗函数。间隔 8 Hz、分辨率 1 Hz谱峰间隔 8 个 bin汉宁窗主瓣宽度 1.44 bin不会造成明显遮蔽。同时汉宁窗旁瓣衰减足够选它。这套参数验证下来100 Hz 到 2000 Hz 范围内的谱峰应该都能清晰呈现8 Hz 间隔的两个峰也能分开。4.3 频谱平均降低方差的有效手段单次 FFT 得到的频谱方差很大谱线看起来毛刺很多。工程上常用的做法是分段平均也就是把长数据切成若干段每段分别做 FFT然后对幅度谱求平均。这里有一个经典权衡段数越多方差越小但每段越短频率分辨率越差。假设你有 10 秒数据采样率 10 kHz总共 100000 点。分 10 段每段 10000 点分辨率 1 Hz平均 10 次分 50 段每段 2000 点分辨率 5 Hz平均 50 次如果你的目标是检测间隔较远的谱峰选第二种如果要分辨靠得很近的谱峰只能牺牲平均次数选第一种。分段时还有一个细节段与段之间通常设置 50% 重叠这样可以在不增加数据总长的前提下增加平均次数同时减少因分段边界造成的信息丢失。这就是 Welch 方法的雏形。4.4 对数谱与线性谱别只会看一种线性幅度谱适合观察谱峰的绝对幅度比如判断某个频率成分是否超标。但线性谱的缺点是动态范围有限弱信号容易被强信号的旁瓣淹没。对数谱dB 谱把幅度取 20 倍对数动态范围一下子拉开。计算公式$$A_{\text{dB}}[k] 20 \log_{10} \left( \frac{A[k]}{A_{\text{ref}}} \right)$$$A_{\text{ref}}$ 通常取 1 或者信号的最大幅度。对数谱适合观察谐波结构、噪声底、弱边带。我自己的习惯是先看线性谱确认主峰位置和幅度再看对数谱找谐波和边带。提示对数谱不能有零值计算前加一个小常数如 $10^{-12}$ 避免 $\log(0)$。5. 常见问题与排查技巧实录5.1 频谱峰位置不对先查频率轴再查信号这是出现频率最高的问题。频谱上的峰位置与理论值对不上通常有三类原因。第一类是频率轴算错了。最常见的是把补零前的点数当成 FFT 点数来算频率轴或者用了双边谱的频率轴去对应单边谱的数据。排查方法很简单输入一个已知频率的正弦波看峰是否落在正确位置。第二类是采样率设置错误。采集卡的实际采样率与代码里写的 $f_s$ 不一致可能是硬件分频、时钟源配置或者驱动层做了重采样。排查方法是采集一个标准信号源输出比如 1 kHz 正弦看频谱峰是否在 1 kHz。第三类是信号本身就不是你以为的那个频率。比如电机转速波动导致基频漂移或者传感器安装方式改变了共振频率。这时候要结合时域波形一起看。现象可能原因排查方法峰位置整体偏移采样率设置错误采集已知频率标准信号验证峰位置成倍数偏差频率轴点数计算错误检查是否用补零点数构建频率轴低频出现大峰直流偏置未去除检查去均值步骤高频出现鬼峰混叠改变采样率观察峰是否移动峰两侧不对称窗函数选择不当换用旁瓣更低的窗5.2 幅度不对归一化系数逐项核对幅度出错的原因比频率出错更隐蔽因为频率轴一眼就能看出对错幅度却需要跟理论值仔细比对。我见过最多的错误是忘记除以 $N$导致幅度大得离谱其次是忘记做单边谱乘 2导致幅度小一半第三是窗函数补偿系数没加幅度小了约一倍。正确的归一化流程是先除以原始数据长度 $N$再除以窗函数补偿系数 $C$最后对非直流分量乘 2。三者的顺序不影响结果但缺一不可。验证方法生成一个幅度为 1 的正弦波做完整流程看输出谱峰是否接近 1。5.3 频谱泄露如何判断和抑制泄露的典型表现是谱峰底部变宽像一座山的裙边拖得很长。如果裙边淹没了旁边的弱信号说明泄露已经影响到分析结果。判断泄露是否严重可以看谱峰两侧的衰减速度。理想情况下远离主瓣后幅度应该快速下降到噪声底。如果下降很慢或者有明显的周期性起伏说明窗函数旁瓣太高或者信号没有整周期截断。抑制泄露的方法有三条选旁瓣更低的窗汉宁换布莱克曼、增加观测时长让信号更接近整周期、如果信号周期已知直接按整周期截断。第三条最彻底但需要知道信号基频适合转速稳定的旋转机械分析。5.4 频率分辨率不够补零救不了你新手最常见的误解是补零能提高分辨率。补零只能让频率轴更密让曲线更平滑不能让两个原本分不开的谱峰分开。真正的解决方法是增加观测时长。假设你要分辨间隔 5 Hz 的两个信号采样率 10 kHz。最少需要观测 0.2 秒即 2000 个点。如果你只有 1000 个点0.1 秒分辨率只有 10 Hz无论怎么补零都分不开。遇到这种情况要么延长采集时间要么用参数估计方法如 MUSIC、ESPRIT做超分辨率分析。后者属于进阶话题适合信噪比高、信号模型明确的场景。5.5 实时频谱分析中的帧长与刷新率权衡做实时频谱显示的时候帧长和刷新率是一对矛盾。帧长越长频率分辨率越高但每帧计算耗时越长刷新率越低。帧长越短刷新快但分辨率粗糙。我的经验值是刷新率保持在 10 到 20 帧每秒比较舒适人眼不会觉得卡顿。基于这个刷新率单帧处理时间要控制在 50 到 100 毫秒以内。对于 10 kHz 采样率、16 位精度、单通道信号这个时间预算足够处理 8192 点左右的 FFT。如果信号采样率很高比如 1 MHz单帧 8192 点只覆盖 8 毫秒分辨率只有 125 Hz很多细节看不到。这时候要么降低刷新率要么用多分辨率分析策略高频段用短帧快速刷新低频段用长帧慢速刷新。6. 进阶话题与个人实践体会6.1 相位谱的价值被严重低估绝大多数人只关注幅度谱把相位谱当成可有可无的附属品。但在某些场景下相位信息才是关键。比如判断两个通道信号的时延用互谱的相位斜率可以精确到采样点以内再比如做模态分析相位关系能区分同频不同振型的信号。相位谱的计算要注意相位是相对于 FFT 起点的如果每次 FFT 的起点没有对齐相位会随机跳变。做相位相关分析时必须保证各次采集的触发时刻一致。6.2 功率谱密度与功率谱的区别功率谱PSD和功率谱密度PSDPower Spectral Density是两个不同的量。功率谱的单位是幅度的平方表示某个频率 bin 内的功率功率谱密度的单位是功率每赫兹表示单位带宽内的功率密度。两者只在频率分辨率归一化上有差异但混用会导致量纲错误。做随机信号分析、噪声测量、振动总级值计算时必须用功率谱密度否则改变 FFT 点数会改变总功率计算结果这是不对的。功率谱密度的归一化系数要除以频率分辨率 $\Delta f$。6.3 我的个人实践清单做了几年频谱分析我总结了一张上电前必查清单每次开始一个新项目都会过一遍传感器带宽是否覆盖关注频率范围抗混叠滤波器截止频率是否低于 $f_s/2$采样率与关注最高频率的比例是否大于 2.5观测时长是否满足频率分辨率要求窗函数类型是否匹配信号特征归一化系数是否逐项核对过是否用已知信号源验证过整条链路对数谱和线性谱是否都看过这八条看起来基础但每一条我都至少踩过一次坑。尤其是最后一条很多隐蔽的谐波和边带只有在对数谱上才看得清楚。6.4 关于学习路径的一点建议如果你刚接触 DFT 和频谱分析我的建议是先不要急着上代码。找一本信号处理教材把 DFT 的定义、性质、帕塞瓦尔定理、循环卷积这几块手推一遍。然后在纸上画一个 8 点 DFT 的蝶形图亲手算一遍输入输出。这些看起来低效的练习会在你后面调参的时候变成直觉。真正开始做工程的时候从正弦波加白噪声这种最简单的信号入手逐步增加复杂度单频加谐波、多频加噪声、调幅信号、调频信号、瞬态冲击。每增加一种信号类型就回头看一下频谱是否符合预期。这个过程大概需要两三个月但走完之后你看到任何频谱图都不会发怵。频谱分析这件事工具和方法都是公开的差距就在对细节的把控上。频率轴多一点少一点、窗函数换一个、归一化差一项结果可能完全不同。这些细节没法从教科书上直接学到只能在反复实践和排查中积累。希望这篇整理能帮你少走几段弯路。