比如说,我有一个例子,信号由三个余弦组成,每个余弦代表4,6和8频带。现在,我用FFT把这个信号扔到频域,在频域上,我切断了6赫兹波段。最后,将信号从频域反演回时域。但是,当我简单地使用numpy.fft.ifft时,我得到了复数数组,这不是进一步分析信号的最佳结果。在进行带通后,我怎样才能反演FFT,从而得到实部和虚部作为一个数字携带的全部信息?我调查了z = sqrt(real^2 + imaginary^2)的事情,但这不是“事情”。
下面我提供一个工作的例子。我会感谢你的帮助。
import numpy as np
from scipy.fftpack import fftfreq
# Define signal.
Fs = 128 # Sampling rate.
Ts = 1 / Fs # Sampling interval.
Time = np.arange(0, 10, Ts) # Time vector.
signal = np.cos(4*np.pi*Time) + np.cos(6*np.pi*Time) + np.cos(8*np.pi*Time)
def spectrum(sig, t):
"""
Represent given signal in frequency domain.
:param sig: signal.
:param t: time scale.
:return:
"""
f = fftfreq(sig.size, d=t[1]-t[0])
y = np.fft.fft(sig)
return f, np.abs(y)
def bandpass(f, sig, min_freq, max_freq):
"""
Bandpass signal in a specified by min_freq and max_freq frequency range.
:param f: frequency.
:param sig: signal.
:param min_freq: minimum frequency.
:param max_freq: maximum frequency.
:return:
"""
return np.where(np.logical_or(f < min_freq, f > max_freq), 0, sig)
freq, spec = spectrum(signal, Time)
signal_filtered = np.fft.ifft(bandpass(freq, spec, 5, 7))
print(signal_filtered)
"""
print(signal_filtered) result:
[ 2.22833798e-15 +0.00000000e+00j 2.13212081e-15 +6.44480810e-16j
1.85209996e-15 +1.23225456e-15j ..., 1.41336488e-15 -1.71179288e-15j
1.85209996e-15 -1.23225456e-15j 2.13212081e-15 -6.44480810e-16j]
"""发布于 2015-04-09 17:34:34
如果你想切断5到7之间的频率,那么你会想把频率保持在
(f < min_freq) | (f > max_freq)这相当于
np.logical_or(f < min_freq, f > max_freq)因此,使用
return np.where(np.logical_or(f < min_freq, f > max_freq), sig, 0)而不是
return np.where(np.logical_or(f < min_freq, f > max_freq), 0, sig)因为np.where的第二个参数包含条件为True时np.where返回的值。
通过这次更改,您的代码将生成
[ 3.00000000 +0.00000000e+00j 2.96514652 +1.24442385e-15j
2.86160515 +2.08976636e-15j ..., 2.69239924 +4.71763845e-15j
2.86160515 +5.88163496e-15j 2.96514652 +6.82134642e-15j]请注意,如果您的信号是真实的,您可以使用rfft对实序列进行离散傅里叶变换,使用irfft进行逆,使用rfftfreq生成频率。
例如,
from __future__ import division
import numpy as np
import scipy.fftpack as fftpack
# Define signal.
Fs = 128 # Sampling rate.
Ts = 1 / Fs # Sampling interval.
Time = np.arange(0, 10, Ts) # Time vector.
signal = np.cos(4*np.pi*Time) + np.cos(6*np.pi*Time) + np.cos(8*np.pi*Time)
def spectrum(sig, t):
"""
Represent given signal in frequency domain.
:param sig: signal.
:param t: time scale.
:return:
"""
f = fftpack.rfftfreq(sig.size, d=t[1]-t[0])
y = fftpack.rfft(sig)
return f, np.abs(y)
def bandpass(f, sig, min_freq, max_freq):
"""
Bandpass signal in a specified by min_freq and max_freq frequency range.
:param f: frequency.
:param sig: signal.
:param min_freq: minimum frequency.
:param max_freq: maximum frequency.
:return:
"""
return np.where(np.logical_or(f < min_freq, f > max_freq), sig, 0)
freq, spec = spectrum(signal, Time)
signal_filtered = fftpack.irfft(bandpass(freq, spec, 5, 7))
print(signal_filtered)收益率
[ 3. 2.96514652 2.86160515 ..., 2.69239924 2.86160515
2.96514652]注意,这里必须使用scipy的fftpack;不要将SciPy的实现与NumPy的实现混为一谈。
发布于 2015-04-10 12:47:04
如果您想要一个严格的实际结果(减去四舍五入误差噪声),则您的IFFT的输入需要是hermician对称的(例如,您需要确保复数组的下半部分是前半部分的复共轭镜像)。看看你最初对真实数据的FFT,你就会看到对称性。
但是看起来你没有过滤负频率,因此向IFFT发送了一个不对称的输入,然后输出一个复杂的结果。
https://stackoverflow.com/questions/29544563
复制相似问题