i=x_n_re[ 1]; x_n_re[ 1]=x_n_re[128]; x_n_re[128]=i;
i=x_n_re[ 2]; x_n_re[ 2]=x_n_re[ 64]; x_n_re[ 64]=i;
i=x_n_re[ 3]; x_n_re[ 3]=x_n_re[192]; x_n_re[192]=i;
i=x_n_re[ 4]; x_n_re[ 4]=x_n_re[ 32]; x_n_re[ 32]=i;
...
i=x_n_re[207]; x_n_re[207]=x_n_re[243]; x_n_re[243]=i;
i=x_n_re[215]; x_n_re[215]=x_n_re[235]; x_n_re[235]=i;
i=x_n_re[223]; x_n_re[223]=x_n_re[251]; x_n_re[251]=i;
i=x_n_re[239]; x_n_re[239]=x_n_re[247]; x_n_re[247]=i;
Radix-2 FFT算法
采样按照位反转方式重新排序后就可进行FFT运算了。本radix-2 FFT应用的固件通过三个主循环执行图2所示的蝶型运算。外循环计数log2(N)级FFT运算。内循环执行每一级的蝶型运算。
FFT
算法的核心部分是执行蝶型运算的一小块代码。程序清单4给出了这一块代码，遗憾的是，它是本应用中唯一“不可移植”的固件。宏MUL_1和MUL_2利用μC的硬件乘法器执行单指令周期乘法运算。这些宏的内容专用于MAXQ2000，可在实际固件中全部看到。
清单4. 用C编写的蝶型运算。
/* (1) Macro MUL_1(A,B,C): C="A"*B (result in Q8.7)*/
/* (2) Macro MUL_2(A,C) : C="A"*last_B (result in Q8.7)*/
MUL_1(cosLUT[tf],x_n_re[b],resultMulReCos);
MUL_2(sinLUT[tf],resultMulReSin);
MUL_1(cosLUT[tf],x_n_im[b],resultMulImCos);
MUL_2(sinLUT[tf],resultMulImSin);
x_n_re[b] = x_n_re[a]-resultMulReCos+resultMulImSin;
x_n_im[b] = x_n_im[a]-resultMulReSin-resultMulImCos;
x_n_re[a] = x_n_re[a]+resultMulReCos-resultMulImSin;
x_n_im[a] = x_n_im[a]+resultMulReSin+resultMulImCos;
复数的极坐标转换
为了便于确定VIN频谱的幅度，我们须要将复数形式的X(k)转换为极坐标形式。实现该转换的固件示于程序清单5。幅度值取代了原始的FFT结果，因为固件不再需要这些数据。
清单5. FFT结果从复数形式转换为极坐标形式。
const unsigned char magnLUT[16][16] =
{
{0x00,0x10,0x20, ... ,0xd0,0xe0,0xf0},
{0x10,0x16,0x23, ... ,0xd0,0xe0,0xf0},
...
{0xe0,0xe0,0xe2, ... ,0xff,0xff,0xff},
{0xf0,0xf0,0xf2, ... ,0xff,0xff,0xff}
};
...
...
/* Compute x_n_re=abs(x_n_re) and x_n_im=abs(x_n_im) */
...
...
x_n_re[0] = magnLUT[x_n_re[0]>>11][0];
for(i=1; ix_n_re[i] = magnLUT[x_n_re[i]>>11][x_n_im[i]>>11];
x_n_re[N_DIV_2] = magnLUT[x_n_re[N_DIV_2]>>11][0];
频谱幅度并非根据式4计算得到，而是通过一个二维LUT查表得到。第一索引为频谱实部的高4位(MSB)，第二索引为频谱虚部的高4位。为得到这些数据，可将带符号的16位数据右移11次。在从频谱的实部和虚部取得索引号前，需首先将它们转换为绝对值。因此，符号位为零。