scipy.signal.

ShortTimeFFT#

class scipy.signal.ShortTimeFFT(win, hop, fs, *, fft_mode='onesided', mfft=None, dual_win=None, scale_to=None, phase_shift=0)[source]#

提供參數化的離散短時傅立葉轉換 (STFT) 及其逆轉換 (ISTFT)。

stft 通過以 hop 增量滑動視窗 (win) 於輸入訊號上,來計算連續的 FFT。它可以用於量化頻譜隨時間的變化。

stft 由複數值矩陣 S[q,p] 表示,其中第 p 列代表一個 FFT,其視窗中心位於時間 t[p] = p * delta_t = p * hop * T,其中 T 是輸入訊號的取樣間隔。第 q 行代表頻率 f[q] = q * delta_f 的值,其中 delta_f = 1 / (mfft * T) 是 FFT 的頻率間隔寬度。

逆 STFT istft 的計算方式是反轉 STFT 的步驟:取 S[q,p] 的第 p 個切片的 IFFT,並將結果與所謂的對偶視窗(參見 dual_win)相乘。將結果平移 p * delta_t,並將結果加到先前的平移結果中以重建訊號。如果只知道對偶視窗且 STFT 是可逆的,則可以使用 from_dual 來實例化此類別。

由於時間 t = 0 的慣例是在輸入訊號的第一個樣本,因此 STFT 值通常具有負時間槽。因此,負索引(如 p_mink_min)並不表示像標準 Python 索引那樣從陣列末端向後計數,而是表示位於 t = 0 的左側。

更多詳細資訊可以在 短時傅立葉轉換 章節的 SciPy 使用者指南 中找到。

請注意,除了 scale_to (使用 scaling) 之外,初始化器的所有參數都具有相同的命名屬性。

參數:
winnp.ndarray

視窗必須是實數或複數值的 1 維陣列。

hopint

視窗在每個步驟中平移的樣本增量。

fsfloat

輸入訊號和視窗的取樣頻率。它與取樣間隔 T 的關係是 T = 1 / fs

fft_mode‘twosided’, ‘centered’, ‘onesided’, ‘onesided2X’

要使用的 FFT 模式(預設為 ‘onesided’)。詳細資訊請參閱屬性 fft_mode

mfft: int | None

使用的 FFT 長度,如果需要零填充 FFT。如果為 None(預設值),則使用視窗 win 的長度。

dual_winnp.ndarray | None

win 的對偶視窗。如果設定為 None,則會在需要時計算。

scale_to‘magnitude’, ‘psd’ | None

如果不是 None(預設值),則會縮放視窗函數,使每個 STFT 列代表 ‘magnitude’ 或功率譜密度 (‘psd’) 頻譜。此參數將屬性 scaling 設定為相同的值。詳細資訊請參閱方法 scale_to

phase_shiftint | None

如果設定,則將線性相位 phase_shift / mfft * f 加到每個頻率 f。預設值 0 確保在第零個切片(其中 t=0 居中)上沒有相位偏移。詳細資訊請參閱屬性 phase_shift

範例

