Main Content

ヒルベルト変換および瞬時周波数

ヒルベルト変換で推定されるのは、単一成分の信号の瞬時周波数のみです。単一成分の信号は、単一の "リッジ" によって時間-周波数平面に表されます。一連の単一成分信号には単一の正弦波とチャープのような信号が含まれます。

1 kHz で 2 秒間サンプリングしたチャープを生成します。初期周波数が 100 Hz で、1 秒後には 200 Hz に増加するようにチャープを指定します。

fs = 1000;
t = 0:1/fs:2-1/fs;
y = chirp(t,100,1,200);

関数 pspectrum に実装された短時間フーリエ変換を使用して、チャープのスペクトログラムを推定します。信号は、各時点の単一のピーク周波数によって明確に表されます。

pspectrum(y,fs,'spectrogram')

Figure contains an axes object. The axes object with title Fres = 10.267 Hz, Tres = 250 ms contains an object of type image.

瞬時周波数を測定するために解析信号を計算して位相を微分します。導関数をスケーリングすることで有効な推定が得られます。

z = hilbert(y);
instfrq = fs/(2*pi)*diff(unwrap(angle(z)));

clf
plot(t(2:end),instfrq)
ylim([0 fs/2])

Figure contains an axes object. The axes object contains an object of type line.

関数 instfreq は、1 ステップで瞬時周波数を計算し、表示します。

instfreq(y,fs,'Method','hilbert')

Figure contains an axes object. The axes object with title Instantaneous Frequency Estimate contains an object of type line.

信号が単一成分ではない場合は、メソッドが失敗します。

1023 Hz で 2 秒間サンプリングされた周波数 60 Hz と 90 Hz の 2 つの正弦波の和を生成します。スペクトログラムを計算してプロットします。各時点で 2 つの成分の存在が示されます。

fs = 1023;
t = 0:1/fs:2-1/fs;
x = sin(2*pi*60*t)+sin(2*pi*90*t);

pspectrum(x,fs,'spectrogram')
yticks([60 90])

Figure contains an axes object. The axes object with title Fres = 10.257 Hz, Tres = 250.2444 ms contains an object of type image.

解析信号を計算し、その位相を微分します。正弦波の周波数を囲む領域でズームインします。解析信号により、正弦波周波数の平均である瞬時周波数が予測されます。

z = hilbert(x);
instfrq = fs/(2*pi)*diff(unwrap(angle(z)));

plot(t(2:end),instfrq)
ylim([60 90])
xlabel('Time (s)')
ylabel('Frequency (Hz)')

Figure contains an axes object. The axes object contains an object of type line.

関数 instfreq は、平均の推定も行います。

instfreq(x,fs,'Method','hilbert')

Figure contains an axes object. The axes object with title Instantaneous Frequency Estimate contains an object of type line.

時間の関数として両方の周波数を推定するには、spectrogram を使用して、パワー スペクトル密度と tfridge を探し、2 つのリッジを追跡します。tfridge で、周波数の変更に対してペナルティを 0.1 に設定します。

[s,f,tt] = pspectrum(x,fs,'spectrogram');

numcomp = 2;
[fridge,~,lr] = tfridge(s,f,0.1,'NumRidges',numcomp);

pspectrum(x,fs,'spectrogram')
hold on
plot3(tt,fridge,abs(s(lr)),'LineWidth',4)
hold off
yticks([60 90])

Figure contains an axes object. The axes object with title Fres = 10.257 Hz, Tres = 250.2444 ms contains 3 objects of type image, line.

参考

|

関連するトピック