对已获取到的音频数据进行复数离散傅里叶分析时，发现俩个在GoldWave观察基本一致的音频数据，经过傅里叶之后的结果却相差甚大
帮忙分析一下，谢谢
音频数据的链接: https://pan.baidu.com/s/1XCUrul6D91tLzbF-Ntc_jA 提取码: nyar
主要运行的函数顺序:
1 int MMITestRealizeFFT(uint16* input, size_t datalen)
2 void cdft(int n, int isgn, double *a, int *ip, double *w)
3 void bitrv2conj(int n, int *ip, double *a)
4 void cftbsub(int n, double *a, double *w)
下面是代码：
uint8 *RecordBuf = malloc(1024 * 2 + 1);
MMITestRealizeFFT((uint16*)RecordBuf, 1024); //RecordBuf 为PCM音频的原始数据
int MMITestRealizeFFT(uint16* input, size_t datalen) //1
{
int ip[64 + 2] = {0};
int lvalue = 0;
int rvalue = 0;
int mvalue = 0;
double *inputfft = malloc((2048 + 1) * sizeof(double));
double ww[2048 * 5 / 4] = {0};
int ii = 0;
ip[0] = 0;
int ret = 0;
static int tmp_se35n = 0;
uint16 outputfft[2048/2 + 1] = {0};
memset(inputfft, 0, (2048 + 1) * sizeof(double));
for(ii= 0; ii< datalen; ii++)
{
inputfft[(ii << 1)] = input[ii];
inputfft[(ii << 1) +1 ] = 0;
}
cdft(datalen<<1, -1, inputfft, ip, ww); // 2 复数离散傅里叶
for(ii= 0; ii< datalen; ii++)
{
outputfft[ii] = (uint16)inputfft[ii << 1];
}
for(ii= 40; ii< 95; ii++) //312-741HZ
{
lvalue += outputfft[ii];
}
lvalue /= 55;
for(ii= 95; ii< 150; ii++) //741-1170HZ
{
mvalue += outputfft[ii];
}
mvalue /= 55;
for(ii= 150; ii< 205; ii++) //1170-1590HZ
{
rvalue += outputfft[ii];
}
rvalue /= 55;
uart_printf("\r\n=========lvalue = %d, mvalue = %d, rvalue = %d===========\n", lvalue,mvalue,rvalue);
free(inputfft);
return CHICK_OK;
}
void cdft(int n, int isgn, double *a, int *ip, double *w)
{
void makewt(int nw, int *ip, double *w);
void bitrv2(int n, int *ip, double *a);
void bitrv2conj(int n, int *ip, double *a);
void cftfsub(int n, double *a, double *w);
void cftbsub(int n, double *a, double *w);
if (n > (ip[0] << 2))
{
makewt(n >> 2, ip, w);
}
if (n > 4)
{
if (isgn >= 0)
{
bitrv2(n, ip + 2, a);
cftfsub(n, a, w);
}
else
{
bitrv2conj(n, ip + 2, a); //3
cftbsub(n, a, w); //4
}
}
else if (n == 4)
{
cftfsub(n, a, w);
}
}
void bitrv2conj(int n, int *ip, double *a)
{
int j, j1, k, k1, l, m, m2;
double xr, xi, yr, yi;
ip[0] = 0;
l = n;
m = 1;
while ((m << 3) < l)
{
l >>= 1;
for (j = 0; j < m; j++)
{
ip[m + j] = ip[j] + l;
}
uart_printf("%d\r ",m);
m <<= 1;
uart_printf("%d\r ", m);
}
m2 = 2 * m;
if ((m << 3) == l)
{
for (k = 0; k < m; k++)
{
for (j = 0; j < k; j++)
{
j1 = 2 * j + ip[k];
k1 = 2 * k + ip[j];
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
j1 += m2;
k1 += 2 * m2;
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
j1 += m2;
k1 -= m2;
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
j1 += m2;
k1 += 2 * m2;
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
}
k1 = 2 * k + ip[k];
a[k1 + 1] = -a[k1 + 1];
j1 = k1 + m2;
k1 = j1 + m2;
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
k1 += m2;
a[k1 + 1] = -a[k1 + 1];
}
}
else
{
a[1] = -a[1];
a[m2 + 1] = -a[m2 + 1];
for (k = 1; k < m; k++)
{
for (j = 0; j < k; j++)
{
j1 = 2 * j + ip[k];
k1 = 2 * k + ip[j];
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
j1 += m2;
k1 += m2;
xr = a[j1];
xi = -a[j1 + 1];
yr = a[k1];
yi = -a[k1 + 1];
a[j1] = yr;
a[j1 + 1] = yi;
a[k1] = xr;
a[k1 + 1] = xi;
}
k1 = 2 * k + ip[k];
a[k1 + 1] = -a[k1 + 1];
a[k1 + m2 + 1] = -a[k1 + m2 + 1];
}
}
}
void cftbsub(int n, double *a, double *w)
{
void cft1st(int n, double *a, double *w);
void cftmdl(int n, int l, double *a, double *w);
int j, j1, j2, j3, j4, j5, j6, j7, l;
double wn4r, x0r, x0i, x1r, x1i, x2r, x2i, x3r, x3i,
y0r, y0i, y1r, y1i, y2r, y2i, y3r, y3i,
y4r, y4i, y5r, y5i, y6r, y6i, y7r, y7i;
l = 2;
if (n > 16)
{
cft1st(n, a, w);
l = 16;
while ((l << 3) < n)
{
cftmdl(n, l, a, w);
l <<= 3;
}
}
if ((l << 2) < n)
{
wn4r = w[2];
for (j = 0; j < l; j += 2)
{
j1 = j + l;
j2 = j1 + l;
j3 = j2 + l;
j4 = j3 + l;
j5 = j4 + l;
j6 = j5 + l;
j7 = j6 + l;
x0r = a[j] + a[j1];
x0i = -a[j + 1] - a[j1 + 1];
x1r = a[j] - a[j1];
x1i = -a[j + 1] + a[j1 + 1];
x2r = a[j2] + a[j3];
x2i = a[j2 + 1] + a[j3 + 1];
x3r = a[j2] - a[j3];
x3i = a[j2 + 1] - a[j3 + 1];
y0r = x0r + x2r;
y0i = x0i - x2i;
y2r = x0r - x2r;
y2i = x0i + x2i;
y1r = x1r - x3i;
y1i = x1i - x3r;
y3r = x1r + x3i;
y3i = x1i + x3r;
x0r = a[j4] + a[j5];
x0i = a[j4 + 1] + a[j5 + 1];
x1r = a[j4] - a[j5];
x1i = a[j4 + 1] - a[j5 + 1];
x2r = a[j6] + a[j7];
x2i = a[j6 + 1] + a[j7 + 1];
x3r = a[j6] - a[j7];
x3i = a[j6 + 1] - a[j7 + 1];
y4r = x0r + x2r;
y4i = x0i + x2i;
y6r = x0r - x2r;
y6i = x0i - x2i;
x0r = x1r - x3i;
x0i = x1i + x3r;
x2r = x1r + x3i;
x2i = x1i - x3r;
y5r = wn4r * (x0r - x0i);
y5i = wn4r * (x0r + x0i);
y7r = wn4r * (x2r - x2i);
y7i = wn4r * (x2r + x2i);
a[j1] = y1r + y5r;
a[j1 + 1] = y1i - y5i;
a[j5] = y1r - y5r;
a[j5 + 1] = y1i + y5i;
a[j3] = y3r - y7i;
a[j3 + 1] = y3i - y7r;
a[j7] = y3r + y7i;
a[j7 + 1] = y3i + y7r;
a[j] = y0r + y4r;
a[j + 1] = y0i - y4i;
a[j4] = y0r - y4r;
a[j4 + 1] = y0i + y4i;
a[j2] = y2r - y6i;
a[j2 + 1] = y2i - y6r;
a[j6] = y2r + y6i;
a[j6 + 1] = y2i + y6r;
}
}
else if ((l << 2) == n)
{
for (j = 0; j < l; j += 2)
{
j1 = j + l;
j2 = j1 + l;
j3 = j2 + l;
x0r = a[j] + a[j1];
x0i = -a[j + 1] - a[j1 + 1];
x1r = a[j] - a[j1];
x1i = -a[j + 1] + a[j1 + 1];
x2r = a[j2] + a[j3];
x2i = a[j2 + 1] + a[j3 + 1];
x3r = a[j2] - a[j3];
x3i = a[j2 + 1] - a[j3 + 1];
a[j] = x0r + x2r;
a[j + 1] = x0i - x2i;
a[j2] = x0r - x2r;
a[j2 + 1] = x0i + x2i;
a[j1] = x1r - x3i;
a[j1 + 1] = x1i - x3r;
a[j3] = x1r + x3i;
a[j3 + 1] = x1i + x3r;
}
}
else
{
for (j = 0; j < l; j += 2)
{
j1 = j + l;
x0r = a[j] - a[j1];
x0i = -a[j + 1] + a[j1 + 1];
a[j] += a[j1];
a[j + 1] = -a[j + 1] - a[j1 + 1];
a[j1] = x0r;
a[j1 + 1] = x0i;
}
}
}
傅里叶变换，是把波形转换成频率
那些频率非常高的波形你用肉眼看能看出来差异，那你也是个电脑
根据提供的代码，你的目标似乎是对音频数据进行傅里叶变换，并计算特定频率范围内的平均值。然而，你提到经过傅里叶变换后的结果与预期的不一致。
首先，傅里叶变换将时域数据转换为频域数据，它将信号分解为一系列频率分量。你的代码似乎是使用库函数来进行快速傅里叶变换（FFT）。然而，代码中有一些问题可能会导致结果与预期不一致。
内存分配问题：在MMITestRealizeFFT函数中，你使用malloc为inputfft和outputfft分配内存。然而，在函数的末尾，你只释放了inputfft的内存，而没有释放outputfft的内存。这可能导致内存泄漏。确保在不再需要这些内存块时进行释放。
傅里叶变换的输入长度：在cdft函数中，你将输入数据长度设置为datalen<<1，但是你的输入数据长度似乎是datalen。这可能会导致访问超出数组边界，从而影响计算结果。
傅里叶变换中的复数数据格式：你的代码使用一个长度为2048+1的双精度浮点数组inputfft来存储复数数据。然而，在傅里叶变换函数cdft中，它使用double *a来接收输入数据，这可能导致类型不匹配的问题。请确保数据类型在整个代码中保持一致。
这些问题可能导致傅里叶变换的结果不正确，进而导致你观察到的结果与预期不一致。你可以尝试修复这些问题并再次运行代码，看看结果是否变得正确。
是精度差异还是错误
不知道你这个问题是否已经解决, 如果还没有解决的话:我可以给您一些排查问题的思路和建议，希望能帮到您。
检查音频数据的采样率和位深度是否一致。如果两个音频数据的采样率或者位深度不同，就算是相同的音频也会有很大差异。您可以使用GoldWave软件查看这些信息，确保两个音频数据是一致的。
检查复数离散傅里叶分析的代码实现是否正确。您可以使用其他的工具(如Python的numpy库)对这两个音频数据进行离散傅里叶变换，比对结果是否一致。如果结果一致，那么就说明您的代码没有问题；如果结果不一致，就需要检查您的代码实现是否存在问题。
检查傅里叶变换后的频谱数据是否有错位。对于傅里叶变换后的频谱数据，我们需要注意两点：频谱的零频分量是否位于频谱数据的中心，还有频谱是否有左右互换的情况。如果您没有注意到这些问题，就可能导致频谱数据的错位，进而导致两个音频数据在频域上的差异。
检查音频数据的前后是否一致。您可以对两个音频数据的前后部分进行比对，看是否存在明显的差异。如果这种差异存在，就可能说明音频数据在处理过程中出现了错误。
尝试对音频数据进行加窗。加窗可以减少频谱泄漏的情况，提高离散傅里叶变换的精度。您可以使用汉明窗或者海明窗等加窗函数，对音频数据进行窗函数处理，然后再进行离散傅里叶变换。这有助于提高分析结果的精度。
如果上述方法都没有解决问题，那么您可能需要对音频数据进行更深入的分析，或者请教其他更专业的人员，以获得更准确的答案。希望能为您提供一些参考，祝您好运。
当对音频数据进行复数离散傅里叶变换（DFT）时，可能会遇到以下一些可能导致结果差异的情况：
数据长度不匹配：确保要进行傅里叶变换的音频数据长度相同。如果长度不一致，将会导致结果差异。
数据类型不匹配：确保要进行傅里叶变换的音频数据具有相同的数据类型。例如，如果一个数据是使用浮点数表示，而另一个是使用整数表示，将会导致结果差异。
数据范围不一致：检查音频数据的范围是否一致。如果一个数据的幅值范围较大，而另一个数据的幅值范围较小，可能导致结果差异。可以尝试将数据进行归一化，使其具有相似的范围。
傅里叶变换参数设置：确保在进行傅里叶变换时，使用相同的参数设置，例如采样率、窗函数类型、窗长等。不同的参数设置可能导致结果差异。
数据处理前的预处理：在进行傅里叶变换之前，对音频数据进行任何预处理或滤波操作时，确保处理步骤相同。不同的预处理步骤可能导致结果差异。
频谱分析的取样点数：傅里叶变换的结果取决于所选的取样点数。请确保在进行频谱分析时，使用相同数量的取样点。
傅里叶变化会放大他们的差异吧
这可能是因为这两个音频数据段中包含的频率成分不同。例如，一个数据段可能包含许多低频成分，而另一个数据段可能包含许多高频成分。在这种情况下，即使两个数据段在时域上看起来非常相似，它们在频域上的差异也会很大。
此外，傅立叶变换的结果还受到信号长度和采样率的影响。如果两个音频数据段的采样率或长度不同，那么它们在傅立叶变换后的结果也会不同。
在进行傅里叶变换前，需要对音频数据进行预处理。如果预处理不充分，例如存在噪声或者截止效应，那么得到的傅里叶变换结果可能会有所偏差，所以检查一下音频数据
检查一下音频数据的采样频率是否一致
如果两段基本一致的音频数据经过傅里叶变换后差异很大，可能是由于以下原因之一：
PCM文件可以以不同的精度进行编码，精度越高，文件的大小就越大，但音频的质量也越高。因此，使用不同的精度来编码声音会导致声音文件的差异变大。例如，16位PCM编码是一种常见的格式，可以提供较高的音频质量，但文件大小较大。而8位PCM编码则会大大减小文件大小，但由于采样精度的减少，音频质量也会相应下降。因此，当使用不同的PCM格式编码声音时，文件大小和音频质量之间会出现差异。
我这里用python尝试了一下：
import wave
import numpy as np
import matplotlib.pyplot as plt
with open('spk_p1Record.pcm', 'rb') as pcmfile:
pcmdata1 = pcmfile.read()
audio1 = np.frombuffer(pcmdata1, np.int8)
print(audio1)
输出：mean(audio1)=-3.8264。 如果要使用int16，同样输出mean(audio1)=-18.2192；如果使用int32，那mean(audio1)=33683879.184。
所以精度不同，会放大不同地方的值，下面是用python的傅里叶变换的图：
import wave
import numpy as np
import matplotlib.pyplot as plt
with open('spk_p1Record.pcm', 'rb') as pcmfile:
pcmdata1 = pcmfile.read()
audio1 = np.frombuffer(pcmdata1, np.int16)
with open('spk_p2Record.pcm', 'rb') as pcmfile:
pcmdata2 = pcmfile.read()
audio2 = np.frombuffer(pcmdata2, np.int16)
audio1_fft = np.abs(np.fft.ifft(np.fft.fftshift(np.fft.fft(audio1))))
audio2_fft = np.abs(np.fft.ifft(np.fft.fftshift(np.fft.fft(audio2))))
tt = audio1_fft - audio2_fft
print(np.mean(tt))
print(np.mean(audio1-audio2))
int8的值分别是：2.1145；-2.0613
int16的值分别是：1064.7189； -17.4793,
int32的值分别是：115403049.13780013；4462099.0158
究其原因还是声音编码的问题！！！
如果有帮助，请采纳，谢谢！
取样频率不同， FFT长度选择不一致，频谱计算方式不一致，频谱图显示参数设置，这几个方面排查
可参考
导致傅里叶变换结果差异很大的可能原因：
数据长度不匹配：输入数据长度为datalen，但是在傅里叶变换中，用datalen<<1表示数据长度。这可能导致数据截断或填充，从而影响傅里叶变换的结果。
傅里叶变换参数的选择：代码中使用的傅里叶变换函数cdft，四个参数分别为datalen<<1、-1、inputfft和ww。参数选择可能不正确，特别是负号的选择对于实数输入数据的处理非常重要。
输入数据的格式：输入数据被视为PCM音频的原始数据。但是，将16位无符号整数转换为双精度浮点数，并在傅里叶变换之前进行了填充。这种数据格式处理可能与实际音频数据的格式不匹配，导致结果差异。
傅里叶变换结果的使用：经过傅里叶变换的结果被转换为uint16类型，然后通过计算加权平均值获得三个范围内的数值。这种数据处理方式可能导致结果差异。
可以试试
确保输入音频数据的格式和采样率与傅里叶变换算法的要求相匹配。
检查傅里叶变换函数cdft的参数设置，确保正确选择参数。可能需要调整参数，以适应实数输入数据的处理。
检查数据的长度匹配。确保正确处理输入数据的长度，并在傅里叶变换中正确设置数据长度参数。
对比代码中的傅里叶变换结果处理和输出部分，确保结果的处理方式符合预期。