python – 使用scipy fft计算信号的自相关性可以得到直接计算的
发布时间:2020-12-20 12:04:54 所属栏目:Python 来源:网络整理
导读:我试图使用自相关是功率谱的逆傅里叶变换的属性来计算信号的自相关.但是,当我使用scipy(或numpy)fft执行此操作并与自相关函数的直接计算进行比较时,我得到了错误的答案,具体而言,fft版本在较大的延迟时间内以较小的负值关闭,这是显然错了. 我的MWE低于输出.
我试图使用自相关是功率谱的逆傅里叶变换的属性来计算信号的自相关.但是,当我使用scipy(或numpy)fft执行此操作并与自相关函数的直接计算进行比较时,我得到了错误的答案,具体而言,fft版本在较大的延迟时间内以较小的负值关闭,这是显然错了.
我的MWE低于输出.我使用fft错了吗? import numpy as np import matplotlib.pyplot as pl from scipy.fftpack import fft,ifft def autocorrelation(x) : xp = (x - np.average(x))/np.std(x) f = fft(xp) p = np.absolute(f)**2 pi = ifft(p) return np.real(pi)[:len(xp)/2]/(len(xp)) def autocorrelation2(x): maxdelay = len(x)/5 N = len(x) mean = np.average(x) var = np.var(x) xp = (x - mean)/np.sqrt(var) autocorrelation = np.zeros(maxdelay) for r in range(maxdelay): for k in range(N-r): autocorrelation[r] += xp[k]*xp[k+r] autocorrelation[r] /= float(N-r) return autocorrelation def autocorrelation3(x): xp = (x - np.mean(x))/np.std(x) result = np.correlate(xp,xp,mode='full') return result[result.size/2:]/len(xp) def main(): t = np.linspace(0,20,1024) x = np.exp(-t**2) pl.plot(t[:200],autocorrelation(x)[:200],label='scipy fft') pl.plot(t[:200],autocorrelation2(x)[:200],label='direct autocorrelation') pl.plot(t[:200],autocorrelation3(x)[:200],label='numpy correlate') pl.legend() pl.show() if __name__=='__main__': main() 解决方法
离散FT假定信号是周期性的.因此,在基于fft的代码中,您正在计算环绕自相关.为了避免这种情况,你必须做一些0-padding形式:
def autocorrelation(x): xp = ifftshift((x - np.average(x))/np.std(x)) n,= xp.shape xp = np.r_[xp[:n//2],np.zeros_like(xp),xp[n//2:]] f = fft(xp) p = np.absolute(f)**2 pi = ifft(p) return np.real(pi)[:n//2]/(np.arange(n//2)[::-1]+n//2) (编辑:李大同) 【声明】本站内容均来自网络,其相关言论仅代表作者个人观点,不代表本站立场。若无意侵犯到您的权利,请及时与联系站长删除相关内容! |