获取与Librosa中的STFT相关的频率

问题描述

使用librosa.stft()计算频谱图时，如何获得相关的频率值？我对生成librosa.display.specshow中的图像不感兴趣，但是我想掌握这些值。

y,sr = librosa.load('../recordings/high_pitch.m4a')
stft = librosa.stft(y,n_fft=256,window=sig.windows.hamming)
spec = np.abs(stft)

spec给了我每个频率的“振幅”或“功率”，但没有给我频率箱本身。我已经看到有一个display.specshow函数将在热图的垂直轴上显示这些频率值，但不会自行返回这些值。

我正在为单个FFT寻找类似于nn.fft.fttfreq()的东西，但是在librosa文档中找不到等效的东西。

解决方法

我想特别指出以下问题和答案：How do I obtain the frequencies of each value in an FFT?。除了查询documentation for the STFT from librosa外，我们知道水平轴是时间轴，而垂直轴是频率。频谱图中的每一列都是时间片的FFT，其中此时间点的中心有一个放置有n_fft=256个分量的窗口。

我们还知道有一个 hop length ，它告诉我们在计算下一个FFT之前需要跳过多少个音频样本。默认情况下，此值为n_fft / 4，因此音频中的每256/4 = 64点，我们将在这个n_fft=256点长的时间点处计算一个新的FFT。如果您想知道每个窗口所处的确切时间点，那就是i / Fs，其中i是音频信号的索引，它是64的倍数。

现在，对于每个FFT窗口，对于真实信号，频谱都是对称的，因此我们仅考虑FFT的正向。这已通过文档进行了验证，其中行数和频率分量的数量为1 + n_fft / 2，其中1为直流分量。既然我们有了这个，请咨询上面从bin编号到相应频率的关系的帖子为i * Fs / n_fft，其中i是bin编号，Fs是采样频率，{{1} }作为FFT窗口中的点数。由于我们仅查看半频谱，而不是n_fft=256从0扩展到i，所以它从0扩展到n_fft，而不是1 + n_fft / 2以外的bin是半频谱的反射版本，因此我们不考虑超过1 + n_fft / 2 Hz的频率分量。

如果您想生成这些频率的NumPy数组，则可以执行以下操作：

Fs / 2

import numpy as np freqs = np.arange(0,1 + n_fft / 2) * Fs / n_fft 是一个将FFT中的bin编号映射到相应频率的数组。作为说明性示例，假设我们的采样频率为16384 Hz，freqs。因此：

n_fft = 256

我们可以看到我们已经生成了In [1]: import numpy as np In [2]: Fs = 16384 In [3]: n_fft = 256 In [4]: np.arange(0,1 + n_fft / 2) * Fs / n_fft Out[4]: array([ 0.,64.,128.,192.,256.,320.,384.,448.,512.,576.,640.,704.,768.,832.,896.,960.,1024.,1088.,1152.,1216.,1280.,1344.,1408.,1472.,1536.,1600.,1664.,1728.,1792.,1856.,1920.,1984.,2048.,2112.,2176.,2240.,2304.,2368.,2432.,2496.,2560.,2624.,2688.,2752.,2816.,2880.,2944.,3008.,3072.,3136.,3200.,3264.,3328.,3392.,3456.,3520.,3584.,3648.,3712.,3776.,3840.,3904.,3968.,4032.,4096.,4160.,4224.,4288.,4352.,4416.,4480.,4544.,4608.,4672.,4736.,4800.,4864.,4928.,4992.,5056.,5120.,5184.,5248.,5312.,5376.,5440.,5504.,5568.,5632.,5696.,5760.,5824.,5888.,5952.,6016.,6080.,6144.,6208.,6272.,6336.,6400.,6464.,6528.,6592.,6656.,6720.,6784.,6848.,6912.,6976.,7040.,7104.,7168.,7232.,7296.,7360.,7424.,7488.,7552.,7616.,7680.,7744.,7808.,7872.,7936.,8000.,8064.,8128.,8192.]) In [5]: freqs = _; len(freqs) Out[5]: 129 元素数组，该数组告诉我们每个对应bin编号的频率。

---

警告

请注意，librosa.display.specshow的默认采样率为22050 Hz，因此，如果您未将采样率（1 + n_fft / 2 = 129）设置为与音频信号相同的采样频率，水平轴将不正确。确保指定sr输入标志以匹配输入音频的采样频率。

,

除了excellent explanation和rayryeng之外，librosa中numpy.fft.fftfreq()的直接等效项是librosa.fft_frequencies()

您可以按以下方式使用它：

y,sr = librosa.load('../recordings/high_pitch.m4a')
Nfft = 256
stft = librosa.stft(y,n_fft=Nfft,window=sig.windows.hamming)
freqs = librosa.fft_frequencies(sr=sr,n_fft=Nfft)

,

您可以按如下方式计算累积能量

samplerate = 48000
Nfft = 8192
freqs = librosa.fft_frequencies(sr=sr,n_fft=Nfft)
plt.loglog(freqs,np.mean(mag**2,axis=1)/(Nfft/2)**2)
plt.xlabel('freq [Hz]')

如果你想对一个频率范围内的能量求和，你可以在频率上使用索引 mag，例如

np.sum(np.mean(mag[(freqs > 1000) & (freqs < 1480),:]**2,axis=1))/(Nfft/2)**2

更一般地，您可以应用过滤器 gain(f)，上面的结果是使用 gain(f) 矩形获得的。

np.sum(np.mean(mag**2,axis=1)*gain(freq))/(Nfft/2)**2

免责声明：我不知道这些比例因子是否适合您。只有形状。