简介《基于傅里叶变换的电力系统谐波分析》是东北电力大学电气工程学院撰写的一份信号与系统课程专题报告面向电气工程专业学生、电力系统初学者及相关工程技术人员聚焦谐波对电网安全稳定运行的影响和傅里叶变换在谐波检测中的核心作用。资源共1个PDF文件压缩包整体452KB内容精炼、结构完整便于直接阅读或打印用于课程汇报。目前已有984人学习下载是同类课程报告中较受欢迎的参考材料。报告在概述傅里叶级数与DFT/FFT原理的基础上重点分析了频谱泄漏、栅栏效应及其加窗抑制方法并结合电缆过热、变压器老化、零线过流、电容器损坏与保护误动等典型谐波危害系统梳理了降低谐波源含量、改善供电环境、安装无源/有源滤波器等治理措施同时给出了在MATLAB中利用FFT函数绘制频谱图的分析思路可作为信号与系统课程设计、谐波专题研究或电力系统谐波问题入门的实用参考资料。1. 谐波分析不是套FFT就完事先看清信号再谈频谱去年帮一个做课程专题报告的同学看代码他的题目就是基于傅里叶变换的电力系统谐波分析。他把电流录波数据丢进fft函数频谱图上 50Hz 附近却堆着一片乱峰5 次谐波标定值 20A画出来只有 17A 多。问题不是 FFT 算错了而是信号本身没处理数据窗没有和电网基波对齐窗函数也没加直流偏置还在。电力系统谐波分析和语音频谱分析不一样基波频率不是固定 50Hz现场录波数据里它可能漂到 49.8Hz采样率稍不匹配谱线就落在两个频率点中间读出来的幅值必然偏小。这篇就按课程报告的实际推进路径讲清楚怎么从录波数据或仿真信号出发把采样率、数据窗、幅值修正、THD 这些环节一步步确认到位让最终的交直流分量、各次谐波幅值和相位经得起追问。适合正在写信号与系统专题报告、或者刚接触电能质量分析想快速跑通一条完整流程的从业者。2. 为什么电力谐波分析绕不开傅里叶变换从周期分解到离散谱2.1 傅里叶级数是谐波分析的起点电力系统正常运行时电压电流是 50Hz 的正弦波但在变频器、整流器、充电桩这类非线性负载接入后电流波形会明显畸变不再是一个干净的正弦。信号与系统课程里讲得很清楚满足狄利克雷条件的周期信号可以分解成直流分量、基波和一系列整数倍频分量的叠加。设基波角频率为 ω02πf0畸变电流可以写为i(t)I0∑_{h1}^{∞}[a_h·cos(hω0t)b_h·sin(hω0t)]工程上更常写成各次谐波单独表示的极坐标形式i(t)∑_{h0}^{∞} I_h·cos(hω0tφ_h)其中 h1 是基波h≥2 是各次谐波。谐波分析要做的就是从这个混合信号里把每个 I_h 和 φ_h 估计出来。傅里叶级数天然适合这个任务因为各次谐波在频域上是离散的谱线50Hz 处是基波250Hz 处是 5 次350Hz 处是 7 次只要频率分辨率足够它们不会叠在一起。这比在时域里用滤波器硬抠某个频带要直观得多。这里要理解一个边界傅里叶级数只适用于周期信号。电网在稳态下的畸变波形近似周期所以课程报告用傅里叶级数展开是合理的但如果是雷击、短路、电弧这类暂态过程波形突变周期假设不成立FFT 给出来的频谱会失真那就得换小波或短时傅里叶。这也是我在报告里一定会写清楚的一句话这次分析针对的是稳态工况不是暂态过程。2.2 从连续积分到 DFT课程报告真正用的是离散算法连续傅里叶变换的定义是积分式数字采集系统拿到的是一串按采样周期 Ts 离散的序列 x(n)没法直接做连续积分。信号与系统课程里对应的方法是 DFS 和 DFT。DFT 的定义式是X(k)∑_{n0}^{N-1} x(n)·e^{−j2πkn/N}这里的 k 是谱线序号k0 对应直流kN/2 对应奈奎斯特频率 fs/2。实际频率和谱线序号的换算关系是f k·fs/N这个公式是整个谐波分析里最关键的换算。fs 是采样率N 是 FFT 点数fs/N 就是频率分辨率 Δf。比如 fs10000Hz、N10000 时Δf1Hz50Hz 正好落在第 50 根谱线上5 次谐波 250Hz 落在第 250 根谱线上。如果 N 选成 1024Δf≈9.77Hz250Hz 对应 k≈25.6谱线对不准能量就会漏到相邻谱线上幅值偏低这在后面专门讲。FFT 并不是一种新变换它只是把 DFT 的计算复杂度从 O(N²) 降到 O(N·logN)。课程报告里不需要自己实现 FFT直接用现成库就好但为了说明原理我习惯让学生先写一个最朴素的 DFT再和库函数对比结果import numpy as np def dft_manual(x): 直接按DFT定义计算仅用于理解原理 N len(x) X np.zeros(N, dtypecomplex) for k in range(N): for n in range(N): X[k] x[n] * np.exp(-2j * np.pi * k * n / N) return X # 用一个短序列验证手写DFT和FFT结果一致 x np.array([1.0, 2.0, 1.0, 0.5, 0.3], dtypefloat) X_manual dft_manual(x) X_fft np.fft.fft(x) print(np.max(np.abs(X_manual - X_fft))) # 约 1e-14可认为一致这段代码的作用不是为了在实际分析中使用而是验证工具链没有方向性错误。DFT 手写实现里有两个循环内层循环是每个时间点 n 对每个频率点 k 求复指数乘积所以复杂度是 N²。np.fft.fft 底层的 FFT 也是算同一个数学对象两者差值在浮点精度范围内说明我们后面的频谱读取可以直接信任库函数的输出。参数上要说明一点上述验证序列 N5这是一个教学演示序列不代表电力谐波场景。实际谐波分析中 N 的来源取决于数据窗长度一般至少覆盖 10 个工频周期也就是 0.2 秒这样每个谐波分量在频谱里有足够的分辨率不会和邻近频率混在一起。如果只取一个周期 20ms频率分辨率是 50Hz5 次谐波 250Hz 和 6 次谐波 300Hz 之间隔了 50Hz刚好能分开但如果存在 5.5 次间谐波就完全看不出来了。2.3 为什么选 FFT而不是滤波器组或其它方法谐波分析的可选方案不止 FFT。常见的还有带通滤波器组、小波变换、瞬时无功功率理论等。课程报告选 FFT最重要的理由是它一次能给出全频谱的所有频率分量幅值和相位都能拿到而且计算量小、实现路径清晰非常适合在 Python 或者 MATLAB 里完成。滤波器组适合嵌入式装置实时监测因为每个通道只计算一个频带硬件实现方便但通道之间的旁瓣泄漏很难做到理想而且相位信息需要额外处理。小波变换在暂态谐波分析里有优势它能同时看到不同时刻的频率变化但是小波基的选取对结果影响很大不同人选的基不一定有可比性在稳态谐波测量中反而没有 FFT 简洁。瞬时无功功率理论更多用于三相电路的实时补偿控制它输出的是时域补偿指令不是直接的各次谐波频谱。FFT 的短板也很明确它对基波频率偏移敏感对非整数次谐波和短时突变信号处理能力弱。这也是我为什么始终强调FFT 不是套上就能出正确结果的算法后面章节的所有参数设置本质上都在弥补这个短板。3. 用 FFT 跑通一条谐波分析流水线采样、加窗与频谱提取3.1 先从仿真数据搭验证环境我每次做谐波分析的第一个建议都是不要直接拿现场录波数据开跑。现场数据没有真值你算出来的 5 次谐波幅值到底对不对没有人能告诉你。课程报告更理想的推进方式是先构造一组已知答案的仿真信号把算法流程跑通确认误差在合理范围内再切换到实际录波数据。这样一旦结果不对能定位是算法问题还是传感器与录波的问题。构造仿真信号时我会模拟一个典型的畸变电流基波 50Hz、100A叠加 5 次谐波 20A、7 次谐波 10A、11 次谐波 4A再给一个幅值远小于谐波的高斯噪声用来模拟测量噪声。import numpy as np fs 10000 # 采样率 10 kHz N 10000 # 采样点数正好 1 秒50 个工频周期 t np.arange(N) / fs f1 50.0 # 基波 5次 7次 11次谐波相位随意设置用于验证相位提取 signal ( 100.0 * np.sin(2 * np.pi * f1 * t) 20.0 * np.sin(2 * np.pi * 5 * f1 * t 0.3) 10.0 * np.sin(2 * np.pi * 7 * f1 * t - 0.8) 4.0 * np.sin(2 * np.pi * 11 * f1 * t 1.2) ) rng np.random.default_rng(42) signal rng.normal(0, 0.5, sizeN) # 模拟测量噪声这里 fs10000、N10000数据窗长 1 秒频率分辨率是 1Hz。基波 50Hz 正好落在第 50 根谱线5 次谐波 250Hz 落在第 250 根谱线等于做了一个“同步采样”避免了频谱泄漏这一最麻烦的误差源。先把最简单的场景跑通再在第五节里故意破坏同步条件演示泄漏这是比较高效的学习路径。噪声幅值 0.5A 相对于 100A 基波是 0.5%不会明显影响谐波幅值读数但可以让代码里的峰值搜索更接近真实情况。随机种子固定成 42是为了结果可以复现课程报告里也应该固定随机种子否则每次运行结果不同答辩时不好解释。3.2 采样率与数据窗怎么定谐波分析里两个最基础但也最容易拍脑袋的参数是采样率和 FFT 点数。采样率决定能看到的最高频率根据奈奎斯特采样定理采样率至少大于最高分析频率的两倍。工程上不会卡在 2 倍通常取 5 到 10 倍以上。如果课程报告的分析范围到 50 次谐波也就是 2500Hz采样率取 10kHz 就足够了2500Hz 对应的谱线还在 5kHz 奈奎斯特频率之内而且留出了几倍余量前端的抗混叠滤波器也更好设计。FFT 点数 N 决定频率分辨率。频率分辨率 Δffs/N 越小相邻谱线之间的距离越近。判断 N 够不够的标准很直接能不能把你想区分的两个频率分量分开。如果信号里有 5 次谐波 250Hz还有 5.5 次间谐波 275Hz它们只差 25Hz那么 Δf 必须小于 25Hz最好在 5Hz 左右。取 fs10000、N2000Δf5Hz正好满足。常见谐波分析仪测 10 个周波再做一次 FFT也就是数据窗长度 0.2 秒本质上也是为了让频率分辨率达到 5Hz。这里给出一个参数对照表方便直接套用场景采样率 fsFFT点数 N频率分辨率 Δf数据窗时长备注演示同步采样10kHz100001Hz1s50Hz和整数倍谐波全部落在整谱线上常用谐波测量10kHz20005Hz0.2s10个工频周期可分辨相邻5Hz间隔关注到50次谐波12.8kHz25605Hz0.2s最高分析2500Hz奈奎斯特6400Hz含间谐波细分10kHz100001Hz1s间谐波分析时希望频率分辨率更高采样率定了之后还要确认前端有没有抗混叠滤波。现在很多采集设备内部自带但低成本采集卡不一定有如果信号里有超过 fs/2 的高频噪声它们会折叠回低频区频谱上凭空多出一堆假谐波。课程报告里没有硬件条件时至少要在文档里写清楚假设采集系统已具备抗混叠措施这一点很关键。3.3 三步预处理去直流、加窗、单边谱修正拿到原始信号后我不会直接丢给 FFT而是先做三步预处理去直流、加窗、取单边谱并做幅值修正。第一步去直流。电流信号如果有直流偏置FFT 的第 0 根谱线会很大而且直流分量可能泄漏到邻近低频谱线影响对低频谐波的判断。处理方式很简单计算整段数据的均值然后每个采样点减去这个均值。注意去直流的前提是直流分量不是我们关心的对象如果分析的目标本身就是直流偏置就不该减要把第 0 根谱线单独读出来。第二步加窗。FFT 处理的是有限长数据相当于加了矩形窗矩形窗主瓣最窄但旁瓣衰减慢一旦基波频率不是处在整谱线上旁瓣能量会污染周围谐波。汉宁窗旁瓣衰减快是谐波分析最常用的折中方案。加了汉宁窗之后主瓣变宽两个频率间隔太近时会分不开所以如果采样严格同步、且只关心整数次谐波用矩形窗也可以如果信号频率可能有偏移就用汉宁窗保底。第三步是单边谱修正。np.fft.fft 输出的是双边谱负频率部分和正频率部分对称。实信号的负频率没有额外信息通常只取正频率一半但幅值要做修正。同时加窗会减小信号能量幅值要按窗函数能量恢复。用汉宁窗时恢复系数不是简单的 2/N而是 N / sum(window)。我一般把这三步写成一个小流程# 1. 去直流 x signal - np.mean(signal) # 2. 加汉宁窗减小频谱泄漏 window np.hanning(N) xw x * window # 3. 单边FFT X np.fft.rfft(xw, nN) freqs np.fft.rfftfreq(N, d1.0 / fs) # 4. 单边谱幅值恢复直流不乘2其余正频率谱线乘2 scale N / window.sum() # 窗函数幅度恢复系数 mag np.abs(X) * 2.0 / window.sum() mag[0] np.abs(X[0]) / window.sum() # 直流分量不乘2 if N % 2 0: mag[-1] np.abs(X[-1]) / window.sum() # 奈奎斯特分量不乘2 # 5. 在幅值谱上找局部峰值提取谐波信息 def find_harmonic_peaks(freqs, mag, min_amp0.5): 最小幅值过滤噪声局部最大判定为谐波谱峰 peaks [] for k in range(1, len(mag) - 1): if mag[k] min_amp and mag[k] mag[k - 1] and mag[k] mag[k 1]: peaks.append((freqs[k], mag[k], k)) return peaks peaks find_harmonic_peaks(freqs, mag, min_amp0.5) for f_peak, amp, k in peaks: h round(f_peak / 50.0) if h 0: print(f检出 {h:2d} 次谐波频率 {f_peak:6.2f} Hz幅值 {amp:6.2f} A)这段代码里最容易出错的是mag[0]的修正。np.fft.rfft 输出的第 0 个元素是直流分量它没有对称的负频率部分所以不需要乘 2。如果 N 是偶数最后一个元素对应奈奎斯特频率同样不乘 2。其余正频率谱线因为双边谱里正负频率各有一半能量取单边时要把另一半加回来所以乘 2。加窗恢复系数N / window.sum()对矩形窗来说等于 1直接写成 2/N 也能用但换成汉宁窗后window.sum()N/2所以恢复系数是 2实际效果相当于幅值乘 4/N这比很多人直接套用的 2/N 正好大了一倍也是“为什么我算出的幅值只有理论值一半”最常见的翻车原因。峰值搜索的min_amp0.5A是一个经验阈值。它要比噪声幅值高同时低于最低关注的谐波幅值。如果现场存在比较大的背景噪声可以先画出幅值谱观察噪声底在什么量级再设定阈值。阈值设大了会漏掉小谐波设小了会把噪声峰当成谐波这是需要反复确认的玄学点之一没有固定公式可套。4. 从频谱到报告指标THD、各次谐波幅值与相位提取4.1 从频谱读到工程上要用的幅值和相位峰值搜索得到了每个谐波对应的频率和幅值但课程报告里还需要把幅值换算成工程上常用的有效值并判断谐波次数。电力系统中电压电流的报表通常用有效值 RMS而 FFT 幅值谱直接给出的是正弦分量的峰值。前面代码里的幅值修正是在把频谱幅度复原为时域正弦波的峰值。如果信号是纯正弦 100AFFT 修正后读出来是 100A对应峰值要得到有效值再除以 √2。工程里常说“电流 100A”通常指有效值但仿真设定时我写的是峰值所以报告里要明确标注单位是峰值还是有效值否则答辩时一定会被问倒。相位提取用np.angle(X[idx], degTrue)就能拿到对应谱线的初相位。这里有一个要注意的地方加窗会改变信号的幅度和相位即使做了幅度修正相位也并不是直接等于原始信号的初相位。具体来说汉宁窗是偶对称窗在同步采样、谱峰正好落在整根谱线上时它不改变该谱线的相位但一旦频率偏移相位误差就会随着偏移量增大。课程报告的谐波相位如果只用于谐波源定位、潮流计算这些后续方向建议报告各次谐波之间的相对相位差也就是以基波相位为参考减去基波相位作为差值。绝对相位受采集时刻的同步影响不同台设备之间没有一致的时间基准对比意义不大。相位提取的代码可以这样写# 找到基波和各次谐波在谱线数组中的索引 target_orders [1, 5, 7, 11] for h in target_orders: f_target h * 50.0 idx int(round(f_target / (fs / N))) # 频率换算到谱线序号 amp mag[idx] if h 1: eff amp / np.sqrt(2.0) # 峰值转有效值仅供参考 else: eff amp / np.sqrt(2.0) phase_deg np.angle(X[idx], degTrue) print(f{h:2d} 次峰值 {amp:6.2f} A有效值 {eff:6.2f} A相位 {phase_deg:7.2f} deg)这里idx int(round(f_target / (fs / N)))是最直接的频率到谱线换算公式它的前提是谐波频率正好落在整根谱线上。前面仿真选了 N10000、Δf1Hz所以 250Hz 对应 idx250完全准确。如果用 N2000、Δf5Hz250Hz 对应 idx50也一样准确。只要保证数据窗长度是工频周期整数倍这个换算就成立否则就要用第六章的插值方法修正。4.2 THD 与各次谐波含有率的计算逻辑衡量波形畸变程度最常用的指标是总谐波畸变率 THD。它的定义是各次谐波有效值的方和根与基波有效值的比值THD√(∑_{h2}^{N_h} I_h²) / I_1谐波含有率 HRU_hI_h/I_1×100%表示某一次谐波相对基波的占比。实际计算时累加范围通常是 2 到 50 次再高次的幅值已经很小而且集肤效应、线路阻抗都会让它们衰减工程仪器一般断到 50 次。课程报告里如果没有硬件仪器的分组滤波更保守的做法是只累加已经识别出的整数次谐波忽略噪声谱线否则把宽带噪声的每个频点都加进 THD数值会明显偏高甚至出现“基波 100A 噪声 1A 但 THD 算成 20%”的离谱结果。代码里我会把 THD 计算限定到一个显式的谐波次数列表fundamentals 100.0 # 由FFT峰值得到的基波幅值 detected_orders [2, 3, 5, 7, 11, 13] harmonic_rms_squares [] for h in detected_orders: idx int(round(h * 50.0 / (fs / N))) amp mag[idx] eff amp / np.sqrt(2.0) harmonic_rms_squares.append(eff ** 2) thd np.sqrt(sum(harmonic_rms_squares)) / (fundamentals / np.sqrt(2.0)) print(f总谐波畸变率 THD {thd * 100:.2f}%)注意代码里基波和所有谐波都换算成有效值再比或者全部用峰值比结果一致但分子分母混用就是错。这个 THD 计算方式只覆盖了检出列表中的谐波次数如果信号里其实还有 9 次谐波但没被列进来结果就不完整。所以更通用的做法是先跑一次峰值搜索再把所有峰值里 h≥2 的分量收集起来参与累加。课程报告中建议两种方法都试一下对比确认 THD 值的差异来源正好可以在报告里体现你对算法边界的理解。仿真信号的理论 THD 可以手算验证只有 5 次 20A、7 次 10A、11 次 4A都相对基波 100A。平方和为 20²10²4²516开根号约 22.7所以 THD 约 22.7%。加噪声后 FFT 读出的值会有微小偏差如果偏差超过 1%就要检查幅值修正或加窗是否正确。4.3 课程报告适合放的三张图、三张表谐波分析报告最忌只有一段结论和一堆数字。评审老师想看到的是从时域到频域的完整证据链。我会建议至少放三张图图一是原始畸变电流波形标注一个工频周期的长度让人直观看到波形不是正弦图二是幅值频谱柱状图横轴频率到 1000Hz 即可覆盖到 19 次谐波或者根据分析范围调整在每个峰上标注谐波次数图三是把数据切成多段窗口每段计算 THD画出 THD 随时间的趋势用来说明分析结论针对的是稳态还是存在波动。三张表对应为表一是算法与数据参数表写清采样率、FFT 点数、频率分辨率、窗函数类型、分析次数范围表二是各次谐波幅值、有效值、含有率、相位的汇总表表三是限值对照表列出各次谐波的含量百分比和对应运行限值标出哪些超标。限值参数我这里给一个示例仅用于课堂演示实际工程要以项目所在电网的最新运行标准为准谐波次数仿真设定峰值(A)FFT测量峰值(A)含有率(%)示例限值(%)1100.0约100.0100.0—520.0约20.020.04.0710.0约10.010.04.0114.0约4.04.02.0130约002.0表里 5 次谐波含有率 20% 远超示例限值 4%这正好能引出结论该场景下存在明显 5 次谐波超标可能来自整流类负载。课程报告要把这个结论和波形畸变对应起来而不是只写“FFT 成功分解了信号”。谐波分析的价值体现在判断来源和评估影响上频率读数本身只是一串数字。5. 谐波分析常见避坑指南频谱泄漏、采样率与幅值修正5.1 谱峰对不准同步采样看似正确频率一漂就翻车现象用 fs10000、N10000 处理 50Hz 基波时一切正常但当现场基波频率变成 49.8Hz同一段代码读出的基波幅值从 100A 掉到 96A5 次谐波从 20A 掉到 18.5A而且频谱图上基波峰旁边出现一排小旁瓣。原因FFT 的频率分辨率 Δffs/N谱线落在 49Hz、50Hz、51Hz 这些离散点上。49.8Hz 不在任何谱线正中能量被分散到邻近谱线这就是频谱泄漏。数据窗长度是 1 秒但窗口内实际不是整周期数首尾相位不连续矩形窗隐含的周期延拓假设被破坏。解决先确认基波频率到底是 50Hz 还是漂移了。常见做法是用锁相环或者过零检测先测出实际基波频率然后调整采样率或数据窗长度使窗口内包含整数个工频周期。更省事的办法是保留汉宁窗因为汉宁窗旁瓣衰减快能让泄漏能量不至于污染到相邻谐波但如果泄漏已经发生仅加窗不能完全恢复幅值需要更进一步做插值修正这在第六章演示。总之把“基波频率是整数”当成默认假设是谐波分析里最容易被反杀的坑。5.2 加窗后的幅值恢复忘掉窗函数恢复系数结果少一半现象用矩形窗跑同步采样信号9 次谐波幅值正确。把窗函数换成汉宁窗后所有谐波幅值变成原来的一半频谱形状没变但数值整体对不上。原因汉宁窗对信号做幅度调制窗内数据的平均能量大约是矩形窗的一半。如果不修正FFT 谱线幅度自然只有真实峰值的一半。我在第三章的代码里用了scale N / window.sum()对汉宁窗来说这个比值是 2等效幅值修正是2/N而矩形窗 window.sum()N比值是 1。所以问题出在把矩形窗的恢复系数套用到汉宁窗上。解决每次换窗函数都要重新计算恢复系数。写代码时不要把2/N写死在程序里而是先算scale N / window.sum()再用mag np.abs(X) * 2.0 / window.sum()。同时要注意直流和奈奎斯特频率分量的器没有乘 2这三处是同一个坑的三种表现做完后用一个已知幅值的正弦信号测试一遍确认读数回得来再进入正式分析。5.3 频率分辨率不足两个分量靠得近FFT 分不开现象信号里有 5 次谐波 250Hz 和 5.5 次间谐波 275Hz幅值分别是 20A 和 2A。用 fs10000、N1000 时频谱上 250Hz 附近只有一个宽峰幅值读到 20.5A但 275Hz 的间谐波完全看不到甚至让 250Hz 的幅值被抬高。原因fs/N10Hz250Hz 和 275Hz 只差 25Hz理论上能分开但 N1000 时数据窗只包含 0.1 秒也就是 5 个工频周期频谱主瓣宽两个峰靠得不够远加上 5 次谐波比间谐波大 10 倍小峰被大峰的旁瓣淹没了。解决提高频率分辨率拉长数据窗到 0.2 秒或 1 秒也就是 N2000 或 10000让 Δf 降到 5Hz 或 1Hz。如果现场只有一段短数据无法拉长窗口可以用零填充增大表观频率分辨率但要注意零填充不能真正提高谱线信息量它只是让谱线插值变细对两个真实分量靠太近导致的合并没有本质帮助。这时候需要考虑 6.2 节的插值方法或直接换成参数化估计方法比如 Prony 分析它不需要整周期数据也能估计频率和衰减因子但课程报告不要求掌握到那么深知道边界即可。5.4 直流分量、奈奎斯特分量两个容易忽略的特殊谱线现象信号带 0.5A 直流偏置FFT 后第 0 根谱线数值巨大低频区域出现一串虚假谱线或者采样率刚好是 2 倍分析频率时最高频率谱线幅值总偏低。原因第 0 根谱线对应直流分量在矩形窗泄漏下会把能量泄露到低频区。奈奎斯特频率分量没有负频率对称部分单边谱修正时如果也乘 2读值就会比实际大一倍反过来如果忘掉它读值又只有一半。解决直流偏置一律先通过减去均值处理。奈奎斯特分量在代码里用if N % 2 0单独处理不乘 2。在课程报告里体现这个细节评审印象分会明显不一样因为它说明你不是只调用了 FFT 函数而是真正理解了单边谱和幅值修正的含义。5.5 噪声阈值设置把噪声峰当谐波THD 虚高现象加噪信号的 THD 计算出来比理论值高好几个百分点查看频谱图发现 17 次、19 次等位置出现了幅值只有 0.3A 的“谐波峰”。这些峰对应频率不是任何典型谐波但峰值搜索把它们纳入了统计。原因高斯噪声的频谱不是一条平直线而是随机起伏某些频点的噪声幅值刚好超过阈值被误识别为谐波。THD 累加时把这些噪声峰当成独立谐波总和就会虚高。解决峰值搜索的min_amp阈值要根据噪声底设定一般取噪声底幅值的 2 到 3 倍。另一个更可靠的方式是只允许整数次谐波参与 THD 统计比如 h2,3,4,...,50把非整数频率点的峰全部忽略。这样即使噪声峰被检出也不会进入指标计算。两者结合后THD 与理论值的偏差可以控制在 0.5% 以内这在课程报告中就足够可信了。6. 短窗下的进阶技巧双谱线插值与基波频率偏移修正前面第四章和第五章都在处理“频率正好落在整谱线上”的理想场景但现实电网中基波频率并不是铁板一块高精度测量时频率偏移会直接影响谐波幅值的准确性。我一般会在课程报告的结尾部分加一个进阶验证用很短的数据窗故意让基波频率偏离整谱线然后用抛物线插值把峰值频率和幅值“救”回来。抛物线插值的思路是当真实峰值落在两根谱线之间时幅值谱上最大的那根谱线和它左右邻居的幅值之间存在一个比例关系可以用二次函数拟合真实峰的位置。设最大谱线索引为 k0对应幅值 y2左右相邻谱线幅值为 y1、y3那么频率偏移量为δ(y3−y1)/(2·(2·y2−y1−y3))修正后的峰值频率为 (k0δ)·fs/N修正后的幅值为 y2 − 0.25·(y1−y3)·δ。这个近似公式对于幅值谱满足局部抛物线形状的情况非常有效汉宁窗下的谱峰比矩形窗更圆润反而更适合抛物线插值。def parabolic_interp(mag, k0): 对幅值谱峰 k0 做抛物线插值返回(频率偏移量, 修正幅值) if k0 0 or k0 len(mag) - 1: return 0.0, mag[k0] y1, y2, y3 mag[k0 - 1], mag[k0], mag[k0 1] denom 2.0 * (2.0 * y2 - y1 - y3) if abs(denom) 1e-12: return 0.0, y2 delta (y3 - y1) / denom amp_corr y2 - 0.25 * (y1 - y3) * delta return delta, amp_corr # 构造基波频率 49.8Hz 的短窗信号 fs 10000 N_short 2000 # 0.2秒10个工频周期Δf5Hz t_short np.arange(N_short) / fs f_off 49.8 signal_off 100.0 * np.sin(2 * np.pi * f_off * t_short) 20.0 * np.sin(2 * np.pi * 5 * f_off * t_short) win np.hanning(N_short) X_off np.fft.rfft(signal_off * win) freq_off np.fft.rfftfreq(N_short, d1.0 / fs) mag_off np.abs(X_off) * 2.0 / win.sum() k0 np.argmax(mag_off[:100]) # 在0-500Hz范围找基波峰 delta, amp_corr parabolic_interp(mag_off, k0) f_est (k0 delta) * fs / N_short print(f基波估计频率 {f_est:.3f} Hz修正前幅值 {mag_off[k0]:.2f} A修正后幅值 {amp_corr:.2f} A)代码输出会显示基波估计频率约 49.8Hz而修正前幅值可能只有 97A 左右修正后回到约 100A。这个对比非常适合放在课程报告的实验验证部分因为它直观说明了一件事FFT 不是不准确而是直接用离散谱线读离散值时不准确插值修正就是把这层误差剥掉。我第一次做谐波分析时就吃过这个亏。某个模拟项目里我用矩形窗直接 FFT基波频率实际是 49.9Hz我没做频率确认就把 50Hz 整数倍谱线的幅值当成基准去算 THD结果基波幅值偏低 3%THD 虚高 0.8%被设备厂商问了一个下午。后来形成习惯每次报告频谱结果前先回答三个问题——采样率是多少、窗口类型是什么、幅值恢复系数用了多少。三个问题答完谐波分析的意外基本就都在控制范围内了。希望这个流程和踩坑记录能帮到你让你少走这一段弯路。本文还有配套的精品资源点击获取