希尔伯特-黄变换说明及程序(标准程序)

希尔伯特黄变换信号处理
希尔伯特黄变换（Hilbert Huang Transform，简称HHT）是一个信
号处理的方法，常常用于分析非线性和非平稳信号。

它是由黄其炎教
授于1996年开发的，因此也叫做黄变换。

HHT的主要目的是将复杂的信号分解成数个瞬时频率相近的固有模态
函数（Intrinsic Mode Function，简称IMF）。

IMF是自然界中任何非线性现象的基本构建块，因此它们的分析在很多领域都非常重要。

HHT算法通常包括以下几个步骤：
1. 将待处理的信号（无论是时域信号还是频域信号）分解成数个组成
部分，即IMF。

2. 对每个IMF进行希尔伯特变换，得到复信号。

3. 计算每个复信号在复平面上的相位角和振幅。

4. 根据每个IMF在时域上的相位角和振幅，重建原信号的相位角和振幅。

5. 最后，将所有IMF的相位角和振幅相加得到原信号的相位角和振幅。

HHT的优点在于它不需要对信号做任何假设或模型。

它可以处理时域
和频域的信号，非常适合于分析非线性和非平稳信号，例如心电图、语音、天气数据和金融数据等。

HHT也有一些缺点，比如计算复杂度比较高，有时候需要选择合适的参数来得到比较准确的结果。

总的来说，希尔伯特黄变换是一个非常有用的信号处理方法，可以帮助我们了解自然界中复杂的现象。

它在科学、工程和医学等领域都得到了广泛应用。

