How to get the FFT of a numpy array to work?

Jessica Flores picture Jessica Flores · Sep 26, 2015 · Viewed 11.9k times · Source

I'm reading a specific column of a csv file as a numpy array. When I try to do the fft of this array I get an array of NaNs. How do I get the fft to work? Here's what I have so far:

#!/usr/bin/env python
from __future__ import division

import numpy as np
from numpy import fft

import matplotlib.pyplot as plt


fileName = '/Users/Name/Documents/file.csv'

#read csv file
df = np.genfromtxt(fileName, dtype = float, delimiter = ',', names = True)

X = df['X'] #get X from file

rate = 1000. #rate of data collection in points per second
Hx = abs(fft.fft(X))
freqX = fft.fftfreq(len(Hx), 1/rate)

plt.plot(freqX,Hx) #plot freqX vs Hx

Answer

Joe Kington picture Joe Kington · Sep 26, 2015

Presumably there are some missing values in your csv file. By default, np.genfromtxt will replace the missing values with NaN.

If there are any NaNs or Infs in an array, the fft will be all NaNs or Infs.

For example:

import numpy as np

x = [0.1, 0.2, np.nan, 0.4, 0.5]
print np.fft.fft(x)

And we'll get:

array([ nan +0.j,  nan+nanj,  nan+nanj,  nan+nanj,  nan+nanj])

However, because an FFT operates on a regularly-spaced series of values, removing the non-finite values from an array is a bit more complex than just dropping them.

pandas has several specialized operations to do this, if you're open to using it (e.g. fillna). However, it's not too difficult to do with "pure" numpy.

First, I'm going to assume that you're working with a continuous series of data because you're taking the FFT of the values. In that case, we'd want to interpolate the NaN values based on the values around them. Linear interpolation (np.interp) may not be ideal in all situations, but it's not a bad default choice:

For example:

import numpy as np

x = np.array([0.1, 0.2, np.nan, 0.4, 0.5])
xi = np.arange(len(x))

mask = np.isfinite(x)
xfiltered = np.interp(xi, xi[mask], x[mask])

And we'll get:

In [18]: xfiltered
Out[18]: array([ 0.1,  0.2,  0.3,  0.4,  0.5])

We can then calculate the FFT normally:

In [19]: np.fft.fft(xfiltered)
Out[19]: 
array([ 1.50+0.j        , -0.25+0.34409548j, -0.25+0.08122992j,
       -0.25-0.08122992j, -0.25-0.34409548j])

...and get a valid result.