Algoritma Invers Short Time Fourier Transform dijelaskan dalam kata-kata

Saya mencoba memahami secara konseptual apa yang terjadi ketika forward dan invers Short Time Fourier Transforms (STFT) diterapkan pada sinyal time-domain diskrit. Saya telah menemukan makalah klasik karya Allen dan Rabiner ( 1977 ), serta artikel Wikipedia ( tautan ). Saya percaya bahwa ada juga artikel bagus lainnya yang dapat ditemukan di sini .

Saya tertarik untuk menghitung transformasi Gabor, yang tidak lebih dari STFT dengan jendela Gaussian.

Inilah yang saya mengerti tentang forward STFT:

1. Sub-urutan dipilih dari sinyal, terdiri dari elemen domain waktu.
2. Sub-urutan dikalikan dengan fungsi jendela menggunakan perkalian titik demi titik dalam domain waktu.
3. Sub-urutan yang dikalikan diambil ke dalam domain frekuensi menggunakan FFT.
4. Dengan memilih sub-urutan yang tumpang tindih berturut-turut, dan mengulangi prosedur di atas, kita mendapatkan matriks dengan m baris dan kolom n . Setiap kolom adalah sub-urutan yang dihitung pada waktu tertentu. Ini dapat digunakan untuk menghitung spektrogram.

Namun, untuk STFT terbalik , makalah berbicara tentang penjumlahan atas bagian analisis yang tumpang tindih. Saya merasa sangat menantang untuk memvisualisasikan apa yang sebenarnya terjadi di sini. Apa yang harus saya lakukan untuk dapat menghitung STFT terbalik (dalam urutan langkah seperti di atas)?

Maju STFT

Saya telah membuat gambar yang menunjukkan apa yang saya pikir sedang terjadi untuk STFT maju. Yang tidak saya mengerti adalah bagaimana cara merakit masing-masing sub-urutan sehingga saya mendapatkan kembali urutan waktu yang asli. Bisakah seseorang memodifikasi gambar ini atau memberikan persamaan yang menunjukkan bagaimana sub-urutan ditambahkan?

Transformasi Terbalik

Inilah yang saya mengerti tentang transformasi terbalik. Setiap jendela berturut-turut diambil kembali ke domain waktu menggunakan IFFT. Kemudian setiap jendela digeser oleh ukuran langkah, dan ditambahkan ke hasil pergeseran sebelumnya. Diagram berikut menunjukkan proses ini. Output terangkum adalah sinyal domain waktu.

Contoh Kode

Kode Matlab berikut menghasilkan sinyal domain waktu sintetis, dan kemudian menguji proses STFT, menunjukkan bahwa kebalikannya adalah ganda dari transformasi maju, dalam kesalahan pembulatan numerik. Awal dan akhir sinyal adalah nol-empuk untuk memastikan bahwa pusat jendela dapat ditempatkan pada elemen pertama dan terakhir dari sinyal domain waktu.

$N + N_0 - 1$$N_0$

% The code computes the STFT (Gabor transform) with step size = 1
% This is most useful when modifications of the signal is required in
% the frequency domain

% The Gabor transform is a STFT with a Gaussian window (w_t in the code)

% written by Nicholas Kinar

% Reference:
% [1] J. B. Allen and L. R. Rabiner, 
% “A unified approach to short-time Fourier analysis and synthesis,” 
% Proceedings of the IEEE, vol. 65, no. 11, pp. 1558 – 1564, Nov. 1977.

% generate the signal
mm = 8192;                  % signal points
t = linspace(0,1,mm);       % time axis

dt = t(2) - t(1);           % timestep t
wSize = 101;                % window size


% generate time-domain test function
% See pg. 156
% J. S. Walker, A Primer on Wavelets and Their Scientific Applications, 
% 2nd ed., Updated and fully rev. Boca Raton: Chapman & Hall/CRC, 2008.
% http://www.uwec.edu/walkerjs/primer/Ch5extract.pdf
term1 = exp(-400 .* (t - 0.2).^2);
term2 = sin(1024 .* pi .* t);
term3 = exp(-400.*(t- 0.5).^2);
term4 = cos(2048 .* pi .* t);
term5 = exp(-400 .* (t-0.7).^2);
term6 = sin(512.*pi.*t) - cos(3072.*pi.*t);
u = term1.*term2  + term3.*term4 + term5.*term6; % time domain signal
u = u';

