将FFT频谱幅度归一化为0dB

Wan*_*ong 8 python fft frequency audacity spectrum

我正在使用FFT从音频文件中提取每个频率分量的幅度.实际上,Audacity中已经有一个名为Plot Spectrum的函数可以帮助解决问题.以3kHz正弦和6kHz正弦组成的音频文件为例,频谱结果如下图所示.你可以看到峰值在3KHz和6kHz,没有额外的频率.

在此输入图像描述

现在我需要实现相同的功能并在Python中绘制类似的结果.我在帮助下接近Audacity结果,rfft但在得到这个结果后仍然有问题需要解决.

在此输入图像描述

  1. 第二张图中振幅的物理意义是什么?
  2. 如何将幅度标准化为0dB,就像Audacity中的那样?
  3. 为什么6kHz以上的频率具有如此高的幅度(≥90)?我可以将这些频率调整到相对较低的水平吗?

相关代码:

import numpy as np
from pylab import plot, show
from scipy.io import wavfile

sample_rate, x = wavfile.read('sine3k6k.wav')
fs = 44100.0

rfft = np.abs(np.fft.rfft(x))
p = 20*np.log10(rfft)
f = np.linspace(0, fs/2, len(p))

plot(f, p)
show()
Run Code Online (Sandbox Code Playgroud)

更新

我将Hanning窗口与整个长度信号相乘(这是正确的吗?)并得到它.裙子的大部分幅度都低于40.

在此输入图像描述

并按照@Mateen Ulhaq的说法将y轴缩放为分贝.结果更接近Audacity one.我可以将低于-90dB的幅度处理得如此之低以至于可以忽略吗?

更新的代码:

fs, x = wavfile.read('input/sine3k6k.wav')
x = x * np.hanning(len(x))

rfft = np.abs(np.fft.rfft(x))
rfft_max = max(rfft)
p = 20*np.log10(rfft/rfft_max)
f = np.linspace(0, fs/2, len(p))
Run Code Online (Sandbox Code Playgroud)

在此输入图像描述


关于赏金

通过上面更新中的代码,我可以用分贝测量频率分量.最高可能值为0dB.但该方法仅适用于特定的音频文件,因为它使用rfft_max此音频.我想像Audacity那样在一个标准规则中测量多个音频文件的频率成分.

我也在Audacity论坛上开始讨论,但我还不清楚如何实现我的目的.

igr*_*nis 5

在对Audacity源代码进行一些反向工程之后,这里有一些答案。首先,他们使用Welch算法估算PSD。简而言之,它将信号拆分为重叠的段,应用一些窗口函数,应用FFT并对结果求平均值。主要是因为存在噪音时,这有助于获得更好的结果。无论如何,在提取必要的参数之后,这里是近似于Audacity频谱图的解决方案:

import numpy as np
from scipy.io import wavfile
from scipy import signal
from matplotlib import pyplot as plt

segment_size = 512

fs, x = wavfile.read('sine3k6k.wav')
x = x / 32768.0  # scale signal to [-1.0 .. 1.0]

noverlap = segment_size / 2
f, Pxx = signal.welch(x,                        # signal
                      fs=fs,                    # sample rate
                      nperseg=segment_size,     # segment size
                      window='hanning',         # window type to use
                      nfft=segment_size,        # num. of samples in FFT
                      detrend=False,            # remove DC part
                      scaling='spectrum',       # return power spectrum [V^2]
                      noverlap=noverlap)        # overlap between segments

# set 0 dB to energy of sine wave with maximum amplitude
ref = (1/np.sqrt(2)**2)   # simply 0.5 ;)
p = 10 * np.log10(Pxx/ref)

fill_to = -150 * (np.ones_like(p))  # anything below -150dB is irrelevant
plt.fill_between(f, p, fill_to )
plt.xlim([f[2], f[-1]])
plt.ylim([-90, 6])
# plt.xscale('log')   # uncomment if you want log scale on x-axis
plt.xlabel('f, Hz')
plt.ylabel('Power spectrum, dB')
plt.grid(True)
plt.show()
Run Code Online (Sandbox Code Playgroud)

有关参数的一些必要说明:

  • wave文件被读取为16位PCM,为了与Audacity兼容,应将其缩放为| A | <1.0
  • segment_size对应Size于Audacity的GUI。
  • 默认窗口类型为“汉宁”,您可以根据需要进行更改。
  • 重叠segment_size/2与Audacity代码相同。
  • 输出窗口的框架遵循Audacity样式。他们扔掉了第一个低频箱,并将所有东西都切掉了-90dB以下

在此处输入图片说明

第二张图片中振幅的物理含义是什么?

它基本上是频率仓中的能量。

如何像Audacity一样将幅度归一化为0dB?

您需要选择一些参考点。分贝图总是与某物相关。当选择最大能量箱作为参考时,您的0db点就是最大能量(显然)。可以将最大振幅的正弦波设置为参考能量。参见ref变量。正弦信号的功率仅是RMS的平方,而要获得RMS,只需将幅度除以sqrt(2)。因此比例因子仅为0.5。请注意,之前的系数log10是10而不是20,这是因为我们处理信号的功率而不是振幅。

我是否可以将-90dB以下的幅度处理得如此之低以至于可以忽略不计?

是的,低于-40dB的任何东西通常被认为是可疏忽的