Documentation

This is machine translation

Translated by Microsoft
Mouseover text to see original. Click the button below to return to the English verison of the page.

Note: This page has been translated by MathWorks. Please click here
To view all translated materals including this page, select Japan from the country navigator on the bottom of this page.

Hilbert Transform and Instantaneous Frequency

The Hilbert transform estimates the instantaneous frequency of a signal for monocomponent signals only. A monocomponent signal is described in the time-frequency plane by a single "ridge." The set of monocomponent signals includes single sinusoids, but also includes signals like chirps.

Generate a chirp sampled at 1 kHz for two seconds. Specify the chirp so its frequency is initially 100 Hz and increases to 200 Hz after one second.

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

Estimate the spectrum of the chirp using the short-time Fourier transform implemented in the spectrogram function. Divide the signal into sections of length 100, windowed with a Hamming window. Specify 80 samples of overlap between adjoining sections and evaluate the spectrum at frequencies.

spectrogram(y,100,80,100,Fs,'yaxis')
view(-27,62)
shading interp

The signal is well described by a single peak frequency at each point in time.

Compute the analytic signal and differentiate its phase to measure the instantaneous frequency. The scaled derivative yields a meaningful estimate.

z = hilbert(y);
instfreq = Fs/(2*pi)*diff(unwrap(angle(z)));

plot(t(2:end),instfreq)
xlabel('Time')
ylabel('Hz')
grid on
title('Instantaneous Frequency')

Repeat the procedure, this time for a signal that is not monocomponent.

Generate a sum of two sinusoids of frequencies 60 Hz and 90 Hz, sampled at 1023 Hz for two seconds. Compute and plot the spectrogram. Use Hamming-windowed sections that have lengths of 200 and sample overlaps of 180. Evaluate the spectrum at frequencies.

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

spectrogram(x,200,180,256,Fs,'yaxis')
ax = gca;
ax.YTick = [60 90];
shading interp

Each time point shows the presence of the two components.

Compute the analytic signal and differentiate its phase. Zoom in to the region enclosing the frequencies of the sinusoids.

z = hilbert(x);
instfreq = Fs/(2*pi)*diff(unwrap(angle(z)));
plot(t(2:end),instfreq)

ylim([60 90])
xlabel('Time')
ylabel('Hz')

The analytic signal predicts an instantaneous frequency that is the average of the sinusoid frequencies.

If you want to estimate both frequencies as functions of time, you can use spectrogram. The ridges are peaks of the power spectral density.

[s,f,tt,p] = spectrogram(x,200,180,256,Fs);

numcomp = 2;
nt = length(tt);

P = zeros(nt,numcomp);
K = zeros(nt,numcomp);

for k = 1:nt
    [pk,lc] = findpeaks(p(:,k),f,'SortStr','descend','NPeaks',numcomp);
    P(k,:) = pk;
    K(k,:) = lc;
end

Plot the results.

spectrogram(x,200,180,256,Fs,'yaxis')
shading interp
view(2)

hold on
for k = 1:numcomp
    plot3(tt,K(:,k),10*log10(P(:,k)),'linewidth',4)
end
hold off
ax = gca;
ax.YTick = [60 90];

See Also

|

Related Topics

Was this topic helpful?