雷德算法

Posted 2020-08-22 prayer521

在实现FFT计算的时候，第一步要做的就是实现倒位序的实现，这里有一种算法，叫做雷德（Rader）算法。

下面假如使用A[I]存的是顺序位序，而B[J]存的是倒位序。I<J的时候需要变序，I>J的时候就不用，不然就白忙活了。

例如   N = 8 的时候，
倒位序 顺序          二进制表示          倒位序         顺序

0  0                                       000          000

4  1                                       100          001

2  2                                       010          010           

6  3                                       110          011

1  4                                       001          100

5  5                                       101          101

3  6                                       011          110

7  7                                       111          111

由上面的表可以看出，按自然顺序排列的二进制数，其下面一个数总是比其上面一个数大1，即下面一个数是上面一个数在最低位加1并向高位进位而得到的。而倒位序二进制数的下面一个数是上面一个数在最高位加1并由高位向低位进位而得到。

I、J都是从0开始，若已知某个倒位序J，要求下一个倒位序数，则应先判断J的最高位是否为0，这可与k=N/2相比较，因为N/2总是等于100..的。如果k>J，则J的最高位为0，只要把该位变为1（J与k=N/2相加即可），就得到下一个倒位序数；如果K<=J，则J的最高位为1，可将最高位变为0（J与k=N/2相减即可）。然后还需判断次高位，这可与k=N\4相比较，若次高位为0，则需将它变为1（加N\4即可）其他位不变，既得到下一个倒位序数；若次高位是1，则需将它也变为0。然后再判断下一位。。。。

下面用C代码实现雷德算法。

#include <stdio.h>

double x[16] = {0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15};

double y[16] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};

int N = 16;

void main()

{

int i,j,k;

double TempR,TempI;

 

int c = 0,m = 0;

 

for(j=0,i=0;i<N-1;i++)    //这里实现了奇偶前后分开排序

{

   

 printf("i = %d \n",i);

 printf("j = %d \n",j);  

 printf("k = %d \n",k);

 if(i<j)                        //如果i<j,即进行变址

 {

  TempR = x[j];

  TempI = y[j];

  x[j]  = x[i];

  y[j]  = y[i];

  x[i]  = TempR;

  y[i]  = TempI;

 }

 k = N/2;                 //求j的下一个倒位序

 while(k<(j+1))        //如果k<=j,表示j的最高位为1 

 {

  j = j-k;                 //把最高位变成0

  k = k/2;               //k/2，比较次高位，依次类推，逐个比较，直到某个位为0

 }

 j = j+k;                //把0改为1

 for(m=0;m<N/2;m++)

 {

  for(c=0;c<2;c++)                                        //一行打印两个数

   printf("        %10.7f+J %10.7f",x[2*m+c],y[2*m+c]);

  printf("\n");

 }

 printf("\n");

}

}