SciPy 1.12 中文文档（三十四）

原文：docs.scipy.org/doc/scipy-1.12.0/index.html

scipy.signal.coherence

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.coherence.html#scipy.signal.coherence

scipy.signal.coherence(x, y, fs=1.0, window='hann', nperseg=None, noverlap=None, nfft=None, detrend='constant', axis=-1)

使用 Welch 方法估计离散时间信号 X 和 Y 的幅度平方相干性估计，Cxy。

Cxy = abs(Pxy)**2/(Pxx*Pyy)，其中 Pxx 和 Pyy 是 X 和 Y 的功率谱密度估计，Pxy 是 X 和 Y 的交叉谱密度估计。

参数：

xarray_like

测量值的时间序列

yarray_like

测量值的时间序列

fsfloat，可选

x 和 y 时间序列的采样频率。默认为 1.0。

windowstr 或者 tuple 或者 array_like，可选

所需使用的窗口。如果 window 是字符串或元组，则传递给 get_window 以生成窗口值，默认情况下为 DFT-even。参见 get_window 获取窗口列表和必需的参数。如果 window 是 array_like，则直接用作窗口，其长度必须为 nperseg。默认为汉宁窗口。

npersegint，可选

每个段的长度。默认为 None，但如果窗口是字符串或元组，则设置为 256，如果窗口是 array_like，则设置为窗口的长度。

noverlap: int, 可选

在分段之间重叠的点数。如果 None，则 noverlap = nperseg // 2。默认为 None。

nfftint，可选

如果需要零填充 FFT，则使用的 FFT 长度。如果 None，则 FFT 长度为 nperseg。默认为 None。

detrendstr 或者 函数 或者 False，可选

指定如何去趋势化每个段。如果 detrend 是一个字符串，则作为 type 参数传递给 detrend 函数。如果它是一个函数，则它接受一个段并返回去趋势化的段。如果 detrend 是 False，则不执行去趋势化。默认为 'constant'。

axisint，可选

计算两个输入信号的相干性的轴；默认为最后一个轴（即 axis=-1）。

返回：

fndarray

样本频率的数组。

Cxyndarray

x 和 y 的幅度平方相干性。

另请参阅

periodogram

简单的，可选修改的周期图

lombscargle

不均匀采样数据的 Lomb-Scargle 周期图

welch

Welch 方法计算的功率谱密度。

csd

Welch 方法的交叉谱密度。

注意事项

适当的重叠量将取决于窗口的选择和您的要求。对于默认的 Hann 窗口，50%的重叠是在准确估计信号功率和不过多计算任何数据之间的合理折衷。更窄的窗口可能需要更大的重叠。

从版本 0.16.0 开始新增。

参考文献

[1]

P. Welch，“用于估计功率谱的快速傅立叶变换的使用：一种基于短期修改周期图平均的方法”，IEEE Trans. Audio Electroacoust. vol. 15, pp. 70-73, 1967 年

[2]

Stoica, Petre 和 Randolph Moses，“信号的频谱分析”，Prentice Hall，2005 年

示例

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

生成两个具有一些共同特征的测试信号。

>>> fs = 10e3
>>> N = 1e5
>>> amp = 20
>>> freq = 1234.0
>>> noise_power = 0.001 * fs / 2
>>> time = np.arange(N) / fs
>>> b, a = signal.butter(2, 0.25, 'low')
>>> x = rng.normal(scale=np.sqrt(noise_power), size=time.shape)
>>> y = signal.lfilter(b, a, x)
>>> x += amp*np.sin(2*np.pi*freq*time)
>>> y += rng.normal(scale=0.1*np.sqrt(noise_power), size=time.shape) 

计算并绘制相干性。

>>> f, Cxy = signal.coherence(x, y, fs, nperseg=1024)
>>> plt.semilogy(f, Cxy)
>>> plt.xlabel('frequency [Hz]')
>>> plt.ylabel('Coherence')
>>> plt.show() 

scipy.signal.spectrogram

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.spectrogram.html#scipy.signal.spectrogram

scipy.signal.spectrogram(x, fs=1.0, window=('tukey', 0.25), nperseg=None, noverlap=None, nfft=None, detrend='constant', return_onesided=True, scaling='density', axis=-1, mode='psd')

使用连续的傅里叶变换计算频谱图。

频谱图可用作可视化非平稳信号频率内容随时间变化的一种方法。

遗留

此函数被视为遗留版本，将不再接收更新。这可能意味着在未来的 SciPy 版本中将被移除。ShortTimeFFT是一个更新的 STFT / ISTFT 实现，具有更多功能，还包括一个spectrogram方法。在SciPy 用户指南的Short-Time Fourier Transform部分中可以找到这些实现之间的比较。

参数：

xarray_like

测量值的时间序列

fsfloat，可选

x时间序列的采样频率。默认为 1.0。

windowstr 或元组或 array_like，可选

期望使用的窗口。如果window是字符串或元组，则会传递给get_window以生成窗口数值，默认情况下为 DFT 偶数。请参阅get_window获取窗口列表和所需参数。如果window是 array_like，则将直接使用作为窗口，并且其长度必须为nperseg。默认为 Tukey 窗口，形状参数为 0.25。

npersegint，可选

每个段的长度。默认为 None，但如果window是字符串或元组，则设置为 256，如果window是 array_like，则设置为窗口的长度。

noverlapint，可选

每个段之间重叠的点数。如果为None，则noverlap = nperseg // 8。默认为None。

nfftint，可选

所使用的 FFT 长度，如果需要零填充 FFT。如果为None，则 FFT 长度为nperseg。默认为None。

detrendstr 或函数或False，可选

指定如何去趋势化每个段。如果detrend是一个字符串，则传递为detrend函数的type参数。如果是一个函数，则接受一个段并返回去趋势化的段。如果detrend为False，则不进行去趋势化。默认为‘constant’。

return_onesidedbool，可选

如果True，返回实数据的单侧频谱。如果False，返回双侧频谱。默认为True，但对于复杂数据，始终返回双侧频谱。

scaling{ ‘density’, ‘spectrum’ }，可选

选择计算功率谱密度（‘density’）或功率谱（‘spectrum’），其中Sxx的单位为 V**2/Hz，如果x以 V 为单位，fs以 Hz 为单位。默认为‘density’。

axisint，可选

计算谱图的轴；默认为最后一个轴（即axis=-1）。

modestr，可选

定义预期的返回值类型。选项有[‘psd’, ‘complex’, ‘magnitude’, ‘angle’, ‘phase’]。‘complex’等同于没有填充或边界扩展的stft的输出。‘magnitude’返回 STFT 的绝对幅度。‘angle’和‘phase’分别返回 STFT 的复角，带有和不带有展开。

