logo标签列表

使用快速傅里叶变换分析音频

pythonaudiosignal-processingfftspectrum
111

我正在尝试用 Python 创建一个图形频谱分析器。

目前，我正在读取 16 位双声道 44,100 Hz 采样率音频流的 1024 字节，并将两个通道的振幅平均值相加。现在我有了一个由 256 个带符号短整数组成的数组。接下来，我想使用类似 numpy 的模块对该数组执行 fft 并使用结果创建图形频谱分析器，首先只会包含 32 条条形图。

我已阅读过快速傅里叶变换（FFT）和离散傅里叶变换（DFT）的维基百科文章，但仍不清楚结果数组代表什么。使用 numpy 对我的数组执行 fft 后，数组看起来像这样：

   [ -3.37260500e+05 +0.00000000e+00j   7.11787022e+05 +1.70667403e+04j
   4.10040193e+05 +3.28653370e+05j   9.90933073e+04 +1.60555003e+05j
   2.28787050e+05 +3.24141951e+05j   2.09781047e+04 +2.31063376e+05j
  -2.15941453e+05 +1.63773851e+05j  -7.07833051e+04 +1.52467334e+05j
  -1.37440802e+05 +6.28107674e+04j  -7.07536614e+03 +5.55634993e+03j
  -4.31009964e+04 -1.74891657e+05j   1.39384348e+05 +1.95956947e+04j
   1.73613033e+05 +1.16883207e+05j   1.15610357e+05 -2.62619884e+04j
  -2.05469722e+05 +1.71343186e+05j  -1.56779748e+04 +1.51258101e+05j
  -2.08639913e+05 +6.07372799e+04j  -2.90623668e+05 -2.79550838e+05j
  -1.68112214e+05 +4.47877871e+04j  -1.21289916e+03 +1.18397979e+05j
  -1.55779104e+05 +5.06852464e+04j   1.95309737e+05 +1.93876325e+04j
  -2.80400414e+05 +6.90079265e+04j   1.25892113e+04 -1.39293422e+05j
   3.10709174e+04 -1.35248953e+05j   1.31003438e+05 +1.90799303e+05j...

我想知道这些数字具体代表什么，以及如何将这些数字转换为32个条形图每个条形的高度的百分比。另外，我是否应该将2个通道平均起来？

4个回答
215
你展示的数组是音频信号的傅里叶变换系数。这些系数可以用于获取音频的频率内容。FFT被定义为复值输入函数，因此你得到的系数将是虚数，即使你的输入都是实数。为了获得每个频率的功率量，你需要计算每个频率的FFT系数的幅值。这不仅仅是系数的实部，你需要计算其实部和虚部平方和的平方根。也就是说，如果你的系数是a + b*j，则它的幅值是sqrt(a^2 + b^2)。
一旦你计算出每个FFT系数的幅值，你需要确定每个FFT系数所属的音频频率。一个N点FFT将在N个等间距频率处给出信号的频率内容，起始频率为0。由于采样频率为44100样本/秒，FFT中点数为256，因此你的频率间隔为44100/256=172 Hz（约）。
您的数组中的第一个系数将是0频率系数。这基本上是所有频率的平均功率水平。其余的系数将从0开始，以172 Hz的倍数递增，直到达到128。在FFT中，您只能测量高达一半样本点的频率。如果您是一个挑战者并且需要知道为什么，请阅读关于Nyquist Frequency和Nyquist-Shannon Sampling Theorem的链接，但基本结果是您的低频将在更高的频率桶中被复制或aliased。因此，频率将从0开始，每个系数增加172 Hz，直到N/2系数，然后减少172 Hz，直到N-1系数。
这应该足够的信息让您开始了。如果您想要比维基百科给出的更容易理解的FFT介绍，可以尝试Understanding Digital Signal Processing: 2nd Ed.。它对我非常有帮助。
所以那些数字就是代表的内容。将每个频率分量的幅度按所有分量幅度的总和进行缩放，即可转换为高度百分比。但这只会给您提供相对频率分布的表示，而不是每个频率的实际功率。您可以尝试按频率分量可能的最大幅度进行缩放，但我不确定是否会显示得很好。找到可行的缩放因子的最快方法是在响亮和柔和的音频信号上进行实验，以找到合适的设置。
最后，如果您想展示整个音频信号的频率内容，则应将两个通道平均在一起。您正在将立体声音频混合成单声道音频，并显示组合频率。如果您想要右侧和左侧频率的两个单独的显示，则需要分别对每个通道执行傅里叶变换。
1
我在网上大多只能找到过于复杂的FFT解释，这是一个很好、简单的解释，说明了采样点数如何影响FFT的结果。谢谢你！ - echolocation
27

虽然这个帖子已经有好几年了，但我觉得非常有帮助。我想为那些正在尝试创建类似东西的人提供我的建议。

至于将数据分成条形图，不应该像Antti建议的那样根据条数平均分配数据。最有用的是将数据分成八度部分，每个八度部分都是前一个频率的两倍。(例如100赫兹是50赫兹的一个八度，而50赫兹又是25赫兹的一个八度)。

根据你想要多少个条形图，你将整个范围分成1/X个八度范围。基于某个中心频率A在条形图上，你可以从以下方法中得到条形图的上下限：

upper limit = A * 2 ^ ( 1 / 2X )
lower limit = A / 2 ^ ( 1 / 2X )

计算下一个相邻中心频率可以使用类似的计算方法：

next lower =  A / 2 ^ ( 1 / X )
next higher = A * 2 ^ ( 1 / X )

接着，您需要对符合这些范围的数据进行平均以得到每个柱形图的振幅。

例如： 我们想将其分成1/3倍频程范围，并从1khz的中心频率开始。

Upper limit = 1000 * 2 ^ ( 1 / ( 2 * 3 ) ) = 1122.5
Lower limit = 1000 / 2 ^ ( 1 / ( 2 * 3 ) ) =  890.9

假设我们有 44100hz 和 1024 个样本（每个数据点之间有43hz），我们应该平均计算从21到26的值。（890.9 / 43 = 20.72约等于21，1122.5 / 43 = 26.10约等于26）

（1/3倍频条可以在 ~40hz 到 ~20khz 之间获得大约30个条。）如你现在所能理解的那样，随着频率升高，我们将会对更多的数字进行平均计算。较低的频率条通常只包含一个或少数几个数据点。而较高的频率条可能是数百个数据点的平均值。原因是 86hz 是 43hz 的八度音......而 10086hz 几乎与 10043hz 听起来一样。

10
您手头拥有的是一个采样时长为256/44100 = 0.00580499秒的样本。这意味着您的频率分辨率为1 / 0.00580499 = 172 Hz。从Python中得到的256个值对应于基本上从86 Hz到255 * 172 + 86 Hz = 43946 Hz的频率。您得到的数字是复数（因此第二个数字末尾带有“j”）。
编辑：修正错误信息
您需要通过计算sqrt（i 2 + j 2 ）来将复数转换为振幅，其中i和j分别是实部和虚部。
如果您想要32个条形图，据我所理解，您应该取四个连续幅度的平均值，以获得256/4 = 32个条形图，就像您想要的那样。
4
请注意，如果c是一个复数，sqrt(c.real2 + c.imag2) == abs(c)。该等式表明，如果您想要计算复数的模长（即距离原点的距离），则可以使用该公式来计算。 - tzot
0

FFT返回N个复数值，您可以计算module=sqrt(real_part^2+imaginary_part^2)。要获取每个频段的值，您必须对频段内所有谐波的模块求和。下面是一个10条柱形谱分析仪的示例。C代码必须被包装以获得pyd Python模块。

float *samples_vett;
float *out_filters_vett;
int Nsamples;
float band_power = 0.0;
float harmonic_amplitude=0.0;
int i, out_index;

out_index=0;


for (i = 0; i < Nsamples / 2 + 1; i++)       
        {
            if (i == 1 || i == 2 || i == 4 || i == 8 || i == 17 || i == 33 || i == 66 || i == 132 || i == 264 || i == 511)
            {
                out_filters_vett[out_index] = band_power; 
                band_power = 0; 
                out_index++;  
            }

            harmonic_amplitude = sqrt(pow(ttfr_out_vett[i].r, 2) + pow(ttfr_out_vett[i].i, 2));
            band_power += harmonic_amplitude;

        }

我使用Python设计并制作了一个整个10 LED灯条谱分析仪。与使用nunmpy库（太大而且没有用，只想得到FFT）不同的是，我创建了一个Python pyd模块（仅27KB）来获取FFT和将整个音频频谱分成频段。

此外，为了读取输出音频，我创建了一个回路WASapi portaudio pyd模块。您可以在图片中看到该项目（块图）10BarsSpectrumAnalyzerWithWASapi.jpg。

我刚刚在我的YouTube频道上添加了一个教程视频：如何设计和制作一个非常聪明的Python谱分析仪10 Led Bar。

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接