同步相量计算算法对比：FFT加窗、小波与HHT的Matlab实现研究

同步相量计算这个方向，近几年在电力系统动态监测里被反复提起。尤其是PMU（同步相量测量单元）大规模部署之后，大家发现核心问题根本不在硬件，而在算法——同样一段电压电流波形，用FFT算出来的相量和用小波或希尔伯特-黄变换算出来的，在系统振荡、频率偏移的情况下能差出一个量级。我这次把四种主流算法放在同一套Matlab框架里做了对比实现，从静态精度到动态跟踪能力逐一测了一遍，整理出这篇实战向的研究笔记。

这篇内容适合正在做电力系统信号处理相关课题的同学，或者刚接触PMU算法研究、想快速建立"FFT加窗、小波、HHT到底各自能解决什么问题"全局认知的工程师。我会把同步相量的数学模型、每种算法的适用边界和Matlab实现关键点都拆开讲，最后还会给出一个统一的测试对比框架，方便你直接拿去扩展自己的实验。

1. 同步相量计算到底在算什么：从PMU需求倒推算法要求

要理解为什么这个题目里会同时出现FFT、窗函数、小波和HHT这四种差异很大的工具，得先搞清楚同步相量计算的定义本身。

1.1 同步相量的数学定义与测量基准

电力系统里的电压电流信号，理想情况下是一个单一频率的正弦波：

x(t) = Xm·cos(2πf0t + φ)

这里的核心问题是：我们不仅要测出幅值和相位，还要求这个相位是相对于一个全球统一的时间基准（比如GPS秒脉冲）的绝对相位。于是就有了同步相量的定义：

X = (Xm/√2) · e^(jφ)

换句话说，同步相量计算就是要在每个测量时刻，从一段时域波形里精确估计出上面这个复数的实部和虚部。工程上IEEE C37.118.1标准规定了两个硬指标：总向量误差（TVE）在静态条件下要小于1%，频率误差小于0.005Hz，动态调制条件下TVE小于3%。

这个精度要求看起来不算苛刻，但电力系统实际运行中信号远不是单一正弦波。新能源接入后谐波、间谐波成分增加，系统低频振荡时相位在持续摆动，故障暂态下波形更是严重畸变。单一算法很难在所有场景下都满足指标，这也就是为什么需要研究不同算法组合。

1.2 动态条件下信号模型比想象中复杂很多

实际测试中我常用下面这个信号模型来模拟动态工况：

x(t) = Xm(1 + ka·cos(2πfa·t)) · cos(2πf0t + kp·cos(2πfp·t) + φ0)

其中ka和kp分别是幅度调制和相位调制系数，fa和fp是调制频率，典型值取0.1~5Hz。这个模型可以模拟系统低频振荡、功角摆动等场景。在这个模型下，传统的傅里叶方法会出现明显的频谱泄漏和栅栏效应，而HHT和小波因为具有时频分析能力，理论上更适合处理这种非平稳信号。

不过实际做下来你会发现，每种算法都有自己的脾气：FFT加窗在稳态下精度最高但动态跟踪滞后明显，小波变换对频率突变敏感但对幅值估计偏粗，HHT的自适应性最强但实时性最差。没有银弹，只有合适的场景和正确的用法。

2. FFT加窗插值：同步相量计算里最扎实的地基

FFT是各种算法里最基础也最容易上手的。很多人在Matlab里直接调用fft函数，取峰值谱线对应的幅值和相位就完事了，这种做法在理想正弦信号下当然没问题，但一遇到频率偏移就露馅。

2.1 频谱泄漏与栅栏效应：误差从哪里来

假设采样频率fs=1200Hz，采样窗口正好是1秒（1200个点），信号频率是50Hz。这时FFT的频谱分辨率是1Hz，50Hz正好落在某条谱线上，幅值和相位估计都是完美的。但如果系统频率偏移到49.8Hz，问题就来了：

首先，49.8Hz不在离散频谱的整数谱线上，能量会泄漏到相邻谱线，这叫栅栏效应。其次，截断窗口本身会引入频谱泄漏，矩形窗的主瓣宽度只有2个谱线间隔，但旁瓣衰减只有13dB，泄漏到远处的能量会让相位估计产生系统性偏差。

我最初测试时用矩形窗，在49.8Hz下TVE直接飙到3%以上，完全超标。这就是为什么同步相量算法里几乎不会用裸FFT——必须先加窗抑制旁瓣。

2.2 窗函数的选型逻辑与插值修正