返回：

fndarray

样本频率的数组。

tndarray

分段时间的数组。

Sxxndarray

x 的谱图。默认情况下，Sxx 的最后一个轴对应于段时间。

另请参阅

periodogram

简单的、可选修改后的周期图

lombscargle

Lomb-Scargle 不规则采样数据的周期图

welch

Welch 方法的功率谱密度。

csd

Welch 方法的交叉谱密度

ShortTimeFFT

提供更多功能的新 STFT/ISTFT 实现，其中还包括一个spectrogram方法。

注释

适当的重叠量取决于窗口的选择和您的需求。与 Welch 方法相反，在计算谱图时，人们可能希望使用较小的重叠（或者根本不重叠），以保持各个段的统计独立性。因此，默认窗口是 Tukey 窗口，每端重叠窗口长度的 1/8。

新版本 0.16.0 中引入。

参考文献

[1]

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

示例

>>> import numpy as np
>>> from scipy import signal
>>> from scipy.fft import fftshift
>>> 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 

计算并绘制谱图。

>>> f, t, Sxx = signal.spectrogram(x, fs)
>>> plt.pcolormesh(t, f, Sxx, shading='gouraud')
>>> plt.ylabel('Frequency [Hz]')
>>> plt.xlabel('Time [sec]')
>>> plt.show() 

注意，如果使用的输出不是单边的话，请使用以下内容：

>>> f, t, Sxx = signal.spectrogram(x, fs, return_onesided=False)
>>> plt.pcolormesh(t, fftshift(f), fftshift(Sxx, axes=0), shading='gouraud')
>>> plt.ylabel('Frequency [Hz]')
>>> plt.xlabel('Time [sec]')
>>> plt.show() 

scipy.signal.lombscargle

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.lombscargle.html#scipy.signal.lombscargle

scipy.signal.lombscargle(x, y, freqs)

计算 Lomb-Scargle 周期图。

Lomb-Scargle 周期图由 Lomb [1]开发，并由 Scargle [2]进一步扩展，用于发现和测试不均匀时间采样中弱周期信号的显著性。

当 normalize 设置为 False（默认值）时，计算得到的周期图未归一化，对于具有振幅 A 的谐波信号，对于足够大的 N，它取值为(A**2) * N/4。

当 normalize 设置为 True 时，计算得到的周期图将通过数据围绕常数参考模型（在零点）的残差进行归一化。

输入数组应为 1-D，并将转换为 float64 类型。

参数：

xarray_like

样本时间。

yarray_like

测量值。

freqsarray_like

输出周期图的角频率。

precenterbool, optional

通过减去均值预置测量值。

normalizebool, optional

计算归一化周期图。

返回：

pgramarray_like

Lomb-Scargle 周期图。

Raises：

ValueError

如果输入数组 x 和 y 的形状不同。

另见

istft

逆短时傅立叶变换

check_COLA

检查是否满足常数重叠加（COLA）约束

welch

Welch 方法的功率谱密度

spectrogram

Welch 方法的谱图

csd

Welch 方法的交叉谱密度

Notes

此子程序使用了由 Townsend 稍作修改的算法来计算周期图[3]，该算法允许在每个频率上仅通过输入数组的一次传递计算周期图。

算法运行时间大致按 O(x * freqs)或 O(N²)缩放，适用于大量样本和频率。

参考文献

[1]

N.R. Lomb，“不等间隔数据的最小二乘频率分析”，《天体物理学和空间科学》，第 39 卷，第 447-462 页，1976 年

[2]

J.D. Scargle，“天文时间序列分析研究 II - 不均匀间隔数据谱分析的统计方面”，《天体物理学期刊》，第 263 卷，第 835-853 页，1982 年

[3]

R.H.D. Townsend，“使用图形处理单元快速计算 Lomb-Scargle 周期图”，《天体物理学期刊增刊》，第 191 卷，第 247-253 页，2010 年

示例

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

首先为信号定义一些输入参数：

>>> A = 2.
>>> w0 = 1.  # rad/sec
>>> nin = 150
>>> nout = 100000 

随机生成样本时间：

>>> x = rng.uniform(0, 10*np.pi, nin) 

绘制所选时间的正弦波：

>>> y = A * np.cos(w0*x) 

定义用于计算周期图的频率数组：

>>> w = np.linspace(0.01, 10, nout) 

计算 Lomb-Scargle 周期图：

>>> import scipy.signal as signal
>>> pgram = signal.lombscargle(x, y, w, normalize=True) 

现在制作输入数据的图表：

>>> fig, (ax_t, ax_w) = plt.subplots(2, 1, constrained_layout=True)
>>> ax_t.plot(x, y, 'b+')
>>> ax_t.set_xlabel('Time [s]') 

然后绘制归一化周期图：

>>> ax_w.plot(w, pgram)
>>> ax_w.set_xlabel('Angular frequency [rad/s]')
>>> ax_w.set_ylabel('Normalized amplitude')
>>> plt.show() 

scipy.signal.vectorstrength

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.vectorstrength.html#scipy.signal.vectorstrength

scipy.signal.vectorstrength(events, period)

确定与给定周期对应的事件的矢量强度。

矢量强度是相位同步的一个度量，表明事件的定时如何与周期信号的单个周期同步。

如果使用多个周期，计算每个的矢量强度。这称为“共振矢量强度”。

参数：

events1D 数组类似

包含事件时间点的时间点数组。

periodfloat 或 array_like

事件应该与之同步的信号周期。周期与 events 单位相同。它也可以是周期数组，此时输出也是相同长度的数组。

返回：

strengthfloat 或 1D 数组

同步的强度。1.0 是完美同步，0.0 是没有同步。如果 period 是一个数组，则这也是一个数组，其中每个元素包含相应周期的矢量强度。

phasefloat 或 array

事件与弧度最强同步的相位。如果 period 是一个数组，则这也是一个数组，其中每个元素包含相应周期的相位。

参考文献：

van Hemmen, JL, Longtin, A 和 Vollmayr, AN。测试共振矢量

strength：听觉系统、电鱼和噪声。混沌 21, 047508 (2011); DOI:10.1063/1.3670512.

van Hemmen, JL。Goldberg、Brown 和 von Mises 后的矢量强度：

生物和数学视角。生物控制。2013 年 8 月;107(4):385-96. DOI:10.1007/s00422-013-0561-7.

van Hemmen, JL 和 Vollmayr, AN。共振矢量强度：发生了什么

当我们改变“探测”频率但保持尖峰时间不变时。生物控制。2013 年 8 月;107(4):491-94。DOI:10.1007/s00422-013-0560-8.