（一）简单的EMD程序function imf = emd(x)% Empiricial Mode Decomposition (Hilbert-Huang Transform)% imf = emd(x)% Func : findpeaksx = transpose(x(:));%转置imf = [];while ~ismonotonic(x) %当x不是单调函数，分解终止条件x1 = x;sd = Inf;%均值%直到x1满足IMF条件，得c1while (sd > 0.1) | ~isimf(x1) %当标准偏差系数sd大于0.1或x1不是固有模态函数时，分量终止条件s1 = getspline(x1);%上包络线s2 = -getspline(-x1);%下包络线x2 = x1-(s1+s2)/2;%此处的x2为文章中的hsd = sum((x1-x2).^2)/sum(x1.^2);x1 = x2;endimf{end+1} = x1;x = x-x1;endimf{end+1} = x;% FUNCTIONSfunction u = ismonotonic(x)%u=0表示x不是单调函数，u=1表示x为单调的u1 = length(findpeaks(x))*length(findpeaks(-x));if u1 > 0, u = 0;else, u = 1; endfunction u = isimf(x)%u=0表示x不是固有模式函数，u=1表示x是固有模式函数N = length(x);u1 = sum(x(1:N-1).*x(2:N) < 0);u2 = length(findpeaks(x))+length(findpeaks(-x));if abs(u1-u2) > 1, u = 0;else, u = 1; endfunction s = getspline(x)%三次样条函数拟合成元数据包络线N = length(x);p = findpeaks(x);s = spline([0 p N+1],[0 x(p) 0],1:N);-------------------------------------------------------------------------------- function n = findpeaks(x)% Find peaks.找到极值% n = findpeaks(x)n = find(diff(diff(x) > 0) < 0);u = find(x(n+1) > x(n));n(u) = n(u)+1;------------------------------------------------------------------------------------------ ---------------------------------------------------------------------------------------- function plot_hht00(x,Ts)% 双边带调幅信号的EMD分解% Plot the HHT.% plot_hht(x,Ts)%% :: Syntax% The array x is the input signal and Ts is the sampling period.% Example on use: [x,Fs] = wavread('Hum.wav');% plot_hht(x(1:6000),1/Fs);% Func : emd% Get HHT.clear all;close all;Ts=0.0005;t=0:Ts:1; % 采样率2000HZ% 调幅信号x=sin(2*pi*t).*sin(40*pi*t);s1 = getspline(x);%上包络线s2 = -getspline(-x);%上包络线x1 = (s1+s2)/2;%此处的x2为文章中的hfigure;plot(t,x);xlabel('Time'), ylabel('Amplitude');title('双边带调幅信号');hold on;plot(t,s1,'-r');plot(t,s2,'-r');plot(t,x1,'g');imf = emd(x);for k = 1:length(imf)b(k) = sum(imf{k}.*imf{k});th = angle(hilbert(imf{k}));d{k} = diff(th)/Ts/(2*pi);end[u,v] = sort(-b);b = 1-b/max(b);% Set time-frequency plots.N = length(x);c = linspace(0,(N-2)*Ts,N-1);%figure;for k = v(1:2)plot(c,d{k},'k.','Color',b([k k k]),'MarkerSize',3); hold on;set(gca,'FontSize',8,'XLim',[0 c(end)],'YLim',[0 50]);%设置x、y轴句柄xlabel('Time'), ylabel('Frequency');title('原信号时频图');end% Set IMF plots.M = length(imf);N = length(x);c = linspace(0,(N-1)*Ts,N);for k1 = 0:4:M-1figurefor k2 = 1:min(4,M-k1), subplot(4,1,k2), plot(c,imf{k1+k2});set(gca,'FontSize',8,'XLim',[0 c(end)]);title('EMD分解结果');endxlabel('Time');end(二)function imf = emd(x)% Empiricial Mode Decomposition (Hilbert-Huang Transform)% imf = emd(x)% Func : findpeaksx = transpose(x(:));%转置imf = [];while ~ismonotonic(x) %当x不是单调函数，分解终止条件x1 = x;sd = Inf;%均值%直到x1满足IMF条件，得c1while (sd > 0.1) | ~isimf(x1) %当标准偏差系数sd大于0.1或x1不是固有模态函数时，分量终止条件s1 = getspline(x1);%上包络线s2 = -getspline(-x1);%下包络线x2 = x1-(s1+s2)/2;%此处的x2为文章中的hsd = sum((x1-x2).^2)/sum(x1.^2);x1 = x2;endimf{end+1} = x1;x = x-x1;endimf{end+1} = x;% FUNCTIONSfunction u = ismonotonic(x)%u=0表示x不是单调函数，u=1表示x为单调的u1 = length(findpeaks(x))*length(findpeaks(-x));if u1 > 0, u = 0;else, u = 1; endfunction u = isimf(x)%u=0表示x不是固有模式函数，u=1表示x是固有模式函数N = length(x);u1 = sum(x(1:N-1).*x(2:N) < 0);u2 = length(findpeaks(x))+length(findpeaks(-x));if abs(u1-u2) > 1, u = 0;else, u = 1; endfunction s = getspline(x)%三次样条函数拟合成元数据包络线N = length(x);p = findpeaks(x);s = spline([0 p N+1],[0 x(p) 0],1:N);--------------------------------------------------------------------------------------------------------------------------------------------------------------- function n = findpeaks(x)% Find peaks.找到极值% n = findpeaks(x)n = find(diff(diff(x) > 0) < 0);u = find(x(n+1) > x(n));n(u) = n(u)+1;------------------------------------------------------------------------------------------ ---------------------------------------------------------------------------------------- function plot_hht00(x,Ts)% 双边带调幅信号的EMD分解% Plot the HHT.% plot_hht(x,Ts)%% :: Syntax% The array x is the input signal and Ts is the sampling period.% Example on use: [x,Fs] = wavread('Hum.wav');% plot_hht(x(1:6000),1/Fs);% Func : emd% Get HHT.clear all;close all;Ts=0.0005;t=0:Ts:1; % 采样率2000HZ% 调幅信号x=sin(2*pi*t).*sin(40*pi*t);s1 = getspline(x);%上包络线s2 = -getspline(-x);%上包络线x1 = (s1+s2)/2;%此处的x2为文章中的hfigure;plot(t,x);xlabel('Time'), ylabel('Amplitude');title('双边带调幅信号');hold on;plot(t,s1,'-r');plot(t,s2,'-r');plot(t,x1,'g');imf = emd(x);for k = 1:length(imf)b(k) = sum(imf{k}.*imf{k});th = angle(hilbert(imf{k}));d{k} = diff(th)/Ts/(2*pi);end[u,v] = sort(-b);b = 1-b/max(b);% Set time-frequency plots.N = length(x);c = linspace(0,(N-2)*Ts,N-1);%figure;for k = v(1:2)plot(c,d{k},'k.','Color',b([k k k]),'MarkerSize',3); hold on;set(gca,'FontSize',8,'XLim',[0 c(end)],'YLim',[0 50]);%设置x、y轴句柄xlabel('Time'), ylabel('Frequency');title('原信号时频图');end% Set IMF plots.M = length(imf);N = length(x);c = linspace(0,(N-1)*Ts,N);for k1 = 0:4:M-1figurefor k2 = 1:min(4,M-k1), subplot(4,1,k2), plot(c,imf{k1+k2}); set(gca,'FontSize',8,'XLim',[0 c(end)]);title('EMD分解结果');endxlabel('Time');endplot_hht00(x,Ts)只画出了时频图，没有三维图。