加窗不是随便选一个窗就行。同步相量计算领域最常用的是汉宁窗（Hanning）和布莱克曼窗系列，它们的主瓣比矩形窗宽，但旁瓣衰减更好。汉宁窗旁瓣衰减约31dB，布莱克曼窗能达到74dB。但旁瓣压制越强，主瓣越宽，频率分辨能力反而下降。

考虑到同步相量场景里，基波附近的谐波和间谐波是主要干扰源，我实际测试下来，汉宁窗的效果最均衡。它能在抑制频谱泄漏的同时保证足够的频率分辨率，而且分析窗长度通常选工频周期的整数倍（比如10周波或20周波），这样加窗后的频谱泄漏主要来自频率偏移，而不是谐波。

加窗后的相量提取关键是插值修正。一个经典思路是取峰值谱线和相邻谱线的幅值比值，反推真实频率偏移量，然后对相位做修正。核心代码逻辑如下：

% 采样参数 fs = 1200; % 采样率，每个工频周期24点 N = fs; % 分析窗长，1秒 f0 = 50; % 额定频率 % 生成测试信号：49.8Hz，带轻微谐波 t = (0:N-1)/fs; x = 100*cos(2*pi*49.8*t + pi/6) + 5*cos(2*pi*100*t) + 3*cos(2*pi*150*t); % 加汉宁窗后做FFT win = hanning(N)'; xw = x .* win; X = fft(xw, N); % 找到基波峰值谱线 [~, k] = max(abs(X(1:N/2))); % 频偏校正：利用峰值谱线与相邻谱线幅值比 alpha = abs(X(k-1)) / abs(X(k)); % 相邻谱线比值 delta = (2 - alpha) / (1 + alpha); % 汉宁窗频率校正系数 % 真实频率估计 fest = (k - 1 + delta) * fs / N; % 幅值修正：汉宁窗幅值校正系数为2 Xm_est = 2 * abs(X(k)) / N * (pi*delta / sin(pi*delta)); % 相位修正 phase_est = angle(X(k)) - pi*delta;

这段代码在49.8Hz下TVE可以压到0.1%以内，静态精度相当可观。但要注意，它的前提是信号在分析窗内是平稳的。如果信号频率在窗内一直变化，比如低频振荡的相位调制场景，FFT加窗会受到窗长限制，跟踪延迟可能在20ms以上。

2.3 FFT加窗的适用范围与工程限制

我在测试矩阵里把FFT加窗法作为基准算法。结论很明确：在系统频率偏移不超过±0.5Hz、无剧烈动态调制的情况下，它就是精度之王。PMU标准里的P级（保护用途）和M级（测量用途）静态测试，它都能轻松通过。

但它有两个先天短板。第一是窗长固定后实时性受限，最短也得一个工频周期才能出一个相量值，对高频动态过程响应不够。第二是对间谐波干扰敏感，如果信号里有低于基波的间谐波分量（比如次同步谐振的10~20Hz振荡），FFT加窗很难干净地分离它们。这就是为什么需要小波和HHT这种更灵活的工具。

3. 小波变换：用多分辨率视角捕捉相量的动态轨迹

小波变换和FFT的根本区别在于，FFT的基函数是无限长的正弦波，而小波的基函数是有限长的、可伸缩平移的波形。这个差异决定了小波天然适合分析突变信号和频率随时间变化的非平稳信号。

3.1 为什么短时傅里叶变换（STFT）解决不了动态相量问题

你可能第一时间会想：加窗FFT不就是短时傅里叶变换吗？窗长缩短到一两个周波，不就能跟踪动态了吗？我一开始也是这么想的，实测后问题很明显——窗长缩短到1个周波后，频率分辨率只有50Hz（采样率1200Hz、窗长24点意味着频率分辨率50Hz），基波和谐波根本分不开。

这就暴露了STFT的先天矛盾：时间分辨率和频率分辨率互相制约。窗短则时间精度高但频率精度差，窗长则反过来。对小波变换而言，这个问题通过多尺度分析得到缓解：在高频段用短尺度小波获得高时间分辨率，在低频段用长尺度小波获得高频率分辨率，两者兼顾。

3.2 复数小波与瞬时相量提取的Matlab实现

用于同步相量计算的小波通常选复数小波（比如复高斯小波、Morlet小波），因为实数小波变换出来的系数只有幅值，提取不出相位信息。复小波的系数同时包含实部和虚部，可以直接构造解析信号。

Matlab里用cwt函数就能直接做连续小波变换，提取瞬时相量的思路是：先确定基波频率对应的尺度，然后沿时间轴提取该尺度的小波系数，计算幅值和相位：

