FFT的C预言实现

/*invert the sequence*/
#include
#include
#define NX 4   /* number of the input figure*/
#define Pi 3.1415926
void main()
{
int M,i,High,J=0,N1,K,NX1,y ;
float xr[NX],T,T1,xi[NX],XR[NX],XI[NX],SUM1,SUM2;
int B,j,P,B1=1;
NX1=NX;
for(i=0;NX1!=0;i++)
NX1=NX1/2;
M=i-1; /*work out the progression*/
/*printf("put the figure M in :\n");
scanf("%d",&M);
for(i=0;iprintf("the figure M is:%d\n",M);
High=NX/2;
/*J=High;*/
/*printf("the current inverse denary figure is :%d\n",J);*/
printf("input the  real part of the figures\n");
for(i=0;i<=NX-1;i++)
{ scanf("%f",&xr[i]);  /*input the real part*/
}
printf("input the image part:\n");
for(i=0;i<=NX-1;i++) xi[i]=0;/* let the image parts be zero*/
/*printf("output the sequence:\n");
for(i=0;i<=NX-1;i++)
printf("%f,",xr[i]); */
printf("\n");
for(i=0;i{
if(i { T=xr[J];xr[J]=xr[i];xr[i]=T;
   }
K=High;
while(J>=K)
{J=J-K;K=K/2;}
J=J+K;
}  /*over*/
printf("output the invert sequence:\n");
for(i=0;i<=NX-1;i++) /* the butterfly algotithm begin*/
printf("%f,",xr[i]);
printf("\n");
for(i=1;i<=M;i++)
{
    B1=1;
    for(j=0;j<=i-1;j++) B1=2*B1;
    printf("tne figure B1 is:%d\n ",B1);B=B1/2;
    for(J=0;J<=B-1;J++)

    { P=NX*J/B1;/*{for(y=0;y        P=P*J;*/
        printf("the figure P is:%d\n",P);

    printf("the valule of J iiiiiiis:%d\n",J);
    printf("cos %f    sin   %f\n",cos(2*Pi*P/NX),sin(2*Pi*P/NX));
   K=J;
   while(K
    { SUM1=cos(2*Pi*P/NX)*xr[K+B]+sin(2*Pi*P/NX)*xi[K+B];
      SUM2=sin(2*Pi*P/NX)*xr[K+B]-cos(2*Pi*P/NX)*xi[K+B];
      XR[K]=xr[K];
      XI[K]=xi[K];
      xr[K]=xr[K]+SUM1;
      printf("the figure %f de %d\n",xr[K],K);
      xi[K]=xi[K]-SUM2;
      printf("the figure %f de %d\n",xi[K],K);
      xr[K+B]=XR[K]-SUM1;
      printf("the figure %f DE %d\n",xr[K+B],K+B);
      xi[K+B]=XI[K]+SUM2;
      printf("the figure %f de %d\n",xi[K+B],K+B);
      K=K+B1;
    }
      }

}
printf("output the sequence after FFT:\n" );
for(i=0;i<=NX-1;i++)
{printf("%f+j %f\n",xr[i],xi[i]); }
getch();
}