目录∙ 1 本质模态函数(IMF)∙ 2 经验模态分解(EMD)∙ 3 结论∙ 4 相关条目∙ 5 参考文献∙ 6 外部链接[编辑]本质模态函数(IMF)任何一个资料，满足下列两个条件即可称作本质模态函数。

⒈局部极大值(local maxima)以及局部极小值(local minima)的数目之和必须与零交越点(zero crossing)的数目相等或是最多只能差1，也就是说一个极值后面必需马上接一个零交越点。

⒉在任何时间点，局部最大值所定义的上包络线(upper envelope)与局部极小值所定义的下包络线，取平均要接近为零。

因此，一个函数若属于IMF，代表其波形局部对称于零平均值。

此类函数类似于弦波（sinusoid-like），但是这些类似于弦波的部分其周期与振幅可以不是固定。

因为，可以直接使用希尔伯特转换，求得有意义的瞬时频率。

[编辑]经验模态分解(EMD)EMD算法流程图建立IMF是为了满足希尔伯特转换对于瞬时频率的限制条件之前置处理，也是一种转换的过程。

我们将IMF来做希尔伯特转换可以得到较良好的特性，不幸的是大部分的资料并不是IMF，而是由许多弦波所合成的一个组合。

如此一来，希尔伯特转换并不能得到正确的瞬时频率，我们便无法准确的分析资料。

为了解决非线性（non-linear）与非稳态（non-stationary）资料在分解成IMF时所遇到的困难，便发展出EMD。

经验模态分解是将讯号分解成IMF的组合。

经验模态分解是借着不断重复的筛选程序来逐步找出IMF。

以讯号为例，筛选程序的流程概述如下:步骤 1 : 找出中的所有局部极大值以及局部极小值，接着利用三次样条(cubic spline)，分别将局部极大值串连成上包络线与局部极小值串连成下包络线。

步骤 2 : 求出上下包络线之平均，得到均值包络线。

步骤 3 : 原始信号与均值包络线相减，得到第一个分量。

步骤 4 : 检查是否符合IMF的条件。

在MATLAB中进行希尔伯特-黄变换（HHT）并得到谱图的步骤如下：
收集数据：首先需要一个信号数据集。

这个数据集可以是地震波、生物信号等非线性和非平稳信号。

经验模态分解（EMD）：使用MATLAB中的EMD函数对数据进行处理，将信号分解为若干固有模态函数（IMF）和一个残差项。

希尔伯特变换：对每个IMF进行希尔伯特变换，得到其解析信号。

计算瞬时频率和瞬时幅度：计算每个解析信号的瞬时频率和瞬时幅度，得到信号在时-频域的表达。