scipy.signal.ShortTimeFFT

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.ShortTimeFFT.html#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)。

stft 通过滑动窗口 (win) 在输入信号上以 hop 增量计算连续的 FFT，可用于量化频谱随时间的变化。

stft 由复数矩阵 S[q,p] 表示，其中第 p 列代表以时间 t[p] = p * delta_t = p * hop * T 居中的窗口的 FFT，其中 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_min或k_min这样的负索引不像标准 Python 索引中的倒数计数从数组末尾开始，而是位于 t = 0 的左侧。

更详细的信息可以在 SciPy 用户指南的短时傅里叶变换部分找到。

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

参数：

winnp.ndarray

窗口必须是一个实数或复数值的一维数组。

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。

示例

以下示例显示了带有变频 (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 

使用的高斯窗口为 50 个样本或 2.5 秒长。参数 mfft=200 在 ShortTimeFFT 中导致频谱过采样 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 的长度，因为在 hop 步骤中 SFT 的长度会增加：

>>> 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

Factor to multiply the STFT values by to scale each frequency slice to a magnitude spectrum.

fac_psd

Factor to multiply the STFT values by to scale each frequency slice to a power spectral density (PSD).

fft_mode

Mode of utilized FFT (‘twosided’, ‘centered’, ‘onesided’ or ‘onesided2X’).

fs

Sampling frequency of input signal and of the window.

hop

Time increment in signal samples for sliding window.

invertible

Check if STFT is invertible.

k_min

The smallest possible signal index of the STFT.

lower_border_end

First signal index and first slice index unaffected by pre-padding.

m_num

Number of samples in window win.

m_num_mid

Center index of window win.

mfft

Length of input for the FFT used - may be larger than window length m_num.

onesided_fft

Return True if a one-sided FFT is used.

p_min

The smallest possible slice index.

phase_shift

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

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, ...])	计算频谱图或交叉谱图。
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)	受后填充影响的第一个信号索引和第一个切片索引。

scipy.signal.stft

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.stft.html#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）。

STFT 可用作量化非平稳信号随时间的频率和相位内容变化的一种方法。

Legacy

此函数被视为传统功能，将不再接收更新。这可能意味着它将在未来的 SciPy 版本中被移除。ShortTimeFFT是一种新的 STFT / ISTFT 实现，具有更多功能。在SciPy 用户指南的教程-STFT部分中可以找到这两种实现的比较。

参数：

x类数组

测量值的时间序列

fs浮点数，可选

x时间序列的采样频率。默认为 1.0。

window字符串或元组或类数组，可选

欲使用的窗口。如果window为字符串或元组，则将其传递给get_window以生成窗口值，默认情况下为 DFT-even。有关窗口和必需参数的列表，请参阅get_window。如果window为类数组，则直接使用它作为窗口，其长度必须为 nperseg。默认为 Hann 窗口。

nperseg整数，可选

每个片段的长度。默认为 256。

noverlap整数，可选

分段之间重叠的点数。如果为None，则noverlap = nperseg // 2。默认为None。当指定时，必须满足 COLA 约束（见下面的说明）。

nfft整数，可选

使用的 FFT 长度，如果需要零填充的 FFT。如果为None，则 FFT 长度为nperseg。默认为None。

detrend字符串或函数或False，可选

指定如何对每个片段进行去趋势化处理。如果detrend为字符串，则将其作为detrend函数的type参数传递。如果它是一个函数，则接受一个片段并返回一个去趋势化的片段。如果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]。

paddedbool, optional

指定输入信号是否在末尾进行零填充，以使信号恰好适合整数个窗段，以便所有信号都包含在输出中。默认为True。如果boundary不是None，并且padded为True（默认情况下是这样），填充将在边界扩展之后进行。

axisint, optional

计算 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() 

比较信号 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 

scipy.signal.istft

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.istft.html#scipy.signal.istft

scipy.signal.istft(Zxx, fs=1.0, window='hann', nperseg=None, noverlap=None, nfft=None, input_onesided=True, boundary=True, time_axis=-1, freq_axis=-2, scaling='spectrum')

执行反短时傅立叶变换（iSTFT）。

传统

此函数被视为遗留版本，将不再接收更新。这也可能意味着它将在未来的 SciPy 版本中删除。ShortTimeFFT是一个新的 STFT / ISTFT 实现，具有更多功能。可以在SciPy 用户指南的短时傅立叶变换部分找到这些实现的比较。

参数：

Zxxarray_like

要重构的信号的 STFT。如果传递的是纯实数组，则将其转换为复杂数据类型。

fsfloat，可选

时间序列的采样频率。默认为 1.0。

windowstr 或 tuple 或 array_like，可选

所需使用的窗口。如果window是字符串或元组，则将其传递给get_window以生成窗口值，默认为 DFT-even。详见get_window获取窗口列表和所需参数。如果window是 array_like，则直接用作窗口，其长度必须为 nperseg。默认为 Hann 窗口。必须与用于生成 STFT 的窗口匹配，以确保忠实反演。

npersegint，可选

数据点数对应于每个 STFT 段。如果每段数据点数为奇数，或者 STFT 通过nfft > nperseg进行填充，则必须指定此参数。如果为None，则其值取决于Zxx和input_onesided的形状。如果input_onesided为 True，则nperseg=2*(Zxx.shape[freq_axis] - 1)。否则，nperseg=Zxx.shape[freq_axis]。默认为None。

noverlapint，可选

点之间重叠的点数。如果为None，则为段长度的一半。默认为None。在指定时，必须满足 COLA 约束（参见下面的注释），并且应与用于生成 STFT 的参数匹配。默认为None。

nfftint，可选

FFT 点数对应于每个 STFT 段。如果 STFT 通过nfft > nperseg进行填充，则必须指定此参数。如果为None，则默认值与nperseg相同，详见上文，但有一例外：如果input_onesided为 True 且nperseg==2*Zxx.shape[freq_axis] - 1，则nfft也取该值。这种情况允许使用nfft=None正确反演奇数长度未填充的 STFT。默认为None。

input_onesidedbool，可选

如果为True，将输入数组解释为单边 FFT，例如由stft返回的return_onesided=True和numpy.fft.rfft。如果为False，将输入解释为双边 FFT。默认为True。

boundarybool, 可选

指定输入信号是否通过向 stft 提供非None boundary 参数来在其边界上扩展。默认为True。

time_axisint, 可选

STFT 的时间段所在位置；默认为最后一轴（即axis=-1）。

freq_axisint, 可选

STFT 的频率轴所在位置；默认为倒数第二轴（即axis=-2）。

scaling: {‘spectrum’, ‘psd’}

