【问题标题】:NumPy Fast Fourier transform (FFT) does not work on sine wave generated in AudacityNumPy 快速傅里叶变换 (FFT) 不适用于 Audacity 中生成的正弦波
【发布时间】:2019-09-22 03:04:25
【问题描述】:

我正在尝试使用 Python 的 NumPy 库进行一些频率分析。我有两个 .wav 文件,它们都包含 440 Hz 正弦波。其中一个是我使用 NumPy 正弦函数生成的,另一个是我在 Audacity 中生成的。 FFT 适用于 Python 生成的 FFT,但对 Audacity 无任何作用。

这是两个文件的链接:

非工作文件:440_audacity.wav

工作文件:440_gen.wav

这是我用来进行傅里叶变换的代码:

import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wave

infile = "440_gen.wav"
rate, data = wave.read(infile)

data = np.array(data)

data_fft = np.fft.fft(data)
frequencies = np.abs(data_fft)

plt.subplot(2,1,1)
plt.plot(data[:800])
plt.title("Original wave: " + infile)

plt.subplot(2,1,2)
plt.plot(frequencies)
plt.title("Fourier transform results")

plt.xlim(0, 1000)

plt.tight_layout()

plt.show()

我有两个 16 位 PCM .wav 文件,一个来自 Audacity,一个使用 NumPy 正弦函数创建。 NumPy 生成的结果如下(正确),峰值为 440Hz:

我用 Audacity 创建的波形虽然看起来相同,但在傅立叶变换上没有给出任何结果:

我承认我在这里不知所措。这两个文件实际上应该包含相同的数据。它们的编码方式相同,波形在上图中显示相同。

这是用于生成工作文件的代码:

import numpy as np
import wave
import struct
import matplotlib.pyplot as plt
from operator import add

freq_one = 440.0
num_samples = 44100
sample_rate = 44100.0
amplitude = 12800

file = "440_gen.wav"

s1 = [np.sin(2 * np.pi * freq_one * x/sample_rate) * amplitude for x in range(num_samples)]

sine_one = np.array(s1)

nframes = num_samples
comptype = "NONE"
compname="not compressed"
nchannels = 1
sampwidth = 2

wav_file = wave.open(file, 'w')
wav_file.setparams((nchannels, sampwidth, int(sample_rate), nframes, comptype, compname))

for s in sine_one:
    wav_file.writeframes(struct.pack('h', int(s)))

【问题讨论】:

  • 有没有机会分享440_gen.wav
  • 在终端尝试用aplay 440_gen.wav录音,看看wav文件是否正确。
  • @Ralph 我正在上传要附加的两个文件。
  • @darkfire613 谢谢。只需 audacity 生成的文件就应该这样做,因为您包含了生成 python 生成的.wav 的代码。此外,“波形在上图中看起来相同”是不正确的。尽管它们的频率相同,但幅度却不同(尽管这无关紧要,我认为这不是问题)
  • @Ralph 我已将 Dropbox 链接添加到这两个文件。

标签: python numpy scipy fft audacity


【解决方案1】:

自从回答了这个问题后,@Konyukh Fyodorov 能够提供更好且合理的解决方案(如下)。


以下内容对我有用,并按预期生成了情节。不幸的是,我无法完全拼凑出它为什么有效,但我分享了这个解决方案,希望它可以帮助其他人实现这一飞跃。

import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wave

infile = "440_gen.wav"
rate, data = wave.read(infile)

data = np.array(data)

# Use first 44100 datapoints in transform
data_fft = np.fft.fft(data[:44100])
frequencies = np.abs(data_fft)

plt.subplot(2,1,1)
plt.plot(data[:800])
plt.title("Original wave: " + infile)

plt.subplot(2,1,2)
plt.plot(frequencies)
plt.title("Fourier transform results")

plt.xlim(0, 1000)

plt.tight_layout()

plt.show()

【讨论】:

  • 谢谢!这也适用于我,使用我在 Audacity 中创建的多个不同的示例文件,并且之前遇到过困难。有必要限制样本数量吗?
  • 太好了,很高兴这对你有用!我不太确定它为什么会起作用;自从我使用傅立叶变换以来已经有很长时间了,我真的不记得它们是如何工作的。我需要考虑一下。
  • “这对我有用”的答案没有说明代码中发生了什么变化,也没有说明它为什么起作用,不帮助任何人学习。这是一个修补代码但不会添加到 Stack Overflow 知识库的答案。
  • @Cris Luengo 尊敬的,我不同意:并非所有问题都以线性方式解决。有时在了解为什么修复之前看到修复可以帮助弥合知识差距并帮助学习。我无法弥合这一差距,但分享了我的补丁,希望它可以帮助其他人实现这一飞跃。您可以选择优先考虑解决方案的“纯度高于实用性”,但不是每个人都这样做。我的功利主义者会很乐意通过一个有效的解决方案而不是等待 2 年,因为我限制性的个人哲学(但这只是我)
【解决方案2】:

让我解释一下为什么您的代码不起作用。以及为什么它适用于[:44100]

首先,你有不同的文件:

440_gen.wav      = 1 sec and 44100  samples (counts)        
440_audacity.wav = 5 sec and 220500 samples (counts)

由于在 FFT 中对于 440_gen.wav,您使用参考点数 N=44100 和采样率 44100,因此您的频率分辨率为 1 Hz(bin 以 1 Hz 为增量跟随)。
因此,在图表上,每个 FFT 样本对应一个等于 1 Hz 的增量。
plt.xlim(0, 1000) 仅对应于 0-1000 Hz 的范围。

但是,对于 FFT 中的 440_audacity.wav,您使用参考点的数量 N=220500 和采样率 44100。您的频率分辨率为 0.2 Hz(bin 以 0.2 Hz 为增量跟随) - 在图表上,每个 FFT 样本对应于 0.2 Hz 增量的频率(最小值-最大值 = +(-) 22500 Hz)。
plt.xlim(0, 1000) 仅对应于 1000x0.2 = 0-200 Hz 的范围。
这就是结果不可见的原因 - 它不在此范围内。

plt.xlim (0, 5000) 将纠正您的情况并将范围扩展到 0-1000 Hz。

jwalton 引入的解决方案[:44100] 实际上只强制 FFT 使用 N = 44100。这重复了计算 440_gen.wav 的情况

更正确的解决您的问题的方法是在代码中使用N (Windows Size) 参数和np.fft.fftfreq() 函数。

下面的示例代码。

我也推荐一篇优秀的文章https://realpython.com/python-scipy-fft/

import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wave

N = 44100 # added

infile = "440_audacity.wav"
rate, data = wave.read(infile)

data = np.array(data)

data_fft = np.fft.fft(data, N)  # added N
frequencies = np.abs(data_fft)
x_freq = np.fft.fftfreq(N, 1/44100)  # added

plt.subplot(2,1,1)
plt.plot(data[:800])
plt.title("Original wave: " + infile)

plt.subplot(2,1,2)
plt.plot(x_freq, frequencies)  # added x_freq 
plt.title("Fourier transform results")

plt.xlim(0, 1000)
plt.tight_layout()
plt.show()

【讨论】:

  • 优秀的答案。欢迎使用 Stack Overflow!
  • 顺便说一句,这也可以在不将信号切割到第一秒的情况下工作,将N设置为信号的长度。
  • 谢谢你这么说,克里斯
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2015-02-10
  • 1970-01-01
  • 1970-01-01
  • 2012-10-14
  • 2011-07-12
  • 2011-04-09
相关资源
最近更新 更多