Hilbert-Huang Transform: matlab 希爾伯特-黃變換: matlab實現


關於Hilbert-Huang的matlab實現,材料匯總,比較雜...感謝所有網絡上的貢獻者們:)

核心:以下代碼計算HHT邊際譜及其對應頻率
工具包要求:G-Rilling EMD Toolbox,TFTB Toolbox
:黃鍔先生課題組開發的工具包(可以在 這里 找到),這里並未用到。

% Empirical mode decomposition, resulting in intrinc mode functions.
% Without parameter 'MAXMODES' the process may be seriously delayed by
% decompose original signals into too many IMFs (not necessary, 9 is 
% enough generally)
imfs = emd(oriSig, 'MAXMODES', 9);
% HHT spectrum: hhtS
[A, f, t] = hhspectrum(imfs);
[hhtS, ~, fCent] = toimage(A, f, t);
% Marginal hilbert spectrum: hhtMS, xf: correspondig frequency
for k = 1 : size(hhtS, 1)
    hhtMS(k) = sum(hhtS(k,  : )) * 1 / fs;
end
xf = fCent(1, :) .* fs;

G-Rilling EMD Toolbox工具包相關配置:了凡春秋

簡單來說,設置好路徑之后輸入 install_emd 即可。


Matlab關於EMD分解后希爾伯特譜分析(相關函數均隸屬 G-Rilling EMD Toolbox)

  • hhspectrum 函數說明(8樓:老老的學生)
