Question

I'm trying to smooth and interpolate some periodic data in python using scipy.fftp. I have managed to take the fft of the data, remove the higher order frequencies above wn (by doing myfft[wn:-wn] = 0) and then reconstruct a "smoothed" version of the data with ifft(myfft). The array created by the ifft has the same number of points as the original data. How can I use that fft to create an array with more points.

x = [i*2*np.pi/360 for i in range(0,360,30)]
data = np.sin(x)
#get fft
myfft = fftp.fft(data)
#kill feqs above wn
myfft[wn:-wn] = 0
#make new series
newdata = fftp.ifft(myfft)

I've also been able to manually recreate the series at the same resolution as demonstrated here Recreating time series data using FFT results without using ifft

but when I tried upping the resolution of the x-values array it didn't give me the right answer either.

Thanks in advance

Niall

Was it helpful?

Solution

What np.fft.fft returns has the DC component at position 0, followed by all positive frequencies, then the Nyquist frequency (only if the number of elements is even), then the negative frequencies in reverse order. So to add more resolution you could add zeros at both sides of the Nyquist frequency:

import numpy as np
import matplotlib.pyplot as plt

y = np.sin(np.linspace(0, 2*np.pi, 32, endpoint=False))

f = np.fft.fft(y)
n = len(f)
f_ = np.concatenate((f[0:(n+1)//2],
                     np.zeros(n//2),
                     [] if n%2 != 0 else f[(n+1)//2:(n+3)//2],
                     np.zeros(n//2),
                     f[(n+3)//2:]))
y_ = np.fft.ifft(f_)
plt.plot(y, 'ro')
plt.plot(y_, 'bo')
plt.show()

enter image description here

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top