scipy.signal.

stft#

scipy.signal.stft(x, fs=1.0, window='hann', nperseg=256, noverlap=None, nfft=None, detrend=False, return_onesided=True, boundary='zeros', padded=True, axis=-1, scaling='spectrum')[源码]#

计算短时傅里叶变换（旧版函数）。

短时傅里叶变换（STFT）可用于量化非平稳信号的频率和相位内容随时间的变化。

旧版

此函数被视为旧版，将不再接收更新。虽然我们目前没有计划将其删除，但建议新代码改用更现代的替代方案。ShortTimeFFT 是一个具有更多功能的更新的 STFT / ISTFT 实现。实现之间的比较可以在 SciPy 用户指南的短时傅里叶变换部分找到。

参数:

x类数组: 测量值的时间序列
fs浮点数，可选: x 时间序列的采样频率。默认为 1.0。
window字符串或元组或类数组，可选: 要使用的窗口。如果 window 是字符串或元组，则将其传递给 get_window 以生成窗口值，这些值默认为 DFT 偶数。有关窗口列表和所需参数，请参阅 get_window。如果 window 是类数组，则将直接用作窗口，其长度必须为 nperseg。默认为 Hann 窗口。
nperseg整数，可选: 每个段的长度。默认为 256。
noverlap整数，可选: 段之间重叠的点数。如果为 None，则 noverlap = nperseg // 2。默认为 None。指定时，必须满足 COLA 约束（见下面的注意事项）。
nfft整数，可选: 如果需要零填充 FFT，则为使用的 FFT 长度。如果为 None，FFT 长度为 nperseg。默认为 None。
detrend字符串或函数或 False，可选: 指定如何去除每个段的趋势。如果 detrend 是字符串，则将其作为 type 参数传递给 detrend 函数。如果它是一个函数，它将接收一个段并返回一个去趋势的段。如果 detrend 为 False，则不进行去趋势。默认为 False。
return_onesided布尔值，可选: 如果为 True，则返回实数数据的一侧频谱。如果为 False，则返回双侧频谱。默认为 True，但对于复数数据，始终返回双侧频谱。
boundary字符串或 None，可选: 指定是否在两端扩展输入信号，以及如何生成新值，以便将第一个加窗段中心对齐到第一个输入点。这有利于在使用窗口函数从零开始时重建第一个输入点。有效选项为 ['even', 'odd', 'constant', 'zeros', None]。默认为 'zeros'，用于零填充扩展。即，对于 nperseg=3，[1, 2, 3, 4] 被扩展为 [0, 1, 2, 3, 4, 0]。
padded布尔值，可选: 指定是否在输入信号末尾进行零填充，以使信号恰好适合整数个窗口段，从而使所有信号都包含在输出中。默认为 True。如果 boundary 不为 None 且 padded 为 True（默认情况），则填充在边界扩展之后发生。
axis整数，可选: 计算 STFT 的轴；默认为最后一个轴（即 axis=-1）。
scaling: {‘spectrum’, ‘psd’}: 默认的 'spectrum' 缩放允许将 Zxx 的每个频率线解释为幅度谱。'psd' 选项将每条线缩放为功率谱密度 - 它允许通过对 abs(Zxx)**2 进行数值积分来计算信号能量。

1.9.0 版新增。

返回:

fndarray: 采样频率数组。
tndarray: 段次数组。
Zxxndarray: x 的 STFT。默认情况下，Zxx 的最后一个轴对应于段次。

另请参阅

istft: 逆短时傅里叶变换
ShortTimeFFT: 提供更多功能的更新的 STFT/ISTFT 实现。
check_COLA: 检查是否满足常量重叠相加（COLA）约束
check_NOLA: 检查是否满足非零重叠相加（NOLA）约束
welch: 基于 Welch 方法的功率谱密度。
spectrogram: 基于 Welch 方法的频谱图。
csd: 基于 Welch 方法的交叉谱密度。
lombscargle: 用于非均匀采样数据的 Lomb-Scargle 周期图

注意事项

为了通过 istft 中的逆短时傅里叶变换实现 STFT 的逆变换，信号加窗必须遵循“非零重叠相加”（NOLA）约束，并且输入信号必须具有完整的加窗覆盖（即 (x.shape[axis] - nperseg) % (nperseg-noverlap) == 0）。可以使用 padded 参数来实现此目的。

给定一个时域信号 \(x[n]\)、一个窗口 \(w[n]\) 和一个跳跃步长 \(H\) = nperseg - noverlap，时间索引 \(t\) 处的加窗帧由下式给出：

\[x_{t}[n]=x[n]w[n-tH]\]

重叠相加（OLA）重建方程由下式给出：

\[x[n]=\frac{\sum_{t}x_{t}[n]w[n-tH]}{\sum_{t}w^{2}[n-tH]}\]

NOLA 约束确保 OLA 重建方程分母中出现的每个归一化项都非零。可以通过 check_NOLA 测试所选的 window、nperseg 和 noverlap 是否满足此约束。

0.19.0 版新增。

参考文献

[1]

Oppenheim, Alan V., Ronald W. Schafer, John R. Buck “Discrete-Time Signal Processing”, Prentice Hall, 1999。

[2]

Daniel W. Griffin, Jae S. Lim “Signal Estimation from Modified Short-Time Fourier Transform”, IEEE 1984, 10.1109/TASSP.1984.1164317

示例

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt
>>> rng = np.random.default_rng()

生成一个测试信号：一个 2 Vrms 正弦波，其频率在 3kHz 附近缓慢调制，并被以 10 kHz 采样且幅度呈指数衰减的白噪声污染。

>>> fs = 10e3
>>> N = 1e5
>>> amp = 2 * np.sqrt(2)
>>> noise_power = 0.01 * fs / 2
>>> time = np.arange(N) / float(fs)
>>> mod = 500*np.cos(2*np.pi*0.25*time)
>>> carrier = amp * np.sin(2*np.pi*3e3*time + mod)
>>> noise = rng.normal(scale=np.sqrt(noise_power),
...                    size=time.shape)
>>> noise *= np.exp(-time/5)
>>> x = carrier + noise

计算并绘制 STFT 的幅度。

>>> f, t, Zxx = signal.stft(x, fs, nperseg=1000)
>>> plt.pcolormesh(t, f, np.abs(Zxx), vmin=0, vmax=amp, shading='gouraud')
>>> plt.title('STFT Magnitude')
>>> plt.ylabel('Frequency [Hz]')
>>> plt.xlabel('Time [sec]')
>>> plt.show()

../../_images/scipy-signal-stft-1_00_00.png

比较信号 x 的能量与其 STFT 的能量

>>> E_x = sum(x**2) / fs  # Energy of x
>>> # Calculate a two-sided STFT with PSD scaling:
>>> f, t, Zxx = signal.stft(x, fs, nperseg=1000, return_onesided=False,
...                         scaling='psd')
>>> # Integrate numerically over abs(Zxx)**2:
>>> df, dt = f[1] - f[0], t[1] - t[0]
>>> E_Zxx = sum(np.sum(Zxx.real**2 + Zxx.imag**2, axis=0) * df) * dt
>>> # The energy is the same, but the numerical errors are quite large:
>>> np.isclose(E_x, E_Zxx, rtol=1e-2)
True