Python中的Parseval定理

Question

我试图抓住Python的fft功能,而我偶然发现的一个奇怪的事情就是Parseval的定理似乎不适用,因为它现在给出了大约50的差异,而它应该是0.

import numpy as np
import matplotlib.pyplot as plt
import scipy.fftpack as fftpack

pi = np.pi

tdata = np.arange(5999.)/300
dt = tdata[1]-tdata[0]

datay = np.sin(pi*tdata)+2*np.sin(pi*2*tdata)
N = len(datay)

fouriery = abs(fftpack.rfft(datay))/N

freqs = fftpack.rfftfreq(len(datay), d=(tdata[1]-tdata[0]))

df = freqs[1] - freqs[0]

parceval = sum(datay**2)*dt - sum(fouriery**2)*df
print parceval

plt.plot(freqs, fouriery, 'b-')
plt.xlim(0,3)
plt.show()

我很确定这是一个规范化因素,但我似乎无法找到它,因为我能找到的关于这个函数的所有信息都是scipy.fftpack.rfft文档.

Answer 1

您的归一化因子来自于尝试将Parseval定理应用于连续信号的傅里叶变换到离散序列.在关于离散傅立叶变换的维基百科文章的侧面板上,讨论了傅里叶变换,傅里叶级数,离散傅里叶变换和狄拉克梳子采样之间的关系.

长话短说,Parseval定理在应用于DFT时,不需要积分,而是求和:a 2*pi通过乘以dt和df你的求和来创建.

还要注意,因为你正在使用scipy.fftpack.rfft,你得到的不是数据的DFT,而只是正数的一半,因为负数将与它对称.如此以来,你只添加一半的数据,再加上0在DC来看,有失踪2去的4*pi发现@unutbu.

在任何情况下,如果datay保持你的序列,你可以验证Parseval定理如下:

fouriery = fftpack.rfft(datay)
N = len(datay)
parseval_1 = np.sum(datay**2)
parseval_2 = (fouriery[0]**2 + 2 * np.sum(fouriery[1:]**2)) / N
print parseval_1 - parseval_2

使用scipy.fftpack.fft或numpy.fft.fft第二次求和不需要采用这种奇怪的形式:

fouriery_1 = fftpack.fft(datay)
fouriery_2 = np.fft.fft(datay)
N = len(datay)
parseval_1 = np.sum(datay**2)
parseval_2_1 = np.sum(np.abs(fouriery_1)**2) / N
parseval_2_2 = np.sum(np.abs(fouriery_2)**2) / N
print parseval_1 - parseval_2_1
print parseval_1 - parseval_2_2