.. note::
    :class: sphx-glr-download-link-note

    Click :ref:`here <sphx_glr_download_intro_scipy_auto_examples_plot_fftpack.py>` to download the full example code
.. rst-class:: sphx-glr-example-title

.. _sphx_glr_intro_scipy_auto_examples_plot_fftpack.py:


=============================================
Plotting and manipulating FFTs for filtering
=============================================

Plot the power of the FFT of a signal and inverse FFT back to reconstruct
a signal.

This example demonstrate :func:`scipy.fftpack.fft`,
:func:`scipy.fftpack.fftfreq` and :func:`scipy.fftpack.ifft`. It
implements a basic filter that is very suboptimal, and should not be
used.




.. code-block:: python


    import numpy as np
    from scipy import fftpack
    from matplotlib import pyplot as plt







Generate the signal
###########################################################



.. code-block:: python


    # Seed the random number generator
    np.random.seed(1234)

    time_step = 0.02
    period = 5.

    time_vec = np.arange(0, 20, time_step)
    sig = (np.sin(2 * np.pi / period * time_vec)
           + 0.5 * np.random.randn(time_vec.size))

    plt.figure(figsize=(6, 5))
    plt.plot(time_vec, sig, label='Original signal')




.. image:: /intro/scipy/auto_examples/images/sphx_glr_plot_fftpack_001.png
    :class: sphx-glr-single-img




Compute and plot the power
###########################################################



.. code-block:: python


    # The FFT of the signal
    sig_fft = fftpack.fft(sig)

    # And the power (sig_fft is of complex dtype)
    power = np.abs(sig_fft)**2

    # The corresponding frequencies
    sample_freq = fftpack.fftfreq(sig.size, d=time_step)

    # Plot the FFT power
    plt.figure(figsize=(6, 5))
    plt.plot(sample_freq, power)
    plt.xlabel('Frequency [Hz]')
    plt.ylabel('plower')

    # Find the peak frequency: we can focus on only the positive frequencies
    pos_mask = np.where(sample_freq > 0)
    freqs = sample_freq[pos_mask]
    peak_freq = freqs[power[pos_mask].argmax()]

    # Check that it does indeed correspond to the frequency that we generate
    # the signal with
    np.allclose(peak_freq, 1./period)

    # An inner plot to show the peak frequency
    axes = plt.axes([0.55, 0.3, 0.3, 0.5])
    plt.title('Peak frequency')
    plt.plot(freqs[:8], power[pos_mask][:8])
    plt.setp(axes, yticks=[])

    # scipy.signal.find_peaks_cwt can also be used for more advanced
    # peak detection




.. image:: /intro/scipy/auto_examples/images/sphx_glr_plot_fftpack_002.png
    :class: sphx-glr-single-img




Remove all the high frequencies
###########################################################

 We now remove all the high frequencies and transform back from
 frequencies to signal.



.. code-block:: python


    high_freq_fft = sig_fft.copy()
    high_freq_fft[np.abs(sample_freq) > peak_freq] = 0
    filtered_sig = fftpack.ifft(high_freq_fft)

    plt.figure(figsize=(6, 5))
    plt.plot(time_vec, sig, label='Original signal')
    plt.plot(time_vec, filtered_sig, linewidth=3, label='Filtered signal')
    plt.xlabel('Time [s]')
    plt.ylabel('Amplitude')

    plt.legend(loc='best')




.. image:: /intro/scipy/auto_examples/images/sphx_glr_plot_fftpack_003.png
    :class: sphx-glr-single-img




**Note** This is actually a bad way of creating a filter: such brutal
cut-off in frequency space does not control distorsion on the signal.

Filters should be created using the scipy filter design code



.. code-block:: python

    plt.show()








**Total running time of the script:** ( 0 minutes  0.110 seconds)


.. _sphx_glr_download_intro_scipy_auto_examples_plot_fftpack.py:


.. only :: html

 .. container:: sphx-glr-footer
    :class: sphx-glr-footer-example



  .. container:: sphx-glr-download

     :download:`Download Python source code: plot_fftpack.py <plot_fftpack.py>`



  .. container:: sphx-glr-download

     :download:`Download Jupyter notebook: plot_fftpack.ipynb <plot_fftpack.ipynb>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.readthedocs.io>`_