第8章 FFT的DSP实现(C55x DSP)

N2 N2 2N ( N 2)L ( N 2)lbN lbN
第8章 FFT的DSP实现
N越大，DIT-FFT运算量就减少得越多，FFT的优越性 就更加突出。例如，当N=256时，直接计算中复数乘法次数 为65 536，FFT算法中复数乘法次数为1024，速度提高倍数 为64。
第8章 FFT的DSP实现

n为偶数

2
nk x(n)WN
n为奇数
x (n ) W
nk N


n 0
N
2 x (2n ) WNnk


n 0
N
2
( x (2n 1) WN2 n 1) k


n 0
N
2
2 x 1 (n ) WNnk

k WN

n 0
N
2
2 x 2 (n ) WNnk
第8章 FFT的DSP实现 因为
2 j 2nk 2 WNnk e N j 2 nk N 2
e
nk WN
2
所以
X( k )

n 0
N
2
nk k x 1 (n)WN WN 2

n 0
N
2
nk x 2 (n)WN
k X1 (k) WN X 2 (k)
(8.6)
2
第8章 FFT的DSP实现
式中，X1(k)与X2(k)分别是x1(n)与x2(n)的N/2点DFT。上 式表明，一个N点DFT可以分解为两个N/2点DFT，此时需要
2. 编译工程并运行程序
在编译工程过程中，如有错误，修改错误，编译成功后， 装载输出文件。
3. 设置断点观察窗口
在程序中的break point处设置断点，共有三个断点，如图8-3 所示。
第8章 FFT的DSP实现
图8-3 在程序中的break point处设置断点
第8章 FFT的DSP实现
4. 设置波形观察窗口 1) 查看输入信号波形
2. FFT流程图 FFT流程图如图8-2所示。
第8章 FFT的DSP实现
图8-2 FFT流程图
第8章 FFT的DSP实现
8.1.5 实验步骤 1. 新建立工程、源文件和命令文件 分别在源文件和命令文件中，输入本章8.1.6小节程序清 单中FFT的C语言源程序和命令文件，并把这两个文件添加 到工程。
2(N/2)2+N/2≈N2/2次复数乘法，N(N/2 -1)+N=N2/2次复数加
法，可见，经过一次分解运算量就减少了接近一半。由于 N/2依然是偶数，故可将两个N/2点的DFT按同样方法分解成 四个N/4点的DFT，四个N/4点的DFT继续分解为八个N/8点 的DFT，如此进行下去，经过L-1次分解后，就把一个N点
第8章 FFT的DSP实现
图8-6 图形属性对话框
第8章 FFT的DSP实现
图8-7 CCS工具的FFT变换结果
第8章 FFT的DSP实现
3) 使用C语言程序，计算FFT变换结果 在主菜单中选择View→Graph→Time/Frequency命令， 出现如图8-8所示的图形属性对话框，按照图8-8所示进行相 应属性修改，修改好后，单击OK确认。将程序运行到第三 个断点处，可以看到如图8-9所示的C语言程序计算FFT变换
试等操作，使用属性窗口设置项目及观察输出波形等。
第8章 FFT的DSP实现
8.1.4 实验原理 1. FFT基本原理 快速傅立叶变换(FFT)是离散傅立叶变换(DFT)的一种高 效运算方法。对于有限长离散数字信号x(n)，0≤n≤N -1，其 离散傅立叶变换公式为
X( k )

