Загрузите гиперболический щебет

Загрузите сигнал, который имеет два гиперболических щебета. Данные производятся на уровне 2 048 Гц. Первый щебет активен между 0,1 и 0,68 секундами, и второй щебет активен между 0,1 и 0,75 секундами. Мгновенная частота (в герц) первого щебета во время $$t$$ $$\frac{15\pi }{{(0.8-t)}^2}/2\pi$$ . Мгновенная частота второго щебета во время $$t$$ $$\frac{5\pi }{{(0.8-t)}^2}/2\pi$$. Постройте сигнал.

load hychirp
plot(t,hychirp)
grid on
title('Signal')
axis tight
xlabel('Time (s)')
ylabel('Amplitude')

Частотно-временной анализ: Преобразование Фурье

Преобразование Фурье (FT) очень хорошо в идентификации частотных составляющих, существующих в сигнале. Однако FT не идентифицирует, когда частотные составляющие происходят.

Постройте спектр величины сигнала. Увеличьте масштаб области между 0 и 200 Гц.

sigLen = numel(hychirp);
fchirp = fft(hychirp);
fr = Fs*(0:1/Fs:1-1/Fs);
plot(fr(1:sigLen/2),abs(fchirp(1:sigLen/2)),'x-')
xlabel('Frequency (Hz)')
ylabel('Amplitude')
axis tight
grid on
xlim([0 200])

Частотно-временной анализ: кратковременное преобразование Фурье

Преобразование Фурье не предоставляет информацию времени. Чтобы определить, когда изменения в частоте происходят, подход кратковременного преобразования Фурье (STFT) сегментирует сигнал на различные фрагменты и выполняет FT на каждом фрагменте. STFT, размещающий рядом в плоскости частоты времени, показывают здесь.

STFT предоставляет некоторую информацию и о синхронизации и о частотах, на которых событие сигнала имеет место. Однако выбор окна (сегмент) размер является ключевым. Для частотно-временного анализа с помощью STFT, выбирая более короткий размер окна помогает получить хорошее разрешение времени за счет разрешения частоты. С другой стороны выбор большего окна помогает получить хорошее разрешение частоты за счет разрешения времени.

Если вы выбираете размер окна, это остается фиксированным для целого анализа. Если можно оценить частотные составляющие, вы ожидаете в своем сигнале, то можно использовать ту информацию, чтобы выбрать размер окна для анализа.

Мгновенные частоты двух щебетов в их начальных моментах времени составляют приблизительно 5 Гц и 15 Гц. Используйте функцию помощника helperPlotSpectrogram построить спектрограмму сигнала с размером окна времени 200 миллисекунд. Исходный код для helperPlotSpectrogram перечислен в приложении. Графики функций помощника мгновенные частоты по спектрограмме как черные сегменты пунктирной линии. Мгновенные частоты разрешены рано в сигнале, но не также позже.

helperPlotSpectrogram(hychirp,t,Fs,200)

Теперь используйте helperPlotSpectrogram построить спектрограмму с размером окна времени 50 миллисекунд. Более высокие частоты, которые возникают позже в сигнале, теперь разрешены, но более низкие частоты в начале сигнала не.

helperPlotSpectrogram(hychirp,t,Fs,50)

Для неустановившихся сигналов как гиперболический щебет, с помощью STFT проблематично. Никакой один размер окна не может разрешить целое содержимое частоты таких сигналов.

Частотно-временной анализ: непрерывный вейвлет преобразовывает

Непрерывный вейвлет преобразовывает (CWT) был создан, чтобы преодолеть проблемы разрешения, свойственные от STFT. Плиточное размещение CWT на плоскости частоты времени показывают здесь.

Плиточное размещение CWT плоскости полезно, потому что много реальных сигналов имеют медленно колеблющееся содержимое, которое происходит по длинным шкалам, в то время как высокочастотные события имеют тенденцию быть резкими или переходными. Однако, если бы было естественно для высокочастотных событий быть длинным в длительности, то затем использование CWT не было бы соответствующим. У вас было бы более плохое разрешение частоты, не получая разрешения времени. Но это - довольно часто не случай. Человеческая слуховая система прокладывает себе путь; у нас есть намного лучшая локализация частоты на более низких частотах и лучшая локализация времени на высоких частотах.

