Hilbert-Huang变换赋能有源电力滤波器：理论、应用与性能优化

Hilbert-Huang变换赋能有源电力滤波器：理论、应用与性能优化一、引言1.1研究背景与意义1.1.1谐波污染问题的严重性随着现代工业和电力电子技术的飞速发展，电力系统中的谐波污染问题日益严重。大量非线性负载，如变频器、整流器、电弧炉等设备的广泛应用，使得电网中的电流和电压波形发生畸变，产生了丰富的谐波成分。据统计，在一些工业发达地区，电网中的谐波含量已经超过了国家标准的限制，严重影响了电能质量和电力系统的安全稳定运行。谐波对电能质量的危害主要体现在以下几个方面。首先，谐波会导致电压畸变，使得电网中的电压波形偏离正弦波，从而影响电气设备的正常运行。例如，谐波会使电动机的转矩脉动增加，降低其效率和使用寿命；会使变压器的铁芯损耗增加，导致其过热，甚至损坏。其次，谐波会增加线路损耗，降低电网的传输效率。由于谐波电流在输电线路中会产生额外的功率损耗，使得电网的传输效率降低，造成能源的浪费。此外，谐波还会对通信系统产生干扰，影响通信质量。谐波电流会在输电线路周围产生电磁场，与通信线路发生电磁耦合，从而干扰通信信号的传输。谐波对设备运行的影响也不容忽视。许多电气设备在设计时是基于正弦波电源的，当电网中存在谐波时，这些设备的性能会受到严重影响。例如，谐波会使电容器发生过电压和过电流，导致其损坏；会使继电器误动作，影响电力系统的保护和控制功能。此外，谐波还会对电子设备产生干扰，如计算机、通信设备等，使其工作不稳定，甚至出现故障。谐波对系统稳定性的威胁同样不可小觑。在电力系统中，谐波可能会引发谐振现象，导致系统的电压和电流急剧增大，甚至引发系统崩溃。此外，谐波还会影响电力系统的动态响应特性，降低系统的稳定性和可靠性。1.1.2有源电力滤波器的重要性有源电力滤波器（APF）作为一种能够动态抑制谐波、全面改善电能质量的电力电子装置，在谐波治理中发挥着关键作用。APF的基本原理是通过实时检测电网中的谐波电流，然后产生与之大小相等、方向相反的补偿电流，注入电网中，从而抵消谐波电流，达到改善电能质量的目的。与传统的无源滤波器相比，APF具有诸多优势。首先，APF能够对变化的谐波进行实时跟踪和补偿，具有良好的动态响应性能。而无源滤波器的补偿特性是固定的，难以适应谐波的变化。其次，APF可以同时补偿多次谐波，而无源滤波器通常只能针对特定次数的谐波进行补偿。此外，APF不受电网阻抗的影响，不会与电网发生谐振，而无源滤波器在某些情况下可能会与电网发生谐振，导致谐波放大。APF的应用范围非常广泛，涵盖了工业、商业和居民等各个领域。在工业领域，APF可以用于谐波污染严重的工厂，如钢铁厂、化工厂、电子厂等，有效地改善电能质量，提高设备的运行效率和可靠性。在商业领域，APF可以用于商场、写字楼等场所，减少谐波对电梯、空调、照明等设备的影响，提高供电质量。在居民领域，APF可以用于住宅小区，降低谐波对家用电器的损害，提高居民的生活质量。1.1.3Hilbert-Huang变换引入的必要性在有源电力滤波器中，准确、实时地检测谐波电流是实现有效谐波补偿的关键。传统的信号处理方法，如快速傅里叶变换（FFT）、小波变换等，在分析电力信号时存在一定的局限性。FFT适用于平稳信号的分析，对于电力系统中存在的大量非线性、非平稳信号，其分析效果不佳，容易出现频谱泄漏和栅栏效应，导致谐波检测精度降低。小波变换虽然在一定程度上能够处理非平稳信号，但由于其基函数的固定性，对于复杂的电力信号，其适应性有限，难以准确地提取信号的特征。Hilbert-Huang变换（HHT）是一种新兴的非线性、非平稳信号处理方法，具有独特的优势。HHT主要由经验模态分解（EMD）和Hilbert谱分析两部分组成。EMD能够将复杂的非平稳信号自适应地分解为一系列固有模态函数（IMF），每个IMF都具有不同的特征尺度和频率成分，能够很好地反映信号的局部特征。然后，对每个IMF进行Hilbert变换，得到信号的Hilbert谱，从而能够在时频域上对信号进行全面、准确的分析。将HHT引入有源电力滤波器中，能够显著提升APF的性能。HHT可以更准确地检测电力信号中的谐波成分，尤其是对于那些频率和幅值随时间变化的谐波，能够实现高精度的检测。通过HHT对谐波电流的准确检测，APF可以生成更精确的补偿电流，从而提高谐波补偿的效果，进一步改善电能质量。此外，HHT还可以用于APF的故障诊断和状态监测，通过对APF运行过程中的信号进行分析，及时发现故障隐患，提高APF的可靠性和稳定性。1.2国内外研究现状在Hilbert-Huang变换理论研究方面，国外学者Huang等人于1998年首次提出了HHT方法，引起了广泛的关注。此后，众多学者对HHT的理论和应用进行了深入研究。Wu和Huang对EMD的特性进行了详细分析，研究了EMD分解过程中的端点效应、模态混叠等问题，并提出了相应的改进方法。Pérez-Cruz等人对HHT进行了全面的综述，总结了HHT在信号处理领域的应用和发展趋势。在国内，也有许多学者对HHT进行了研究。例如，赵春晖等人对HHT在电力系统故障诊断中的应用进行了研究，提出了基于HHT的故障特征提取方法，取得了较好的效果。在有源电力滤波器技术发展方面，国外的研究起步较早，技术相对成熟。ABB、西门子等公司在APF的研发和生产方面处于领先地位，其产品具有高性能、高可靠性等特点。国内对APF的研究始于20世纪90年代，经过多年的发展，取得了显著的成果。一些高校和科研机构，如清华大学、浙江大学等，在APF的控制策略、拓扑结构等方面进行了深入研究，提出了许多新的方法和技术。同时，国内也涌现出了一批优秀的APF生产企业，如荣信电力电子股份有限公司、上海思源电气股份有限公司等，其产品在国内市场占据了一定的份额。在Hilbert-Huang变换在有源电力滤波器中应用的研究方面，国内外学者也进行了大量的工作。国外学者Kashif等人将HHT应用于APF的谐波检测中，通过实验验证了HHT在谐波检测中的有效性和准确性。国内学者李鹏等人提出了一种基于HHT的APF控制策略，该策略能够根据电网谐波的变化实时调整APF的补偿电流，提高了APF的补偿性能。然而，当前的研究仍存在一些不足之处。在HHT理论方面，EMD分解过程中的端点效应和模态混叠问题尚未得到彻底解决，这在一定程度上影响了HHT的分析精度。在APF技术方面，如何进一步提高APF的效率、降低成本、增强可靠性，仍然是需要解决的问题。在HHT在APF中应用方面，如何将HHT与APF的控制策略更好地结合，实现APF的智能化控制，还有待进一步研究。1.3研究目标与内容本研究的目标是通过引入Hilbert-Huang变换，提升有源电力滤波器的性能，实现更精准的谐波检测与补偿，有效改善电能质量。具体研究内容包括以下几个方面：Hilbert-Huang变换原理分析：深入研究HHT的基本原理，包括EMD分解和Hilbert谱分析的过程和特点。分析EMD分解过程中可能出现的端点效应和模态混叠问题，探讨相应的改进方法，以提高HHT的分析精度。有源电力滤波器系统建模：建立APF的数学模型，包括主电路拓扑结构和控制策略的模型。分析APF的工作原理和性能特点，为后续的研究提供理论基础。基于Hilbert-Huang变换的APF算法改进：将HHT应用于APF的谐波检测中，提出基于HHT的谐波检测算法。结合APF的控制策略，对算法进行优化，提高APF的动态响应性能和补偿精度。实验验证：搭建APF实验平台，对基于HHT的APF算法进行实验验证。通过实验数据的分析，评估算法的性能，验证研究的有效性和可行性。1.4研究方法与技术路线本研究采用理论分析、仿真模拟与实验验证相结合的研究方法。理论分析：对Hilbert-Huang变换和有源电力滤波器的相关理论进行深入研究，分析其工作原理和性能特点，为后续的研究提供理论基础。仿真模拟：利用MATLAB/Simulink等仿真软件，建立基于HHT的APF仿真模型，对算法进行仿真分析。通过仿真，优化算法参数，验证算法的可行性和有效性。实验验证：搭建APF实验平台，进行实验研究。通过实验数据的采集和分析，评估算法的实际性能，进一步改进和完善算法。技术路线图如下：首先，进行文献调研，了解国内外在Hilbert-Huang变换和有源电力滤波器领域的研究现状和发展趋势。然后，对HHT的原理进行深入研究，分析其在处理电力信号时的优势和存在的问题，并提出相应的改进方法。接着，建立APF的数学模型，分析其工作原理和性能特点。在此基础上，将HHT应用于APF的谐波检测中，提出基于HHT的APF算法，并进行仿真分析。根据仿真结果，优化算法参数，进一步提高算法的性能。最后，搭建APF实验平台，进行实验验证，通过实验数据的分析，评估算法的实际性能，总结研究成果，提出未来的研究方向。二、Hilbert-Huang变换理论基础2.1Hilbert-Huang变换概述Hilbert-Huang变换（Hilbert-HuangTransform，HHT）是由美国国家航空航天局（NASA）的NordenE.Huang等人于1998年提出的一种全新的非线性、非平稳信号分析方法。在传统的信号处理领域，傅里叶变换长期占据主导地位，它基于信号是线性和平稳的假设，通过将信号分解为不同频率的正弦和余弦波的叠加，来分析信号的频率成分。然而，在实际应用中，许多信号，如电力系统中的谐波信号、生物医学中的心电信号、地震工程中的地震波信号等，往往具有非线性和非平稳的特性，傅里叶变换在处理这些信号时存在很大的局限性。例如，傅里叶变换无法准确地反映信号的局部特征，对于非平稳信号的分析结果可能会产生严重的偏差。为了克服传统信号处理方法的局限性，学者们不断探索新的方法。小波变换便是其中一种重要的改进，它通过选择合适的小波基函数，能够在一定程度上处理非平稳信号，对信号进行多分辨率分析，提供了信号在不同时间尺度上的信息。但是，小波变换的基函数是固定的，在面对复杂多变的信号时，其适应性仍然有限。HHT的出现为非线性、非平稳信号的处理提供了新的思路和方法。它打破了传统方法对信号线性和平稳性的假设，能够更加准确地分析信号的时频特性。HHT主要由经验模态分解（EmpiricalModeDecomposition，EMD）和Hilbert谱分析（HilbertSpectrumAnalysis，HSA）两部分组成。EMD是HHT的核心部分，它能够将复杂的非平稳信号自适应地分解为一系列固有模态函数（IntrinsicModeFunctions，IMF）。每个IMF分量都代表了信号在不同时间尺度上的局部特征，并且满足一定的条件，如在整个数据段内，极值点的数量与过零点的数量相等或至多相差1；信号由局部极大值确定的包络线和由局部极小值确定的包络线均值为零。这种自适应的分解方式使得EMD能够根据信号自身的特点进行分解，不需要预先设定基函数，与传统的傅里叶变换和小波变换相比，具有更强的自适应性。在得到IMF分量后，对每个IMF进行Hilbert变换，从而得到信号的Hilbert谱。Hilbert谱能够全面地展示信号的时间-频率-能量分布，使得我们可以直观地观察到信号在不同时刻的频率组成和能量变化情况。通过Hilbert谱分析，可以清晰地获取信号的瞬时频率和瞬时幅值等重要信息，为信号的深入分析提供了有力的工具。自HHT提出以来，在众多领域得到了广泛的应用和深入的研究。在电力系统中，它被用于谐波分析、故障诊断、电能质量检测等方面；在生物医学领域，可用于心电信号分析、脑电信号处理等；在地震工程中，能够对地震波进行分析，帮助研究地震的特性和传播规律。随着研究的不断深入，HHT的理论和应用也在不断发展和完善，为解决各种实际问题提供了更有效的方法和手段。2.2经验模态分解（EMD）2.2.1EMD基本原理经验模态分解（EMD）的核心思想是将复杂的非平稳信号分解为若干个固有模态函数（IMF）之和，每个IMF分量都代表了信号中不同的特征尺度和振动模式，能够反映信号的局部特性。其分解过程基于以下假设：任何复杂信号都是由若干个不同频率的本征模态函数组成，这些本征模态函数可以是线性的，也可以是非线性的。IMF必须满足两个条件：一是在整个数据序列中，极值点（极大值点和极小值点）的数量与过零点的数量相等，或最多相差不能多于一个；二是在任一时间点上，信号由局部极大值定义的上包络和由局部极小值定义的下包络的平均值为零。第一个条件类似于传统平稳高斯过程关于窄带的定义，第二个条件则是将传统的全局限定变为局部限定，可有效去除由于波形不对称而造成的瞬时频率的波动，使得IMF的瞬时频率具有实际物理意义。EMD的分解过程是一个反复筛选的过程，具体步骤如下：对于给定的原始信号x(t)，首先找出信号的所有局部极大值点和极小值点。利用三次样条插值函数，分别对极大值点和极小值点进行插值，得到信号的上包络线e_{max}(t)和下包络线e_{min}(t)。计算上下包络线的平均值m_1(t)=\frac{e_{max}(t)+e_{min}(t)}{2}。将原始信号x(t)减去平均包络m_1(t)，得到一个新的信号h_1(t)=x(t)-m_1(t)。判断h_1(t)是否满足IMF的两个条件。如果满足，则h_1(t)就是第一个IMF分量，记为c_1(t)=h_1(t)；如果不满足，则将h_1(t)作为新的原始信号，重复步骤1至4，直到得到满足IMF条件的分量。将第一个IMF分量c_1(t)从原始信号x(t)中分离出来，得到残余信号r_1(t)=x(t)-c_1(t)。把残余信号r_1(t)当作新的原始信号，重复上述步骤，依次得到第二个IMF分量c_2(t)、第三个IMF分量c_3(t)……直到残余信号r_n(t)变为单调函数或只存在一个极点为止。此时，原始信号x(t)可以表示为x(t)=\sum_{i=1}^{n}c_i(t)+r_n(t)，其中c_i(t)为第i个IMF分量，r_n(t)为最终的残余分量。通过这种方式，EMD能够将复杂信号中不同尺度的波动分量逐级分解出来，每个IMF分量都包含了信号在特定时间尺度上的主要特征，为后续的信号分析提供了便利。2.2.2EMD算法实现EMD算法的具体实现步骤如下：初始化：输入原始信号x(t)，设置分解的终止条件，如残余信号的长度小于某个阈值或者残余信号为单调函数等。寻找极值点：采用合适的算法，如峰值检测算法，找出信号x(t)的所有局部极大值点和极小值点。这一步骤的关键在于准确地识别出信号的极值点，以确保后续包络线的拟合精度。包络线拟合：利用三次样条插值函数，对极大值点和极小值点分别进行插值，得到信号的上包络线e_{max}(t)和下包络线e_{min}(t)。三次样条插值函数能够在保证曲线光滑的前提下，较好地拟合极值点，从而得到准确的包络线。计算均值：计算上下包络线的平均值m(t)=\frac{e_{max}(t)+e_{min}(t)}{2}。计算差值：将原始信号x(t)减去平均包络m(t)，得到差值信号h(t)=x(t)-m(t)。IMF判断：判断差值信号h(t)是否满足IMF的两个条件。如果满足，则h(t)就是一个IMF分量；如果不满足，则将h(t)作为新的原始信号，返回步骤2，继续进行筛选。在判断过程中，需要严格按照IMF的条件进行检验，以确保得到的IMF分量的准确性。分离IMF：将得到的IMF分量c(t)从原始信号x(t)中分离出来，得到残余信号r(t)=x(t)-c(t)。重复分解：将残余信号r(t)当作新的原始信号，重复步骤2至7，直到满足分解的终止条件。在实际编程实现中，可以使用各种编程语言，如MATLAB、Python等。以MATLAB为例，实现EMD算法的核心代码如下：function[imf,res]=emd(x)%x为输入的原始信号%imf为分解得到的IMF分量%res为最终的残余分量N=length(x);imf=[];res=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endend%x为输入的原始信号%imf为分解得到的IMF分量%res为最终的残余分量N=length(x);imf=[];res=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endend%imf为分解得到的IMF分量%res为最终的残余分量N=length(x);imf=[];res=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endend%res为最终的残余分量N=length(x);imf=[];res=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendN=length(x);imf=[];res=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendimf=[];res=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendres=x;whiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendwhiletrueh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendh=res;whiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendwhiletrue%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endend%寻找极值点max_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));m=(upper_envelope+lower_envelope)/2;ifmean(abs(m))>1e-6flag=false;return;endflag=true;endfunctionflag=isTerminationCondition(res)%判断残余信号是否满足终止条件%这里简单地判断残余信号是否为单调函数iflength(res)<10||all(diff(res)>=0)||all(diff(res)<=0)flag=true;elseflag=false;endendmax_index=findpeaks(h);min_index=findpeaks(-h);%拟合包络线upper_envelope=spline(max_index,h(max_index),1:N);lower_envelope=spline(min_index,-h(min_index),1:N);%计算均值m=(upper_envelope+lower_envelope)/2;%计算差值h_new=h-m;%判断是否满足IMF条件ifisIMF(h_new)break;endh=h_new;end%分离IMFimf=[imf;h];res=res-h;%判断是否满足终止条件ifisTerminationCondition(res)break;endendendfunctionflag=isIMF(h)%判断信号h是否为IMF%条件1：极值点个数与过零点个数相等或最多相差1zero_crossings=find(diff(sign(h))~=0);num_zero_crossings=length(zero_crossings);max_index=findpeaks(h);min_index=findpeaks(-h);num_extrema=length(max_index)+length(min_index);ifabs(num_zero_crossings-num_extrema)>1flag=false;return;end%条件2：上下包络均值为零max_index=findpeaks(h);min_index=findpeaks(-h);upper_envelope=spline(max_index,h(max_index),1:length(h));lower_envelope=spline(min_index,-h(min_index),1:length(h));