默认的'spectrum'缩放允许解释Zxx的每个频率线为幅度谱。'psd'选项将每行缩放到功率谱密度 - 允许通过数值积分计算信号的能量 abs(Zxx)**2。

返回:

tndarray

输出数据数组的时间。

xndarray

Zxx的逆短时傅立叶变换。

另请参阅

stft

短时傅立叶变换

ShortTimeFFT

更多功能的新 STFT/ISTFT 实现。

check_COLA

检查是否满足 Constant OverLap Add (COLA)约束

check_NOLA

检查是否满足 Nonzero Overlap Add (NOLA)约束

注意事项

为了通过istft反转 STFT 以进行反 STFT，信号窗必须遵守“非零重叠添加”（NOLA）约束：

[\sum_{t}w^{2}[n-tH] \ne 0]

这确保了出现在重叠添加重建方程分母中的归一化因子

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

不为零。使用check_NOLA函数可以检查 NOLA 约束。

已修改的 STFT（通过掩蔽或其他方式）不能保证与确切可实现信号对应。该函数通过最小二乘估计算法实现了 iSTFT，该算法详细说明见[2]，其生成的信号最小化了返回信号的 STFT 和修改后 STFT 之间的均方误差。

版本 0.19.0 中的新功能。

参考文献

[1]

Oppenheim, Alan V., Ronald W. Schafer, John R. Buck “离散时间信号处理”，Prentice Hall，1999 年。

[2]

Daniel W. Griffin, Jae S. Lim “从修改后的短时傅里叶变换估计信号”, 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 的 50Hz 正弦波，受 1024 Hz 采样的 0.001 V**2/Hz 白噪声的影响。

>>> fs = 1024
>>> N = 10*fs
>>> nperseg = 512
>>> amp = 2 * np.sqrt(2)
>>> noise_power = 0.001 * fs / 2
>>> time = np.arange(N) / float(fs)
>>> carrier = amp * np.sin(2*np.pi*50*time)
>>> noise = rng.normal(scale=np.sqrt(noise_power),
...                    size=time.shape)
>>> x = carrier + noise 

计算 STFT，并绘制其幅度

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

将幅度为载波幅度的 10%或更少的分量置零，然后通过逆 STFT 转换回时间序列

>>> Zxx = np.where(np.abs(Zxx) >= amp/10, Zxx, 0)
>>> _, xrec = signal.istft(Zxx, fs) 

将清理后的信号与原始和真实的载波信号进行比较。

>>> plt.figure()
>>> plt.plot(time, x, time, xrec, time, carrier)
>>> plt.xlim([2, 2.1])
>>> plt.xlabel('Time [sec]')
>>> plt.ylabel('Signal')
>>> plt.legend(['Carrier + Noise', 'Filtered via STFT', 'True Carrier'])
>>> plt.show() 

注意，清理后的信号并不像原始信号那样突然开始，因为某些瞬态的系数也被移除了：

>>> plt.figure()
>>> plt.plot(time, x, time, xrec, time, carrier)
>>> plt.xlim([0, 0.1])
>>> plt.xlabel('Time [sec]')
>>> plt.ylabel('Signal')
>>> plt.legend(['Carrier + Noise', 'Filtered via STFT', 'True Carrier'])
>>> plt.show() 

scipy.signal.check_COLA

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.check_COLA.html#scipy.signal.check_COLA

scipy.signal.check_COLA(window, nperseg, noverlap, tol=1e-10)

检查常量重叠添加（COLA）约束是否满足。

参数：

window字符串或元组或 array_like

所需使用的窗口。如果 window 是字符串或元组，则将其传递给 get_window 以生成窗口值，默认情况下为 DFT-even。有关窗口和所需参数的列表，请参见 get_window。如果 window 是 array_like，则将其直接用作窗口，其长度必须为 nperseg。

nperseg整数

每个片段的长度。

noverlap整数

段之间重叠的点数。

tol浮点数，可选

每个频段加权和与中位数频段和的允许方差。

返回：

verdict布尔值

True 如果选择的组合在 tol 范围内满足 COLA，否则 False

请参见

check_NOLA

检查是否满足非零重叠添加（NOLA）约束

stft

短时傅里叶变换

istft

逆短时傅里叶变换

注释

为了通过逆短时傅里叶变换中的逆 STFT 实现 STFT 的反演，在 istft 中，只需确保信号窗口符合“常数重叠添加”（COLA）的约束即可。这确保了输入数据中的每个点都被等权重，从而避免混叠，并允许完全重建。

满足 COLA 的一些窗口示例：

• 重叠为 0、1/2、2/3、3/4 等的矩形窗口

• Bartlett 窗口在 1/2、3/4、5/6 等重叠时

• Hann 窗口在 1/2、2/3、3/4 等重叠时

• 任何 Blackman 家族窗口的 2/3 重叠

• 任何具有 noverlap = nperseg-1 的窗口

在[2]中可以找到其他窗口的非常全面的列表，在“幅度平坦度”为单位时满足 COLA 条件。

从版本 0.19.0 开始新增。

参考文献

[1]

Julius O. Smith III，《Spectral Audio Signal Processing》，W3K Publishing，2011 年，ISBN 978-0-9745607-3-1。

[2]

G. Heinzel, A. Ruediger and R. Schilling，《Spectrum and spectral density estimation by the Discrete Fourier transform (DFT)，including a comprehensive list of window functions and some new at-top windows》，2002 年，hdl.handle.net/11858/00-001M-0000-0013-557A-5

示例

>>> from scipy import signal 

确认 75%（3/4）重叠的矩形窗口的 COLA 条件：

>>> signal.check_COLA(signal.windows.boxcar(100), 100, 75)
True 

对于 25%（1/4）重叠，COLA 不成立：

>>> signal.check_COLA(signal.windows.boxcar(100), 100, 25)
False 

“对称”Hann 窗口（用于滤波器设计）不满足 COLA：

>>> signal.check_COLA(signal.windows.hann(120, sym=True), 120, 60)
False 

“周期性”或“DFT-even”的 Hann 窗口（用于 FFT 分析）在 1/2、2/3、3/4 等重叠情况下是 COLA 的：

>>> signal.check_COLA(signal.windows.hann(120, sym=False), 120, 60)
True 

>>> signal.check_COLA(signal.windows.hann(120, sym=False), 120, 80)
True 

>>> signal.check_COLA(signal.windows.hann(120, sym=False), 120, 90)
True 

scipy.signal.check_NOLA

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.check_NOLA.html#scipy.signal.check_NOLA

scipy.signal.check_NOLA(window, nperseg, noverlap, tol=1e-10)

检查是否满足非零重叠添加（NOLA）约束。

参数：

