将Rosseta代码FFT实现到VBA Excel中
我试图将FFT Rosetta代码实现到VBA excel中。 我无法完全重build与Rosetta代码完全相同的输出数据。 起初,我认为这是types转换不匹配,但通过1.1缩放input值也将输出值缩放1.1。 我认为可能存在问题的唯一方面是我如何在代码中实现引用的数组,而不是Rosetta写入的。 我发现Rosetta通过在它的奇recursion调用中写出+ step和buf + step来移位数组的地址。 我不了解VBA中的类似构造,所以我简单地传递了一个额外的recursion参数移位,它跟踪地址在传递到下一个recursion调用时的移位。
我的class次执行有什么问题,还是其他的?
Rosetta声称它的FFT是内存,但是他们做了额外的数据拷贝并将其存储在out []中。 我将不胜感激解释为什么这是如此。
FFT中的Rosetta代码
void _fft(cplx buf[], cplx out[], int n, int step) { if (step < n) { _fft(out, buf, n, step * 2); _fft(out + step, buf + step, n, step * 2); for (int i = 0; i < n; i += 2 * step) { cplx t = cexp(-I * PI * i / n) * out[i + step]; buf[i / 2] = out[i] + t; buf[(i + n)/2] = out[i] - t; } } } void fft(cplx buf[], int n) { cplx out[n]; for (int i = 0; i < n; i++) out[i] = buf[i]; _fft(buf, out, n, 1); }
我在VBA中尝试过FFT
Private Sub rec_fft(ByRef buf, ByRef out, ByVal n, ByVal step, ByVal shift) If step < n Then Dim ii As Long Dim pi As Double pi = 3.14159265358979 + 3.23846264338328E-15 Dim c1(1 To 2) As Long Dim c2(1 To 2) As Double Call rec_fft(out, buf, n, step * 2, shift) Call rec_fft(out, buf, n, step * 2, shift + step) For i = 1 To n / (2 * step) ii = 2 * step * (i - 1) c1(1) = Cos(-1 * pi * CDbl(ii) / CDbl(n)) c1(2) = Sin(-1 * pi * CDbl(ii) / CDbl(n)) c2(1) = c1(1) * out(shift + 1 + ii + step, 1) - c1(2) * out(shift + 1 + ii + step, 2) c2(2) = c1(1) * out(shift + 1 + ii + step, 2) + c1(2) * out(shift + 1 + ii + step, 1) buf(shift + 1 + ii / 2, 1) = out(shift + 1 + ii, 1) + c2(1) buf(shift + 1 + ii / 2, 2) = out(shift + 1 + ii, 2) + c2(2) buf(shift + 1 + (ii + n) / 2, 1) = out(shift + 1 + ii, 1) - c2(1) buf(shift + 1 + (ii + n) / 2, 2) = out(shift + 1 + ii, 2) - c2(2) Next i End If End Sub Private Sub fft(ByRef buf, ByVal n) Dim out() As Double ReDim out(1 To n, 1 To 2) For i = 1 To n out(i, 1) = buf(i, 1) out(i, 2) = buf(i, 2) Next i Call rec_fft(buf, out, n, 1, 0) End Sub
输出比较
Input Data: 1 1 1 1 0 0 0 0 Rosetta FFT : 4 (1, -2.41421) 0 (1, -0.414214) 0 (1, 0.414214) 0 (1, 2.41421) My FFt : 4 (1, -3) 0 (1, -1) 0 (1, 1) 0 (1, 3)
您已将c1声明为Long
:
"Dim c1(1 To 2) As Long"
将其更改为Double
,它的工作原理是:
"Dim c1(1 To 2) As Double"