Python 看似错误频率的频率过滤

Python Frequency filtering with seemingly wrong frequencies

下面的脚本通过切割所有大于 6 的频率来过滤频率。 然而,不是使用看似正确的函数 rfftfreq,而是使用 fftfreq

据我了解rfftfreq应该和rfft一起使用。为什么此代码可以使用 fftfreqrfft

import numpy as np
from scipy.fftpack import rfft, irfft, fftfreq

time   = np.linspace(0,10,2000)
signal = np.cos(5*np.pi*time) + np.cos(7*np.pi*time)

W = fftfreq(signal.size, d=time[1]-time[0])
f_signal = rfft(signal)

# If our original signal time was in seconds, this is now in Hz    
cut_f_signal = f_signal.copy()
cut_f_signal[(W<6)] = 0

cut_signal = irfft(cut_f_signal)

背景:rfft 给出了一个数组,在单独的条目中对具有真实和想象的傅立叶模式进行排序。如[R(0), R(1), I(1), ... R(N/2),I(N/2)]R(n)I(n)分别是傅立叶模的实部和虚部。 (假设一个偶数条目数组)

因此,rfftfreq 产生与该数组对应的频率数组,例如(假设输入数组为偶数,采样间距为 1):

[0, 1/n, 1/n, 2/n, 2/n ... n/(2n), n/(2n)]

但是此代码适用于 fftfreq,其中函数的输出为

[0, 1/n, ... n/(2n), -n/(2n), ..., -1/n]

显然,fftfreqrfft 一起使用时会导致错误的结果,因为频率和 FFT 区间不匹配。

你 mis-specified 原始信号中的频率。

一个正弦波根据这个等式(来自Wikipedia)进行参数化:

signal = np.cos(5*np.pi*time) + np.cos(7*np.pi*time) 的定义中缺少因子 2。因此,实际频率为

5*pi*t = 2*pi*t * f
f = (5*pi*t) / (2*pi*t) = 5 / 2

7*pi*t = 2*pi*t * f
f = (7*pi*t) / (2*pi*t) = 7 / 2

换句话说,这两个频率是你认为的一半。具有讽刺意味的是,这就是为什么它似乎适用于 fftfreq 而不是 rfftfreq。前者覆盖了两倍的频率范围(正频率和负频率),因此补偿了缺失的因子 2。

这是正确的代码:

signal = np.cos(5 * 2*np.pi * time) + np.cos(7 * 2*np.pi * time)
W = rfftfreq(signal.size, d=time[1]-time[0])