numpy.fft.fft() implementation in Python

Question

I am trying to obtain the spectrum of periodic signals using the fft function. Then plot the magnitude and phase of the transform. Magnitude plots are OK, but phase plots are completely unexpected. For example I used the functions Sin³(t) and Cos³(t). The code I have used is:

import matplotlib.pyplot as plt
import numpy.fft as nf
import math
import numpy as np
pi = math.pi

N=512

# Sin³(x)
t=np.linspace(-4*pi,4*pi,N+1);t=t[:-1]
y=(np.sin(t))**3
Y=nf.fftshift(nf.fft(y))/N
w=np.linspace(-64,64,N+1);w=w[:-1]

plt.figure("0")
plt.subplot(2,1,1)
plt.plot(w,abs(Y),'ro',lw=2)
plt.xlim((-4,4))
plt.ylabel(r"$|Y|$",size=16)
plt.title("Spectrum of sin\u00B3(t)")
plt.grid(True)
plt.subplot(2,1,2)
ii=np.where(abs(Y)>1e-3)
plt.plot(w[ii],np.angle(Y[ii]),'go',lw=2)
plt.xlim((-4,4))
plt.ylabel(r"Phase of $Y$",size=16)
plt.xlabel(r"$\omega$",size=16)
plt.grid(True)

# Cos³(x)

t=np.linspace(-4*pi,4*pi,N+1);t=t[:-1]
y=(np.cos(t))**3
Y=nf.fftshift(nf.fft(y))/N
w=np.linspace(-64,64,N+1);w=w[:-1]

plt.figure("1")
plt.subplot(2,1,1)
plt.plot(w,abs(Y),'ro',lw=2)
plt.xlim((-4,4))
plt.ylabel(r"$|Y|$",size=16)
plt.title("Spectrum of cos\u00B3(t)")
plt.grid(True)
plt.subplot(2,1,2)
ii=np.where(abs(Y)>1e-3)
plt.plot(w[ii],(np.angle(Y[ii])),'go',lw=2)
plt.xlim((-4,4))
plt.ylabel(r"Phase of $Y$",size=16)
plt.xlabel(r"$\omega$",size=16)
plt.grid(True)

The plots obtained are:

i) For Sin³(t) - Magnitude and Phase of Spectrum of Sin³(t)

ii) For Cos³(t) - Magnitude and Phase of Spectrum of Cos³(t)

As You can see in the links above, The magnitude of Both functions are fine. The phase of Spectrum of Sin³(t) is correct as expected.

Since Cos³(t) is real and even, and expanding in terms of complex exponentials reveal that the Phase of the Coefficients is 0. But the plotted graph shows entirely different answer (see 2): at w=-3 the phase is near 5 radians. What is the mistake I have done and what is the correct way to implement the FFT.

jwalton · Accepted Answer · 2020-04-19 17:34:56Z

3

np.fft is behaving as expected; it's your plot which is causing the confusion. A call to plt.tight_layout() should help clear things up for you:

If you take a careful look at the y-axis of the phase for cos³(t) you'll see that the values all have a 1e-17 prefactor. So at w = -3 the phase is not near 5 radians, it's actually near 5e-17 radians (which is for all intents and purposes zero).

I tidied up your code whilst I was looking into this:

import matplotlib.pyplot as plt
import numpy.fft as nf
import numpy as np

plt.rcParams['axes.grid'] = True

pi = np.pi
N = 512

t = np.linspace(-4*pi, 4*pi, N+1)[:-1]


def fft_func(x, func, N):
    y = func(x) ** 3
    Y = nf.fftshift(nf.fft(y)) / N
    w = np.linspace(-64, 64, N+1)[:-1]
    return y, Y, w


def plot_mag_phase(w, Y, func_name):
    fig, ax = plt.subplots(2, 1)

    ii = np.where(abs(Y) > 1e-3)

    ax[0].plot(w, abs(Y), 'ro', lw=2)
    ax[1].plot(w[ii], np.angle(Y[ii]), 'go', lw=2)

    ax[0].set_xlim(-4, 4)
    ax[1].set_xlim(-4, 4)

    ax[0].set_ylabel("$|Y|$", size=16)
    ax[1].set_ylabel("Phase of $Y$", size=16)

    ax[0].set_title("Spectrum of {}\u00B3(t)".format(func_name))
    ax[1].set_xlabel(r'$\omega$', size=16)

    fig.tight_layout()

# sin x ^ 3
y, Y, w = fft_func(t, np.sin, N)
plot_mag_phase(w, Y, 'sin')

# cos x ^ 3
y, Y, w = fft_func(t, np.cos, N)
plot_mag_phase(w, Y, 'cos')

edited Apr 19, 2020 at 17:34

answered Apr 19, 2020 at 17:08

jwalton

5,7561 gold badge22 silver badges45 bronze badges

Sign up to request clarification or add additional context in comments.

4 Comments

Grentouce Over a year ago

This is my first post in overflow , your answer gave me additional insights in plotting,and formatting, thanks a lot !

jwalton Over a year ago

You're welcome, I'm glad this was useful :) When you're coding you should try to avoid repeating the same / similar bits of code over and over. Instead, try to capture the logic you're using and implement it as some function or method.

Grentouce Over a year ago

Idk y but what is the relation between w and t , i.e. why w ranges from -64 to 64 ? Thanks.

jwalton Over a year ago

I'm not sure why the range of w is -64 to 64, I only copied it from your code. Given that you set the x-limits of the plot to -4 -> 4, you could drop the w altogether and use plot_mag_phase(t, Y, 'sin') for the same effect.

Collectives™ on Stack Overflow

numpy.fft.fft() implementation in Python

1 Answer 1

4 Comments

Your Answer

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

4 Comments

Your Answer

Sign up or log in

Post as a guest

Related