保存在文件中的双样本的 FFTW
FFTW of double samples saved in file
我正在尝试使用 FFTW 库计算 53k 双样本的 FFT,并在此基础上猜测信号的基频是多少。样本由sndfile库在wav输入文件的基础上生成(程序加载wav文件,生成double数据的样本并保存到文本文件)。每个样本的数量级范围从 -0.0009 到 +0.0009。使用下面介绍的函数计算 FFT 后,我收到 in[0][0] = -2.142。这是一种不正确的数据检索方法,或者输入文件不正确。我究竟做错了什么?传递给 FFTW 函数的数据是错误的,文件中存储的数据不正确还是我编写的函数不正确?缓冲区是一个包含样本的浮点数组。
static fftw_complex* calculateFourier(float* buffer, const unsigned int bufferLen)
{
unsigned int i;
//const unsigned short windowSize = 1024;
fftw_plan result;
fftw_complex *out, *in;
in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * bufferLen);
out = (fftw_complex*) fftw_malloc( sizeof(fftw_complex) * bufferLen);
double* doubleBuffer = castToDouble(buffer, bufferLen);
for(i = 0; i < bufferLen; i++)
{
in[i][0] = doubleBuffer[i];
in[i][1] = 0.000000000000;
}
result = fftw_plan_dft_1d(bufferLen, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
//result = fftw_plan_dft_r2c_1d(fileSize, castToDouble(buffer, bufferLen), out, FFTW_MEASURE);
fftw_execute(result);
fftw_destroy_plan(result);
return out;
}
static double* castToDouble(float* buffer, const unsigned int bufferLen)
{
unsigned int i;
double* doubleBuffer = (double*)malloc( sizeof(double) * bufferLen);
for(i = 0; i < bufferLen; i++)
doubleBuffer[i] = (double)buffer[i];
return doubleBuffer;
}
我认为问题在于您没有对原始数据应用 window function。因此,傅里叶变换具有很强的频率分量,对应于样本缓冲区的长度。
尝试通过替换此行来应用 Hanning window:
in[i][0] = doubleBuffer[i];
有了这个:
in[i][0] = doubleBuffer[i] * 0.5 * (1.0 - cos(2.0 * M_PI * i / (bufferLen - 1)));
(显然这个window函数如果需要重复使用可以缓存。)
好好看看下面这段代码有什么问题:
#include <stdio.h>
void prn (float *a, int n);
int main (void)
{
float floatarray[] = { 0.0, 5.0, 10.0, 15.0, 20.0, 25.0 };
prn (floatarray, sizeof floatarray/sizeof *floatarray);
return 0;
}
void prn (float *a, int n)
{
int i;
double *d = (double *)a;
for (i = 0; i < n; i++)
printf (" a[%d] : %01.0lf\n", i, d[i]);
}
运行它,发生了什么?为什么?
什么是sizeof
一个float
?什么是sizeof
一个double
?当演员 double *d = (double *)a;
发生时会发生什么? d
的每个元素需要多少字节?那里有多少?你的
double* doubleBuffer = castToDouble(buffer, bufferLen);
导致完全相同的问题。
逐个元素Cast/Assignment需要
无论您使用分配有malloc
(或calloc
)的新double
数组,都必须逐个元素 cast/assignment 完成数组的转换。这是相同的代码,但使用可变长度数组进行了一些添加以避免动态分配内存:
#include <stdio.h>
void castfloatdouble (double *d, float *f, int n);
void prn (float *a, int n);
void prndouble (double *d, int n);
int main (void)
{
float floatarray[] = { 0.0, 5.0, 10.0, 15.0, 20.0, 25.0 };
int n = (int)(sizeof floatarray/sizeof *floatarray);
double doublearray[n]; /* memset is a good idea to zero a VLA */
prn (floatarray, n); /* incorrect results due to typesize mismatch */
putchar ('\n');
castfloatdouble (doublearray, floatarray, n);
prndouble (doublearray, n); /* correct results after deep cast/assign */
putchar ('\n');
return 0;
}
void castfloatdouble (double *d, float *f, int n)
{
int i;
for (i = 0; i < n; i++)
d[i] = (double)f[i];
}
void prn (float *a, int n)
{
int i;
double *d = (double *)a;
for (i = 0; i < n; i++)
printf (" a[%d] : %01.0lf\n", i, d[i]);
}
void prndouble (double *d, int n)
{
int i;
for (i = 0; i < n; i++)
printf (" a[%d] : %01.0lf\n", i, d[i]);
}
输出显示不正确和正确的结果
$ ./bin/castfloatdouble2
a[0] : 2048
a[1] : 16777220
a[2] : 805306499
a[3] : 0
a[4] : 0
a[5] : 0
a[0] : 0
a[1] : 5
a[2] : 10
a[3] : 15
a[4] : 20
a[5] : 25
我正在尝试使用 FFTW 库计算 53k 双样本的 FFT,并在此基础上猜测信号的基频是多少。样本由sndfile库在wav输入文件的基础上生成(程序加载wav文件,生成double数据的样本并保存到文本文件)。每个样本的数量级范围从 -0.0009 到 +0.0009。使用下面介绍的函数计算 FFT 后,我收到 in[0][0] = -2.142。这是一种不正确的数据检索方法,或者输入文件不正确。我究竟做错了什么?传递给 FFTW 函数的数据是错误的,文件中存储的数据不正确还是我编写的函数不正确?缓冲区是一个包含样本的浮点数组。
static fftw_complex* calculateFourier(float* buffer, const unsigned int bufferLen)
{
unsigned int i;
//const unsigned short windowSize = 1024;
fftw_plan result;
fftw_complex *out, *in;
in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * bufferLen);
out = (fftw_complex*) fftw_malloc( sizeof(fftw_complex) * bufferLen);
double* doubleBuffer = castToDouble(buffer, bufferLen);
for(i = 0; i < bufferLen; i++)
{
in[i][0] = doubleBuffer[i];
in[i][1] = 0.000000000000;
}
result = fftw_plan_dft_1d(bufferLen, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
//result = fftw_plan_dft_r2c_1d(fileSize, castToDouble(buffer, bufferLen), out, FFTW_MEASURE);
fftw_execute(result);
fftw_destroy_plan(result);
return out;
}
static double* castToDouble(float* buffer, const unsigned int bufferLen)
{
unsigned int i;
double* doubleBuffer = (double*)malloc( sizeof(double) * bufferLen);
for(i = 0; i < bufferLen; i++)
doubleBuffer[i] = (double)buffer[i];
return doubleBuffer;
}
我认为问题在于您没有对原始数据应用 window function。因此,傅里叶变换具有很强的频率分量,对应于样本缓冲区的长度。
尝试通过替换此行来应用 Hanning window:
in[i][0] = doubleBuffer[i];
有了这个:
in[i][0] = doubleBuffer[i] * 0.5 * (1.0 - cos(2.0 * M_PI * i / (bufferLen - 1)));
(显然这个window函数如果需要重复使用可以缓存。)
好好看看下面这段代码有什么问题:
#include <stdio.h>
void prn (float *a, int n);
int main (void)
{
float floatarray[] = { 0.0, 5.0, 10.0, 15.0, 20.0, 25.0 };
prn (floatarray, sizeof floatarray/sizeof *floatarray);
return 0;
}
void prn (float *a, int n)
{
int i;
double *d = (double *)a;
for (i = 0; i < n; i++)
printf (" a[%d] : %01.0lf\n", i, d[i]);
}
运行它,发生了什么?为什么?
什么是sizeof
一个float
?什么是sizeof
一个double
?当演员 double *d = (double *)a;
发生时会发生什么? d
的每个元素需要多少字节?那里有多少?你的
double* doubleBuffer = castToDouble(buffer, bufferLen);
导致完全相同的问题。
逐个元素Cast/Assignment需要
无论您使用分配有malloc
(或calloc
)的新double
数组,都必须逐个元素 cast/assignment 完成数组的转换。这是相同的代码,但使用可变长度数组进行了一些添加以避免动态分配内存:
#include <stdio.h>
void castfloatdouble (double *d, float *f, int n);
void prn (float *a, int n);
void prndouble (double *d, int n);
int main (void)
{
float floatarray[] = { 0.0, 5.0, 10.0, 15.0, 20.0, 25.0 };
int n = (int)(sizeof floatarray/sizeof *floatarray);
double doublearray[n]; /* memset is a good idea to zero a VLA */
prn (floatarray, n); /* incorrect results due to typesize mismatch */
putchar ('\n');
castfloatdouble (doublearray, floatarray, n);
prndouble (doublearray, n); /* correct results after deep cast/assign */
putchar ('\n');
return 0;
}
void castfloatdouble (double *d, float *f, int n)
{
int i;
for (i = 0; i < n; i++)
d[i] = (double)f[i];
}
void prn (float *a, int n)
{
int i;
double *d = (double *)a;
for (i = 0; i < n; i++)
printf (" a[%d] : %01.0lf\n", i, d[i]);
}
void prndouble (double *d, int n)
{
int i;
for (i = 0; i < n; i++)
printf (" a[%d] : %01.0lf\n", i, d[i]);
}
输出显示不正确和正确的结果
$ ./bin/castfloatdouble2
a[0] : 2048
a[1] : 16777220
a[2] : 805306499
a[3] : 0
a[4] : 0
a[5] : 0
a[0] : 0
a[1] : 5
a[2] : 10
a[3] : 15
a[4] : 20
a[5] : 25