Обнаружьте близко расположенные синусоиды

Рассмотрите синусоиду, f(x)=ej2πνx, оконный с окном Gaussian, g(t)=e-πt2. Кратковременное преобразование

Vgf(t,η)=ej2πνt-e-π(x-t)2e-j2π(x-t)(η-ν)dx=e-π(η-ν)2ej2πνt.

Когда просматривается как функция частоты, преобразование комбинирует константу (вовремя) колебание в ν с Гауссовым затуханием далеко от него. synchrosqueezing оценка мгновенной частоты,

Ωgf(t,η)=1j2πe-π(η-ν)2tej2πνte-π(η-ν)2ej2πνt=ν,

равняется значению, полученному при помощи стандартного определения, (2π)-1dаргументf(x)/dx. Для суперпозиции синусоид,

f(x)=k=1KAkej2πνkx,

кратковременное преобразование становится

Vgf(t,η)=k=1KAke-π(η-νk)2ej2πνkt.

Создайте 1 024 выборки сигнала, состоящего из двух синусоид. Одна синусоида имеет нормированную частоту ω0=π/5 рад/выборка. Другая синусоида имеет три раза частоту и три раза амплитуду.

N = 1024;
n = 0:N-1;

w0 = pi/5;
x = exp(1j*w0*n)+3*exp(1j*3*w0*n);

Вычислите кратковременное преобразование Фурье сигнала. Используйте окно Gaussian с 256 выборками с α=20, 255 выборок перекрытия между смежными разделами и точки ДПФ 10:24. Постройте абсолютное значение преобразования.

Nw = 256;
nfft = 1024;
alpha = 20;

[s,w,t] = spectrogram(x,gausswin(Nw,alpha),Nw-1,nfft,'centered');

surf(t,w/pi,abs(s),'EdgeColor','none')
view(2)
axis tight
xlabel('Samples')
ylabel('Normalized Frequency (\times\pi rad/sample)')

Synchrosqueezed преобразование Фурье приводит к более резкой, лучшей локализованной оценке спектра.

[ss,sw,st] = fsst(x,[],gausswin(Nw,alpha));

fsst(x,'yaxis')

Синусоиды отображаются как постоянные колебания в ожидаемых значениях частоты. Чтобы видеть, что затухание далеко от гребней является Гауссовым, постройте мгновенное значение преобразования и наложите два экземпляра Гауссова. Выразите Гауссово амплитудное и стандартное отклонение в терминах α и длина окна. Вспомните, что стандартное отклонение окна частотного диапазона является обратной величиной стандартного отклонения окна временного интервала.

rstdev = (Nw-1)/(2*alpha);
amp = rstdev*sqrt(2*pi);

instransf = abs(s(:,128));

plot(w/pi,instransf)
hold on
plot(w/pi,[1 3]*amp.*exp(-rstdev^2/2*(w-[1 3]*w0).^2),'--')
hold off
xlabel('Normalized Frequency (\times\pi rad/sample)')
lg = legend('Transform','First sinusoid','Second sinusoid');
lg.Location = 'best';

Synchrosqueezed преобразование Фурье концентрирует энергетическое содержимое сигнала на предполагаемых мгновенных частотах.

stem(sw/pi,abs(ss(:,128)))
xlabel('Normalized Frequency (\times\pi rad/sample)')
title('Synchrosqueezed Transform')

synchrosqueezed оценки мгновенной частоты допустимы, только если синусоиды разделяются больше, чем 2Δ, где

Δ=1σ2журнал2

для окна Gaussian и σ стандартное отклонение.

Повторите предыдущее вычисление, но теперь укажите, что вторая синусоида имеет нормированную частоту ω0+1.9Δ рад/выборка.

D = sqrt(2*log(2))/rstdev;

w1 = w0+1.9*D;

x = exp(1j*w0*n)+3*exp(1j*w1*n);

[s,w,t] = spectrogram(x,gausswin(Nw,alpha),Nw-1,nfft,'centered');
instransf = abs(s(:,20));

plot(w/pi,instransf)
hold on
plot(w/pi,[1 3]*amp.*exp(-rstdev^2/2*(w-[w0 w1]).^2),'--')
hold off
xlabel('Normalized Frequency (\times\pi rad/sample)')
lg = legend('Transform','First sinusoid','Second sinusoid');
lg.Location = 'best';

Synchrosqueezed преобразование Фурье не может разрешить синусоиды хорошо потому что |ω1-ω0|<2Δ. Спектральные оценки значительно уменьшаются в значении.

[ss,sw,st] = fsst(x,[],gausswin(Nw,alpha));

stem(sw/pi,abs(ss(:,128)))
xlabel('Normalized Frequency (\times\pi rad/sample)')
title('Synchrosqueezed Transform')

Смотрите также

| | |

Похожие темы

Для просмотра документации необходимо авторизоваться на сайте