快速傅立叶变换(FFT)的实现
从傅里叶变换到快速傅里叶变换的基本实现方法

从傅里叶变换到快速傅里叶变换的基本实现方法(原创实用版4篇)目录(篇1)I.傅里叶变换的概念和意义1.傅里叶变换是一种将时域信号转换为频域信号的数学变换方法2.在信号处理、图像处理、通信等领域有着广泛的应用II.快速傅里叶变换(FFT)的基本原理1.傅里叶变换的乘法运算导致计算效率低下2.快速傅里叶变换利用了周期函数的周期性性质,将乘法运算转化为加法运算3.FFT的基本算法思想:基于递归的方式,将大的傅里叶变换问题分解为更小的子问题III.FFT的具体实现方法1.迭代实现方法:主要用于离散傅里叶变换(DFT)的实现2.迭代实现方法的优化:使用蝶形图表示FFT的运算过程,便于理解和计算3.直接实现方法:对于特定的离散序列,可以直接计算其FFT结果,不需要进行迭代正文(篇1)一、傅里叶变换的概念和意义傅里叶变换是一种将时域信号转换为频域信号的数学变换方法。
它可以将一个时域信号表示为不同频率的正弦和余弦函数的线性组合,使得信号的频域分析变得更加方便。
傅里叶变换在信号处理、图像处理、通信等领域有着广泛的应用。
二、快速傅里叶变换(FFT)的基本原理傅里叶变换的乘法运算导致计算效率低下,快速傅里叶变换(FFT)利用了周期函数的周期性性质,将乘法运算转化为加法运算。
FFT的基本算法思想是:基于递归的方式,将大的傅里叶变换问题分解为更小的子问题。
FFT算法可以分为迭代实现方法和直接实现方法,其中迭代实现方法主要用于离散傅里叶变换(DFT)的实现。
三、FFT的具体实现方法1.迭代实现方法:迭代实现方法的主要思想是将大的傅里叶变换问题分解为更小的子问题,通过递归的方式逐步求解。
迭代实现方法可以使用蝶形图表示FFT的运算过程,便于理解和计算。
2.迭代实现方法的优化:迭代实现方法的优化主要是为了减少计算量,例如使用树形结构来存储中间结果,减少重复计算。
3.直接实现方法:对于特定的离散序列,可以直接计算其FFT结果,不需要进行迭代。
快速傅里叶变换fft的c程序代码实现

快速傅里叶变换fft的c程序代码实现标题:一种高效实现快速傅里叶变换(FFT)的C语言程序代码导言:快速傅里叶变换(Fast Fourier Transform,FFT)是一种在信号处理、图像处理、通信系统等领域广泛应用的重要算法。
它通过将输入信号从时域转换到频域,实现了对信号的频谱分析和频率成分提取。
在实际应用中,为了获得高效的FFT计算过程,我们需要使用合适的算法和优化技巧,并将其转化为高质量的C语言代码。
本文将介绍一种基于Cooley-Tukey算法的快速傅里叶变换的C语言程序代码实现。
我们将从原理开始详细讲解FFT算法,然后逐步引入代码实现的步骤,并进行相关优化。
我们将总结整个实现过程,并分享一些个人对FFT算法的理解和观点。
一、快速傅里叶变换(FFT)的原理(1)傅里叶级数与离散傅里叶变换傅里叶级数是将一个周期函数分解为一系列正弦和余弦函数的和的方法。
然而,实际数字信号往往是离散的。
我们需要离散傅里叶变换(Discrete Fourier Transform,DFT)来对离散信号进行频谱分析。
(2)DFT的定义及其计算复杂度离散傅里叶变换通过对离散信号的变换矩阵进行乘法运算,得到其频谱表示。
然而,直接使用定义式计算DFT的时间复杂度为O(N^2),其中N为信号长度,这对于大规模信号计算是不可接受的。
(3)引入快速傅里叶变换 (FFT)Cooley-Tukey算法是一种最常用的FFT算法,通过将DFT分解为多个较小规模的DFT计算来降低计算复杂度。
FFT的时间复杂度为O(NlogN),大大提高了计算效率。
二、快速傅里叶变换(FFT)的C语言实现(1)算法流程和数据结构设计以一维FFT为例,我们需要定义合适的数据结构来表示复数和存储输入输出信号,同时设计实现FFT的主要流程。
(2)递归实现方法递归实现是最直观的FFT实现方法,基于Cooley-Tukey算法的思想。
将输入信号分为偶数和奇数两部分,然后递归计算它们的FFT。
快速傅里叶变换FFT的C语言实现及应用

