在不使用ifft的情况下使用FFT结果重新创建时间序列数据

推荐答案调用fft(wolfer)时，您告诉转换假定一个等于数据长度的基本周期.要重建数据，您必须使用相同基本周期= 2*pi/N的基函数.同样，您的时间索引xs必须在原始信号的时间采样范围内.When you called fft(wolfer), you told the transform to assume a fundamental period equal to the length of the data. To reconstruct the data, you have to use basis functions of the same fundamental period = 2*pi/N. By the same token, your time index xs has to range over the time samples of the original signal.另一个错误是忘记进行完整的复数乘法.更容易将其视为Y[omega]*exp(1j*n*omega/N).Another mistake was in forgetting to do to the full complex multiplication. It's easier to think of this as Y[omega]*exp(1j*n*omega/N).这是固定代码.请注意，为了避免与sqrt(-1)混淆，我将i重命名为ctr，并且将n重命名为N以遵循通常的信号处理约定，即对样本使用小写字母，对总样本长度使用大写字母. .我还导入了__future__ division以避免混淆整数除法.Here's the fixed code. Note I renamed i to ctr to avoid confusion with sqrt(-1), and n to N to follow the usual signal processing convention of using the lower case for a sample, and the upper case for total sample length. I also imported __future__ division to avoid confusion about integer division. 忘了早点添加:请注意，SciPy的fft在累加后不会除以N.在使用Y[n]之前，我没有对此进行划分；如果您想找回相同的数字，而不只是看到相同的形状，就应该这样做.forgot to add earlier: Note that SciPy's fft doesn't divide by N after accumulating. I didn't divide this out before using Y[n]; you should if you want to get back the same numbers, rather than just seeing the same shape.最后，请注意，我正在对整个频率系数范围进行求和.当我绘制np.abs(Y)时，它看起来在较高频率中有明显的值，至少直到样本70左右.我认为，通过对整个范围求和，查看正确的结果，然后回算系数并查看会发生什么，将更容易理解结果.And finally, note that I am summing over the full range of frequency coefficients. When I plotted np.abs(Y), it looked like there were significant values in the upper frequencies, at least until sample 70 or so. I figured it would be easier to understand the result by summing over the full range, seeing the correct result, then paring back coefficients and seeing what happens.from __future__ import divisionimport numpy as npfrom scipy import *from matplotlib import pyplot as gpltfrom scipy import fftpackdef f(Y,x, N): total = 0 for ctr in range(len(Y)): total += Y[ctr] * (np.cos(x*ctr*2*np.pi/N) + 1j*np.sin(x*ctr*2*np.pi/N)) return real(total)tempdata = np.loadtxt("sunspots.dat")year=tempdata[:,0]wolfer=tempdata[:,1]Y=fft(wolfer)N=len(Y)print(N)xs = range(N)gplt.plot(xs, [f(Y, x, N) for x in xs])gplt.show() 这篇关于在不使用ifft的情况下使用FFT结果重新创建时间序列数据的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持！ 上岸，阿里云！