windowstr 或 tuple 或 array_like

要使用的期望窗口。如果window是字符串或元组，则将其传递给get_window以生成窗口值，默认情况下为 DFT-even。有关窗口和所需参数的列表，请参见get_window。如果window是 array_like，则将其直接用作窗口，其长度必须为 nperseg。

npersegint

每个段的长度。

noverlapint

分段之间重叠的点数。

tolfloat, 可选

一个频段的加权和与中位数频段的加权和的允许方差。

返回：

verdictbool

如果所选组合符合tol内的 NOLA 约束条件，则返回True，否则返回False

参见

check_COLA

检查是否满足恒定重叠添加（COLA）约束

stft

短时傅立叶变换

istft

逆短时傅立叶变换

注释

为了通过istft中的逆 STFT 启用 STFT 的反演，信号窗必须遵守“非零重叠添加”（NOLA）约束：

[\sum_{t}w^{2}[n-tH] \ne 0]

对于所有的(n)，其中(w)是窗口函数，(t)是帧索引，(H)是跨步大小（(H) = nperseg - noverlap）。

这确保了重叠添加反演方程中的归一化因子不为零。只有非常异常的窗口才会不满足 NOLA 约束。

1.2.0 版中的新功能。

参考文献

[1]

Julius O. Smith III，“音频信号谱分析”，W3K Publishing，2011 年，ISBN 978-0-9745607-3-1。

[2]

G. Heinzel, A. Ruediger and R. Schilling, “离散傅立叶变换（DFT）估计的频谱和谱密度，包括详细的窗函数列表和一些新的顶部窗口”，2002 年，hdl.handle.net/11858/00-001M-0000-0013-557A-5

示例

>>> import numpy as np
>>> from scipy import signal 

确认 75%（3/4）重叠的矩形窗口的 NOLA 条件：

>>> signal.check_NOLA(signal.windows.boxcar(100), 100, 75)
True 

对于 25%（1/4）重叠，NOLA 也成立：

>>> signal.check_NOLA(signal.windows.boxcar(100), 100, 25)
True 

“对称”Hann 窗（用于滤波器设计）也满足 NOLA：

>>> signal.check_NOLA(signal.windows.hann(120, sym=True), 120, 60)
True 

只要有重叠，就需要非常异常的窗口才能不满足 NOLA：

>>> w = np.ones(64, dtype="float")
>>> w[::2] = 0
>>> signal.check_NOLA(w, 64, 32)
False 

如果重叠不足，末端带有零的窗口将无法工作：

>>> signal.check_NOLA(signal.windows.hann(64), 64, 0)
False
>>> signal.check_NOLA(signal.windows.hann(64), 64, 1)
False
>>> signal.check_NOLA(signal.windows.hann(64), 64, 2)
True 

scipy.signal.czt

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/czt-function.html#scipy.signal.czt

scipy.signal.czt(x, m=None, w=None, a=1 + 0j, *, axis=-1)

计算 Z 平面中螺旋周围的频率响应。

参数：

x：array

要变换的信号。

m：int，可选

所需输出点的数量。默认为输入数据的长度。

w：complex，可选

在每个步骤中点之间的比率。这必须精确，否则累积误差将使输出序列的尾部退化。默认为整个单位圆周围均匀分布的点。

a：complex，可选

复平面中的起始点。默认为 1+0j。

axis：int，可选

计算 FFT 的轴。如果未给出，则使用最后一个轴。

返回：

out：ndarray

一个与 x 相同尺寸的数组，但是变换轴的长度设置为 m。

另见

CZT

创建可调用的啁啾 z 变换函数的类。

zoom_fft

部分 FFT 计算的便捷函数。

注释

默认值选取为 signal.czt(x) 等同于 fft.fft(x)，如果 m > len(x)，则 signal.czt(x, m) 等同于 fft.fft(x, m)。

如果需要重复变换，请使用 CZT 来构建一个专门的变换函数，可以在不重新计算常量的情况下重复使用。

一个示例应用是在系统识别中，重复评估系统的 Z 变换的小片段，以精炼估计极点的真实位置。

参考文献

[1]

Steve Alan Shilling，《啁啾 z 变换及其应用研究》，第 20 页（1970）krex.k-state.edu/dspace/bitstream/handle/2097/7844/LD2668R41972S43.pdf

示例

生成一个正弦波：

>>> import numpy as np
>>> f1, f2, fs = 8, 10, 200  # Hz
>>> t = np.linspace(0, 1, fs, endpoint=False)
>>> x = np.sin(2*np.pi*t*f2)
>>> import matplotlib.pyplot as plt
>>> plt.plot(t, x)
>>> plt.axis([0, 1, -1.1, 1.1])
>>> plt.show() 

其离散傅里叶变换的能量全集中在单一频率箱中：

>>> from scipy.fft import rfft, rfftfreq
>>> from scipy.signal import czt, czt_points
>>> plt.plot(rfftfreq(fs, 1/fs), abs(rfft(x)))
>>> plt.margins(0, 0.1)
>>> plt.show() 

然而，如果正弦波是对数衰减的：

>>> x = np.exp(-t*f1) * np.sin(2*np.pi*t*f2)
>>> plt.plot(t, x)
>>> plt.axis([0, 1, -1.1, 1.1])
>>> plt.show() 

DFT 将具有频谱泄漏：

>>> plt.plot(rfftfreq(fs, 1/fs), abs(rfft(x)))
>>> plt.margins(0, 0.1)
>>> plt.show() 

尽管 DFT 总是在单位圆周围采样 Z 变换，啁啾 z 变换允许我们沿任何对数螺旋（例如半径小于单位的圆）采样 Z 变换：

>>> M = fs // 2  # Just positive frequencies, like rfft
>>> a = np.exp(-f1/fs)  # Starting point of the circle, radius < 1
>>> w = np.exp(-1j*np.pi/M)  # "Step size" of circle
>>> points = czt_points(M + 1, w, a)  # M + 1 to include Nyquist
>>> plt.plot(points.real, points.imag, '.')
>>> plt.gca().add_patch(plt.Circle((0,0), radius=1, fill=False, alpha=.3))
>>> plt.axis('equal'); plt.axis([-1.05, 1.05, -0.05, 1.05])
>>> plt.show() 

使用正确的半径，这将转换衰减正弦波（以及具有相同衰减率的其他波形），而不会出现频谱泄漏：

>>> z_vals = czt(x, M + 1, w, a)  # Include Nyquist for comparison to rfft
>>> freqs = np.angle(points)*fs/(2*np.pi)  # angle = omega, radius = sigma
>>> plt.plot(freqs, abs(z_vals))
>>> plt.margins(0, 0.1)
>>> plt.show() 