figure;
plot(u)

Nmid = (wSize - 1) / 2 + 1;    % midway point in the window
hN = Nmid - 1;                 % number on each side of center point       


% stores the output of the Gabor transform in the frequency domain
% each column is the FFT output
Umat = zeros(wSize, mm);     


% generate the Gaussian window 
% [1] Y. Wang, Seismic inverse Q filtering. Blackwell Pub., 2008.
% pg. 123.
T = dt * hN;                    % half-width
sp = linspace(dt, T, hN); 
targ = [-sp(end:-1:1) 0 sp];    % this is t - tau
term1 = -((2 .* targ) ./ T).^2;
term2 = exp(term1);
term3 = 2 / (T * sqrt(pi));
w_t = term3 .* term2;
wt_sum = sum ( w_t ); % sum of the wavelet


% sliding window code
% NOTE that the beginning and end of the sequence
% are padded with zeros 
for Ntau = 1:mm

    % case #1: pad the beginning with zeros
    if( Ntau <= Nmid )
        diff = Nmid - Ntau;
        u_sub = [zeros(diff,1); u(1:hN+Ntau)];
    end

    % case #2: simply extract the window in the middle
    if (Ntau < mm-hN+1 && Ntau > Nmid)
        u_sub = u(Ntau-hN:Ntau+hN);
    end

    % case #3: less than the end
    if(Ntau >= mm-hN+1)
        diff = mm - Ntau;
        adiff = hN - diff;
        u_sub = [ u(Ntau-hN:Ntau+diff);  zeros(adiff,1)]; 
    end   

    % windowed trace segment
    % multiplication in time domain with
    % Gaussian window  function
    u_tau_omega = u_sub .* w_t';

    % segment in Fourier domain
    % NOTE that this must be padded to prevent
    % circular convolution if some sort of multiplication
    % occurs in the frequency domain
    U = fft( u_tau_omega );

    % make an assignment to each trace
    % in the output matrix
    Umat(:,Ntau) = U;

end

% By here, Umat contains the STFT (Gabor transform)

