我的FFT算法有什么问题? [英] what’ wrong with My FFT algorithm?
问题描述
大家好，最近我遇到了有关FFT算法的问题:我试图参考教科书在C语言中对该算法进行编程，但它不能作为matlab的fft函数使用.
首先，我在"Complex.h"中定义了用于复杂操作的数据类型.
Hello,everyone,recently I had met a problem about FFT algorithm:I tried to program the algorithm in C with reference to the textbook,but it didn''t work as the matlab'' fft function.
Firstly,I defined a data type for complex operation in "Complex.h"
typedef struct
{
double x_real;
double x_imag;
}Complex;
// the multiplication for complex x1 and x2
Complex CMul(Complex x1,Complex x2)
{
Complex result;
result.x_real=x1.x_real*x2.x_real-x1.x_imag*x2.x_imag;
result.x_imag=x1.x_real*x2.x_imag+x1.x_imag*x2.x_real;
return result;
}
//the addition for the complex x1 and x2
Complex CAdd(Complex x1,Complex x2)
{
Complex result;
result.x_real=x1.x_real+x2.x_real;
result.x_imag=x1.x_imag+x2.x_imag;
return result;
}
//the subtract for the complex x1 and x2
Complex CSub(Complex x1,Complex x2)
{
Complex result;
result.x_real=x1.x_real-x2.x_real;
result.x_imag=x1.x_imag-x2.x_imag;
return result;
}
然后我开始在"fft.cpp"中对FFT进行编程
Then I began to program the FFT in"fft.cpp"
#include <stdio.h>
#include <math.h>
#include "Complex.h"
#define pi 3.1415926
//bit-reverse operation for the input
void RSort(Complex x[],int length )
{
int LH,J,N1,I,K;
Complex T;
LH=length/2;
J=LH;
N1=length-2;
for (I=1;I<=N1;I++)
{
if (I<J)
{
T=x[I];
x[I]=x[J];
x[J]=T;
}
K=LH;
while(J>=K)
{
J=J-K;
K=K/2;
}
J=J+K;
}
}
/*
void FFT(Complex x[],int N ,int M)---compute the FFT of the origin signal
x[]--the input of origin and output of its FFT
N--the length of x[];
M--the power of N;
*/
void FFT(Complex x[],int N ,int M)
{
/*
L--the level of the Butterfly operation
B--the distance between two signal in every Butterfly operation
wn[]---the rotation factor
p--the exponent of rotation factor
*/
int L,B,J,p,k,i;
Complex wn[4];
//to get the rotation factor
for (i=0;i<N/2;i++)
{
wn[i].x_real=cos(2*pi*i/N);
wn[i].x_imag=-sin(2*pi*i/N);
}
//bit-reversal for the input signal
RSort(x,N);
//the main part of the FFT algorithm
for (L=1;L<=M;L++)
{
B=int(pow(2,L-1));
for (J=0;J<=B-1;J++)
{
p=int(J*pow(2,M-L));
for (k=J;k<=N-1;)
{
x[k]=CAdd(x[k],CMul(wn[p],x[k+B]));
x[k+B]=CSub(x[k],CMul(wn[p],x[k+B]));
k=k+int(pow(2,L));
}
}
}
}
void main()
{
Complex *x;
int i;
for (i=0;i<=15;i++)
{
x[i].x_real=i;
x[i].x_imag=0;
}
printf("the input signal：");
for (i=0;i<=15;i++)
{
printf("%lf ",x[i].x_real);
}
printf("\n");
FFT(x,16,4);
printf("the result after FFT：");
for (i=0;i<=15;i++)
{
printf("%lf i%lf ",x[i].x_real,x[i].x_imag);
}
printf("\n");
}
为方便起见，我将输入信号定义为信号x [] = [0,1,2,3,4,5,6,7].执行程序后的结果输出为:
For convenience,I defined the input signal to be a signal x[]=[0,1,2,3,4,5,6,7].The output of the result after executing program is :
输入信号:0.000000 1.000000 2.000000 3.000000 4.000000 5.000000 6.
000000 7.000000
FFT后的结果:28.000000 i0.000000 -1.414213 i-4.828427 4.000000 i-6.00
0000 -0.707107 i-0.707107 12.000000 i0.000000 0.000000 i-2.000000 4.000000 i
0.000000 0.000000 i0.000000
按任意键继续
the input signal：0.000000 1.000000 2.000000 3.000000 4.000000 5.000000 6.
000000 7.000000
the result after FFT：28.000000 i0.000000 -1.414213 i-4.828427 4.000000 i-6.00
0000 -0.707107 i-0.707107 12.000000 i0.000000 0.000000 i-2.000000 4.000000 i
0.000000 0.000000 i0.000000
Press any key to continue
而且在matlab中使用FFT函数的结果是:
And the result using FFT function in matlab is：
>> x = [0 1 2 3 4 5 6 7];
>> y = fft(x)
y =
28.0000 -4.0000 + 9.6569i -4.0000 + 4.0000i -4.0000 + 1.6569i -4.0000 -4.0000-1.6569i -4.0000-4.0000i -4.0000-9.6569i
>> x=[0 1 2 3 4 5 6 7];
>> y=fft(x)
y =
28.0000 -4.0000 + 9.6569i -4.0000 + 4.0000i -4.0000 + 1.6569i -4.0000 -4.0000 - 1.6569i -4.0000 - 4.0000i -4.0000 - 9.6569i
我非常期待您的意见，非常感谢！^ _ ^^ _ ^ ...
I''m looking forward to your opinions,Thanks very much!^_^^_^...
推荐答案
请注意:您的问题并不是真的要求专家意见或类似的东西；这确实是为您做大量工作的要求:测试，调试，检测错误并尝试修复它们.也许有人会为您找到一些东西，但是对于像我这样的许多人来说，从头开始实现该算法会容易得多.
所有这些实际上都超出了CodeProject快速问答的格式.答案.
所以我不知道还有什么帮助...
至少使用C ++复杂库，而不是基本的实现. FFT有许多可以使用的实现方式.阅读一些文章，例如:
如何实现FFT算法 [ http://www.librow.com/articles/article-10 [
Just a note: your question is not really a request for expert advice or something like that; this is really a request to do considerable amount of job for you: testing, debugging, detecting bugs and trying to fix them. Maybe someone will find something for you, but for many like myself it would be much easier to implement the algorithm from scratch.
All of this really goes beyond to format of CodeProject Quick Questions & Answers.
So I don''t know how else to help…
At least use C++ complex library instead you rudimentary implementation. FFT has many implementations you could use. Read some articles, such as:
How to implement the FFT algorithm[^],
http://www.librow.com/articles/article-10[^].
You can find a lot more.
—SA
作为逻辑错误的错误发生在FFT()函数中:
The errors which are logical errors occur in the function of FFT():
x[k]=CAdd(x[k],CMul(wn[p],x[k+B]));
x[k+B]=CSub(x[k],CMul(wn[p],x[k+B]));
这些代码应修改为(在定义Complex的这些变量T1和T2之前):
These codes should be modified to(before these variable T1 and T2 of Complex should be defined):
T2=CMul(wn[p],x[k+B]);
T1=x[k];
x[k]=CAdd(T1,T2);
x[k+B]=CSub(T1,T2);
然后，它工作正常！对此我感到很兴奋.现在，我首先要感谢SAKryukov的帮助！:-) :-)我已经开始对自己充满信心了！再次感谢.
该程序无法有效运行，仍然需要改进！让我们继续吧！
And then,it works correctly!I felt excited about that.Now I wanted to give thanks to SAKryukov for his help firstly!:-):-)I had began to gain strong confidence in myself! Thanks again.
This program didn''t work efficiently and it still needed to be improved well!Let''s go for it!
这篇关于我的FFT算法有什么问题?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
