tmporaries
tmporaries

Reputation: 1553

Phase spectrum with python and FFT

I'm trying to calculate a phase spectrum of sinusoid. The following code generates 1Hz sinusoid with zero initial phase.

import numpy
from numpy import pi, sin, arange
from pylab import plot, show, xlabel, ylabel, xlim, grid

sampling_rate = 500
sampling_time = 1 / sampling_rate
length = 1 # in seconds
n = sampling_rate * length # number of points

time = arange(0, n * sampling_time, sampling_time)
# Generate sinusoid: frequency=1Hz, phase=0
signal = sin(2 * pi * time)

fft = numpy.fft.fft(signal)
fft_phase = numpy.angle(fft)
fft_freq = numpy.arange(n) * sampling_rate / n

plot(fft_freq, fft_phase)
ylabel("FFT Angle")
xlabel("Frequency (Hz)")
xlim(left=0, right=5)
grid(True)
show()

But result doesn't match my expectations. It has non-zero phase of 1 Hz component:

enter image description here

It shows incorrect phase of 1 Hz harmonic. What's wrong with the code (or approach)?

Upvotes: 1

Views: 4777

Answers (1)

Cris Luengo
Cris Luengo

Reputation: 60761

When the magnitude is zero, then the phase is given by numerical imprecision.

If you display the values computed by fft you’ll see that the values you expect to be 0 are actually in the order of 1e-16 or something like that. This is numerical imprecision caused by rounding in the floating-point computations.

The solution is to compute both the magnitude and the phase, and ignore the phase component if the magnitude is too small.

Upvotes: 1

Related Questions