% Notice how the Fourier transform is symmetrical 
% (we only need the first N/2+1
% points, but I've plotted the full transform here
figure;
imagesc( (abs(Umat)).^2 )


% now let's try to get back the original signal from the transformed
% signal

% use IFFT on matrix along the cols
us = zeros(wSize,mm);
for i = 1:mm 
    us(:,i) = ifft(Umat(:,i));
end

figure;
imagesc( us );

% create a vector that is the same size as the original signal,
% but allows for the zero padding at the beginning and the end of the time
% domain sequence
Nuu = hN + mm + hN;
uu = zeros(1, Nuu);

% add each one of the windows to each other, progressively shifting the
% sequence forward 
cc = 1; 
for i = 1:mm
   uu(cc:cc+wSize-1) = us(:,i) + uu(cc:cc+wSize-1)';
   cc = cc + 1;
end

% trim the beginning and end of uu 
% NOTE that this could probably be done in a more efficient manner
% but it is easiest to do here

% Divide by the sum of the window 
% see Equation 4.4 of paper by Allen and Rabiner (1977)
% We don't need to divide by L, the FFT transform size since 
% Matlab has already taken care of it 
uu2 = uu(hN+1:end-hN) ./ (wt_sum); 

figure;
plot(uu2)

% Compare the differences bewteen the original and the reconstructed
% signals.  There will be some small difference due to round-off error
% since floating point numbers are not exact
dd = u - uu2';

figure;
plot(dd);

Pasangan transformasi STFT dapat ditandai dengan 4 parameter yang berbeda:

1. Ukuran FFT (N)
2. Ukuran langkah (m)
3. Jendela analisis (ukuran N)
4. Jendela sintesis (ukuran N)

Prosesnya adalah sebagai berikut:

1. Ambil N (ukuran fft) sampel dari lokasi input saat ini
2. Terapkan jendela analisis
3. Lakukan FFT
4. Lakukan apa pun yang ingin Anda lakukan di domain frekuensi
5. FFT terbalik
6. Terapkan jendela sintesis
7. Tambahkan ke output di lokasi output saat ini
8. Majukan input dan output lokasi dengan sampel M (ukuran langkah)

Algoritma add yang tumpang tindih adalah contoh yang baik untuk itu. Dalam hal ini ukuran langkah adalah N, ukuran FFT adalah 2 * N, jendela analisis adalah persegi panjang dengan N yang diikuti oleh N nol dan jendela sintesis hanyalah semua yang.

Ada banyak pilihan lain untuk itu dan dalam kondisi tertentu transfer maju / terbalik sepenuhnya merekonstruksi (yaitu Anda bisa mendapatkan sinyal asli kembali).

Poin kunci di sini adalah bahwa setiap sampel keluaran biasanya menerima kontribusi aditif dari lebih dari satu FFT terbalik. Output perlu diakumulasikan pada beberapa frame. Jumlah frame yang berkontribusi hanya diberikan oleh ukuran FFT dibagi dengan ukuran step (dibulatkan ke atas, jika perlu).

Tujuh tahun setelah pertanyaan ini pertama kali diajukan, saya mengalami kebingungan serupa dengan @Nicholas Kinar. Di sini saya ingin memberikan beberapa ide dan penjelasan perseptual "tidak resmi" dan "benar tidak sepenuhnya terjamin".

Judul pernyataan berikut ini dilebih-lebihkan untuk kejelasan yang lebih baik.

1. Proses maju STFT tidak benar-benar dimaksudkan untuk mempertahankan sinyal asli.
   • Saat menggunakan STFT dengan jendela non-sepele (Tidak semua-yang), sinyal input ke FFT adalah versi miring / membentang dari fragmen sinyal asli.
   • Ini bagus untuk ekstraksi fitur, di mana data yang tidak berguna / redundan disaring. Seperti dalam deteksi suku kata, tidak semua data temporal diperlukan untuk mendeteksi beberapa nada tertentu dalam sebuah pidato.
   • Puncak dalam vektor jendela mewakili minoritas posisi dalam sinyal audio di mana algoritma harus memperhatikan.
2. Jadi hasil mentah dari STFT terbalik bisa menjadi sesuatu yang mungkin tidak kita harapkan secara intuitif.
   • Seharusnya berupa fragmen sinyal berjendela yang menyerupai fitur STFT.
3. Untuk mendapatkan fragmen sinyal un-windowed yang asli, seseorang dapat menerapkan invers-window ke output mentah ifft.
   • Sangat mudah untuk merancang fungsi pemetaan yang dapat membatalkan efek jendela hann / hamming.
4. Jendela sintesis kemudian dilibatkan untuk menangani fragmentasi temporal yang tumpang tindih
   • Karena fragmen sinyal yang tidak berjendela yang asli dapat dilihat telah diperoleh, "bobot transisi" apa pun dapat digunakan untuk menginterpolasi bagian yang tumpang tindih.
5. Jika Anda ingin mempertimbangkan bahwa pidato pidato berjendela mungkin kurang menghargai sinyal lemah tetapi memuja sinyal kuat itu, maka mungkin ada cara untuk merancang jendela sintesis yang sesuai.
6. Juga, algoritma pembuatan jendela sintesis sintesis dapat diberikan dengan menerapkan prinsip-prinsip berikut:
   • Berat lebih tinggi posisi di jendela sintesis jika nilai jendela analisis untuk posisi ini tinggi, dibandingkan dengan fragmen lain yang tumpang tindih dengan posisi ini.
   • Bobot lebih rendah posisi di jendela sintesis jika nilai jendela analisis untuk posisi ini rendah, dan fragmen yang tumpang tindih lainnya lebih menghargai posisi ini dengan nilai jendela analisis yang lebih besar.