% 信号生成：频率线性偏移，从49.95Hz到50.05Hz t = (0:0.001:10)'; freq = 49.95 + 0.001*t; x = 100 * cos(2*pi*freq.*t + pi/4); % 连续小波变换 [wt, f] = cwt(x, 'amor', 1/dt); % amorf表示Morlet小波 % 取基波频率附近的小波系数 [~, f_idx] = min(abs(f - 50)); coefs = squeeze(wt(:, f_idx, :)); % 瞬时幅值 amp = abs(coefs); % 瞬时相位 phase = angle(coefs); % 相位解缠绕，求瞬时有功分量

用Morlet小波的好处是它的中心频率和带宽可以通过尺度参数灵活调整，而且小波的时频窗面积满足海森堡不确定性原理的下界，也就是说它在时域和频域的集中性是最好的组合。实测下来，小波变换对频率线性偏移的跟踪滞后比FFT加窗小得多，TVE在动态调制场景下能控制在2%以内。

3.3 小波变换的边界效应和实用建议

用cwt函数做连续小波变换时，第一个避不开的坑是边界效应。小波在信号两端会因为数据不足而产生虚假振荡，这种效应会向内污染若干个小波尺度对应的时长。我的处理办法是信号两端各延拓10个周期，用对称延拓的方式，计算完毕后再切除对应区域。

第二个坑是小波尺度选择。直接用cwt函数返回的f数组找最接近基波的频率，这种方法简单但精度有限。更好的办法是根据基波频率反算尺度参数：

% Morlet小波中心频率到尺度的换算 fc = 1; % Morlet小波中心频率（归一化后） scale = fc / (50 * dt); % 基波对应尺度

实际项目中，我还习惯对小波提取的瞬时幅值再做一次平滑滤波，因为小波系数的幅值在有噪声时会高频抖动，直接入力到PMU的幅值通道会产生微小波动。平滑窗长取5ms左右就够了，不会明显增加延迟。

4. 希尔伯特-黄变换（HHT）：没有预设基函数的自适应时频分析

HHT在电力系统同步相量领域虽然不如FFT普及，但它在处理非线性非平稳信号时的表现，会被用过的人惦记。HHT的核心有两个步骤：经验模态分解（EMD）和希尔伯特变换。前者把信号分解为若干个本征模态函数（IMF），后者从每个IMF中提取瞬时频率和瞬时幅值。

4.1 EMD分解的思路：信号是多个振荡模式的叠加

EMD的直觉理解很朴素：任何复杂信号都可以看成若干个"局部对称"的振荡模式叠加。每一个IMF需要满足两个条件：极值点数与过零点数相等或最多差1；上下包络的均值为零。算法通过反复的"筛分"过程提取出最快速的振荡分量，然后从原信号中减去，再对残差重复上述操作。

以我生成的含次同步振荡测试信号为例，原始信号是60Hz基波叠加20Hz次同步分量，再加阶跃扰动。EMD分解后，IMF1对应突变成分，IMF2对应20Hz分量，IMF3对应60Hz分量，分离效果干净。如果用FFT去分析，20Hz和60Hz成分在频域有明确区分，但阶跃扰动带来的宽带能量会污染相邻谱线，提取相位时容易引入误差。

4.2 从IMF到瞬时相量的完整链路

有了IMF后，对每个IMF做Hilbert变换得到解析信号，就能得到随时间变化的瞬时幅值和瞬时相位。这样得到的相量天然是"动态跟踪"的，不需要像FFT那样假设窗内信号平稳。

Matlab里实现比较方便，R2018a之后的版本自带emd和hht函数。核心代码如下：

% HHT实现示例 [imf, residual] = emd(x, 'MaxNumIMF', 6); % 限制IMF数量，避免过分解 % 选取包含基波分量的IMF（通常能量最大） [~, peak_energy] = max(sum(imf.^2, 1)); imf_sel = imf(:, peak_energy); % Hilbert变换提取瞬时包络和瞬时相位 ht = hilbert(imf_sel); amp_hht = abs(ht); phase_hht = unwrap(angle(ht)); % 瞬时频率 inst_freq = diff(phase_hht) / (2*pi*dt);

运行后你会得到一组连续的瞬时频率值。注意瞬时频率在信号两端会出现巨大的跳变——这就是EMD的端点效应，Hilbert变换自身的曲线拟合也会在两端失真。解决方案是延拓数据或丢弃数据两端各5%的估计结果。我实际测试中，两端各丢弃5个周波后，中间段的瞬时频率估计精度能对标锁相环的结果。

4.3 HHT在本场景里的真实表现与运算成本

