从 FFT 访问处理后的值

Posted 2023-02-25

【中文标题】从 FFT 访问处理后的值

【英文标题】：Accessing processed values from FFT

【发布时间】：2016-09-07 00:39:03

【问题描述】：

我正在尝试使用Lomont FFT 来返回复数以使用 c# 构建频谱图/频谱密度图。

我无法理解如何从类中返回值。 这是我迄今为止编写的代码，似乎正在运行。

    int read = 0;
    Double[] data;
    byte[] buffer = new byte[1024];

    FileStream wave = new FileStream(args[0], FileMode.Open, FileAccess.Read);  
    read = wave.Read(buffer, 0, 44);                                            
    read = wave.Read(buffer, 0, 1024);                                          
    data = new Double[read];                                                    


    for (int i = 0; i < read; i+=2)
    
        data[i] = BitConverter.ToInt16(buffer, i) / 32768.0;
        Console.WriteLine(data[i]);
    
    LomontFFT LFFT = new LomontFFT();
    LFFT.FFT(data, true);

我不清楚的是，如何将 Lomont FFT 实现中的值返回/访问回我的应用程序（控制台）？

作为 c# 开发的新手，我想我可能缺少了解如何从 Lomont 类的实例中检索已处理值的基本方面，或者甚至可能不正确地调用它。

Console.WriteLine(LFFT.A); // Returns 0
Console.WriteLine(LFFT.B); // Returns 1

我一直在寻找代码 sn-p 或如何执行此操作的解释，但到目前为止，我没有找到任何我理解或解释我所面临问题的特定方面的内容。任何指导将不胜感激。

上面代码中提到的data 数组中的结果子集可以在下面找到，根据我目前的理解，它似乎是有效的：

0.00531005859375
0.0238037109375
0.041473388671875
0.0576171875
0.07183837890625
0.083465576171875
0.092193603515625
0.097625732421875
0.099639892578125
0.098114013671875
0.0931396484375
0.0848388671875
0.07354736328125
0.05963134765625
0.043609619140625
0.026031494140625
0.007476806640625
-0.011260986328125
-0.0296630859375
-0.047027587890625
-0.062713623046875
-0.076141357421875
-0.086883544921875
-0.09454345703125
-0.098785400390625
-0.0994873046875
-0.0966796875
-0.090362548828125
-0.080810546875
-0.06842041015625
-0.05352783203125
-0.036712646484375
-0.0185546875

我实际上想做什么？（视角）

我希望将波形文件加载到控制台应用程序中，并以 jpg/png 格式返回频谱图/频谱密度图/图像以供进一步处理。

我正在阅读的波形文件是单声道格式

---

更新 1

根据使用的 FFT，我收到的结果略有不同。

使用 RealFFT

    for (int i = 0; i < read; i+=2)
    
        data[i] = BitConverter.ToInt16(buffer, i) / 32768.0;
        //Console.WriteLine(data[i]);
    

    LomontFFT LFFT = new LomontFFT();
    LFFT.RealFFT(data, true);

    for (int i = 0; i < buffer.Length / 2; i++)
    
        System.Console.WriteLine("0",
          Math.Sqrt(data[2 * i] * data[2 * i] + data[2 * i + 1] * data[2 * i + 1]));
    

RealFFT 的部分结果

0.314566983321381
0.625242818210924
0.30314888696868
0.118468857708093
0.0587697011760449
0.0369034115568654
0.0265842582236275
0.0207195964060356
0.0169601273233317
0.0143745438577886
0.012528799609089
0.0111831275153128
0.0102313284519146
0.00960198279358434
0.00920236001619566

