最新のリリースでは、このページがまだ翻訳されていません。 このページの最新版は英語でご覧になれます。

hht

ヒルベルト・ファン変換

構文

hs = hht(imf,fs)

[hs,f,t] = hht(imf,fs)

[hs,f,t,imfinsf,imfinse] = hht(___)

[___] = hht(___,Name,Value)

hht(___)

hht(___,freqlocation)

説明

hs = hht(imf,fs) は、固有モード関数 imf とオプションのサンプリング周波数 fs によって指定される信号のヒルベルト スペクトル hs を返します。hs は、スペクトル成分が経時変化する信号の混成から成る信号の解析に役立ちます。hht を使用して信号のヒルベルトスペクトル解析を実行し、局所化された特徴を特定します。

[hs,f,t] = hht(imf,fs) は、hs に加えて周波数ベクトル f と時間ベクトル t を返します。

[hs,f,t,imfinsf,imfinse] = hht(___) は、信号の診断のために、固有モード関数の瞬時周波数 imfinsf と瞬間エネルギー imfinse も返します。

[___] = hht(___,Name,Value) では、1 つ以上の Name,Value のペア引数で指定された追加オプションを使用して、ヒルベルト スペクトル パラメーターを推定します。

出力引数なしで hht(___) を使用すると、現在の Figure ウィンドウにヒルベルト スペクトルがプロットされます。この構文は、これより前の構文の任意の入力引数で使用できます。

hht(___,freqlocation) では、周波数軸の場所を指定するためにオプションの freqlocation 引数を使用してヒルベルト スペクトルをプロットします。既定の設定では、周波数は Y 軸上に表されます。

例

経験的モード分解の実行と信号のヒルベルト スペクトルの可視化

この例では、周波数にはっきりした変化が含まれる正弦波で構成される非定常連続信号について考えます。ジャックハンマーの振動、または花火の音は、非定常連続信号の例です。

サンプリング周波数 fs で非定常信号のデータを読み込み、混合正弦波信号を可視化します。

load('sinusoidalSignalExampleData.mat','X','fs');
t = (0:length(X)-1)/fs;
plot(t,X);
xlabel('Time(s)');

混合信号に異なる振幅と周波数値をもつ正弦波が含まれることを観察します。

ヒルベルト スペクトル プロットを作成するには、信号の IMF が必要です。経験的モード分解を実行して、固有モード関数と信号の残差を計算します。信号が滑らかではないため、'pchip' を Interpolation メソッドとして指定します。

[imf,residual,info] = emd(X,'Interpolation','pchip');

Current IMF  |  #Sift Iter  |  Relative Tol  |  Stop Criterion Hit  
      1      |        2     |     0.026352   |  SiftMaxRelativeTolerance
      2      |        2     |    0.0039573   |  SiftMaxRelativeTolerance
      3      |        1     |     0.024838   |  SiftMaxRelativeTolerance
      4      |        2     |      0.05929   |  SiftMaxRelativeTolerance
      5      |        2     |      0.11317   |  SiftMaxRelativeTolerance
      6      |        2     |      0.12599   |  SiftMaxRelativeTolerance
      7      |        2     |      0.13802   |  SiftMaxRelativeTolerance
      8      |        3     |      0.15937   |  SiftMaxRelativeTolerance
      9      |        2     |      0.15923   |  SiftMaxRelativeTolerance
The decomposition stopped because the number of extrema of the residual signal is less than 'MaxNumExtrema'.

コマンド ウィンドウ内に生成されたテーブルは、生成された各 IMF のシフト反復の数、相対許容誤差、およびシフト停止基準を示します。この情報は info にも含まれます。Display を 0 と指定することによってテーブルを非表示にできます。

経験的モード分解を使用して取得した imf 成分を使用してヒルベルト スペクトル プロットを作成します。

hht(imf,fs);

周波数対時間のプロットは、IMF の各点における瞬間エネルギーを示す垂直方向のカラー バーがあるスパース プロットです。このプロットは、元の混合信号から分解された各成分の瞬時周波数スペクトルを表します。プロットから、1 秒での周波数に明瞭な変化がある 3 つの IMF が観察されます。

信号のヒルベルト スペクトル パラメーターの計算