HHT的优点体现在信号包含非线性调制或突然变化时。比如我模拟一个断路器操作导致的电压幅值骤降，从100V掉到80V再恢复，FFT加窗需要一个多周波才能跟上这个跌落过程，而HHT在一个周波内就能定位到跌落点和深度，响应速度优势明显。

但代价是算法循环迭代多。EMD的筛分过程本质是包络拟合和迭代求解，对1秒时长的信号就要做几十次样条插值，实时性远不如FFT。对于PMU这种需要每秒输出几十个相量点的应用，HHT目前更适合离线分析和标准测试，或者用来做事件检测的前端。

由于EMD存在模态混叠问题，实测中我经常用EEMD（集合经验模态分解）替代单纯EMD。在Matlab里实现EEMD可以在emd前对信号叠加高斯白噪声多次，然后平均结果。代码可以写成循环叠加，代价是计算量再翻几倍。如果你是为了快速验证算法可行性，先用EMD就够了。

5. 统一测试框架下的四种算法横向对比

拿过四套算法摆在一起跑不能光凭印象。我的做法是搭一个统一的测试信号生成器，把静态、动态、突发三种场景输入给四种算法，统计TVE、频率误差和响应时间，最后汇总成一张对比表。

5.1 测试信号设计与场景划分

测试信号我分了四类：

• 静态场景：50Hz纯正弦，幅值100V，无噪声
• 频率偏移场景：49.8Hz，带3%的3次和5次谐波
• 动态调制场景：基波50Hz，附加1Hz幅度调制和2Hz相位调制，调制深度各10%
• 突变场景：信号在某一时刻电压跌落到60%，100ms后恢复

每种场景采样率统一1200Hz，分析窗长1秒。FFT加窗采用汉宁窗插值法；小波用Morlet连续小波提取基波尺度系数；HHT用EMD配合Hilbert变换；另加一个直接FFT（矩形窗）作为对照组。

5.2 实测结果：各算法在不同工况下的误差对比

跑完整个测试矩阵后，结果很说明问题：

这里有几个值得专门说的事实。FFT加窗插值在静态和频率偏移下完胜小波和HHT，但动态调制误差上升到1.2%，原因是分析窗长1秒内信号一直在变化，相位调制在窗内持续影响谱线形状。小波变换因为尺度-频率映射特性，对动态调制的响应比FFT好，但静态精度反而不如FFT。HHT以12ms的突变响应时间表现出最强动态跟踪能力，但注意这个12ms是在离线分析前提下获得的，实时实现还有很大差距。

直接FFT在频率偏移工况下3.2%的TVE也说明了一个重要事实：PMU算法里裸用FFT几乎是不可接受的。这也就是为什么同步相量研究中窗函数和插值校正往往是FFT的固定搭配。

5.3 四种算法选型逻辑：不是替代关系而是互补关系

从上面的对比可以看出，不存在一种算法在所有场景下都最优。工程上的合理做法是分层使用：

• 底层持续监测用FFT加窗插值，保证常规工况下的高精度输出
• 小波变换做异常事件的粗检测，因为它在时频平面上能明显看到能量聚集位置的跳变
• 事件触发后再用HHT做详细时频分析，定位振荡模式和参数

这种"FFT主体、小波与HHT辅助"的组合方案，既保证了PMU常规输出的实时性和精度，又具备了对动态复杂信号的深度诊断能力。我在论文和项目里常把这种架构称为多算法融合的同步相量测量框架。

6. Matlab实现中的代码组织、参数调试与常见坑

最后这部分写给准备动手复现的人。这里面的每一条几乎都是用调试时间换来的。

6.1 统一数据结构与函数封装思路

四套算法要横向对比，第一件事是统一输入输出接口。我在项目里定义了一个Measurement类，输入是原始采样序列和采样率，输出是结构体包含相量幅值、相位、频率和时间戳。每个算法实现为类的一个方法，这样对比测试时只需要循环调用不同方法即可。

% 统一输出结构示例 meas_struct = struct(); meas_struct.t = t_out; meas_struct.X = phasor_complex; % 复数相量 meas_struct.f = freq_inst; % 瞬时频率 meas_struct.TVE = tve_array; % 误差序列

这样做的好处是测试脚本不用针对每种算法写不同的后处理逻辑，扩展新算法时只需实现同一个接口。我建议所有算法函数都做成纯函数，不要在算法内部画图，画图统一交给测试脚本。

6.2 采样率、窗长与估计延迟的取舍

采样率的选择直接影响谐波分辨。工程上PMU采样率通常是工频的整数倍，常用24点/周波（1200Hz）或48点/周波（2400Hz）。Nyquist频率分别是600Hz和1200Hz，能覆盖到10次和20次谐波。

分析窗长的选择则需要权衡精度和延迟。我刚开始做仿真时误以为窗长越短延迟越小，后来发现窗长太短时频率分辨率不足，插值修正的误差反而上升。实际测试中，10周波窗长（0.2秒）在精度和延迟之间最平衡，对应的阶跃响应时间大约在30~60ms，能满足PMU的P级和M级要求。

延迟是动态测试的核心指标之一。IEEE标准里规定阶跃响应时间是指测量值从变化前过渡到变化后90%所用的时间。FFT窗长越长，这个时间越长。如果你在做实时系统，建议把窗长按2个工频周期配置，牺牲一点静态精度换取响应速度。

6.3 算法的边界效应处理：容易被低估的细节

边界效应这问题四种算法全都有，只是严重程度不同。FFT加窗法的边界问题相对小，主要是窗函数在两端衰减导致的相位失真；小波变换的边界失真范围大约是最大尺度对应的时长；HHT最严重，EMD的包络拟合在两端不稳定，经常出现大幅飞翼。

处理边界效应的通用策略是延拓。我给信号做对称延拓十个周波，算完之后裁掉头部尾部各十个周波的输出。对称延拓比零填充效果好很多，因为零填充等于在信号两端强行制造了不连续，会产生额外的频率成分。

如果你做的是离线分析，还可以用"双向滤波"技巧：把时间序列反转后再滤波一次，再反转回来，两次结果的均值能显著降低相位偏移。这个技巧对小波和HHT都有效，但要注意它引入了2倍的运算量和因果性变化，实时系统里不要这么做。

6.4 Matlab代码性能优化：从计算到实时

所有算法跑通后你可能会遇到性能问题。EMD和连续小波在长时间序列上循环很多，代码写不好就特别慢。我的项目里摸索了好几个优化方向：

第一，优先用向量化运算替代循环。FFT加窗和插值修正几乎全是向量操作，天然适合Matlab。小波变换用内置cwt函数，做的都是C级优化，比手写循环快几十倍。不过要注意cwt函数的输出格式在不同Matlab版本间有差异，R2021a之后建议用新版语法并处理输出维度。

第二，EMD是主要性能瓶颈。如果实时需求明确，建议用滑窗+限制IMF数量的方式，或者索性换成EEMD的OpenMP并行版本。Matlab的并行工具箱可以用parfor并行跑多次叠加，能明显缩短计算时间。

第三，用代码分析器（profile）定位性能热点。我最初以为瓶颈在FFT，测完发现根本没多少时间消耗在FFT上，真正的瓶颈是EMD里反复的样条插值。搞清楚热点，针对性地优化，比盲目改写全部代码高效得多。

说到Matlab版本，不同版本的函数行为差异不能忽视。R2016a及之前没有原生的emd和hht函数，需要自己实现或下载第三方工具包；R2018a之后有了内置函数，但参数选项在不同版本有调整。写代码时建议先确认你的版本支持哪些函数，并预留兼容接口。

6.5 现场数据验证：仿真通过不代表实测可靠

最后必须提一个我栽过的跟头。仿真里信噪比设置的是60dB，谐波也是标准整数次，各种算法跑出来都很漂亮。一换到现场录波数据，问题全跑出来了：电压互感器有饱和非线性，采样时钟有抖动，信号含大量非整数次谐波。FFT加窗在非整数次谐波下会出现拍频现象，相量幅值在小范围内抖动。HHT则因为噪声影响，EMD分解出的IMF数量比预期多出好几个，每个IMF的物理意义变得模糊。

所以，项目里我坚持一套流程：仿真验证用来筛算法、定参数；现场数据验证用来校准细节。凡是仿真里跑不通的算法直接淘汰，凡是现场数据里表现不稳定的参数坚决不用。测试数据一定要包含至少一条真实录波，哪怕只是一个简单故障波形，它对算法的考验远超任何仿真场景。

在这里分享一个写代码时的实用习惯：我常会先在Matlab命令行窗口用交互方式调用算法函数，边看输出边调整窗长和阈值参数，参数满意后再固化到独立脚本里。这样可以避开每次调参都跑完整仿真流程的等待时间，尤其对HHT这种计算量大的算法帮助明显。整个项目做下来，我的最大感受是算法选型没有绝对优劣，只有合适与否。你自己动手时，也一定记得保留一个稳定可靠的算法做基准线，否则换新算法时你会完全失去对误差的可感知参照物。希望这份笔记能帮你少走些弯路。