scipy.signal.zoom_fft

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.zoom_fft.html#scipy.signal.zoom_fft

scipy.signal.zoom_fft(x, fn, m=None, *, fs=2, endpoint=False, axis=-1)

仅计算范围在 fn 中的频率的 x 的 DFT。

参数：

x：数组

要变换的信号。

fn：类似数组

长度为 2 的序列 [f1, f2] 给出频率范围，或者一个标量，其中假设范围为 [0, fn]。

m：整数，可选

要评估的点数。默认为 x 的长度。

fs：浮点数，可选

采样频率。例如，如果 fs=10 表示 10 kHz，那么 f1 和 f2 也应该以 kHz 表示。默认采样频率为 2，因此 f1 和 f2 应在 [0, 1] 范围内以保持变换低于奈奎斯特频率。

endpoint：布尔值，可选

如果为 True，f2 是最后一个样本。否则，不包括它。默认为 False。

axis：整数，可选

计算 FFT 的轴。如果未给出，则使用最后一个轴。

返回：

out：ndarray

转换后的信号。傅里叶变换将在点 f1, f1+df, f1+2df, …, f2 处计算，其中 df=(f2-f1)/m。

另请参阅

ZoomFFT

创建一个可调用的部分 FFT 函数的类。

注意事项

默认选择这样，使得 signal.zoom_fft(x, 2) 等价于 fft.fft(x)，如果 m > len(x)，那么 signal.zoom_fft(x, 2, m) 等价于 fft.fft(x, m)。

要绘制结果变换的幅度图，请使用：

plot(linspace(f1, f2, m, endpoint=False), abs(zoom_fft(x, [f1, f2], m))) 

如果需要重复变换，请使用 ZoomFFT 构建一个专门的变换函数，可以在不重新计算常数的情况下重复使用。

示例

要绘制变换结果，请使用类似以下的方法：

>>> import numpy as np
>>> from scipy.signal import zoom_fft
>>> t = np.linspace(0, 1, 1021)
>>> x = np.cos(2*np.pi*15*t) + np.sin(2*np.pi*17*t)
>>> f1, f2 = 5, 27
>>> X = zoom_fft(x, [f1, f2], len(x), fs=1021)
>>> f = np.linspace(f1, f2, len(x))
>>> import matplotlib.pyplot as plt
>>> plt.plot(f, 20*np.log10(np.abs(X)))
>>> plt.show() 

scipy.signal.CZT

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.CZT.html#scipy.signal.CZT

class scipy.signal.CZT(n, m=None, w=None, a=1 + 0j)

创建一个可调用的啁啾变换函数。

转换以计算螺旋周围的频率响应。此类对象是可调用的，可以在其输入上计算啁啾变换。此对象预先计算给定变换中使用的恒定啁啾。

参数：

nint

信号的大小。

mint，可选

所需的输出点数。默认为n。

wcomplex，可选

每个步骤中点之间的比例。这必须是精确的，否则累积误差将降低输出序列的尾部。默认为整个单位圆周围均匀分布的点。

acomplex，可选

复平面中的起始点。默认为 1+0j。

返回：

fCZT

用于在x上计算啁啾变换的可调用对象f(x, axis=-1)。

另请参见

czt

用于快速计算 CZT 的便利函数。

ZoomFFT

创建可调用的部分 FFT 函数的类。

注意事项

默认值选为使f(x)等同于fft.fft(x)，如果m > len(x)，则使f(x, m)等同于fft.fft(x, m)。

如果w不位于单位圆上，则变换将围绕指数增长半径的螺旋进行。无论如何，角度将线性增加。

对于位于单位圆上的变换，当使用ZoomFFT时，精度更高，因为w中的任何数值误差在长数据长度上累积，偏离单位圆。

与等效零填充 FFT 相比，啁啾变换可能更快。尝试使用您自己的数组大小进行测试。

然而，啁啾变换的精度明显低于等效的零填充 FFT。

由于此 CZT 使用 Bluestein 算法实现，因此可以在 O(N log N)时间内计算大素数长度的傅里叶变换，而不是直接 DFT 计算所需的 O(N**2)时间。(scipy.fft也使用 Bluestein 的算法。)

（“啁啾变换”名称来自 Bluestein 算法中使用的啁啾。它不像其他带有“啁啾”名称的变换那样将信号分解为啁啾。）

参考文献

[1]

Leo I. Bluestein，“离散傅里叶变换计算的线性滤波方法”，东北电子研究与工程会议记录 10，218-219（1968）。

[2]

Rabiner、Schafer 和 Rader，“啁啾变换算法及其应用”，贝尔系统技术杂志 48，1249-1292（1969）。

示例

计算多个素数长度 FFT：

>>> from scipy.signal import CZT
>>> import numpy as np
>>> a = np.random.rand(7)
>>> b = np.random.rand(7)
>>> c = np.random.rand(7)
>>> czt_7 = CZT(n=7)
>>> A = czt_7(a)
>>> B = czt_7(b)
>>> C = czt_7(c) 

显示计算 FFT 的点：

>>> czt_7.points()
array([ 1.00000000+0.j        ,  0.62348980+0.78183148j,
 -0.22252093+0.97492791j, -0.90096887+0.43388374j,
 -0.90096887-0.43388374j, -0.22252093-0.97492791j,
 0.62348980-0.78183148j])
>>> import matplotlib.pyplot as plt
>>> plt.plot(czt_7.points().real, czt_7.points().imag, 'o')
>>> plt.gca().add_patch(plt.Circle((0,0), radius=1, fill=False, alpha=.3))
>>> plt.axis('equal')
>>> plt.show() 

方法

__call__(x, *[, axis])	计算信号的奇异变换。
points()	返回进行奇异变换的点的位置。

scipy.signal.ZoomFFT

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.ZoomFFT.html#scipy.signal.ZoomFFT

class scipy.signal.ZoomFFT(n, fn, m=None, *, fs=2, endpoint=False)

创建一个可调用的变焦 FFT 变换函数。

这是圆周单位周围等间距频率的啁啾变换（CZT）的特化，用于比计算整个 FFT 并截断更有效地计算 FFT 的一部分。

参数：

n整数

信号的大小。

fn类似数组

长度为 2 的序列[f1, f2]表示频率范围，或者标量，假定范围[0, fn]。

m整数，可选

评估点数。默认为n。

fs浮点数，可选

采样频率。例如，如果fs=10表示 10 kHz，则f1和f2也应以 kHz 为单位。默认的采样频率为 2，因此f1和f2的范围应在[0, 1]之间，以使变换保持在奈奎斯特频率以下。

endpoint布尔值，可选