この例では、周波数にはっきりした変化が含まれる正弦波で構成される非定常連続信号について考えます。ジャックハンマーの振動、または花火の音は、非定常連続信号の例です。

サンプリング周波数 fs で非定常信号のデータを読み込み、混合正弦波信号を可視化します。

load('sinusoidalSignalExampleData.mat','X','fs');
t = (0:length(X)-1)/fs;
plot(t,X);
xlabel('Time(s)');

混合信号に異なる振幅と周波数値をもつ正弦波が含まれることを観察します。

ヒルベルト スペクトル パラメーターを計算するには、信号の IMF が必要です。経験的モード分解を実行して、固有モード関数と信号の残差を計算します。信号が滑らかではないため、'pchip' を Interpolation メソッドとして指定します。

[imf,residual,info] = emd(X,'Interpolation','pchip');

Current IMF  |  #Sift Iter  |  Relative Tol  |  Stop Criterion Hit  
      1      |        2     |     0.026352   |  SiftMaxRelativeTolerance
      2      |        2     |    0.0039573   |  SiftMaxRelativeTolerance
      3      |        1     |     0.024838   |  SiftMaxRelativeTolerance
      4      |        2     |      0.05929   |  SiftMaxRelativeTolerance
      5      |        2     |      0.11317   |  SiftMaxRelativeTolerance
      6      |        2     |      0.12599   |  SiftMaxRelativeTolerance
      7      |        2     |      0.13802   |  SiftMaxRelativeTolerance
      8      |        3     |      0.15937   |  SiftMaxRelativeTolerance
      9      |        2     |      0.15923   |  SiftMaxRelativeTolerance
The decomposition stopped because the number of extrema of the residual signal is less than 'MaxNumExtrema'.

コマンド ウィンドウ内に生成されたテーブルは、生成された各 IMF のシフト反復の数、相対許容誤差、およびシフト停止基準を示します。この情報は info にも含まれます。Display を 0 と指定することによってテーブルを非表示にできます。

ヒルベルト スペクトル パラメーターである、ヒルベルト スペクトル hs、周波数ベクトル fs、時間ベクトル t、瞬時周波数 imfinsf, および瞬間エネルギー imfinse を計算します。

[hs,f,t,imfinsf,imfinse] = hht(imf,fs);

計算したヒルベルト スペクトル パラメーターを時間-周波数解析と信号の診断に使用します。

入力引数

imf — 固有モード関数

配列 | timetable

固有モード関数。配列または timetable として指定します。imf は、包絡線がゼロに関して対称である極値またはゼロクロッシングの数以下の任意の信号です。emd は、複雑な信号をヒルベルト スペクトル解析を実行するために必要な固有モード関数の有限数に分解および単純化するために使用されます。

hht は、imf の各列を固有モード関数として処理します。imf の計算の詳細は、emd を参照してください。

fs — サンプリング周波数

2π (既定値) | 正のスカラー

サンプリング周波数。正のスカラーとして指定します。サンプリング周波数またはサンプリング レートは単位時間あたりのサンプル数です。

fs が指定されない場合、2π の正規化された周波数がヒルベルト スペクトルの計算に使用されます。imf が timetable として指定されない場合、これからサンプリング時間が推論されます。

freqlocation — プロット上の周波数軸の位置

'yaxis' (既定値) | 'xaxis'

プロット上の周波数軸の位置。'yaxis' または 'xaxis' として指定します。プロットの Y 軸または X 軸上に周波数データを表示するには、それぞれ 'yaxis' または 'xaxis' として freqlocation を指定します。

名前と値のペアの引数

オプションの引数 Name,Value のコンマ区切りペアを指定します。Name は引数名で、Value は対応する値です。Name は引用符で囲まなければなりません。Name1,Value1,...,NameN,ValueN のように、複数の名前と値のペアの引数を、任意の順番で指定できます。

例: 'FrequencyResolution',1

'FrequencyLimits' — ヒルベルト スペクトルを計算する周波数限界

[0,fs/2] (既定値) | 1x2 正の整数配列

ヒルベルト スペクトルを計算する周波数限界。'FrequencyLimits' と 1x2 の正の整数配列で構成されるコンマ区切りのペアとして指定します。FrequencyLimits は Hz で指定します。