汇总所有IMF的时-频域表达：对每个IMF的时-频域表达进行求和，得到原始信号的时-频域表达。

这就是希尔伯特-黄变换得到的信号谱。

Python希尔伯特黄变换（Python Hilbert-Huang Transform，简称HHT）是一种复杂非线性信号分析方法，结合了希尔伯特变换和黄变换的优势，能够有效地对非线性和非平稳信号进行时频谱分析。

本文将从HHT的原理、基本步骤和Python实现方法三个方面进行介绍。

一、HHT的原理1.希尔伯特变换希尔伯特变换是一种将实数信号转换为解析信号的数学方法，通过对原信号进行傅立叶变换得到频谱信息，再对频谱信息进行一定的处理得到解析频谱，从而实现信号的解析表示。

希尔伯特变换的核心是求出原信号的解析函数，即原信号的复数形式，其中实部是原信号本身，虚部是原信号的希尔伯特变换。

希尔伯特变换在信号处理领域有着广泛的应用，能够提取信号的瞬时特征，对非平稳信号进行时频分析具有很高的效果。

2.黄变换黄变换是一种局部线性和非线性信号分解方法，可以将非线性和非平稳信号分解成若干个固有模态函数（Intrinsic Mode Function，简称IMF）的线性组合。

黄变换首先对原信号进行极值点的提取，然后通过极值点之间的插值得到包络线，再将原信号减去包络线得到一维信号，并对得到的一维信号进行数据挑选和插值，最终得到IMF。

多次重复以上步骤，直到原信号能够被分解为若干个IMF，再通过IMF的线性组合得到原信号的近似表示。

3.HHT的结合HHT将希尔伯特变换和黄变换结合在一起，利用希尔伯特变换提取信号的瞬时特征，再通过黄变换将信号分解成若干个IMF，从而能够更准确地描述信号的时频特性。

HHT的优势在于能够适用于非线性和非平稳信号，对信号的局部特征具有很好的描述能力，因此在振动信号分析、生物医学信号处理等领域有着广泛的应用。

二、HHT的基本步骤1.信号分解HHT首先对原信号进行希尔伯特变换，得到信号的瞬时频率特征，然后通过黄变换将信号分解成若干个IMF。

2.IMF的提取针对得到的IMF，需要对每个IMF进行较为严格的判别，确定其是否符合IMF的特征：极值点交替出现、包络线对称、局部频率单调。

希尔伯特黄变换及其应用希尔伯特黄变换及其应用希尔伯特黄变换（Hilbert-Huang Transform，HHT）是一种用于分析非线性和非平稳信号的方法，它由黄其森（Norden E. Huang）和希尔伯特（Hilbert）共同提出。

该方法通过将信号分解为一组固有模态函数（Intrinsic Mode Functions，IMF）来提取信号中的模式和趋势。

本文将介绍希尔伯特黄变换的应用，并详细讲解其中的几个应用领域。

应用一：信号处理•希尔伯特黄变换可以用于音频信号处理，通过提取信号的固有模态函数，可以分离出音频信号中的主要频率成分，从而实现去噪、降噪等处理。

•在图像处理中，希尔伯特黄变换可以用于边缘检测和纹理分析。

通过提取图像的固有模态函数，可以分离出图像中的纹理信息和边缘信息，从而实现图像增强和分割等操作。

应用二：地震学•地震学中的信号分析是一项重要的任务，希尔伯特黄变换可以用于地震信号的分析和处理。

通过将地震信号分解为固有模态函数，可以提取出地震信号中的地震波的时频特征，从而实现地震信号的分类和识别。

•希尔伯特黄变换还可以用于地震信号的时频谱分析，通过将地震信号分解为固有模态函数，并对每个分量进行傅里叶变换，可以得到地震信号的时频谱图，从而更好地理解地震信号的时频特性。

