scipy.signal.

ShortTimeFFT

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

提供参数化的离散短时傅里叶变换（STFT）及其逆变换（ISTFT）。

通过以 hop 增量在输入信号上滑动窗口（win），STFT 计算顺序的 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 实例化此类。

默认情况下，使用所谓的规范双窗口。它是所有可能的双窗口中能量最小的窗口。from_win_equals_dual 和 closest_STFT_dual_window 提供了利用替代双窗口的方法。请注意，win 也始终是 dual_win 的双窗口。

由于约定将时间 t = 0 设为输入信号的第一个采样点，STFT 值通常具有负时间槽。因此，像 p_min 或 k_min 这样的负索引并不表示像标准 Python 索引那样从数组末尾倒数，而是表示在 t = 0 的左侧。

更多详细信息请参见 SciPy 用户指南 的短时傅里叶变换 部分。

请注意，初始化器的所有参数，除了 scale_to（它使用 scaling）之外，都具有同名的属性。

参数:

winnp.ndarray

窗口必须是实值或复值的 1d 数组。

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 列代表“幅度”或功率谱密度（“PSD”）频谱。此参数将属性 scaling 设置为相同的值。详见方法 scale_to。

phase_shiftint | None

如果设置，则向每个频率 f 添加线性相位 phase_shift / mfft * f。默认值 0 确保在第零个切片（其中 t=0 居中）上没有相位偏移。有关更多详细信息，请参见属性 phase_shift。

属性:

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

如果设置，则向频率为 f 的每个 FFT 切片添加线性相位 phase_shift / mfft * f。

scaling

应用于窗口函数的归一化（“magnitude”、“psd”、“unitary”或 None）。

win

作为实值或复值 1d 数组的窗口函数。

方法

extent(n[, axes_seq, center_bins])	返回最小和最大时频值。
from_dual(dual_win, hop, fs, *[, fft_mode, ...])	仅通过提供双窗口实例化 ShortTimeFFT。
from_win_equals_dual(desired_win, hop, fs, *)	创建窗口及其双窗口在缩放因子上相等的实例。
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 的“幅度”或“psd”缩放。
spectrogram(x[, y, detr, p0, p1, k_offset, ...])	计算频谱图或交叉频谱图。
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)	受后填充影响的第一个信号索引和第一个切片索引。

备注

典型的 STFT 应用是创建各种类型的时频图，通常归纳为“频谱图”这一术语。请注意，该术语也明确指 STFT 的绝对平方 [11]，如 spectrogram 中所示。

STFT 也可用于滤波和滤波器组，如 [12] 中所述。

参考文献

[11]

Karlheinz Gröchenig：《时频分析基础》，Birkhäuser Boston 2001，10.1007/978-1-4612-0003-1

[12]

Julius O. Smith III，《频谱音频信号处理》，在线书籍，2011，https://www.dsprelated.com/freebooks/sasp/

示例

在浏览器中尝试！

以下示例展示了一个具有变化的频率 \(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 会产生超出 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()

使用 istft 重建信号很简单，但请注意，应该指定 x1 的长度，因为 STFT 长度以 hop 步长增加。

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

可以计算信号部分的 STFT

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

组合顺序 STFT 部分时，需要考虑重叠。

>>> p0_ub = SFT.upper_border_begin(N2)[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

也可以计算信号部分的 istft

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

返回在新标签页中打开