问Python意外的FFT结果
EN

Stack Overflow用户

提问于 2018-07-19 10:42:01

回答 1查看 2.8K关注 0票数 2

我正在分析什么是呼吸波形，由3种不同的形状构成(数据来源于磁共振成像，因此使用了多次回波时间；如果您想要一些快速的背景信息，请参阅这里 )。以下是为某些上下文绘制的两个波形的几个片段：

对于每个波形，我试图执行一个DFT，以发现主导频率或呼吸频率。

我的问题是，当我绘制我所执行的FFT时，我得到了两件事之一，取决于我使用的FFT库。此外，和都不代表我所期望的。我知道数据并不总是我们想要的样子，但是我的数据中显然有波形，所以我希望离散傅里叶变换在某个合理的地方产生一个频率峰值。供参考，我预计约80至130赫兹。

我的数据存储在一个pandas数据帧中，每个回波时间的数据都存储在一个单独的系列中。我正在应用所选择的FFT函数(请参阅下面的代码)，并为每个函数接收不同的结果。

SciPy (fftpack)

import pandas as pd
import scipy.fftpack

# temporary copy to maintain data structure
lead_pts_fft_df = lead_pts_df.copy()

# apply a discrete fast Fourier transform to each data series in the data frame
lead_pts_fft_df.magnitude = lead_pts_df.magnitude.apply(scipy.fftpack.fft)
lead_pts_fft_df

NumPy：

import pandas as pd
import numpy as np 

# temporary copy to maintain data structure
lead_pts_fft_df = lead_pts_df.copy()

# apply a discrete fast Fourier transform to each data series in the data frame
lead_pts_fft_df.magnitude = lead_pts_df.magnitude.apply(np.fft.fft)
lead_pts_fft_df

其余的相关守则：

ECHO_TIMES = [0.080, 0.200, 0.400] # milliseconds

f_s = 1 / (0.006) # 0.006 = time between samples
freq = np.linspace(0, 29556, 29556) * (f_s / 29556) # (29556 = length of data)

# generate subplots
fig, axes = plt.subplots(3, 1, figsize=(18, 16))

for idx in range(len(ECHO_TIMES)):
    axes[idx].plot(freq, lead_pts_fft_df.magnitude[idx].real, # real part
                   freq, lead_pts_fft_df.magnitude[idx].imag, # imaginary part

    axes[idx].legend() # apply legend to each set of axes

# show the plot
plt.show()

后DFT (SciPy fftpack)：

后DFT (NumPy)

这里是用于创建这些图的数据集(在Dropbox上)的链接，而这里是指向木星笔记本的链接。

编辑：

我使用了张贴的建议，并采取了大小(绝对值)的数据，并绘制了一个对数y轴。新的结果公布如下。看来我的信号里有一些环绕。我没有使用正确的频率刻度吗？更新后的代码和情节如下。

# generate subplots
fig, axes = plt.subplots(3, 1, figsize=(18, 16))

for idx in range(len(ECHO_TIMES)):
    axes[idx].plot(freq[1::], np.log(np.abs(lead_pts_fft_df.magnitude[idx][1::])),
                   label=lead_pts_df.index[idx], # apply labels
                   color=plot_colors[idx]) # colors
    axes[idx].legend() # apply legend to each set of axes

# show the plot
plt.show()

python

numpy

scipy

fft

dft

电商音视频处理解决方案

为电商客户场景提供短视频处理、音频加密处理等场景方案，提高客户效率和用户体验

回答 1

Stack Overflow用户

回答已采纳

发布于 2018-08-08 10:15:18

我已经解决了我的问题。

鲁安戈对此链接非常有帮助，它帮助我确定如何正确地缩放我的频率箱，并适当地绘制DFT。

另外:我忘了考虑信号中的偏移量。它不仅在DFT中引起0 Hz的巨峰问题，而且对变换后的信号中的大部分噪声也起着一定的作用。我利用scipy.signal.detrend()消除了这个问题，得到了一个非常合适的DFT。

# import DFT and signal packages from SciPy
import scipy.fftpack
import scipy.signal

# temporary copy to maintain data structure; original data frame is NOT changed due to ".copy()"
lead_pts_fft_df = lead_pts_df.copy()

# apply a discrete fast Fourier transform to each data series in the data frame AND detrend the signal
lead_pts_fft_df.magnitude = lead_pts_fft_df.magnitude.apply(scipy.signal.detrend)
lead_pts_fft_df.magnitude = lead_pts_fft_df.magnitude.apply(np.fft.fft)
lead_pts_fft_df

相应地安排频率箱：

num_projections = 29556
N = num_projections
T = (6 * 3) / 1000 # 6*3 b/c of the nature of my signal: 1 pt for each waveform collected every third acquisition
xf = np.linspace(0.0, 1.0 / (2.0 * T), num_projections / 2)

然后情节：

# generate subplots
fig, axes = plt.subplots(3, 1, figsize=(18, 16))

for idx in range(len(ECHO_TIMES)):
    axes[idx].plot(xf, 2.0/num_projections * np.abs(lead_pts_fft_df.magnitude[idx][:num_projections//2]),
                   label=lead_pts_df.index[idx], # apply labels
                   color=plot_colors[idx]) # colors
    axes[idx].legend() # apply legend to each set of axes

# show the plot
plt.show()