应用三：医学工程•在医学工程中，希尔伯特黄变换可以用于生物信号的分析和处理，如心电图（ECG）和脑电图（EEG）等。

通过将生物信号分解为固有模态函数，可以提取出信号中的重要特征，如心跳频率、脑电波的频率等，从而实现疾病的诊断和监测。

•希尔伯特黄变换还可以用于生物信号的时频谱分析，通过将生物信号分解为固有模态函数，并对每个分量进行傅里叶变换，可以得到信号的时频谱图，从而更好地分析信号的时频特性。

应用四：金融市场•在金融市场中，希尔伯特黄变换可以用于股票价格的分析和预测。

通过将股票价格分解为固有模态函数，可以提取出股票价格的趋势和周期成分，从而更好地预测股票价格的走势。

希尔伯特黄变换和经验模态分解下载温馨提示：该文档是我店铺精心编制而成，希望大家下载以后，能够帮助大家解决实际的问题。

文档下载后可定制随意修改，请根据实际需要进行相应的调整和使用，谢谢!并且，本店铺为大家提供各种各样类型的实用资料，如教育随笔、日记赏析、句子摘抄、古诗大全、经典美文、话题作文、工作总结、词语解析、文案摘录、其他资料等等，如想了解不同资料格式和写法，敬请关注!Download tips: This document is carefully compiled by the editor. I hope that after you download them, they can help you solve practical problems. The document can be customized and modified after downloading, please adjust and use it according to actual needs, thank you!In addition, our shop provides you with various types of practical materials, such as educational essays, diary appreciation, sentence excerpts, ancient poems, classic articles, topic composition, work summary, word parsing, copy excerpts, other materials and so on, want to know different data formats and writing methods, please pay attention!希尔伯特黄变换和经验模态分解：理论与应用导言希尔伯特黄变换（Hilbert-Huang Transform，HHT）和经验模态分解（Empirical Mode Decomposition，EMD）是近年来在信号处理领域备受关注的两大方法。

HHT-希尔伯特·黄变换1998年，Norden E. Huang等人提出了经验模态分解方法，并引入了Hilbert谱的概念和Hilbert谱分析的方法，美国国家航空和宇航局（NASA）将这一方法命名为Hilbert-Huang Transform，简称HHT，即希尔伯特-黄变换。

HHT主要内容包含两部分，第一部分为经验模态分解（Empirical Mode Decomposition，简称EMD），它是由Huang提出的；第二部分为Hilbert谱分析（Hilbert Spectrum Analysis，简称HAS）。

简单说来，HHT处理非平稳信号的基本过程是：首先利用EMD方法将给定的信号分解为若干固有模态函数（以Intrinsic Mode Function或IMF表示，也称作本征模态函数），这些IMF是满足一定条件的分量；然后，对每一个IMF进行Hilbert变换，得到相应的Hilbert谱，即将每个IMF表示在联合的时频域中；最后，汇总所有IMF的Hilbert谱就会得到原始信号的Hilbert谱。

与传统的信号或数据处理方法相比，HHT具有如下特点：(1)HHT能分析非线性非平稳信号。

传统的数据处理方法，如傅立叶变换只能处理线性非平稳的信号，小波变换虽然在理论上能处理非线性非平稳信号，但在实际算法实现中却只能处理线性非平稳信号。

历史上还出现过不少信号处理方法，然而它们不是受线性束缚，就是受平稳性束缚，并不能完全意义上处理非线性非平稳信号。

HHT则不同于这些传统方法，它彻底摆脱了线性和平稳性束缚，其适用于分析非线性非平稳信号。

(2)HHT具有完全自适应性。

HHT能够自适应产生“基”，即由“筛选”过程产生的IMF。

这点不同于傅立叶变换和小波变换。

傅立叶变换的基是三角函数，小波变换的基是满足“可容性条件”的小波基，小波基也是预先选定的。

在实际工程中，如何选择小波基不是一件容易的事，选择不同的小波基可能产生不同的处理结果。