如果为 True，则f2为最后一个样本。否则，不包括在内。默认为 False。

返回：

fZoomFFT

可调用对象f(x, axis=-1)用于计算x上的变焦 FFT。

另请参阅

zoom_fft

用于计算变焦 FFT 的便捷函数。

注：

默认设置使得f(x, 2)等同于fft.fft(x)，如果m > len(x)，那么f(x, 2, m)等同于fft.fft(x, m)。

采样频率是信号x中样本之间的时间步长的倒数。单位圆对应从 0 到采样频率的频率。默认的采样频率为 2，因此f1和f2的值应在范围[0, 1)内，以保持变换在奈奎斯特频率以下。

请记住，变焦 FFT 只能插值现有 FFT 的点。它无法帮助解决两个分开的附近频率。只能通过增加采集时间来增加频率分辨率。

这些函数使用 Bluestein 算法实现（就像scipy.fft一样）。[2]

参考文献

[1]

Steve Alan Shilling，“啁啾变换及其应用研究”，第 29 页（1970 年）krex.k-state.edu/dspace/bitstream/handle/2097/7844/LD2668R41972S43.pdf

[2]

Leo I. Bluestein，“离散傅立叶变换的线性滤波方法”，东北电子研究与工程会议记录第 10 卷，218-219 页（1968 年）。

示例

要绘制变换结果，请使用类似以下的内容：

>>> import numpy as np
>>> from scipy.signal import ZoomFFT
>>> t = np.linspace(0, 1, 1021)
>>> x = np.cos(2*np.pi*15*t) + np.sin(2*np.pi*17*t)
>>> f1, f2 = 5, 27
>>> transform = ZoomFFT(len(x), [f1, f2], len(x), fs=1021)
>>> X = transform(x)
>>> f = np.linspace(f1, f2, len(x))
>>> import matplotlib.pyplot as plt
>>> plt.plot(f, 20*np.log10(np.abs(X)))
>>> plt.show() 

方法

__call__ (x, *[, axis])	计算信号的奇异变换。
points()	返回进行奇异变换计算的点。

scipy.signal.czt_points

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.czt_points.html#scipy.signal.czt_points

scipy.signal.czt_points(m, w=None, a=1 + 0j)

返回进行啁啾变换的点。

参数：

mint

所需点的数量。

w复数，可选

每个步骤中点的比率。默认为均匀分布在整个单位圆周围的点。

a复数，可选

复平面中的起始点。默认为 1+0j。

返回：

outndarray

Z 平面中的点，CZT 在调用时以复数m、w和a作为参数进行 z 变换采样。

另请参阅：

CZT

创建一个可调用的啁啾变换函数的类。

czt

用于快速计算 CZT 的便捷函数。

示例

绘制 16 点 FFT 的点：

>>> import numpy as np
>>> from scipy.signal import czt_points
>>> points = czt_points(16)
>>> import matplotlib.pyplot as plt
>>> plt.plot(points.real, points.imag, 'o')
>>> plt.gca().add_patch(plt.Circle((0,0), radius=1, fill=False, alpha=.3))
>>> plt.axis('equal')
>>> plt.show() 

和一个穿过单位圆的 91 点对数螺旋：

>>> m, w, a = 91, 0.995*np.exp(-1j*np.pi*.05), 0.8*np.exp(1j*np.pi/6)
>>> points = czt_points(m, w, a)
>>> plt.plot(points.real, points.imag, 'o')
>>> plt.gca().add_patch(plt.Circle((0,0), radius=1, fill=False, alpha=.3))
>>> plt.axis('equal')
>>> plt.show() 

稀疏矩阵（scipy.sparse）

原文：docs.scipy.org/doc/scipy-1.12.0/reference/sparse.html

SciPy 二维稀疏数组包，用于数值数据。

注意

此软件包正在切换到与 NumPy 数组兼容的数组接口，而不再使用旧的矩阵接口。我们建议您对所有新工作使用数组对象（bsr_array, coo_array 等）。

在使用数组接口时，请注意：

• x * y 现在不再执行矩阵乘法，而是执行元素级乘法（与 NumPy 数组类似）。为了使代码同时适用于数组和矩阵，使用 x @ y 来进行矩阵乘法。

• 诸如 sum 的操作，原先生成密集矩阵，现在生成数组，其乘法行为类似但有所不同。

• 稀疏数组目前必须是二维的。这也意味着这些对象上的所有 切片 操作必须产生二维结果，否则将导致错误。这将在未来版本中解决。

构造实用程序（eye, kron, random, diags 等）尚未移植完成，但可以将其结果封装成数组：

A = csr_array(eye(3)) 

内容

稀疏数组类

bsr_array(arg1[, shape, dtype, copy, blocksize])	块稀疏行（Block Sparse Row）格式的稀疏数组。
coo_array(arg1[, shape, dtype, copy])	COOrdinate 格式的稀疏数组。
csc_array(arg1[, shape, dtype, copy])	压缩稀疏列（Compressed Sparse Column）数组。
csr_array(arg1[, shape, dtype, copy])	压缩稀疏行（Compressed Sparse Row）数组。
dia_array(arg1[, shape, dtype, copy])	带有对角线存储的稀疏数组。
dok_array(arg1[, shape, dtype, copy])	基于键的字典（Dictionary Of Keys）稀疏数组。
lil_array(arg1[, shape, dtype, copy])	基于行的列表（LIst of Lists）稀疏数组。
sparray()	该类为所有稀疏数组提供基类。

稀疏矩阵类：

bsr_matrix(arg1[, shape, dtype, copy, blocksize])	块稀疏行格式的稀疏矩阵。
coo_matrix(arg1[, shape, dtype, copy])	COO 格式的稀疏矩阵。
csc_matrix(arg1[, shape, dtype, copy])	压缩稀疏列矩阵。
csr_matrix(arg1[, shape, dtype, copy])	压缩稀疏行矩阵。
dia_matrix(arg1[, shape, dtype, copy])	带有对角线存储的稀疏矩阵。
dok_matrix(arg1[, shape, dtype, copy])	基于键的字典稀疏矩阵。
lil_matrix(arg1[, shape, dtype, copy])	基于行的链表稀疏矩阵。
spmatrix()	该类为所有稀疏矩阵类提供基类。

函数：

构建稀疏数组：