快速傅里叶变换FFT的C语言实现及应用快速傅里叶变换(Fast Fourier Transform, FFT)是一种快速计算离散傅里叶变换(Discrete Fourier Transform, DFT)的算法。
它能够在较短时间内计算出巨大数据集的傅里叶变换,广泛应用于信号处理、图像处理、通信等领域。
C语言是一种广泛应用于嵌入式系统和科学计算的编程语言,拥有高效、灵活和可移植等特点。
下面将介绍FFT的C语言实现及应用。
首先,我们需要了解离散傅里叶变换的定义。
离散傅里叶变换将一组离散的时域信号转换成一组对应的频域信号,可以表示为以下公式:X(k) = Σ[ x(n) * W^(kn) ]其中,X(k)是频域信号,x(n)是时域信号,W是单位复数旋转因子,其计算公式为W=e^(-j*2π/N),其中j是虚数单位,N是信号的长度。
实现FFT算法的关键在于计算旋转因子的值,一种常用的计算方法是采用蝶形算法(butterfly algorithm)。
蝶形算法将DFT分解成多个子问题的求解,通过递归调用可以快速计算出结果。
以下是一种基于蝶形算法的FFT实现的示例代码:```c#include <stdio.h>#include <math.h>typedef structfloat real;float imag;if (N <= 1)return;}for (int i = 0; i < N/2; i++)even[i] = signal[2*i];odd[i] = signal[2*i + 1];}fft(even, N/2);fft(odd, N/2);for (int k = 0; k < N/2; k++)signal[k].real = even[k].real + temp.real;signal[k].imag = even[k].imag + temp.imag;signal[k + N/2].real = even[k].real - temp.real; signal[k + N/2].imag = even[k].imag - temp.imag; }int maiint N = sizeof(signal) / sizeof(signal[0]);fft(signal, N);printf("频域信号:\n");for (int i = 0; i < N; i++)printf("%f + %fi\n", signal[i].real, signal[i].imag);}return 0;```以上代码实现了一个简单的4点FFT算法,输入时域信号为{1,0,1,0},输出为对应的频域信号。
(完整word版)基于DSP的快速傅立叶变换(FFT)的实现(汇编语言)

