/* VC 下运行通过 2007。4。26 把radix4中某些层for->while等效 14:51*//*****************fft programe*********************/void rad4(double x[],double y[],double TAB[],int i,int n2,int adr1,int adr2,int adr3){ int i1,i2,i3; double co1,co2,co3,si1,si2,si3; double r1,r2,r3,r4,s1,s2,s3,s4;     i1=i+n2;i2=i1+n2;i3=i2+n2;                       /* 数据位置 */  co1=TAB[2*adr1];co2=TAB[2*adr2];co3=TAB[2*adr3]; /* 旋转因子 */ si1=TAB[2*adr1+1];si2=TAB[2*adr2+1];si3=TAB[2*adr3+1];  r1=x[i]+x[i2];r3=x[i]-x[i2]; s1=y[i]+y[i2];s3=y[i]-y[i2]; r2=x[i1]+x[i3];r4=x[i1]-x[i3]; s2=y[i1]+y[i3];s4=y[i1]-y[i3];  x[i]=r1+r2; y[i]=s1+s2;   r2=r1-r2; r1=r3-s4;r3=r3+s4; s2=s1-s2; s1=s3+r4;s3=s3-r4; x[i1]=co1*r3+si1*s3;y[i1]=co1*s3-si1*r3; x[i2]=co2*r2+si2*s2;y[i2]=co2*s2-si2*r2; x[i3]=co3*r1+si3*s1;y[i3]=co3*s1-si3*r1;} void switchadr(double x[],double y[],int n) // 变换后,如果想顺序输出就调用这个{ double tempx,tempy; int i,j,k; for(i=0,j=0;i<=n-1;i++) {  if(i<j)  {   tempx=x[j];  tempy=y[j]; x[j]=x[i];y[j]=y[i]; x[i]=tempx;  y[i]=tempy;  }  k=n/4;     while(3*k<(j+1)&&(j>0) )  {   j=j-3*k; k=k/4;  }  j=j+k; }} void radix4(double x[],double y[],double TAB[],int n) // 核心变换代码{ int k,e,n1,n2; int a,b,c;    n2=n; e=1;    /* 1 -- n 有效 */ k=n; while(k>1) /*内部循环没有使用K值，只是用来控制循环次数而已,log4\n次*/  {  int j;  n1=n2; n2=n2/4;   for(j=1;j<=n2;j++)  {   int i;   a=(j-1)*e;b=a+a; c=a+b;   a &= (n-1);  b &= (n-1); c &= (n-1); /* 旋转因子周期 N */             i = (j-1);    do {    rad4(x,y,TAB,i,n2,a,b,c);    i+=n1;   }while(i<n);  }  e <<= 2; k >>= 2; }}

评论