numpy.fft.fft是如何工作的?

7
我现在正在试图理解numpy中的fft函数。为此,我测试了以下假设:
我有两个函数:f(x)= x ^ 2g(x)= f'(x)= 2 * x。根据傅里叶变换定律和wolfram alpha,应该是G(w)= 2pi * i * F(w)(前置因子可能会有所变化,但只应有一个常数因子)。在python中实现时,我写入:
import numpy as np
def x2(x):
    return x*x
def nx(x):
    return 2*x

a = np.linspace(-3, 3, 16)
a1 = x2(a)
a2 = nx(a)

b1 = np.fft.fft(a1)
b2 = np.fft.fft(a2)

c = b1/b2

现在我期望得到一个几乎恒定的值c,但实际上却得到了:

array([  1.02081592e+16+0.j        ,   1.32769987e-16-1.0054679j ,
         4.90653893e-17-0.48284271j,  -1.28214041e-16-0.29932115j,
        -1.21430643e-16-0.2j       ,   5.63664751e-16-0.13363573j,
        -5.92271642e-17-0.08284271j,  -4.21346622e-16-0.03978247j,
        -5.55111512e-16-0.j        ,  -5.04781597e-16+0.03978247j,
        -6.29288619e-17+0.08284271j,   8.39500693e-16+0.13363573j,
        -1.21430643e-16+0.2j       ,  -0.00000000e+00+0.29932115j,
        -0.00000000e+00+0.48284271j,   1.32769987e-16+1.0054679j ])

我的错误在哪里,我该如何使用fft?


2
这里错误的部分是 G(w) = 2 pi i F(w)。应该是 G(w) = 2 pi i w F(w) - Dietrich Epp
@DietrichEpp:这些数组中的w是什么,它也是虚数吗? - arc_lupus
@arc_lupus w是对应于x的频率空间等价物。 在这种情况下,您在原始空间中有x ^ 2,它本身没有特定的频率(至少不是以易于理解的方式),因此您将得到多个不同的峰来表示x ^ 2。例如,如果f(x)= sin(x),则您将在1处获得一个漂亮简单的delta函数,以捕获sin(x)中的简单频率。 - Mike Williamson
1个回答

6
你所提供的属性适用于连续傅里叶变换(CFT)。FFT计算的是离散傅里叶变换(DFT),它与CFT有关但并不完全等同。
确实,在某些条件下,DFT与CFT成比例:即对于足够采样的函数,其在采样限制之外为零(例如,请参见本书的附录E)。
上述函数均不满足这两个条件,因此DFT与CFT不成比例,您的数值结果反映了这一点。
这里有一些代码,使用适当采样的带限函数通过FFT确认您感兴趣的关系:
import numpy as np

def f(x):
    return np.exp(-x ** 2)
def fprime(x):
    return -2 * x * f(x)

a = np.linspace(-10, 10, 100)
a1 = f(a)
a2 = fprime(a)

b1 = np.fft.fft(a1)
b2 = np.fft.fft(a2)
omega = 2 * np.pi * np.fft.fftfreq(len(a), a[1] - a[0])

np.allclose(b1 * 1j * omega, b2)
# True

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接