傅里叶或时域积分

Question

我在努力理解一个信号的数值积分问题。基本上，我有一个信号，我想集成或执行，并反消去，作为时间的函数（集成拾音器线圈获得磁场）。我试过两种不同的方法，原则上应该是一致的，但它们不是。我使用的代码如下。注意，代码中的信号y之前已经使用巴特沃斯滤波（类似于这里的http://wiki.scipy.org/Cookbook/ButterworthBandpass）进行了高通滤波。信号和时间基准可以在这里下载（https://www.dropbox.com/s/fi5z38sae6j5410/trial.npz?dl=0）

import scipy as sp
from scipy import integrate
from scipy import fftpack
data = np.load('trial.npz')
y = data['arr_1'] # this is the signal 
t = data['arr_0']
# integration using pfft
bI = sp.fftpack.diff(y-y.mean(),order=-1)
bI2= sp.integrate.cumtrapz(y-y.mean(),x=t)

现在这两个信号（除了可以得到的最终不同的线性趋势）是不同的，或者更好地说，它们在相同的振荡时间下非常相似，但是在这两个信号之间有一个大约30倍的因子，在这个意义上bI2比bI低30倍（大约）。顺便说一句，我已经减去了这两个信号的平均值，以确保它们是零均值信号，并在IDL中进行积分（在等效的cumsumtrapz和fourier域中）得到了与bI2兼容的值。任何线索都很受欢迎

Answer 1

很难知道scipy.fftpack.diff()在引擎盖下做什么。在

为了解决你的问题，我挖掘了一个旧的频域积分函数，这是我不久前写的。值得指出的是，在实践中，人们通常希望对某些参数的控制比scipy.fftpack.diff()给你的更多一些控制。例如，myintf()函数的f_lo和f_hi参数允许您对输入进行频带限制，以排除可能有噪声的非常低或非常高的频率。特别是在集成过程中，噪声低频会“爆炸”，并压倒信号。您可能还需要在时间序列的开始和结束处使用一个窗口来阻止频谱泄漏。在

我计算了bI2和一个结果，bI3，使用以下代码与intf()集成一次（为了简单起见，我假设了一个平均采样率）：

import intf
from scipy import integrate
data = np.load(path)
y = data['arr_1']
t = data['arr_0']

bI2= sp.integrate.cumtrapz(y-y.mean(),x=t)
bI3 = intf.intf(y-y.mean(), fs=500458, f_lo=1, winlen=1e-2, times=1)

我画了bI2和bI3：

尽管bI2存在明显的分段线性趋势，但这两个时间序列具有相同的数量级和大致相同的形状。我知道这并不能解释scipy函数中发生了什么，但至少这表明它不是频域方法的问题。在

intf的代码完整粘贴在下面。在

^{pr2}$