% [A,f,tt] = HHSPECTRUM(imf,t,l,aff)
% Input:
%- imf : matrix with one IMF per row    % 將emd分解得到的IMF代入就可以,就是你的程序中寫的c變量,不用加最后一行的趨勢項
%- t   : time instants    % 瞬時時間或持續時間 ??(寫[1:信號長度]就可以,真實的時間可以根據采樣率轉換
%- l   : estimation parameter for instfreq    % 瞬時頻率的估計參數 ??(寫1就可以,決定瞬時頻率估計時的邊界從第幾個點開始
%- aff : if 1, displays the computation evolution    % 顯示計算進程選項,不想顯示寫0就可以
%
% Output:
%   - A   : amplitudes of IMF rows
%   - f   : instantaneous frequencies
%   - tt  : truncated time instants    % 截止時間 ??(截斷時間,返回的是瞬時頻率對應的時間,要比原來信號的時間按短,由前面的l值決定)
%
% calls:
%   - hilbert  : computes the analytic signal
%   - instfreq : computes the instantaneous frequency    % 瞬時頻率
%
% Example:
[A,f,tt] = hhspectrum(c(1:end-1,:),[1:n],1,0);
[im,tt,ff] = toimage(A,f,tt,512);
disp_hhs(im,tt,[],fs);
  • 9樓:老老的學生

關於時頻圖的概念,我認為是與諸如小波,gabor等聯合時頻分析方法聯系在一起的。

小波,gabor等具有多尺度分析的概念,得到的時頻分布是一個二維的矩陣(橫軸時間,縱軸頻率,可以用不同的顏色(光譜圖spectrogram)或瀑布圖表示不同的幅度)。

對HHT來講,我覺得氣時頻圖的概念是有些不同的。EMD分解的作用就是把復雜的信號分界為簡單的單分量的信號,使其可以應用瞬時頻率的概念,hilbert變換的目的就是分析出瞬時頻率。所以HHT在每一時刻得到的只有一個值,而不是像小波之類的得到一系列的值(多尺度分析)。所以我們從其時頻分布圖上看到的是一條線,而不是一幅圖。

  • 27樓:songzy41
    HHT的重要意義是通過EMD分解后得到IMP,用IMP做hilbert,對每個固有模態函數才有瞬時信息。原信號是有許多固有模態函數之和,所有原信號沒有瞬時信息可言。這是HHT的基本出發點。

基於G-Rilling EMD Toolbox的Hilbert譜計算示例代碼:nike815

%示例程序
N=1000;
n=1:N;
fs=1000;
t=n/fs;
fx=50;
fy=150;
x=2*sin(2*pi*fx*t);
y=sin(2*pi*fy*t);
z=x+y;
data=z;
imf=emd(data);  % 對輸入信號進行EMD分解    
[A,f,t]=hhspectrum(imf);    % 對IMF分量求取瞬時頻率與振幅:A:是每個IMF的振幅向量,f:每個IMF對應的瞬時頻率,t:時間序列號
[E,t,Cenf]=toimage(A,f);    % 將每個IMF信號合成求取Hilbert譜,E:對應的振幅值,Cenf:每個網格對應的中心頻率。這里橫軸為時間,縱軸為頻率。即時頻圖(用顏色表示第三維值的大小)和三維圖(三維坐標系:時間,中心頻率,振幅)        
cemd_visu(data,1:length(data),imf); % 顯示每個IMF分量及殘余信號
disp_hhs(E);    % 希爾伯特譜
% 畫出邊際譜
colormap(flipud(gray)); % 黑白顯示
%N=length(Cenf);    % 設置頻率點數。完全從理論公式出發。網格化后中心頻率很重要,大家從連續數據變為離散的角度去思考,相信應該很容易理解
for k=1:size(E,1)
    bjp(k)=sum(E(k,:))*1/fs;
end
figure(3);
plot(Cenf(1,:)*fs,bjp);  % 作邊際譜圖   進行求取Hilbert譜時頻率已經被抽樣成具有一定窗長的離散頻率,所以此時的頻率軸已經是中心頻率
xlabel('頻率 / Hz');
ylabel('幅值');

% 繪制瞬時包絡圖和瞬時頻率圖
figure;
subplot(221),plot(t/N,A(1,:));xlabel('時間 t/s');ylabel('幅值');title('imf1分量瞬時包絡');
subplot(222),plot(t/N,f(1,:)*fs);xlabel('時間 t/s');ylabel('頻率');title('imf1分量瞬時頻率');
subplot(223),plot(t/N,A(2,:));xlabel('時間 t/s');ylabel('幅值');title('imf2分量瞬時包絡');
subplot(224),plot(t/N,f(2,:)*fs);xlabel('時間 t/s');ylabel('頻率');title('imf2分量瞬時頻率');

Hilbert變換-Hilbert譜、Hilbert邊際譜和包絡譜:於青民

Hilbert譜:信號的希爾伯特變換后做fft,表示信號幅值在整個頻率段上隨時間和頻率的變化規律;
Hilbert邊際譜:對hilbert譜做積分,表示信號幅值在整個頻率段上隨頻率的變化情況,它相當於傅里葉譜,但比傅里葉譜具有更高的頻率分辨率。Hilbert邊際譜是通過對Hilbert譜(在時間上)積分得到的;
Hilbert包絡譜:希爾伯特變換后做包絡后再fft,不同於Hilbert譜和Hilbert邊際譜,是直接對信號進行Hilbert變換后構造解析函數,然后依據解析函數求模值,求的模值即為包絡,然后對信號包絡進行FFT后得到的即為Hilbert包絡譜。


學習筆記:rogen

邊際譜從統計意義上表示了整組數據每個頻率點的積累幅值分布,而傅立葉譜的某點幅值表示在整個信號里有一個含有此頻率的三角函數組分,而且幅值越大只是說明在整個數據段上,局部存在的可能性越大。
再看得到的圖形,FFT 表示的是整個數據中,能量在一個頻率上分布的可能性地描述,而邊際譜表示在在每一個頻率上幅值的積累,如果想知道具體時間那么就看HHT譜,這個時間-幅值-頻率的三維譜。說到瞬時頻率,傅立葉變換不強調局部性,而是強調全局性。咱們的HHT才提出一個唯一的瞬時頻率的定義。因此拿瞬時頻率來衡量傅立葉變換也是不公平的。
對於邊際譜,就是hilbert譜對時間的積分,從積分的角度來講,就相當於對任意一階頻率把所有的時間上的幅值都加起來了,這就反映這階頻率在所有時間上的幅值積累,而對於他們對於頻率所對應的幅值,描述為該頻率在整個時間上出現的可能,我個人認為是既然出現了,那就是存在的,不能說是一種可能,只能說他出現的次數的多少,經過累積以后就變成了他的幅值,而對於fourier來講,他只能說明在用正余弦函數擬合這個信號的時候需要這一階頻率,幅值大,則說明他對擬合這個信號的貢獻大,而不是一種出現的可能性的度量。
能量譜,可以理解為對邊際譜的平方,這個只是具有能量的形式,而具體是否能代表能量,需要更進一步的探討,還有一個瞬時能量,這些名詞確實很誘人,但究竟如何,還需要大家的努力。


另一個比較詳細的HHT資料,附matlab代碼:cwjy

  • 邊際譜與傅里葉譜的比較
    意義不同:邊際譜從統計意義上表征了整組數據每個頻率點的累積幅值分布,而傅里葉頻譜的某一點頻率上的幅值表示在整個信號里有一個含有此頻率的三角函數組分。
    作用不同:邊際譜可以處理非平穩信號,如果信號中存在某一頻率的能量出現,就表示一定有該頻率的振動波出現,也就是說,邊際譜能比較准確地反映信號的實際頻率成分。而傅里葉變換只能處理平穩信號。
  • HHT & Hilbert變換
    Hilbert變換只是單純地求信號的瞬時振幅,頻率和相位,有可能出現沒有意義的負頻率;HHT變換先將信號進行EMD分解,得到的是各個不同尺度的分量,對每一個分量進行Hilbert變換后得到的是有實際意義的瞬時頻率。


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM