执行 FFT 后如何将复数转换回“正常”数

Posted 2023-02-25

【中文标题】执行 FFT 后如何将复数转换回“正常”数

【英文标题】：How to convert complex numbers back into "normal" numbers after performing FFT

【发布时间】：2015-04-09 16:56:41

【问题描述】：

说，我有一个示例信号由三个余弦组成，每个余弦代表 4、6 和 8 个频带。现在，我使用 FFT 将此信号放入频域，并在频域中切断不需要的 6 Hz 频带。最后，我想将频域的信号反相回到时域。但是当我简单地使用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]
"""

【参考方案1】：

如果你想去掉 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 的第二个参数包含np.where 在条件为True 时返回的值。

通过这一更改，您的代码生成

[ 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 的实现混在一起。

【参考方案2】：

如果您想要一个严格真实的结果（减去舍入误差噪声），您的 IFFT 输入需要是赫米斯对称的（例如，您需要确保复数数组的后半部分是前半部分的复共轭镜像)。查看真实数据的初始 FFT，您会看到对称性。

但看起来您没有过滤负频率，因此向 IFFT 发送了一个不对称输入，然后输出一个复杂的结果。