diags_array(diagonals, /, *[, offsets, ...])	从对角线构造稀疏数组。
eye_array(m[, n, k, dtype, format])	稀疏数组格式中的单位矩阵
random_array(shape, *[, density, format, ...])	返回一个 0, 1) 范围内均匀随机数的稀疏数组
[block_array(blocks, *[, format, dtype])	从稀疏子块构建稀疏数组

构建稀疏矩阵：

eye(m[, n, k, dtype, format])	对角线上有 1 的稀疏矩阵
identity(n[, dtype, format])	稀疏格式中的单位矩阵
diags(diagonals[, offsets, shape, format, dtype])	从对角线构造稀疏矩阵。
spdiags(data, diags[, m, n, format])	从对角线返回稀疏矩阵。
bmat(blocks[, format, dtype])	从稀疏子块构建稀疏数组或矩阵
random(m, n[, density, format, dtype, ...])	生成给定形状和密度的稀疏矩阵，值为随机分布。
rand(m, n[, density, format, dtype, ...])	生成给定形状和密度的稀疏矩阵，值均匀分布。

从更小的结构（数组或矩阵）构建更大的结构

kron(A, B[, format])	稀疏矩阵 A 和 B 的 Kronecker 乘积
kronsum(A, B[, format])	方阵稀疏矩阵 A 和 B 的 Kronecker 和
block_diag(mats[, format, dtype])	从提供的矩阵构建块对角稀疏矩阵或数组
tril(A[, k, format])	返回稀疏数组或矩阵的下三角部分
triu(A[, k, format])	返回稀疏数组或矩阵的上三角部分
hstack(blocks[, format, dtype])	水平堆叠稀疏矩阵（按列）
vstack(blocks[, format, dtype])	垂直堆叠稀疏数组（按行）

保存和加载稀疏矩阵：

save_npz(file, matrix[, compressed])	使用.npz格式将稀疏矩阵或数组保存到文件中。
load_npz(file)	使用.npz格式从文件加载稀疏数组/矩阵。

稀疏工具：

find(A)	返回矩阵非零元素的索引和值

辨识稀疏数组：

• 使用 isinstance(A, sp.sparse.sparray) 检查是否为数组或矩阵。

• 使用 A.format == ‘csr’ 来检查稀疏格式

辨识稀疏矩阵：

issparse(x)	x 是否为稀疏数组或稀疏矩阵类型？
isspmatrix(x)	x 是否为稀疏矩阵类型？
isspmatrix_csc(x)	x 是否为 csc_matrix 类型？
isspmatrix_csr(x)	x 是否为 csr_matrix 类型？
isspmatrix_bsr(x)	x 是否为 bsr_matrix 类型？
isspmatrix_lil(x)	x 是否为 lil_matrix 类型？
isspmatrix_dok(x)	x 是否为 dok_array 类型？
isspmatrix_coo(x)	x 是否为 coo_matrix 类型？
isspmatrix_dia(x)	x 是否为 dia_matrix 类型？

子模块

csgraph	压缩稀疏图例程 (scipy.sparse.csgraph)
linalg	稀疏线性代数 (scipy.sparse.linalg)

异常情况

SparseEfficiencyWarning
SparseWarning

使用信息

有七种可用的稀疏数组类型：

1. csc_array: 压缩稀疏列格式
3. csr_array: 压缩稀疏行格式
5. bsr_array: 块稀疏行格式
7. lil_array: 列表列表格式
9. dok_array: 键字典格式
11. coo_array: COO 格式（即 IJV，三元组格式）
13. dia_array: 对角线格式

要高效构造数组，请使用dok_array或者lil_array。lil_array类支持基本切片和与 NumPy 数组类似语法的花式索引。正如下文所示，COO 格式也可用于高效构造数组。尽管它们与 NumPy 数组相似，强烈不建议直接在这些数组上使用 NumPy 函数，因为 NumPy 可能无法正确转换它们以进行计算，导致意外（和错误）的结果。如果确实要在这些数组上应用 NumPy 函数，请首先检查 SciPy 是否有适用于给定稀疏数组类的自己的实现，或者在应用方法之前将稀疏数组转换为 NumPy 数组（例如，使用类的toarray方法）。

要执行诸如乘法或求逆之类的操作，首先将数组转换为 CSC 或 CSR 格式。lil_array格式是基于行的，因此转换为 CSR 是有效的，而转换为 CSC 则不太有效。

在 CSR、CSC 和 COO 格式之间的所有转换都是高效的、线性时间的操作。

矩阵向量乘积

要在稀疏数组和向量之间进行向量乘积，简单地使用数组的dot方法，如其文档字符串中所述：

>>> import numpy as np
>>> from scipy.sparse import csr_array
>>> A = csr_array([[1, 2, 0], [0, 0, 3], [4, 0, 5]])
>>> v = np.array([1, 0, -1])
>>> A.dot(v)
array([ 1, -3, -1], dtype=int64) 

警告

从 NumPy 1.7 版本开始，np.dot不知道稀疏数组，因此使用它将导致意外的结果或错误。应该首先获得相应的密集数组：

>>> np.dot(A.toarray(), v)
array([ 1, -3, -1], dtype=int64) 

但这样一来，所有的性能优势都会丧失。

CSR 格式特别适合快速矩阵向量乘积。

示例 1

构造一个 1000x1000 的lil_array并给它添加一些值：

>>> from scipy.sparse import lil_array
>>> from scipy.sparse.linalg import spsolve
>>> from numpy.linalg import solve, norm
>>> from numpy.random import rand 

>>> A = lil_array((1000, 1000))
>>> A[0, :100] = rand(100)
>>> A[1, 100:200] = A[0, :100]
>>> A.setdiag(rand(1000)) 

现在将其转换为 CSR 格式并解决 A x = b 得到 x：

>>> A = A.tocsr()
>>> b = rand(1000)
>>> x = spsolve(A, b) 

将其转换为密集数组并求解，并检查结果是否相同：

>>> x_ = solve(A.toarray(), b) 

现在我们可以计算误差的范数：

>>> err = norm(x-x_)
>>> err < 1e-10
True 

应该很小 :)

示例 2

在 COO 格式中构造一个数组：

>>> from scipy import sparse
>>> from numpy import array
>>> I = array([0,3,1,0])
>>> J = array([0,3,1,2])
>>> V = array([4,5,7,9])
>>> A = sparse.coo_array((V,(I,J)),shape=(4,4)) 

注意索引不需要排序。

在转换为 CSR 或 CSC 时，重复的(i,j)条目将被求和。

>>> I = array([0,0,1,3,1,0,0])
>>> J = array([0,2,1,3,1,0,0])
>>> V = array([1,1,1,1,1,1,1])
>>> B = sparse.coo_array((V,(I,J)),shape=(4,4)).tocsr() 

这对于构造有限元刚度和质量矩阵非常有用。

进一步细节

CSR 列索引不一定排序。同样适用于 CSC 行索引。当需要排序索引时，请使用.sorted_indices()和.sort_indices()方法（例如，将数据传递给其他库时）。