Создайте 1024 выборки сигнала, который состоит из суммы синусоид, встроенных в белый гауссов шум. Нормализованные частоты синусоид составляют $$\mathrm{}2\delta /5\mathrm{}$$рад/образец и $$\mathrm{}4\lambda /5\mathrm{}$$рад/образец. Синусоида более высокой частоты в 3 раза превышает амплитуду другой синусоиды.

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

w0 = 2*pi/5;
x = sin(w0*n)+3*sin(2*w0*n);

Вычислите синхронизированное преобразование Фурье сигнала. Постройте график результата.

[s,w,n] = fsst(x);

mesh(n,w/pi,abs(s))

axis tight
view(2)
colorbar

Вычислите обычное кратковременное преобразование Фурье сигнала для сравнения. Использовать значения по умолчанию spectrogram. Постройте график результата.

[s,w,n] = spectrogram(x);
 
surf(n,w/pi,abs(s),'EdgeColor','none')

axis tight
view(2)
colorbar

Создайте сигнал, состоящий из двух чирпов. Дискретизируют сигнал при частоте 3 кГц в течение одной секунды. Первая чирпа имеет начальную частоту 400 Гц и достигает 800 Гц в конце выборки. Вторая чирпа начинается с 500 Гц и достигает 1000 Гц в конце. Вторая чирп имеет вдвое большую амплитуду, чем первая чирп.

fs = 3000;
t = 0:1/fs:1-1/fs;

x1 = chirp(t,400,t(end),800);
x2 = 2*chirp(t,500,t(end),1000);

Вычислите и постройте график синхронизированного преобразования Фурье сигнала.

fsst(x1+x2,fs,'yaxis')

Сравните синхронизированное преобразование с кратковременным преобразованием Фурье (STFT). Вычислите STFT с помощью spectrogram функция. Укажите параметры по умолчанию, используемые fsst:

[stft,f,t] = spectrogram(x1+x2,kaiser(256,10),255,256,fs);

Постройте график абсолютного значения STFT.

mesh(t,f,abs(stft))

xlabel('Time (s)') 
ylabel('Frequency (Hz)')
title('Short-Time Fourier Transform')
axis tight
view(2)

Вычислите и отобразите синхронизированное преобразование Фурье квадратичного чирпа, которое начинается со 100 Гц и пересекает 200 Гц при t = 1. Укажите частоту дискретизации 1 кГц. Выразить время выборки как duration скаляр.

fs = 1000;
t = 0:1/fs:2;
ts = duration(0,0,1/fs);

x = chirp(t,100,1,200,'quadratic');

fsst(x,ts,'yaxis')

title('Quadratic Chirp')

Алгоритм синхроквизирования работает в предположении, что частота сигнала изменяется медленно. Таким образом, спектр лучше концентрируется в ранние времена, когда скорость изменения частоты меньше.

Вычислите и отобразите синхронизированное преобразование Фурье линейной частотной частоты, которая начинается при постоянном токе и пересекает 150 Гц при t = 1. Используйте 256-образное окно Хэмминга.

x = chirp(t,0,1,150);

fsst(x,ts,hamming(256),'yaxis')

title('Linear Chirp')

Вычислите и отобразите синхронизированное преобразование Фурье логарифмической чирпы. Чирп отбирается при 1 кГц, начинается при 20 Гц и пересекает 60 Гц при t = 1. Используйте 256-образное окно Кайзера с β = 20.

x = chirp(t,20,1,60,'logarithmic');

[s,f,t] = fsst(x,fs,kaiser(256,20));

clf
mesh(t,f,(abs(s)))

title('Logarithmic Chirp') 
xlabel('Time (s)')
ylabel('Frequency (Hz)')
view(2)

Используйте логарифмическую шкалу для частотной оси. Преобразование становится прямой линией.

ax = gca;
ax.YScale = 'log';
axis tight

Загрузите речевой сигнал, дискретизированный при $${}_{\mathrm{Fs}}\mathrm{=}\mathrm{}\mathrm{}\mathrm{}$$7418 Гц. Файл содержит запись женского голоса, говорящего слово «MATLAB ®».

load mtlb

% To hear, type sound(mtlb,Fs)

Вычислите синхронизированное преобразование сигнала. Используйте окно Ханна длиной 256. Отображение времени по оси X и частоты по оси Y.

fsst(mtlb,Fs,hann(256),'yaxis')

Использовать ifsst для инвертирования преобразования. Сравните исходные и восстановленные сигналы.

sst = fsst(mtlb,Fs,hann(256));

xrc = ifsst(sst,hann(256));

plot((0:length(mtlb)-1)/Fs,[mtlb xrc xrc-mtlb])
legend('Original','Reconstructed','Difference')

% To hear, type sound(xrc-mtlb,Fs)