以下範例顯示了頻率變化的正弦波 \(f_i(t)\) 的 STFT 幅度(在圖中以紅色虛線標記)

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from scipy.signal import ShortTimeFFT
>>> from scipy.signal.windows import gaussian
...
>>> T_x, N = 1 / 20, 1000  # 20 Hz sampling rate for 50 s signal
>>> t_x = np.arange(N) * T_x  # time indexes for signal
>>> f_i = 1 * np.arctan((t_x - t_x[N // 2]) / 2) + 5  # varying frequency
>>> x = np.sin(2*np.pi*np.cumsum(f_i)*T_x) # the signal

使用的 Gaussian 視窗為 50 個樣本或 2.5 秒長。ShortTimeFFT 中的參數 mfft=200 會使頻譜被過採樣 4 倍

>>> g_std = 8  # standard deviation for Gaussian window in samples
>>> w = gaussian(50, std=g_std, sym=True)  # symmetric Gaussian window
>>> SFT = ShortTimeFFT(w, hop=10, fs=1/T_x, mfft=200, scale_to='magnitude')
>>> Sx = SFT.stft(x)  # perform the STFT

在圖中,訊號 x 的時間範圍以垂直虛線標記。請注意,SFT 產生Values超出 x 的時間範圍。左側和右側的陰影區域表示邊界效應,這是由於該區域中的視窗切片未完全在 x 的時間範圍內所導致的

>>> fig1, ax1 = plt.subplots(figsize=(6., 4.))  # enlarge plot a bit
>>> t_lo, t_hi = SFT.extent(N)[:2]  # time range of plot
>>> ax1.set_title(rf"STFT ({SFT.m_num*SFT.T:g}$\,s$ Gaussian window, " +
...               rf"$\sigma_t={g_std*SFT.T}\,$s)")
>>> ax1.set(xlabel=f"Time $t$ in seconds ({SFT.p_num(N)} slices, " +
...                rf"$\Delta t = {SFT.delta_t:g}\,$s)",
...         ylabel=f"Freq. $f$ in Hz ({SFT.f_pts} bins, " +
...                rf"$\Delta f = {SFT.delta_f:g}\,$Hz)",
...         xlim=(t_lo, t_hi))
...
>>> im1 = ax1.imshow(abs(Sx), origin='lower', aspect='auto',
...                  extent=SFT.extent(N), cmap='viridis')
>>> ax1.plot(t_x, f_i, 'r--', alpha=.5, label='$f_i(t)$')
>>> fig1.colorbar(im1, label="Magnitude $|S_x(t, f)|$")
...
>>> # Shade areas where window slices stick out to the side:
>>> for t0_, t1_ in [(t_lo, SFT.lower_border_end[0] * SFT.T),
...                  (SFT.upper_border_begin(N)[0] * SFT.T, t_hi)]:
...     ax1.axvspan(t0_, t1_, color='w', linewidth=0, alpha=.2)
>>> for t_ in [0, N * SFT.T]:  # mark signal borders with vertical line:
...     ax1.axvline(t_, color='y', linestyle='--', alpha=0.5)
>>> ax1.legend()
>>> fig1.tight_layout()
>>> plt.show()
../../_images/scipy-signal-ShortTimeFFT-1_00_00.png

使用 istft 重建訊號很簡單,但請注意,應指定 x1 的長度,因為 SFT 長度會以 hop 步長增加

>>> SFT.invertible  # check if invertible
True
>>> x1 = SFT.istft(Sx, k1=N)
>>> np.allclose(x, x1)
True

可以計算訊號部分的 SFT

>>> p_q = SFT.nearest_k_p(N // 2)
>>> Sx0 = SFT.stft(x[:p_q])
>>> Sx1 = SFT.stft(x[p_q:])

當組裝連續的 STFT 部分時,需要考慮重疊

>>> p0_ub = SFT.upper_border_begin(p_q)[1] - SFT.p_min
>>> p1_le = SFT.lower_border_end[1] - SFT.p_min
>>> Sx01 = np.hstack((Sx0[:, :p0_ub],
...                   Sx0[:, p0_ub:] + Sx1[:, :p1_le],
...                   Sx1[:, p1_le:]))
>>> np.allclose(Sx01, Sx)  # Compare with SFT of complete signal
True

也可以計算訊號部分的 itsft

>>> y_p = SFT.istft(Sx, N//3, N//2)
>>> np.allclose(y_p, x[N//3:N//2])
True
屬性:
T

輸入訊號和視窗的取樣間隔。

delta_f

STFT 的頻率間隔寬度。

delta_t

STFT 的時間增量。

dual_win

正規對偶視窗。

f

STFT 的頻率值。

f_pts

沿頻率軸的點數。

fac_magnitude

將 STFT 值乘以的因子,以將每個頻率切片縮放為幅度譜。

fac_psd

將 STFT 值乘以的因子,以將每個頻率切片縮放為功率譜密度 (PSD)。

fft_mode

所用 FFT 的模式(‘twosided’、‘centered’、‘onesided’ 或 ‘onesided2X’)。

fs

輸入訊號和視窗的取樣頻率。

hop

滑動視窗的訊號樣本中的時間增量。

invertible

檢查 STFT 是否可逆。

k_min

STFT 的最小可能訊號索引。

lower_border_end

不受預先填充影響的第一個訊號索引和第一個切片索引。

m_num

視窗 win 中的樣本數。

m_num_mid

視窗 win 的中心索引。

mfft

所用 FFT 的輸入長度 - 可能大於視窗長度 m_num

onesided_fft

如果使用單邊 FFT,則返回 True。

p_min

最小可能的切片索引。

phase_shift

如果設定,則將線性相位 phase_shift / mfft * f 加到每個頻率 f 的 FFT 切片。

scaling

應用於視窗函數的正規化(‘magnitude’、‘psd’ 或 None)。

win

作為實數或複數值 1 維陣列的視窗函數。

方法

extent(n[, axes_seq, center_bins])

返回最小和最大值時間-頻率值。

from_dual(dual_win, hop, fs, *[, fft_mode, ...])

僅通過提供對偶視窗來實例化 ShortTimeFFT

from_window(win_param, fs, nperseg, noverlap, *)

通過使用 get_window 來實例化 ShortTimeFFT

istft(S[, k0, k1, f_axis, t_axis])

逆短時傅立葉轉換。

k_max(n)

訊號結束後第一個未被時間切片觸及的樣本索引。

nearest_k_p(k[, left])

返回最接近的樣本索引 k_p,對於該索引 t[k_p] == t[p] 成立。

p_max(n)

對於 n 個樣本輸入,第一個非重疊上部時間切片的索引。

p_num(n)

具有 n 個樣本的輸入訊號的時間切片數。

p_range(n[, p0, p1])

確定並驗證切片索引範圍。

scale_to(scaling)

縮放視窗以獲得 STFT 的 ‘magnitude’ 或 ‘psd’ 縮放。

spectrogram(x[, y, detr, p0, p1, k_offset, ...])

計算 spectrogram 或 cross-spectrogram。

stft(x[, p0, p1, k_offset, padding, axis])

執行短時傅立葉轉換。

stft_detrend(x, detr[, p0, p1, k_offset, ...])

短時傅立葉轉換,並預先從每個區段中減去趨勢。

t(n[, p0, p1, k_offset])

具有 n 個樣本的輸入訊號的 STFT 時間。

upper_border_begin(n)

受後填充影響的第一個訊號索引和第一個切片索引。