.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "api/gallery/signal_processing/spectral/plot_03_gabor_stft.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code or to run this example in your browser via JupyterLite. .. rst-class:: sphx-glr-example-title .. _sphx_glr_api_gallery_signal_processing_spectral_plot_03_gabor_stft.py: Gabor's time-frequency atoms (1946) =================================== A Fourier spectrum says which frequencies a signal contains but not when. Dennis Gabor (1946) described signals instead by "atoms": Gaussian-windowed sinusoids, each occupying a cell of the time-frequency plane. Sliding a window along the signal and taking a DFT of each frame gives the short-time Fourier transform (:func:`~engineeringkit.signal_processing.stft`). The cells cannot be made arbitrarily small. For RMS widths :math:`\sigma_t` and :math:`\sigma_f` of a window's energy in time and frequency, :math:`\sigma_t\sigma_f \ge 1/(4\pi)`, and only the Gaussian reaches the bound. .. GENERATED FROM PYTHON SOURCE LINES 17-23 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from scipy import signal from engineeringkit.signal_processing import stft .. GENERATED FROM PYTHON SOURCE LINES 24-31 The uncertainty product of a window ----------------------------------- The RMS widths are computed from the window's energy :math:`|w|^2` and from its finely sampled spectrum :math:`|W|^2`. The Gaussian meets the bound, the smooth Hann window comes within a few percent, and the abrupt edges of the rectangular window spread its energy far out in frequency. .. GENERATED FROM PYTHON SOURCE LINES 31-55 .. code-block:: Python def spread(window, fs): """RMS duration [s], RMS bandwidth [Hz] and their product.""" n = window.size t = (np.arange(n) - (n - 1) / 2) / fs e_t = window**2 / np.sum(window**2) sigma_t = np.sqrt(np.sum(e_t * t**2)) W = np.fft.fftshift(np.fft.fft(window, 64 * n)) f = np.fft.fftshift(np.fft.fftfreq(64 * n, 1 / fs)) e_f = abs(W) ** 2 / np.sum(abs(W) ** 2) sigma_f = np.sqrt(np.sum(e_f * f**2)) return sigma_t, sigma_f, sigma_t * sigma_f fs, n = 1000.0, 256 bound = 1 / (4 * np.pi) products = {} for name, w in {"gaussian": signal.get_window(("gaussian", 32), n), "hann": signal.get_window("hann", n), "boxcar": np.ones(n)}.items(): st, sf, products[name] = spread(w, fs) print(f"{name:8s}: sigma_t = {1e3 * st:5.1f} ms, sigma_f = {sf:6.1f} Hz, product = {products[name] / bound:5.3f} x 1/(4 pi)") assert np.isclose(products["gaussian"], bound, rtol=1e-3) assert products["hann"] > 1.01 * bound and products["boxcar"] > 10 * bound .. rst-class:: sphx-glr-script-out .. code-block:: none gaussian: sigma_t = 22.6 ms, sigma_f = 3.5 Hz, product = 1.000 x 1/(4 pi) hann : sigma_t = 36.2 ms, sigma_f = 2.3 Hz, product = 1.026 x 1/(4 pi) boxcar : sigma_t = 73.9 ms, sigma_f = 16.6 Hz, product = 15.381 x 1/(4 pi) .. GENERATED FROM PYTHON SOURCE LINES 56-62 Short windows see time, long windows see frequency -------------------------------------------------- A chirp sweeping from 50 to 400 Hz, plus a click at 1.5 s. A 32-sample (32 ms) Gaussian window locates the click but blurs the chirp's frequency; a 512-sample window does the opposite. .. GENERATED FROM PYTHON SOURCE LINES 62-80 .. code-block:: Python t = np.arange(int(2 * fs)) / fs x = signal.chirp(t, 50, t[-1], 400) x[int(1.5 * fs)] += 20.0 fig, axes = plt.subplots(1, 2, figsize=(10, 3.8), sharey=True) for ax, length in zip(axes, (32, 512)): result = stft(x, fs, length, hop=8, window=("gaussian", length / 6)) result.plot(ax=ax) info = result.summary() ax.set_title(f"{length}-sample window: $\\Delta f$ = {info['frequency_resolution']:.0f} Hz") # The ridge of the chirp follows its instantaneous frequency. middle = np.argmin(abs(result.time - 1.0)) ridge = result.frequency[np.argmax(abs(result.values[:, middle]))] print(f"{length:3d}-sample window: ridge at t = 1 s is {ridge:5.1f} Hz (chirp: 225 Hz), {info['frames']} frames") assert abs(ridge - 225) <= info["frequency_resolution"] fig.tight_layout() plt.show() .. image-sg:: /api/gallery/signal_processing/spectral/images/sphx_glr_plot_03_gabor_stft_001.png :alt: 32-sample window: $\Delta f$ = 31 Hz, 512-sample window: $\Delta f$ = 2 Hz :srcset: /api/gallery/signal_processing/spectral/images/sphx_glr_plot_03_gabor_stft_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none 32-sample window: ridge at t = 1 s is 218.8 Hz (chirp: 225 Hz), 247 frames 512-sample window: ridge at t = 1 s is 224.6 Hz (chirp: 225 Hz), 187 frames .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.072 seconds) .. _sphx_glr_download_api_gallery_signal_processing_spectral_plot_03_gabor_stft.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: lite-badge .. image:: images/jupyterlite_badge_logo.svg :target: ../../../../lite/lab/index.html?path=api/gallery/signal_processing/spectral/plot_03_gabor_stft.ipynb :alt: Launch JupyterLite :width: 150 px .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_03_gabor_stft.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_03_gabor_stft.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_03_gabor_stft.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_