%短時(shí)傅里葉變換STFT
%依據(jù)FFT手動(dòng)實(shí)現(xiàn)STFT
clear
clc
close?all
Fs?=?1000;????????????%?Sampling?frequency
T?=?1/Fs;?????????????%?Sampling?period
L?=?1000;?????????????%?Length?of?signal
t?=?（0:L-1）*T;????????%?Time?vector
S?=?20*cos（100*2*pi*t）+40*cos（50*2*pi*t）;%0.2-0.7*cos（2*pi*50*t+20/180*pi）?+?0.2*cos（2*pi*100*t+70/180*pi）?;

%%?所需參數(shù)
%主要包含:信號(hào)分割長(zhǎng)度（默認(rèn)分割8個(gè)窗口），海明窗口，重疊率，N點(diǎn)采樣
%默認(rèn)設(shè)置：
%?nsc=floor（L/4.5）;%海明窗的長(zhǎng)度
%?nov=floor（nsc/2）;%重疊率
%?nff=max（2562^nextpow2（nsc））;%N點(diǎn)采樣長(zhǎng)度
%也可手動(dòng)設(shè)置
nsc=100;%海明窗的長(zhǎng)度即每個(gè)窗口的長(zhǎng)度
nov=0;%重疊率
nff=max（2562^nextpow2（nsc））;%N點(diǎn)采樣長(zhǎng)度
%%?手動(dòng)實(shí)現(xiàn)STFT
h=hamming（nsc?‘periodic‘）;%計(jì)算海明窗的數(shù)值，給窗口內(nèi)的信號(hào)加權(quán)重
coln?=?1+fix（（L-nsc）/（nsc-nov））;%信號(hào)被分成了多少個(gè)片段
%如果nfft為偶數(shù)，則S的行數(shù)為（nfft/2+1），如果nfft為奇數(shù)，則行數(shù)為（nfft+1）/2
%因?yàn)閙atlab的FFT結(jié)果是對(duì)稱的，只需要一半
rown=nff/2+1;
STFT_X=zeros（rowncoln）;%初始化最終結(jié)果
%對(duì)每個(gè)片段進(jìn)行fft變換
index=1;%當(dāng)前片段第一個(gè)信號(hào)位置在原始信號(hào)中的索引
for?i=1:coln
????%提取當(dāng)前片段信號(hào)值并用海明窗進(jìn)行加權(quán)
????temp_S=S（index:index+nsc-1）.*h‘;
????%進(jìn)行N點(diǎn)FFT變換
????temp_X=fft（temp_Snff）;
????%取一半
????STFT_X（:i）=temp_X（1:rown）‘;
????%將索引后移
????index=index+（nsc-nov）;
end

%%?matlab自帶函數(shù)
spectrogram（Shamming（nsc?‘periodic‘）novnffFs）;
title（‘spectro