Main Content

hilbert

ヒルベルト変換を使用した離散時間解析信号

ページ内をすべて折りたたむ

構文

x = hilbert(xr)
x = hilbert(xr,n)

説明

例

x = hilbert(xr) は、実数データ シーケンス、xr からの解析信号 x を返します。xr が行列の場合、hilbert では、各列に対応する解析信号を検出します。

例

x = hilbert(xr,n) では、n 点の高速フーリエ変換 (FFT) を使用してヒルベルト変換が計算されます。入力データは、長さ n になるように、ゼロが付加されるか切り捨てられます。

例

すべて折りたたむ

ライブ スクリプトを開く

シーケンスを定義し、hilbert を使用してその解析信号を計算します。

xr = [1 2 3 4];
x = hilbert(xr)
x = 1×4 complex

   1.0000 + 1.0000i   2.0000 - 1.0000i   3.0000 - 1.0000i   4.0000 + 1.0000i

x の虚数部は、xr のヒルベルト変換であり、実数部は xr 自体です。

imx = imag(x)
imx = 1×4

     1    -1    -1     1


rex = real(x)
rex = 1×4

     1     2     3     4

x の離散フーリエ変換 (DFT) の後半はゼロです(この例では、変換の最後の半分がちょうどその最後の要素です)。fft(x) の DC 要素とナイキスト要素は純粋な実数です。

dft = fft(x)
dft = 1×4 complex

  10.0000 + 0.0000i  -4.0000 + 4.0000i  -2.0000 + 0.0000i   0.0000 + 0.0000i
ライブ スクリプトを開く

関数 hilbert は、有限データ ブロックに正確に一致する解析信号を検出します。また、有限インパルス応答 (FIR) ヒルベルト変換フィルターを使用して解析信号を生成し、虚数部に対する近似を計算することもできます。

周波数が 203 Hz、721 Hz および 1001 Hz である 3 つの正弦波で構成されるシーケンスを生成します。シーケンスは 10 kHz で約 1 秒間サンプリングされています。関数 hilbert を使用して解析信号を計算します。これを 0.01 秒と 0.03 秒の間でプロットします。

fs = 1e4;
t = 0:1/fs:1; 

x = 2.5 + cos(2*pi*203*t) + sin(2*pi*721*t) + cos(2*pi*1001*t);

y = hilbert(x);

plot(t,real(y),t,imag(y))
xlim([0.01 0.03])
legend('real','imaginary')
title('hilbert Function')
xlabel('Time (s)')

元のシーケンスおよび解析信号のパワー スペクトル密度のウェルチ推定を計算します。ハミング ウィンドウを適用した長さ 256 のオーバーラップのないセクションにシーケンスを分割します。負の周波数では解析信号にパワーがないことを確認します。

