【发布时间】:2013-10-09 16:13:01
【问题描述】:
我一直在这个 DFT 上撞墙。它应该打印出: 8,0,0,0,0,0,0,0 但我得到 8 然后是非常非常小的数字。这些是舍入错误吗?有什么我可以做的吗?我的 Radix2 FFT 给出了正确的结果,看起来愚蠢的 DFT 也不能工作。
我从复数开始,所以我知道有一点遗漏,我试图将其剥离以说明问题。
#include <cstdlib>
#include <math.h>
#include <iostream>
#include <complex>
#include <cassert>
#define SIZE 8
#define M_PI 3.14159265358979323846
void fft(const double src[], double dst[], const unsigned int n)
{
for(int i=0; i < SIZE; i++)
{
const double ph = -(2*M_PI) / n;
const int gid = i;
double res = 0.0f;
for (int k = 0; k < n; k++) {
double t = src[k];
const double val = ph * k * gid;
double cs = cos(val);
double sn = sin(val);
res += ((t * cs) - (t * sn));
int a = 1;
}
dst[i] = res;
std::cout << dst[i] << std::endl;
}
}
int main(void)
{
double array1[SIZE];
double array2[SIZE];
for(int i=0; i < SIZE; i++){
array1[i] = 1;
array2[i] = 0;
}
fft(array1, array2, SIZE);
return 666;
}
【问题讨论】:
-
小数字有多小?
-
@PaulR 在 1e-16 的顺序上
-
这听起来很正确 - 16 或 17 位小数大约是双精度的限制。请注意,您应该真正使用
中的 M_PI 而不是自己滚动,但我认为在这种情况下并没有太大区别。 -
是的,这真的很烦人。我可以使用 std::cout.setf ( std::ios::fixed, std::ios::floatfield ) 修复输出,但这只能修复显示。在大型数据集上,累积误差实际上加起来是有意义的。