Main Content

希尔伯特变换与瞬时频率

打开实时脚本

希尔伯特变换仅可估计单分量信号的瞬时频率。单分量信号在时频平面中用单一“脊”来描述。单分量信号包括单一正弦波信号和 chirp 等信号。

生成以 1 kHz 采样的时长为两秒的啁啾信号。指定啁啾信号的最初频率为 100 Hz，一秒后增加到 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, xlabel Time (s), ylabel Frequency (Hz) 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 函数只需一步即可计算并显示瞬时频率。

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

Figure contains an axes object. The axes object with title Instantaneous Frequency Estimate, xlabel Time (s), ylabel Frequency (Hz) contains an object of type line.

当信号不是单分量时，该方法会失败。

生成频率为 60 Hz 和 90 Hz 的两个正弦波的总和，以 1023 Hz 采样两秒。计算并绘制频谱图。在每个时间点都显示存在两个分量。

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, xlabel Time (s), ylabel Frequency (Hz) 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 with xlabel Time (s), ylabel Frequency (Hz) 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, xlabel Time (s), ylabel Frequency (Hz) contains an object of type line.

要采用时间的函数来估算这两个频率，请使用 spectrogram 求功率频谱密度，使用 tfridge 跟踪两个脊。在 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, xlabel Time (s), ylabel Frequency (Hz) contains 3 objects of type image, line.

另请参阅

|

相关主题