离散傅里叶变换DFT
发布时间
阅读量:
阅读量
离散傅里叶变换DFT
DFT和IDFT
时域中长度为N的序列x[n]的离散傅里叶变换(DFT)和逆变换(IDFT)
为了方便记忆,常用一个中间变量
这样DFT和IDFT变成下面容易记忆的形式:
进一步,正反变换都可以看成是两个矩阵的乘法
观察到对 W 矩阵进行求逆操作转换为 W^{-1} 的过程只需将每个元素取倒数,并除以 即可完成。无需直接计算其逆矩阵,在理论上这种转换过程非常高效。
numpy实现
代码参考自这里。
波形显示函数
def show(ori_func, ft, sampling_period = 5):
n = len(ori_func)
interval = sampling_period / float(n)
print interval
# 绘制原始函数
plt.subplot(2, 1, 1)
plt.plot(np.arange(0, sampling_period, interval), ori_func, 'black')
plt.xlabel('Time'), plt.ylabel('Amplitude')
# 绘制变换后的函数
plt.subplot(2,1,2)
frequency = np.arange(n / 2) / (n * interval)
nfft = abs(ft[range(int(n / 2))] / n )
plt.plot(frequency, nfft, 'red')
plt.xlabel('Freq (Hz)'), plt.ylabel('Amp. Spectrum')
plt.show()
多个谐波叠加波形
time = np.arange(0, 5, 0.005)
x1 = np.sin(2 * np.pi * 1 * time) # 频率为1
x2 = np.sin(2 * np.pi * 20 * time) # 频率为20
x3 = np.sin(2 * np.pi * 60 * time) # 频率为60
x = x1 + x2 + x3 # 叠加
y = np.fft.fft(x)
show(x, y)

验证numpy.fft.fft的计算原理
# 傅里叶变换DFT
x = np.random.random(500)
N = len(x)
n = np.arange(N)
k = n.reshape((N, 1))
W = np.exp(-2j * np.pi * k * n / N) # 500*500 矩阵
y = np.dot(W, x)
print np.allclose(y, np.fft.fft(x)) # 应该是 True
验证numpy.fft.ifft的计算原理
# 傅里叶逆变换IDFT
W_inv = np.exp(2j * np.pi * k * n / N) / N
x = np.dot(y, W_inv)
print np.allclose(x, np.fft.ifft(y)) # 应该是 True
参考资料
博客
博客
API
书籍
全部评论 (0)
还没有任何评论哟~
