傅里叶变换频谱图代码实操：采样率、频率分辨率与频谱泄漏详解

简介《傅里叶变换频谱图你必须知道的》视频的配套代码包面向正在学习信号处理、音频分析与傅里叶变换的读者尤其适合希望把理论公式落到实际频谱图中的初学者。压缩包内共2个文件包含1个MATLAB脚本和1个MATLAB数据文件整体仅128KB下载和运行都很轻量。脚本基于MATLAB内置fft函数演示如何对时间域信号进行快速傅里叶变换并绘制幅度谱内容覆盖示例信号构造、频率轴生成、单边幅度谱计算等关键步骤可直观看到5Hz与10Hz频率分量的强度对比。数据文件提供实际采集的音频或图像信号读者可直接load后复用同一套流程体验从离散采样到频域可视化的完整过程已有2125人学习下载。结合视频讲解操作既能理解傅里叶变换、DFT与逆变换的关系也能在音频编辑、降噪等真实场景中迁移应用。 最近有朋友私信我说看了《傅里叶变换频谱图你必须知道的》那期视频代码下载了但不知道从哪一步开始跑数据文件打开了也是一脸懵。这期内容本来就是把傅里叶变换频谱图从概念到实操串了一遍视频里为了节奏砍掉了很多细节配套的代码和数据反而是最值得反复琢磨的部分。所以今天这篇文章就把视频里用的代码、数据、还有没来得及展开的参数逻辑一次性说清楚不管你是刚接触频谱分析还是已经写过几个FFT脚本都应该能从里面拿到点实在东西。这篇文章不是单纯贴代码而是把代码为什么这么写、数据为什么这么组织、参数为什么这么选都拆开讲。你跟着走一遍能搞清楚采样率、频率分辨率、频谱泄漏这些概念在代码里到底怎么体现看完可以自己动手改参数、换数据做出属于你自己的频谱图。1. 视频里的核心内容频谱图到底在“看”什么1.1 视频讲了什么那期视频的核心是讲清楚傅里叶变换频谱图怎么读、怎么用。很多人一开始接触FFT代码跑通了图也画出来了但横轴纵轴代表什么、峰值说明什么、噪声底在哪里完全没概念。视频解决的问题就是把这些“必须知道”的东西讲透。配套的代码和数据是为了让你能亲手复现视频里的每一张图。代码不是摆设数据也不是随便造出来的每一处都对应视频里一个知识点。比如频谱图为什么左边有个大峰值、右边却是一堆毛刺这背后是采样率选择的问题再比如同一个信号为什么视频里画的频谱图和你自己跑出来的不一样很可能是补零和窗函数的区别。1.2 代码和数据的设计思路视频里用的信号本身不复杂是一个多频率叠加的合成信号加上了一点点噪声。这么做的好处是你明知道这个信号“应该”长什么样再看频谱图就能直观理解FFT输出的每一部分。等你理解了再把真实数据套进来就能举一反三。配套数据文件里包含了原始时域信号和预处理后的频谱数据。代码部分分成了三个脚本数据生成、频谱分析、可视化出图。这样拆开是故意的方便你单独改某一部分比如想换一个信号只动数据生成脚本就行不必把分析和绘图全部重写。2. 动手之前必须搞懂的三个参数2.1 采样率决定频谱图能看多“宽”采样率是所有频谱分析的基石。奈奎斯特采样定理告诉我们采样率必须大于信号最高频率的两倍否则高频成分会折叠到低频区间产生混叠。视频里用的采样率是1000Hz也就是说频谱图能有效显示的最高频率是500Hz。# 信号生成脚本中的采样率设置 fs 1000 # 采样率 1000Hz奈奎斯特频率 500Hz t np.arange(0, 1, 1/fs) # 生成1秒时长的时间轴这里为什么选1000Hz而不是800Hz或者2000Hz因为视频里要展示的三个频率成分分别是50Hz、120Hz和300Hz1000Hz采样率轻轻松松满足了奈奎斯特定理的要求同时计算量又很小单次FFT几毫秒就出结果。如果你处理的信号最高频率是5kHz那采样率至少要10kHz这是硬性要求。2.2 频率分辨率决定频谱图能看多“细”频率分辨率等于采样率除以FFT点数。如果你的采样率是1000HzFFT点数也是1000那频率分辨率就是1Hz。这意味着你能分辨的最小频率差是1Hz两个频率如果只差0.5Hz在频谱图上会叠成一个峰。提高频率分辨率的办法有两个一是增加采样时长二是做补零。但补零并不增加真实信息只是让频谱曲线更平滑。很多人误以为补零能提高分辨率这是个经典的坑。视频里有张对比图专门讲这个配套代码里也能复现。# 频谱分析脚本中的FFT参数 n len(signal) # 1000个点 freq np.fft.rfftfreq(n, d1/fs) # 频率轴 fft_vals np.fft.rfft(signal) # 单边FFT2.3 幅值归一到dB让大峰小峰都能看清线性幅值下50Hz的大峰可能会把300Hz的小峰压到看不见。视频里展示的频谱图都用了dB刻度本质是对幅值取对数压缩动态范围。代码里就一行# 幅值转dB参考幅值设为1 db_vals 20 * np.log10(np.abs(fft_vals) / n 1e-12)这个1e-12是防止取对数时出现负无穷。实际处理信号时可能还会遇到DC分量特别大的情况这时通常会用减去均值或者加高通滤波的方式处理视频里没有深入讲但代码里的数据预处理环节预留了这个口子。3. 配套代码与数据的完整解读3.1 代码文件结构视频配套的代码一共三个脚本逻辑非常清晰文件功能对应视频章节generate_signal.py生成多频叠加信号并保存为CSV什么是频谱图fft_analysis.py读取CSV做FFT分析如何计算频谱plot_spectrum.py可视化出图频谱图怎么读这么拆的好处是你可以像搭积木一样组合使用。比如你已经有了自己的数据只要格式和CSV对得上完全可以直接跳过第一步把数据丢给fft_analysis.py去分析。3.2 数据说明CSV文件是三列结构时间、原始信号、含噪信号。视频里大部分频谱图都是用含噪信号画的这样更贴近真实场景。数据生成的核心代码逻辑如下# 多频叠加信号 signal (np.sin(2 * np.pi * 50 * t) 0.8 * np.sin(2 * np.pi * 120 * t) 0.5 * np.sin(2 * np.pi * 300 * t)) # 加一点高斯白噪声模拟真实采集 noisy_signal signal 0.2 * np.random.randn(len(t))这里三个频率的幅值故意设置成不一样的50Hz是1.0120Hz是0.8300Hz是0.5。这么做是为了让你看频谱图时能直观地通过峰高对比出信号的幅值比例这也是频谱图“定量读法”的一个练习。3.3 依赖库和版本兼容代码依赖numpy、pandas、matplotlib三个库。版本方面没有特别苛刻的要求numpy只要不是特别老的版本都能跑。如果你用的是Anaconda环境这三个库通常已经装好了。推荐用conda建一个隔离环境安装依赖避免和系统环境冲突conda create -n fft_demo python3.10 conda activate fft_demo pip install numpy pandas matplotlib4. 完整实操从原始信号到可发布的频谱图4.1 信号生成与保存第一步是运行数据生成脚本把时域信号落盘。这里面有几个容易踩的坑一是时间轴用np.arange还是np.linspace有的资料推荐linspace但我这里用arange是为了保证采样点数严格等于采样率乘以时长二是随机噪声种子代码里用了np.random.seed不设置seed的话每次跑出来的数据都不一样后面对比效果就不好复现。# 保存数据 df pd.DataFrame({time: t, signal: signal, noisy_signal: noisy_signal}) df.to_csv(fft_demo_data.csv, indexFalse)跑完之后你会得到1000行数据用Excel打开也能看不过1000行看趋势还行真要分析还是得靠代码。4.2 频谱分析核心逻辑频谱分析脚本做了三件事读CSV、做FFT、把复数频谱转成幅值谱。需要注意rfft和fft的选择。信号是实数序列时fft输出的后半部分是前半部分的共轭对称没有额外信息所以用rfft只取前半部分计算量直接减半。# 完整频谱分析流程 def analyze_spectrum(file_path, fs1000): df pd.read_csv(file_path) signal df[noisy_signal].values # 去直流分量避免0Hz处出现超大峰值 signal_detrend signal - np.mean(signal) # 加汉宁窗减少频谱泄漏 window np.hanning(len(signal_detrend)) signal_windowed signal_detrend * window # FFT计算 n len(signal_windowed) fft_vals np.fft.rfft(signal_windowed) freq np.fft.rfftfreq(n, d1/fs) # 幅值校正加窗后要补偿窗函数的能量损失 amplitude 2.0 * np.abs(fft_vals) / n amplitude / np.mean(window) # 窗函数幅值补偿 return freq, amplitude这里有几个细节值得展开解释一下。首先是去直流如果没有这一步0Hz处会有一个超级大的峰值其他频率的小峰值全都看不见。然后是加窗视频里特别演示了不加窗和加窗的对比不加窗时频谱图会有明显的“裙边”这是频谱泄漏加汉宁窗后旁瓣大幅降低主瓣变宽了一点点但整体清晰度提升明显。4.3 可视化出图可视化脚本里最核心的是画频谱图的方式。视频里展示了两种线性坐标和dB坐标。两种图表达的信息侧重点不同线性坐标适合看“哪些频率存在”dB坐标适合看“频率成分之间的相对强弱”。# 频谱可视化 plt.figure(figsize(12, 5)) # 子图1线性坐标 plt.subplot(1, 2, 1) plt.plot(freq, amplitude, linewidth0.8) plt.xlabel(Frequency (Hz)) plt.ylabel(Amplitude) plt.xlim(0, 500) # 子图2dB坐标 plt.subplot(1, 2, 2) plt.plot(freq, 20*np.log10(amplitude 1e-12), linewidth0.8) plt.xlabel(Frequency (Hz)) plt.ylabel(Magnitude (dB)) plt.xlim(0, 500) plt.ylim(-80, 0) plt.tight_layout() plt.savefig(spectrum_compare.png, dpi150)输出结果你应该能清楚看到三个峰值分别落在50Hz、120Hz、300Hz附近。dB坐标下50Hz峰大约在-12dB120Hz在-15dB左右300Hz在-20dB左右。这和你设置的幅值比例是吻合的因为20倍log10后幅值减半大约对应-6dB的落差。4.4 真实数据的替换套路学会看合成信号后下一步就是拿真实数据试试。视频评论区有不少人问能不能用这个代码处理自己的数据。答案是完全可以前提是数据得是(时间, 幅值)的两列结构比如传感器采集的振动信号、音频片段、心电信号等。举个例子如果你想分析一段录音的频谱用sounddevice库录一段10秒的音频采样率44100Hz保存成numpy数组然后直接用频谱分析脚本替换数据源就行。唯一要改的是fs参数和信号长度相关的逻辑其他都不用动。# 替换数据源示例读取语音文件 import soundfile as sf data, fs sf.read(speech.wav) # 读取音频 if data.ndim 1: data data.mean(axis1) # 多声道转单声道 signal data[:fs*5] # 取前5秒 # 后续直接复用analyze_spectrum逻辑5. 常见问题与排查技巧实录5.1 频谱图“毛刺”太多峰值不明显这是加窗之前最常遇到的问题。如果信号不是整周期截断FFT默认加的矩形窗会产生严重的频谱泄漏表现为主峰周围出现很多旁瓣噪声底抬高。解决办法就是加窗。汉宁窗是通用选择旁瓣衰减快主瓣宽度也还能接受。如果频率分辨率要求很高可以选布莱克曼窗但主瓣会再宽一点。实际处理中先加窗、再FFT频谱会干净很多。5.2 峰值位置跟预期不一致这个坑我遇到过好几次。比如你明明造了一个300Hz的信号频谱图上峰值却落在305Hz。原因通常是频率分辨率不够300Hz落在两个频率点中间能量被分散到附近多个点上了。解决办法有两个方向一是增加采样时长从1秒加到2秒频率分辨率就从1Hz变成了0.5Hz峰值能更贴近300Hz二是用rfftfreq时检查频率轴生成是否正确。另外还有一个容易忽略的问题如果你的信号长度不是2的整数次幂有些老代码会直接报错建议用numpy自带的FFT它对非2次幂的长度兼容性要好得多。5.3 加窗后的幅值变了很多人第一次加窗后发现峰值幅值比原来小了以为是bug。这不是bug是窗函数本身的特性。汉宁窗中间是1、两端是0信号整体能量被削弱了所以幅值会变小。解决办法是在加窗后做幅值补偿除以窗函数的平均值这样幅值就能恢复到接近真实值。这个补偿逻辑在之前的代码里已经体现了。5.4 数据量太大FFT跑不动如果一个信号是10分钟长的音频采样率44100Hz总点数超过2600万。直接做FFT内存开销巨大而且频谱图把2600万个点画出来也没意义。实用做法是分段处理每段一两秒做一次FFT把多次结果平均起来这就是经典的Welch方法scipy.signal.welch里直接封装好了。from scipy import signal as scipy_signal f, Pxx scipy_signal.welch(signal, fsfs, nperseg2048, windowhann, scalingspectrum)5.5 CSV文件读取报错如果你用Excel打开过CSV又保存了一次可能会发现Python读出来的数据不对了。这是因为Excel在部分区域设置下会把CSV里的逗号分隔符替换成分号。如果pd.read_csv报错或者读出来的数据只有一列用记事本打开CSV看看分隔符是什么然后给read_csv传对应的sep参数就行。6. 一些想单独拎出来说的实操心得视频里很多结论是靠代码验证出来的但我实际做的时候还是踩了不少坑。比如最开始做频谱图时峰值旁边总有一圈一圈的波纹查了资料才知道是指数窗和矩形窗的区别。后来我把视频里所有图和代码对应起来用不同窗函数对比着看才算真正理解了频谱泄漏这件事。再分享一个比较实用的小技巧当你想观察频谱图上某个不起眼的小峰时可以在plot的时候加上plt.yscale(log)或者用dB刻度小峰马上就能看出来。另外如果你的原始信号里有明显的直流偏置画图前先减均值不然0Hz处的大峰可能会把y轴范围撑爆其他细节全是平的。还有一个容易被忽略的点是时域信号的长度。FFT处理很“吃”长度长度太短的话频率分辨率不够峰值会糊在一起长度太长又会影响计算速度和内存。做实际项目时建议先估算信号的最低频率成分确保一个FFT块里至少包含几个周期否则低频成分根本看不出来。这个项目虽然是以视频配套资料的形式发布的但代码和数据本身完全可以当做一个微型频谱分析工具箱来用。不管是做音频处理、振动分析还是学习信号处理都可以拿这套代码做底子按需改一改。以后如果我遇到有意思的信号处理案例也会继续把代码和数据更新到这个仓库里大家可以持续关注。本文还有配套的精品资源点击获取