使用 FFT

    for (int i = 0; i < read; i+=2)
    
        data[i] = BitConverter.ToInt16(buffer, i) / 32768.0;
        //Console.WriteLine(data[i]);
    

    double[] bufferB = new double[2 * data.Length];

    for (int i = 0; i < data.Length; i++)
    
        bufferB[2 * i] = data[i];
        bufferB[2 * i + 1] = 0;
    

    LomontFFT LFFT = new LomontFFT();
    LFFT.FFT(bufferB, true);

    for (int i = 0; i < bufferB.Length / 2; i++)
    
        System.Console.WriteLine("0",
          Math.Sqrt(bufferB[2 * i] * bufferB[2 * i] + bufferB[2 * i + 1] * bufferB[2 * i + 1]));
    

FFT 的部分结果：

0.31456698332138
0.625242818210923
0.303148886968679
0.118468857708092
0.0587697011760447
0.0369034115568653
0.0265842582236274
0.0207195964060355
0.0169601273233317
0.0143745438577886
0.012528799609089
0.0111831275153127
0.0102313284519146
0.00960198279358439
0.00920236001619564

【参考方案1】：

查看LomontFFT.FFT 文档：

计算数据的正向或反向傅里叶变换，使用 包含复值数据的数据作为交替实数和 虚构部分。长度必须是 2 的幂。数据是 原地修改。

这告诉我们一些事情。首先，该函数需要复值数据，而您的数据是真实的。对此的快速解决方法是创建另一个两倍大小的缓冲区并将所有虚部设置为 0：

double[] buffer = new double[2*data.Length];
for (int i=0; i<data.Length; i++)

  buffer[2*i] = data[i];
  buffer[2*i+1] = 0;

文档还告诉我们，计算已就地完成。这意味着在对FFT 的调用返回后，输入数组将替换为计算结果。因此，您可以使用以下命令打印光谱：

LomontFFT LFFT = new LomontFFT();
LFFT.FFT(buffer, true);

for (int i = 0; i < buffer.Length/2; i++)

  System.Console.WriteLine("0",
    Math.Sqrt(buffer[2*i]*buffer[2*i]+buffer[2*i+1]*buffer[2*i+1]));

请注意，由于您的输入 data 是实值，您也可以使用 LomontFFT.RealFFT。在这种情况下，给定一个稍微不同的打包规则，您将获得 FFT 结果：

LomontFFT LFFT = new LomontFFT();
LFFT.RealFFT(data, true);

System.Console.WriteLine("0", Math.Abs(data[0]);
for (int i = 1; i < data.Length/2; i++)

  System.Console.WriteLine("0",
    Math.Sqrt(data[2*i]*data[2*i]+data[2*i+1]*data[2*i+1]));

System.Console.WriteLine("0", Math.Abs(data[1]);

这将为您提供频谱的非冗余下半部分（与提供整个频谱的LomontFFT.FFT 不同）。此外，相对于LomontFFT.FFT，可以预期在double 精度（大约1e-16 乘以频谱峰值）数量级上的数值差异。

【讨论】：

感谢您的解释，您所指出的对我来说很有意义。根据我使用的 FFT 算法（RealFFT 还是 FFT），我收到的结果似乎略有不同。我不确定这实际上有多重要？ OP 中的更新 1 对此进行了描述。 @SleuthEye

@SleutheEye，所以在这里澄清一下，如果我只能期望来自 RealFFT 的 nyquist 和以下结果，为什么我们在 for 循环中将 data.Length 除以 2 data.Length/2？这是由于从 RealFFT 返回的虚数以及将data.Length 减半用作过滤高于奈奎斯特的所有项目的一种手段吗？

data.Length/2 for 循环是由于复数由虚部和实部组成（即每个复数有 2 个 double 值）。 RealFFT 仅返回频谱的下半部分这一事实意味着 data.Length == buffer.Length/2 （对应的等效 buffer 与 FFT 调用一起使用）。因此，如果您有 1024 个实数值，您将有 data.Length==1024、buffer.Length==2048 并且您将执行 RealFFT sn-p 的循环体 511 次（我们跳过索引 0，其中包含 2 个 0Hz 的纯实分量和奈奎斯特频率）。