meFFT >>> 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 The utilized Gaussian window is 50 samples or 2.5 s long. The parameter ``mfft=200`` in `ShortTimeFFT` causes the spectrum to be oversampled by a factor of 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 In the plot, the time extent of the signal `x` is marked by vertical dashed lines. Note that the SFT produces values outside the time range of `x`. The shaded areas on the left and the right indicate border effects caused by the window slices in that area not fully being inside time range of `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() Reconstructing the signal with the `~ShortTimeFFT.istft` is straightforward, but note that the length of `x1` should be specified, since the SFT length increases in `hop` steps: >>> SFT.invertible # check if invertible True >>> x1 = SFT.istft(Sx, k1=N) >>> np.allclose(x, x1) True It is possible to calculate the SFT of signal parts: >>> N2 = SFT.nearest_k_p(N // 2) >>> Sx0 = SFT.stft(x[:N2]) >>> Sx1 = SFT.stft(x[N2:]) When assembling sequential STFT parts together, the overlap needs to be considered: >>> 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 It is also possible to calculate the `itsft` for signal parts: >>> y_p = SFT.istft(Sx, N//3, N//2) >>> np.allclose(y_p, x[N//3:N//2]) True Ú