'FrequencyResolution' — 周波数範囲を離散化する周波数分解能

(f_high-f_low)/100 (既定値) | 正のスカラー

周波数限界を離散化する周波数分解能。'FrequencyResolution' と正のスカラーとで構成されるコンマ区切りのペアとして指定します。

FrequencyResolution は Hz で指定します。FrequencyResolution が指定されない場合、(f_high-f_low)/100 の値は FrequencyLimits から推定されます。ここで、f_high は FrequencyLimits の上限、f_low は下限です。

'MinThreshold' — ヒルベルト スペクトルの最小しきい値

-inf (既定値) | スカラー

ヒルベルト スペクトルの最小しきい値。'MinThreshold' とスカラーで構成されるコンマ区切りのペアとして指定します。

MinThreshold は、$$10{\log }_{10}\left(hs\right)$$ の対応する要素が MinThreshold より小さい場合に hs の要素を 0 に設定します。

出力引数

hs — 信号のヒルベルト スペクトル

スパース配列

信号のヒルベルト スペクトル。スパース配列として返されます。hs を使用して時間-周波数を解析し、信号の局所的な特徴を特定します。

f — 信号の周波数ベクトル

配列 | timetable

信号の周波数ベクトル。配列または timetable として返されます。hht は、周波数ベクトル f と時間ベクトル t を使用してヒルベルト スペクトル プロットを作成します。

f は以下として返されます。

• imf を配列として指定した場合は配列。

• imf を一様にサンプリングされた timetable として指定した場合は timetable。

数学的には、f は次のように表されます。

f = f_low:FrequencyResolution:f_high

t — 信号の時間ベクトル

配列 | timetable

信号の時間ベクトル。配列または timetable として返されます。hht は、時間ベクトル t と周波数ベクトル f を使用してヒルベルト スペクトル プロットを作成します。

t は以下として返されます。

• imf を配列として指定した場合は配列。

• imf を一様にサンプリングされた timetable として指定した場合は timetable。

imfinsf — 各 imf の瞬時周波数

配列 | timetable

各 imf の瞬時周波数。配列または timetable として返されます。

insf は imf と同じ列数であり、以下として返されます。

• imf を配列として指定した場合は配列。

• imf を一様にサンプリングされた timetable として指定した場合は timetable。

imfinse — 各 imf の瞬間エネルギー

配列 | timetable

各 imf の瞬間エネルギー。配列または timetable として返されます。

inse は imf と同じ列数であり、以下として返されます。

• imf を配列として指定した場合は配列。

• imf を一様にサンプリングされた timetable として指定した場合は timetable。

アルゴリズム

ヒルベルト スペクトル パラメーターは次の方法で計算されます。

1. 複雑な信号 X_c(t) は、経験的モード分解を使用して N 個の固有モード関数に分解されます。

2. 瞬時周波数 ω_i(t) と対応するエネルギー a_i²(t) は次のように計算されます。

   $$X_c\left(t\right)=H\left\{{\displaystyle \sum _{i=1}^N{c}_i\left(t\right)}\right\}={\displaystyle \sum _{i=1}^N{a}_i\left(t\right)e^{{\displaystyle {\int }_0^t{\omega }_i\left(t\right)dt}}}.$$

3. hht は、N個の時系列 $$\left[t,{\omega }_i\left(t\right),a_i^2\left(t\right)\right]$$ を離散化された時間-周波数平面に重ね合わせることによってヒルベルト スペクトル プロットを作成します。ここで、$$imfinsf=\left[{\omega }_i\left(t\right)\right]$$ および $$imfinse=\left[a_i^2\left(t\right)\right]$$ です。

参照

[1] Huang, Norden E., Shen, Samuel S. P. "Hilbert-Huang transform and its applications." World Scientific. 2014. Vol. 16.

[2] Huang, Norden E., Wu, Zhaohuaau, Long, Steven R., Arnold, Kenneth C., Chen, Xianyaoau & Blank, Karin. "On instantaneous frequency." Advances in Adaptive Data Analysis. World Scientific. 1.02 (2009): 177-229.

参考

emd | hilbert

トピック

• 解析信号とヒルベルト変換
• ヒルベルト変換および瞬時周波数

R2018a で導入