快速傅立叶变换(FFT )的实现一、实验目的1.了解FFT 的原理及算法;2.了解DSP 中FFT 的设计及编程方法;3.熟悉FFT 的调试方法;二、实验原理FFT 是一种高效实现离散付立叶变换的算法,把信号从时域变换到频域,在频域分析处理信息。
对于长度为N 的有限长序列x (n ),它的离散傅里叶变换为:(2/)j N nk N W e π-=,称为旋转因子,或蝶形因子。
在x (n )为复数序列的情况下,计算X (k ):对某个k 值,需要N 次复数乘法、(N -1)次复数加法;对所有N 个k 值,需要2N 次复数乘法和N (N -1)次复数加法。
对于N 相当大时(如1024)来说,直接计算它的DFT 所作的计算量是很大的,FFT 的基本思想在于: 利用2()j nk N N W e π-=的周期性即:k N k N N W W +=对称性:/2k k N N N W W +=-将原有的N 点序列分成两个较短的序列,这些序列的DFT 可以很简单的组合起来得到原序列的DFT 。
按时间抽取的FFT ——DIT FFT 信号流图如图5.1所示:图5.1 时间抽取的FFT —DIT FFT 信号流图FFT 算法主要分为以下四步。
第一步 输入数据的组合和位倒序∑=-=10)()(N n nk N W n x k X把输入序列作位倒序是为了在整个运算最后的输出中得到的序列是自然顺序。
第二步 实现N 点复数FFT第一级蝶形运算;第二级蝶形运算;第三级至log2N 级蝶形运算;FFT 运算中的旋转因子N W 是一个复数,可表示:为了实现旋转因子N W 的运算,在存储空间分别建立正弦表和余弦表,每个表对应从0度到180度,采用循环寻址来对正弦表和余弦表进行寻址。
第三步 功率谱的计算X (k )是由实部()R X k 和虚部()I X k 组成的复数:()()()R I X k X k jX k =+;计算功率谱时只需将FFT 变换好的数据,按照实部()R X k 和虚部()I X k 求它们的平方和,然后对平方和进行开平方运算。
编程实现快速傅里叶变换(fft)

一、概述傅里叶变换是信号处理和数据压缩中常用的数学工具,它可以将时域信号转换为频域信号,从而便于分析和处理。
而快速傅里叶变换(FFT)则是一种高效的计算傅里叶变换的方法,可以大大提高计算效率,广泛应用于信号处理、图像处理、通信系统等领域。
二、傅里叶变换原理傅里叶变换的基本思想是将一个时域信号分解为不同频率的正弦和余弦函数的叠加,从而得到该信号的频谱图。
具体来说,对于一个连续信号x(t),它的傅里叶变换X(ω)定义为:X(ω) = ∫[0,∞]x(t)e^(-jωt)dt其中,ω为频率变量,X(ω)表示在频率ω处的信号能量。
而对于离散信号x[n],它的傅里叶变换X[k]则定义为:X[k] = ∑[n=0,N-1]x[n]e^(-j2πkn/N)其中,N为信号的采样点数,k为频率域的序号。
上述公式称为离散傅里叶变换(DFT),计算复杂度为O(N^2)。
而快速傅里叶变换则通过巧妙的算法设计,将计算复杂度降低到O(NlogN)。
三、快速傅里叶变换算法概述快速傅里叶变换的算法最早由Cooley和Tukey在1965年提出,它的基本思想是将一个长度为N的DFT分解为两个长度为N/2的DFT的组合,通过递归地分解和合并,最终实现对整个信号的快速计算。
下面我们来介绍一种常用的快速傅里叶变换算法:递归式分治法。
四、递归式分治法递归式分治法是一种高效的计算DFT的方法,它的基本思想是将长度为N的DFT分解为两个长度为N/2的DFT,并通过递归地调用自身,最终实现对整个信号的傅里叶变换。
具体来说,假设有一个长度为N的信号x[n],对其进行快速傅里叶变换的过程可以分为如下几个步骤:1. 将长度为N的信号x[n]分为长度为N/2的偶数序号和奇数序号的两个子信号x_even[n]和x_odd[n];2. 对子信号x_even[n]和x_odd[n]分别进行快速傅里叶变换,得到它们的频域表示X_even[k]和X_odd[k];3. 结合X_even[k]和X_odd[k],计算原信号的频域表示X[k]。
快速傅里叶变换FFT的FPGA设计与实现--电科1704 郭衡

快速傅里叶变换FFT的FPGA设计与实现学生姓名郭衡班级电科1704学号17419002064指导教师谭会生成绩2020年5 月20 日快速傅里叶变换FFT 的设计与实现一、研究项目概述非周期性连续时间信号x(t)的傅里叶变换可以表示为:=)(ϖX dt tj et x ⎰∞∞--1)(ϖ,式中计算出来的是信号x(t)的连续频谱。
但是,在实际的控制系统中能够式中计算出来的是信号x(t)的连续频谱。
但是,在实际的控制系统中能够算信号x(t)的频谱。
有限长离散信号x(n),n=0,1,…,N-1的DFT 定义为:∑-=-=-==102,1.....10)()(N n Nj N knNeW N k W n x K X π、、。
可以看出,DFT 需要计算大约N2次乘法和N2次加法。
当N 较大时,这个计算量是很大的。
利用WN 的对称性和周期性,将N 点DFT 分解为两个N /2点的DFT ,这样两个N /2点DFT 总的计算量只是原来的一半,即(N /2)2+(N /2)2=N2/2,这样可以继续分解下去,将N /2再分解为N /4点DFT 等。
对于N=2m 点的DFT 都可以分解为2点的DFT ,这样其计算量可以减少为(N /2)log2N 次乘法和Nlog2N 次加法。
图1为FFT 与DFT-所需运算量与计算点数的关系曲线。
由图可以明显看出FFT 算法的优越性。
图1 FFT 与DFT 所需乘法次数比较X[1]将x(n)分解为偶数与奇数的两个序列之和,即x(n)=x1(n)+x2(n)。
x1(n)和x2(n)的长度都是N /2,x1(n)是偶数序列,x2(n)是奇数序列,则∑∑=--=-=+2)12(1202)1.....,0()(2)(1)(N n kn N N n km N N k W n x W n x K X所以)1...,0()(2)(1)(1222120-=+=∑∑-=-=N k W n x W W n x K X N n km N k N km N Nn由于kmN N jkm Njkm NW eeW2/2/2222===--ππ,则)1.....,0)((2)(1)(2)(1)(122/1202/-=+=+=∑∑-=-=N k k X W k X W n x W W n x K X kN N n km N k N Nn kn N其中X1(k)和X2(k)分别为x1(n)和x2(n)的N /2点DFT 。
快速傅里叶变换FFT算法源码经典

快速傅里叶变换FFT算法源码经典以下是一个经典的快速傅里叶变换(FFT)算法的源码,包含详细的注释解释每个步骤的作用。
```pythonimport cmath#递归实现快速傅里叶变换def fft(x):N = len(x)#基本情况:如果输入向量只有一个元素,则直接返回该向量if N <= 1:return x#递归步骤:#将输入向量分成两半even = fft(x[0::2]) # 偶数索引的元素odd = fft(x[1::2]) # 奇数索引的元素T = [cmath.exp(-2j * cmath.pi * k / N) * odd[k] for k in range(N // 2)]#组合结果return [even[k] + T[k] for k in range(N // 2)] + \[even[k] - T[k] for k in range(N // 2)]#逆傅里叶变换def ifft(X):N = len(X)#将输入向量取共轭X_conj = [x.conjugate( for x in X]#应用快速傅里叶变换x_conj = fft(X_conj)#将结果取共轭并归一化return [(x.conjugate( / N).real for x in x_conj]#示例测试if __name__ == "__main__":x=[1,2,3,4]X = fft(x)print("快速傅里叶变换结果:", X)print("逆傅里叶变换恢复原始向量:", ifft(X))```这个源码实现了一个经典的快速傅里叶变换(FFT)算法。
首先,`fft`函数实现了递归的快速傅里叶变换,接收一个输入向量`x`作为参数,返回傅里叶变换后的结果`X`。
如果输入向量只有一个元素,则直接返回。
否则,将输入向量分成两半,分别对偶数索引和奇数索引的元素递归应用FFT。
[C++]频谱图中FFT快速傅里叶变换C++实现
![[C++]频谱图中FFT快速傅里叶变换C++实现](https://img.taocdn.com/s3/m/f1be4c19773231126edb6f1aff00bed5b9f373d0.png)
[C++]频谱图中FFT快速傅⾥叶变换C++实现在项⽬中,需要画波形频谱图,因此进⾏查找,不是很懂相关知识,下列代码主要是针对这篇⽂章。
//快速傅⾥叶变换/*⼊⼝参数:inv: =1,傅⾥叶变换; =-1,逆傅⾥叶变换N:输⼊的点数,为偶数,⼀般为2的幂次级,2,4,8,16...k: 满⾜N=2^k(k>0),实质上k是N个采样数据可以分解为偶次幂和奇次幂的次数real[]: inv=1时,存放N个采样数据的实部,inv=-1时,存放傅⾥叶变换的N个实部imag[]: inv=1时,存放N个采样数据的虚部,inv=-1时,存放傅⾥叶变换的N个虚部出⼝参数:real[]: inv=1时,返回傅⾥叶变换的实部,inv=-1时,返回逆傅⾥叶变换的实部imag[]: inv=1时,返回傅⾥叶变换的虚部,inv=-1时,返回逆傅⾥叶变换的虚部*/void FFT::dealFFT(double real[], double imag[], double dSp[], int N, int k, int inv){int i, j, k1, k2, m, step, factor_step;double temp_real, temp_imag, factor_real, factor_imag;if (inv != 1 && inv != -1)return;//double *real = new double[N];//double *imag = new double[N];//倒序j = 0;for (i = 0; i < N; i++){if (j>i){temp_real = real[j];real[j] = real[i];real[i] = temp_real;temp_imag = imag[j];imag[j] = imag[i];imag[i] = temp_imag;}m = N / 2;while (j >= m&&m != 0){j -= m;m >>= 1;}j += m;}//蝶形运算for (i = 0; i < k; i++){step = 1 << (i + 1);factor_step = N >> (i + 1); //旋转因数变化速度//初始化旋转因⼦factor_real = 1.0;factor_imag = 0.0;for (j = 0; j < step / 2; j++){for (k1 = j; k1 < N; k1 += step){k2 = k1 + step / 2; //蝶形运算的两个输⼊/* temp_real = real[k1] + real[k2] * factor_real - imag[k2] * factor_imag;temp_imag = imag[k1] + real[k2] * factor_imag + imag[k2] * factor_real;real[k2] = real[k1] - (real[k2] * factor_real - imag[k2] * factor_imag);imag[k2] = imag[k1] - (real[k2] * factor_imag + imag[k2] * factor_real);real[k1] = temp_real;imag[k1] = temp_imag;*///上⾯⽅法虽然直⽩,但效率太低,稍微改变结构如下:temp_real = real[k2] * factor_real - imag[k2] * factor_imag;temp_imag = real[k2] * factor_imag + imag[k2] * factor_real;real[k2] = real[k1] - temp_real;imag[k2] = imag[k1] - temp_imag;real[k1] = real[k1] + temp_real;imag[k1] = imag[k1] + temp_imag;}factor_real = inv*cos(-2 * PI*(j + 1)*factor_step / N);factor_imag = inv*sin(-2 * PI*(j + 1)*factor_step / N);}}if (inv == -1){for (i = 0; i <= N - 1; i++){real[i] = real[i] / N;imag[i] = imag[i] / N;}}for (i = 0; i<N;i++){dSp[i] = sqrt(real[i] * real[i] + imag[i] * imag[i]);}}⼀般好像需要进⾏下转换,即后半部分和前半部分置换,即1234变成3412.void FFT::FFTShift(double dp[], int len){for (int i = 0; i < len / 2; i++){double tmp = dp[i];dp[i] = dp[i + len / 2];dp[i + len / 2] = tmp;}}此时得到的应该是实部和虚部解出来的频谱图的Y轴电压值,⼀般频谱图Y轴为dB,因此需要进⾏转换void FFT::getFFT(double dRe[], double dIm[], double dSp[], int len, int nBits, double dWorkingImpedance){dealFFT(dRe, dIm, dSp, len, nBits, 1);FFTShift(dSp,len); //此时得到的应该是实部和虚部解出来的频谱图的Y轴电压值,还需要转换////dBW = 10lg(电压^2/阻抗);dBm =dBW+30,注意电压单位是Vfor (int i = 0; i<len; i++){dSp[i] = 10 * log10(dSp[i] * dSp[i] / dWorkingImpedance)+30;}}getFFT()输出之后的dp才是要的频谱图Y轴值,频谱图X轴的坐标得到通过以下⽅式://X轴精确度,采样频率/数据个数 = 步长m_DeltaX_S = m_dataPara.nSampleFrequency / nDataNumOfPage_S;data_SX[i / 2] = m_dataPara.nCenterFrequency + count*m_DeltaX_S - m_dataPara.nWorkingBandWidth/2;//中⼼频率+当前点*步长-带宽/2在项⽬中,实际代码如下:int count = 0;for (int i = 0; i < nDataNumOfPage_S * 2; i++){if (i % 2 == 0)data_SQ[i / 2] = data_S[i] * m_DeltaY_S;elsedata_SI[i / 2] = data_S[i] * m_DeltaY_S;if (i % 2 == 0){count++;data_SX[i / 2] = m_dataPara.nCenterFrequency + count*m_DeltaX_S - m_dataPara.nWorkingBandWidth/2;}}m_dataPara.nWorkingImpedance = 50;FFT fft;int nBits = log10(nDataNumOfPage_S) / log10(2);//因为参数需要是2的N次⽅fft.getFFT(data_SQ, data_SI, data_SS, nDataNumOfPage_S, nBits, m_dataPara.nWorkingImpedance); LoadData_S(data_SX, data_SS, nDataNumOfPage_S);。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
电子科技大学通信与信息工程学院标准实验报告(实验)课程名称DSP设计与实践电子科技大学教务处制表电 子 科 技 大 学实 验 报 告学生姓名: 学 号 指导教师:实验地点: 实验时间: 一、实验室名称: 科B341 二、实验项目名称:快速傅立叶变换(FFT )的实现 三、实验学时:4 四、实验原理:基—2按时间抽取FFT 算法对于有限长离散数字信号{x[n]},0 ≤ n ≤ N-1,其离散谱{x[k]}可以由离散付氏变换(DFT )求得。
DFT 的定义为可以方便的把它改写为如下形式:不难看出,W N 是周期性的,且周期为N ,即W N 的周期性是DFT 的关键性质之一。
为了强调起见,常用表达式W N 取代W 以便明确其周期是N 。
由DFT 的定义可以看出,在x[n]为复数序列的情况下,完全直接运算N 点DFT 需要(N-1)2次复数乘法和N (N-1)次加法。
因此,对于一些相当大的N 值(如1024)来说,直接计算它的DFT 所作的计算量是很大的。
FFT 的基本思想在于,将原有的N 点序列序列分成两个较短的序列,这些序列的DFT 可以很简单的组合起来得到原序列的DFT 。
例如,若N 为偶数,将原有的N 点序列分成两个(N/2)点序列,那么计算N 点DFT 将只需要约[(N/2)2 ·2]=N 2/2次复数乘法。
即比直接计算少作一半乘法。
因子(N/2)2表示直接计算(N/2)点DFT 所需要的乘法次数,而乘数2代表必须完成两个DFT 。
()1,...,1,0][)2(10-==--=∑N k en x k X nk Nj N n π()1,...,1,0][10-==∑-=N k W n x k X nkNN n ...2,1,0,))((±±==++l m W W nk NlN k mN n N上述处理方法可以反复使用,即(N/2)点的DFT 计算也可以化成两个(N/4)点的DFT (假定N/2为偶数),从而又少作一半的乘法。
这样一级一级的划分下去一直到最后就划分成两点的FFT 运算的情况。
比如,一个N = 8点的FFT 运算按照这种方法来计算FFT 可以用下面的流程图来表示:x(0)x(1)x(2)x(3)x(4)x(5)x(6)x(7)X(7)X(6)X(5)X(4)X(3)X(2)X(1)X(0)关于蝶形结运算的具体原理及其推导可以参照讲义,在此就不再赘述。
实数FFT 运算对于离散傅立叶变换(DFT )的数字计算,FFT 是一种有效的方法。
一般假定输入序列是复数。
当实际输入是实数时,利用对称性质可以使计算DFT 非常有效。
一个优化的实数FFT 算法是一个组合以后的算法。
原始的2N 个点的实输入序列组合成一个N 点的复序列,之后对复序列进行N 点的FFT 运算,最后再由N 点的复数输出拆散成2N 点的复数序列,这2N 点的复数序列与原始的2N 点的实数输入序列的DFT 输出一致。
使用这种方法,在组合输入和拆散输出的操作中,FFT 运算量减半。
这样利用实数FFT 算法来计算实输入序列的DFT 的速度几乎是一般复FFT 算法的两倍。
本实验就用这种方法实现了一个256点实数FFT (2N = 256)运算。
⒈ 实数FFT 运算序列的存储分配如何利用有限的DSP 系统资源,合理的安排好算法使用的存储器是一个比较重要的问题。
参见FFT 实验程序的CMD 文件:MEMORY {PAGE 0: IPROG: origin = 0x3080, len = 0x1F80 VECT: origin = 0x3000, len = 0x80 EPROG: origin = 0x38000, len = 0x8000PAGE 1: USERREGS: origin = 0x60, len = 0x1c BIOSREGS: origin = 0x7c, len = 0x4 IDA TA: origin = 0x80, len = 0xB80 EDATA: origin = 0xC00, len = 0x1400}SECTIONS{.vectors: {} > VECT PAGE 0.sysregs: {} > BIOSREGS PAGE 1.trcinit: {} > IPROG PAGE 0.gblinit: {} > IPROG PAGE 0.bios: {} > IPROG PAGE 0frt: {} > IPROG PAGE 0.text: {} > IPROG PAGE 0.cinit: {} > IPROG PAGE 0.pinit: {} > IPROG PAGE 0.sysinit: {} > IPROG PAGE 0.data {} > EDA TA PAGE 1.bss: {} > IDA TA PAGE 1.far: {} > IDA TA PAGE 1.const: {} > IDATA PAGE 1.switch: {} > IDATA PAGE 1.sysmem: {} > IDATA PAGE 1.cio: {} > IDA TA PAGE 1.MEM$obj: {} > IDA TA PAGE 1.sysheap: {} > IDATA PAGE 1}从上面的连接定位CMD文件可以了解到,程序代码安排在0x3000开始的存储器中。
其中0x3000-0x3080存放中断向量表。
FFT程序使用的正弦表、余弦表数据(.data段)安排在0xc00开始的地方。
变量(.bss段定义)存放在0x80开始的地址中。
另外,本256点实数FFT程序的输入数据缓冲为0x2300-0x23ff,FFT后功率谱的计算结果存放在0x2200-0x22ff中。
⒉基二实数FFT运算的算法该算法主要分为四步:第一步,输入数据的组合和位倒序把输入序列作位倒序,是为了在整个运算最后的输出中得到的序列是自然顺序。
首先,原始的输入的2N = 256个点的实数序列复制放到标记有“d_input_addr”的相邻单元,当成N = 128点的复数序列d[n]。
奇数地址是d[n]的实部,偶数地址是d[n]的虚部。
这个过程叫做组合(n是从0到无穷,指示时间的变量,N是常量)。
然后,复数序列经过位倒序,存储在数据处理缓冲器中,标记为“fft-data”。
①如图2,输入实数序列为a[n],n=0,1,2,3,…,255。
分离a[n]成两个序列,如图3所示。
原始的输入序列是从地址0x2300到0x23FF,其余的从0x2200到0x22FF的是经过位倒序之后的组合序列:n=0,1,2,3, (127)②d[n]表示复合FFT的输入,r[n]表示实部,i[n]表示虚部:d[n]=r[n]+j i[n]按位倒序的方式存储d[n]到数据处理缓冲中,如图2。
图2* 编程技巧:在用C54x进行位倒序组合时,使用位倒序寻址方式可以大大提高程序执行的速度和使用存储器的效率。
在这种寻址方式中,AR0存放的整数N是FFT点数的一半,一个辅助寄存器指向一数据存放的单元。
当使用位倒序寻址把AR0加到辅助寄存器中时,地址以位倒序的方式产生,即进位是从左到右,而不是从右到左。
例如,当AR0 = 0x0060,AR2 = 0x0040时,通过指令:MAR AR2+0B我们就可以得到AR2位倒序寻址后的值为0x0010。
下面是0x0060(1100000)与0x0040(1000000)以位倒序方式相加的过程:1 1 0 0 0 0 01 0 0 0 0 0 0+0 0 1 0 0 0 0实现256点数据位倒序存储的具体程序段如下:bit_rev:STM #d_input_addr,ORIGINAL_INPUT ;在AR3(ORIGINAL_INPUT)中;放入输入地址STM #fft_data,DATA_PROC_BUF ;在AR7(DATA_PROC_BUF)中;放入处理后输出的地址MVMM DA TA_PROC_BUF,REORDERED_DA TA ;AR2(REORDERED_DA TA);中装入第一个位倒序数据指针STM #K_FFT_SIZE-1,BRCSTM #K_FFT_SIZE,AR0 ;AR0的值是输入数据数目的一半=128RPTB bit_rev_endMVDD *ORIGINAL_INPUT+,*REORDERED_DATA+ ;将原始输入缓冲中的数据;放入到位倒序缓冲中去之;后输入缓冲(AR3)指针加1;位倒序缓冲(AR2)指针也加;一MVDD *ORIGINAL_INPUT-,*REORDERED_DA TA+ ;将原始输入缓冲中的数据;放入到位倒序缓冲中去之;后输入缓冲(AR3)指针减一;位倒序缓冲(AR2)指针加一;以保证位倒序寻址正确MAR *ORIGINAL_INPUT+0B ;按位倒序寻址方式修改AR3 bit_rev_end:注意,在上面的程序中。
输入缓冲指针AR3(即ORIGINAL_INPUT)在操作时先加一再减一,是因为我们把输入数据相邻的两个字看成一个复数,在用寄存器间接寻址移动了一个复数(两个字的数据)之后,对AR3进行位倒序寻址之前要把AR3的值恢复到这个复数的首字的地址,这样才能保证位倒序寻址的正确。
第二步,N点复数FFT在数据处理缓冲器里进行N点复数FFT运算。
由于在FFT运算中要用到旋转因子W N,它是一个复数。
我们把它分为正弦和余弦部分,用Q15格式将它们存储在两个分离的表中。
每个表中有128项,对应从0度到180度。
因为采用循环寻址来对表寻址,128 = 27 < 28,因此每张表排队的开始地址就必须是8个LSB位为0的地址。
参照前面图1 DES 系统的存储区分配,所以我们必须把正弦表第一项“sine_table”放在0x0D00的位置,余弦表第一项“cos_table ”放在0x0E00的位置。
① 根据公式利用蝶形结对d[n]进行N=128点复数FFT 运算,其中所需的正弦值和余弦值分别以Q15的格式存储于内存区以0x0D00开始的正弦表和以0x0E00开始的余弦表中。
我们把128点的复数FFT 分为七级来算,第一级是计算两点的FFT 蝶形结,第二级是计算四点的FFT 蝶形结,然后是八点、十六点、三十二点六十四点、一百二十八点的蝶形结计算。
最后所得的结果表示为: D[k] = F{d[n]} = R[k] + j I[k]其中,R[k]、I[k]分别是D[k]的实部和虚部。
图3② FFT 完成以后,结果序列D[k]就存储到数据处理缓冲器的上半部分,如图3()1,...,1,0][10-==∑-=N k W n d k D nkNN n)2sin()2cos()2(nk Nj nk N eWnk Nj nk Nπππ-==-所示,下半部分任然保留原始的输入序列a[n],这半部分将在第三步中被改写。