Постройте scalogram CWT. scalogram является абсолютным значением CWT, построенного как функция времени и частоты. График использует логарифмическую ось частоты, потому что частоты в CWT являются логарифмическими. Присутствие двух гиперболических щебетов в сигнале ясно из scalogram. С CWT можно точно оценить мгновенные частоты в длительности сигнала, не волнуясь о выборе длины сегмента.

cwt(hychirp,Fs)

Белая пунктирная линия отмечает то, что известно как конус влияния. Конус влияния показывает области в scalogram, потенциально затронутом граничными эффектами. Для получения дополнительной информации смотрите Граничные эффекты и Конус Влияния.

Чтобы получить смысл того, как быстро величина коэффициентов вейвлета растет, используйте функцию помощника helperPlotScalogram3d построить scalogram как 3-D поверхность. Исходный код для helperPlotScalogram3d перечислен в приложении.

helperPlotScalogram3d(hychirp,Fs)

Используйте функцию помощника helperPlotScalogram построить scalogram сигнала и мгновенных частот. Исходный код для helperPlotScalogram перечислен в приложении. Мгновенные частоты выравниваются хорошо с функциями scalogram.

helperPlotScalogram(hychirp,Fs)

Приложение – функции помощника

helperPlotSpectrogram

function helperPlotSpectrogram(sig,t,Fs,timeres)
% This function is only intended to support this wavelet example.
% It may change or be removed in a future release.

[px,fx,tx] = pspectrum(sig,Fs,'spectrogram','TimeResolution',timeres/1000);
hp = pcolor(tx,fx,20*log10(abs(px)));
hp.EdgeAlpha = 0;
ylims = hp.Parent.YLim;
yticks = hp.Parent.YTick;
cl = colorbar;
cl.Label.String = 'Power (dB)';
axis tight
hold on
title(['Time Resolution: ',num2str(timeres),' ms'])
xlabel('Time (s)')
ylabel('Hz');
dt  = 1/Fs;
idxbegin = round(0.1/dt);
idxend1 = round(0.68/dt);
idxend2 = round(0.75/dt);
instfreq1 = abs((15*pi)./(0.8-t).^2)./(2*pi);
instfreq2 = abs((5*pi)./(0.8-t).^2)./(2*pi);
plot(t(idxbegin:idxend1),(instfreq1(idxbegin:idxend1)),'k--');
hold on;
plot(t(idxbegin:idxend2),(instfreq2(idxbegin:idxend2)),'k--');
ylim(ylims);
hp.Parent.YTick = yticks;
hp.Parent.YTickLabels = yticks;
hold off
end

helperPlotScalogram

function helperPlotScalogram(sig,Fs)
% This function is only intended to support this wavelet example.
% It may change or be removed in a future release.
[cfs,f] = cwt(sig,Fs);

sigLen = numel(sig);
t = (0:sigLen-1)/Fs;

hp = pcolor(t,log2(f),abs(cfs));
hp.EdgeAlpha = 0;
ylims = hp.Parent.YLim;
yticks = hp.Parent.YTick;
cl = colorbar;
cl.Label.String = 'magnitude';
axis tight
hold on
title('Scalogram and Instantaneous Frequencies')
xlabel('Seconds');
ylabel('Hz');
dt  = 1/2048;
idxbegin = round(0.1/dt);
idxend1 = round(0.68/dt);
idxend2 = round(0.75/dt);
instfreq1 = abs((15*pi)./(0.8-t).^2)./(2*pi);
instfreq2 = abs((5*pi)./(0.8-t).^2)./(2*pi);
plot(t(idxbegin:idxend1),log2(instfreq1(idxbegin:idxend1)),'k--');
hold on;
plot(t(idxbegin:idxend2),log2(instfreq2(idxbegin:idxend2)),'k--');
ylim(ylims);
hp.Parent.YTick = yticks;
hp.Parent.YTickLabels = 2.^yticks;
end

helperPlotScalogram3d

function helperPlotScalogram3d(sig,Fs)
% This function is only intended to support this wavelet example.
% It may change or be removed in a future release.
figure
[cfs,f] = cwt(sig,Fs);

sigLen = numel(sig);
t = (0:sigLen-1)/Fs;
surface(t,f,abs(cfs));
xlabel('Time (s)')
ylabel('Frequency (Hz)')
zlabel('Magnitude')
title('Scalogram In 3-D')
set(gca,'yscale','log')
shading interp
view([-40 30])
end