n 0
N 1
第8章 FFT的DSP实现
第8章 FFT的DSP实现
8.1 利用标准C语言实现FFT实验
8.2 利用DSP库函数实现FFT实验
第8章 FFT的DSP实现
8.1 利用标准C语言实现FFT实验
8.1.1 实验目的 ■了解FFT基本原理。 ■学习c55x DSP中用标准C实现FFT的设计和编程思想。
第8章 FFT的DSP实现
法(DIT-FFT)，另一类是按频率抽取的FFT算法(DIF-FFT)。
这里只介绍DIT-FFT，至于DIF-FFT原理类似。 设有限长序列x(n)的长度为N=2L，L为整数，显然N为 偶数，若不满足该条件，加零补足。此时定义两个x(n)的偶 数项和奇数项序列，长度均为N/2，即：
第8章 FFT的DSP实现
第8章 FFT的DSP实现
init_fft_tab(); input_data(); for ( i=0;i<sample_l;i++ ) { fWaveR[i]=input[i]; //break point
fWaveI[i]=0.0f;
w[i]=0.0f; }
fft(fWaveR,fWaveI);
k N WN WN k
(8.2)
2
W WN
k N
kN
(8.3)
第8章 FFT的DSP实现
利用周期性与对称性，一方面可以在DFT的运算过程中 合并某些项，另一方面可以把长序列的DFT分解成若干个短
序列的DFT。由于DFT运算量与变换长度的平方成正比，因
此分解成短序列的DFT可以大大减少运算量。 常用的FFT算法有两大类，一类是按时间抽取的FFT算
第8章 FFT的DSP实现
{
x0=x1=x2=x3=x4=x5=x6=0;
x0=i&0x01; x1=(i/2)&0x01; x2=(i/4)&0x01; x3=(i/8)&0x01; x4=(i/16)&0x01; x5=(i/32)&0x01; x6=(i/64)&0x01;x7=(i/128)&0x01;
8.1.2 实验设备 ■PC兼容机一台；操作系统为Windows 2000(或 Windows NT、Windows 98、Windows XP)。
■计算机安装CCS 5000或CCS v 3.1。 8.1.3 实验要求 使用CCS集成仿真环境，完成建立工程、源文件、命令
文件，保存和添加文件到工程，进行编译、链接、运行和调
nk x (n )WN
k=0,1,2,…,N-1
(8.1)
第8章 FFT的DSP实现
式中
WN e
j 2 N
称为旋转因子或蝶形因子。由上式可以看
出，计算N点DFT的所有X(k)值需要N2次复数乘法和N(N -1)
次复数加法。当N较大时，N(N-1)≈N2，计算X(k)的运算量
几乎与N2成正比，因此直接计算DFT的运算量很大，即使采 用计算机也很难做到实时处理，必须加以改进。 FFT主要利用DFT旋转因子的周期性与对称性来减少运 算量，周期性和对称性分别为
图8-5 输入信号波形
第8章 FFT的DSP实现
2) 使用CCS提供的工具，观察输入信号FFT变换结果 在主菜单中选择View→Graph→Time/Frequency命令，出现 如图8-6所示的图形属性对话框，按照图8-6所示进行相应属 性修改，修改好后，单击OK确认。可以看到使用CCS提供 的工具，对输入信号进行FFT变换的结果如图8-7所示。
在主菜单中选择View→Graph→Time/Frequency命令，
出现如图8-4所示的图形属性对话框，按照图8-4所示进行相 应属性修改，修改好后，单击OK确认。将程序运行到第二 个断点处，得到输入信号波形，如图8-5所示。
第8章 FFT的DSP实现
图8-4 图形属性对话框
第8章 FFT的DSP实现
DFT分解成N/2个2点的DFT，至此分解结束。图8-1所示为
一个完整的N=8点DIT-FFT运算蝶形图。
第8章 FFT的DSP实现
图8-1 N=8点DIT-FFT运算蝶形图
第8章 FFT的DSP实现
DIT-FFT运算量为：当N=2L时，经过L -1次分解，整 个运算过程有L级蝶形，每一级蝶形有N/2个蝶形运算，每
第8章 FFT的DSP实现
float fWaveR[sample_l],fWaveI[sample_l],w[sample_l]; float sin_tab[sample_l];
float cos_tab[sample_l];
void init_fft_tab(); void input_data(); void fft(float dataR[sample_l],float dataI[sample_l]); void main() { int i;
xx=x0*128+x1*64+x2*32+x3*16+x4*8+x5*4+x6*2+x7;
dataI[xx]=dataR[i]; }
第8章 FFT的DSP实现
for ( i=0;i<sample_l;i++ ) {
dataR[i]=dataI[i]; dataI[i]=0;
} /************** following code FFT ******************************************/ for ( L=1;L<=8;L++ )
第8章 FFT的DSP实现
8.1.6 FFT的标准C语言参考程序清单 1. FFT.c文件
#include <math.h>
#define pi 3.1415926 #define sample_l 256 #define signal_1_f 60; #define signal_2_f 200; #define signal_sample_f 512; int input[sample_l];
while(1);
//break point
//break point
第8章 FFT的DSP实现
} void init_fft_tab() { float wt1; float wt2;
int i;
for (i=0;i<sample_l;i++) {
第8章 FFT的DSP实现
wt1=2*pi*i*signal_1_f; wt1=wt1/signal_sample_f;
x 1 ( n ) x ( 2n )
(8.4) N n 0,1,2, , 1 2
N n 0,1,2, , 1 (8.5) x 2 (n ) x (2n 1) 2
第8章 FFT的DSP实现
则x(n)的N点序列DFT可写为
X( k )

n 0
N 1
nk x (n ) WN
一个蝶形运算有一次复数乘法和两次复数加法，所以整个运
算过程的运算量如下：
复数乘法
N N mF L 1 lbN 2 2
F L N 2 N lbN 2
复数加法
第8章 FFT的DSP实现
直接计算DFT需要N2次复数乘法，N(N-1)次复数加法， 直接计算DFT与DIT-FFT复数乘法的运算量之比为：
结果。
第8章 FFT的DSP实现
图8-8 图形属性对话框
wenku.baidu.com
第8章 FFT的DSP实现
图8-9 C语言程序计算FFT变换结果
第8章 FFT的DSP实现
对使用CCS工具的FFT变换结果和使用C语言程序计算 的FFT变换结果进行比较分析，可以看到两个结果是一致的。 读者可以修改两个信号的频率，重新编译、链接、装载 输出文件和运行程序，比较使用CCS工具的FFT变换结果和 使用C语言程序计算的FFT变换结果，看两个结果是否 一致？
wt2=2*pi*i*signal_2_f;
wt2=wt2/signal_sample_f; input[i]=(cos(wt1)+cos(wt2))/2*32768; } }
第8章 FFT的DSP实现
void input_data() {
int i;
for(i=0;i<sample_l;i++) { sin_tab[i]=sin(2*pi*i/sample_l); cos_tab[i]=cos(2*pi*i/sample_l); } }
{
/* for(1) */
第8章 FFT的DSP实现
b=1; i=L-1; while ( i>0 ) { b=b*2; i--; } /* b= 2^(L-1) */
第8章 FFT的DSP实现
void fft(float dataR[sample_l],float dataI[sample_l]) {
int x0,x1,x2,x3,x4,x5,x6,x7,xx;
int i,j,k,b,p,L; float TR,TI,temp; /********** following code invert sequence ***********************************/ for ( i=0;i<sample_l;i++ )