ドキュメンテーション

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

hht

ヒルベルト・ファン変換

説明

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) が必要です。経験的モード分解を実行して、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' を内挿法として指定します。

[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 にも含まれます。Display0 と指定することによってテーブルを非表示にできます。

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

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

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

入力引数

すべて折りたたむ

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

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

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

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

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

名前と値のペアの引数

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

例: 'FrequencyResolution',1

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

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

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

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

MinThreshold は、10log10(hs) の対応する要素が MinThreshold より小さい場合に hs の要素を 0 に設定します。

出力引数

すべて折りたたむ

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

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

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

f = f_low:FrequencyResolution:f_high

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

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

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

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

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

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

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

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

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

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

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

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

アルゴリズム

ヒルベルト・ファン変換は、非定常および非線形データの時間-周波数解析の実行に役立ちます。ヒルベルト・ファン手順は以下のステップから構成されます。

  1. emd はデータセット x を固有モード関数の有限個に分解します。

  2. 各固有モード関数 xi の場合、関数 hht は次のとおりです。

    1. hilbert を使用して解析信号 zi(t)=xi(t)+jH{xi(t)} を計算します。ここで、H{xi}xi のヒルベルト変換です。

    2. zizi(t)=ai(t)ejθi(t) として表します。ここで、ai(t) は瞬間的な振幅で、θi(t) は瞬間的な位相です。

    3. 瞬間エネルギー |ai(t)|2、と瞬時周波数 ωi(t)dθi(t)/dt を計算します。サンプルレートを指定した場合、hht は、ωi(t) を周波数 (Hz) に変換します。

    4. imfinse の瞬間エネルギーと imfinsf の瞬時周波数を出力します。

  3. 出力引数なしで呼び出した場合、hht は、信号のエネルギーを、振幅に比例する色をもつ時間と周波数の関数としてプロットします。

参照

[1] Huang, Norden E., and Samuel S. P. Shen. Hilbert-Huang Transform and Its Applications. Singapore: World Scientific, 2014.

[2] Huang, Norden E., Zhaohua Wu, Steven R. Long, Kenneth C. Arnold, Xianyao Chen, and Karin Blank. "On Instantaneous Frequency." Advances in Adaptive Data Analysis. Vol. 1, No. 2, 2009, pp. 177–229.

R2018a で導入