pwelch([x;y].',256,0,[],fs,'centered')
legend('Original','hilbert')

関数 designfilt を使用して 60 次のヒルベルト変換器 FIR フィルターを設計します。400 Hz の遷移幅を指定します。フィルターの周波数応答を可視化します。

fo = 60;

d = designfilt('hilbertfir','FilterOrder',fo, ...
       'TransitionWidth',400,'SampleRate',fs); 

freqz(d,1024,fs)

正弦波シーケンスをフィルター処理し、解析信号の虚数部を近似します。

hb = filter(d,x);

フィルターの群遅延 grd はフィルターの次数の 1/2 と等価です。この遅延を補正します。虚数部の最初の grd サンプルと実数部の最後の grd サンプルおよび時間ベクトルを削除します。0.01 秒と 0.03 秒の間の結果をプロットします。

grd = fo/2;
   
y2 = x(1:end-grd) + 1j*hb(grd+1:end);
t2 = t(1:end-grd);

plot(t2,real(y2),t2,imag(y2))
xlim([0.01 0.03])
legend('real','imaginary')
title('FIR Filter')
xlabel('Time (s)')

近似の解析信号のパワー スペクトル密度 (PSD) を推定し、hilbert の結果と比較します。

pwelch([y;[y2 zeros(1,grd)]].',256,0,[],fs,'centered')
legend('hilbert','FIR Filter')

入力引数

すべて折りたたむ

入力信号。実数値のベクトルまたは行列で指定します。xr が複素数の場合、hilbert では虚数部が無視されます。

例: sin(2*pi*(0:15)/16) は、正弦波の 1 周期を指定します。

例: sin(2*pi*(0:15)'./[16 8]) は、2 チャネルの正弦波信号を指定します。

データ型: single | double

DFT の長さ。正の整数スカラーで指定します。

データ型: single | double

出力引数

すべて折りたたむ

解析信号。ベクトルまたは行列として返されます。

詳細

すべて折りたたむ

解析信号

hilbert は、実数データ シーケンスから、複素螺旋シーケンス ("解析信号" とも呼ばれる) を返します。

解析信号 x = xr + jxi は、実数部 xr と虚数部 xi から成り、実数部は元のデータで、虚数部はそのヒルベルト変換です。虚数部は、元の実数シーケンスを 90°位相シフトしたものです。したがって、正弦値は余弦値に変換され、逆に余弦値は正弦値に変換されます。ヒルベルト変換された級数は、元のシーケンスと同じ振幅と周波数成分をもちます。元の位相に対応した位相情報を含みます。

ヒルベルト変換は、振幅や周波数などの時系列の瞬間的な特性を計算するのに役立ちます。瞬間的な振幅は、複素ヒルベルト変換の振幅であり、瞬間的な周波数は、瞬間的な位相角の時間変化率です。純粋な正弦波の場合、瞬間的な振幅と瞬間的な周波数は一定値です。ただし、瞬間的な位相は、局所的な位相角が単一のサイクルで線形に変化する様子を反映して、ノコギリ形です。混合正弦波の場合、属性は多くても 2～3 点の短期的すなわち局所的な平均です。例については、ヒルベルト変換および瞬時周波数を参照してください。

参考文献[1]には、時系列のスペクトル密度の対数のヒルベルト変換を実行して、最小位相を復元する Kolmogorov の方法が記述されています。このツールボックスの関数 rceps はこの復元を行います。

アルゴリズム

シーケンス xr の解析信号は "片側フーリエ変換" です。つまり、変換は負の周波数に対してはゼロになります。解析信号を近似するため、hilbert では、入力シーケンスの FFT が計算され、負の周波数に対応する FFT 係数がゼロに置き換えられて、その結果の逆 FFT が計算されます。

hilbert は、以下の 4 つの手順のアルゴリズムを使用します。

  1. 入力シーケンスの FFT が計算され、その結果がベクトル x に格納されます。

  2. 要素 h(i) が以下の値をもつ、ベクトル h が作成されます。

    • 1 (i = 1, (n/2)+1)

    • 2 (i = 2, 3, … , (n/2))

    • 0 (i = (n/2)+2, … , n)

  3. x と h の要素単位の積が計算されます。

  4. 手順 3 で得られたシーケンスの逆 FFT が計算され、結果から最初の n 個の要素が返されます。

このアルゴリズムは、まず[2]で紹介されました。この手法は、入力信号 x が有限データ ブロックであることを前提としています。この前提により、x のスペクトル冗長性をこの関数で厳密に排除することが可能になります。FIR フィルター処理に基づいた手法は解析信号の近似しかできませんが、データに連続して作用するという利点があります。FIR フィルターを使用して計算されるヒルベルト変換の他の例については、単側波帯振幅変調を参照してください。

参照

[1] Claerbout, Jon F. Fundamentals of Geophysical Data Processing with Applications to Petroleum Prospecting. Oxford, UK: Blackwell, 1985.

[2] Marple, S. L. “Computing the Discrete-Time Analytic Signal via FFT.” IEEE® Transactions on Signal Processing. Vol. 47, 1999, pp. 2600–2603.

[3] Oppenheim, Alan V., Ronald W. Schafer, and John R. Buck. Discrete-Time Signal Processing. 2nd Ed. Upper Saddle River, NJ: Prentice Hall, 1999.

拡張機能

C/C++ コード生成
MATLAB® Coder™ を使用して C および C++ コードを生成します。

GPU コード生成
GPU Coder™ を使用して NVIDIA® GPU のための CUDA® コードを生成します。

バージョン履歴

R2006a より前に導入

すべて展開する

参考

| |