TI库实现点或更大点数FFT.docx
- 文档编号:6452143
- 上传时间:2023-01-06
- 格式:DOCX
- 页数:12
- 大小:65.62KB
TI库实现点或更大点数FFT.docx
《TI库实现点或更大点数FFT.docx》由会员分享,可在线阅读,更多相关《TI库实现点或更大点数FFT.docx(12页珍藏版)》请在冰豆网上搜索。
TI库实现点或更大点数FFT
在DSP运算中,经常需要把输入时域信号在频域进行处理之后,再还原为时域信号,这样就需要进行FFT和IFFT运算:
x(n)-> FFT ->X(f)->频域处理->Y(f)->IFFT->y(n)
而一般的DSP芯片只支持整数运算,也就是说只能进行定点小数计算。
N点FFT计算出0…N-1,N个复数:
0,A,N/2,A*,A为(N/2-1)个复数,A*为A的共轭复数。
FFT的公式为:
N
X(k)= sum x(n)*exp(-j*2*pi*(k-1)*(n-1)/N)、1<= k<= N.
n=1
IFFT的公式为:
N
x(n)=(1/N)sum X(k)*exp(j*2*pi*(k-1)*(n-1)/N)、1<= n<= N.
k=1
假设我们对ADC转换器转换的数字信号进行FFT运算,若输入数据为16bit的短整型数,我们可以把它看作Q15的从-1到1之间的小数。
根据FFT的公式我们可以知道,FFT变换之后的结果将超出这个范围。
例如在matlab中输入fft(sin([1:
8]*0.5)),可以看到结果:
2.8597,-0.8019-3.0216i,0.4312-0.8301i,0.5638-0.3251i,0.5895,0.5638+0.3251i,0.4312+0.8301i,-0.8019+3.0216i
实际上,FFT变换之后的数据的范围在-N到N之间,N为FFT的点数。
为了正确地表达-N到N之间的数值,输出数据的Q值将变小,例如若N=1024,输入数据为Q15的话,那么输出数据则必须为Q5才能够确保结果不会溢出。
这样的结果将丢失很多信息,以至于IFFT无法还原为原来的数据。
如下图所示:
Q15->1024FFT->Q5->1024IFFT->Q5
这样,经过FFT和IFFT变换之后,数据从Q15变成了Q5,丢失了10bit的信息。
为了使得定点FFT和IFFT之后能够还原为原来的数据,必须使用32bit定点FFT,输入数据虽然是16bit,在进行FFT的时候将将它转换为32bit运算,这样输入的数据为Q31,1024点FFT的输出数据为Q21,满足16bit的精度要求。
这样整个计算流程变为:
Q15->Q31->1024FFT32bit->Q21->1024IFFT32bit->Q21->Q15
在TI的DSPLIB中的FFT和IFFT都提供了两种选择:
SCALE和NOSCALE。
仍以1024点计算为例,它们的意思分别为:
Q15->NOSCALE FFT ->Q15(结果可能溢出)
Q15->SCALE FFT ->Q5(结果不会溢出,但是精度降低)
Q15->NOSCALEIFFT->Q15(结果可能溢出)
Q15->SCALEIFFT ->Q5(结果不会溢出,但是精度降低)
32位定点FFT,IFFT与此类似,不再重复。
因此,用DSPLIB进行FFT和IFFT计算时,注意必须采用32bit精度,而且FFT变换时采用SCALE,而IFFT变换时采用NOSCALE。
关于TIFFT库的使用
2010-01-2512:
53星期一
这篇文章是应一些找我讨论DSP的同学所写,贴在这里大家一起学习。
曾有不少论坛上的同学(包括DSP算法讨论群里的一些同学)问过我关于TI的FFT库的使用,这里我将我使用过的一些经验说一下。
TI的这个FFT库在计算速度、计算精度以及数据存储等方面是做了不少优化的,比如数据存储,若作N点的FFT,供查表用的旋转因子必须有N/2点的正弦值与N/2点的余弦值,这个库将其压缩成3N/4点的正弦值,因此就节省了N/4点的存储空间;另外计算N点实数FFT时,一般简单的做法是将N点实数的虚部全化为0来处理,而这个库则把N点实数数据打包成N/2点复数数据来处理,在计算速度和存储空间都有很大改进。
之前我浏览helloDSP论坛的帖子,有很多人发出疑问:
为什么我计算出的mag值全为零?
这样的帖子真不少见;TI的官方工程师论坛不少老外也发问:
WhydidIgetallzeros?
事实上我想主要原因是输入数据格式不对。
我认为使用这个库主要注意一下两点:
1.数据输入输出的Q格式;
2.存储空间分配。
下面以32位实数FFT为例来说明。
注意到文档的40页有如下说明:
1.在函数voidcalc(RFFT32_handle)有如下一句:
NotethattheinputandoutputdataareinQ31format.;
2.在函数voidmag(FFT128R_handle)有如下一句:
NotethatthemagnitudeoutputisstoredinQ30format.
因为28x系列DSP是定点处理器,而FFT计算涉及到不少浮点计算,TI使用Q格式来解决这个问题(Q格式说明可参考sprc087_IQmath)。
事实上输入数据采用Q31格式能在避免计算溢出前提下获得最好的计算精度。
对AD采样的数据进行FFT计算,定义计算缓冲区数组:
longipcb[N+2];
因为AD结果寄存器是12位的,在数据左对齐的情况下直接左移15位即可:
ipcb[ConversionCount]=((unsignedlong)AdcRegs.ADCRESULT0)<<15;
另外我记得有份文档提到在某些存储器下是右对齐,此时则需左移19位,大家在使用时注意这个问题。
当采样完成后,按照FFT库的文档上的说明或者仿照文档附带的例程,进行相应函数调用来实现自己的FFT计算,比如按计算点数来定义各个变量数组,是否加窗,是否求解幅值平方值等。
在进行存储器分配时,文档上要求(128点实数FFT为例):
FFTipcbALIGN(256):
{}>L0L1RAMPAGE1
FFTmag>L0L1RAMPAGE1
FFTtf>NVMEMPAGE0/*Nonvolatilememory*/
.econst>NVMEMPAGE0/*Nonvolatilememory*/
注意两点:
FFT计算缓冲区FFTipcb需在page1上连续分配2N个位置(以ALIGN来指定),FFTtf(旋转因子)需放在Nonvolatilememory的page0内(事实上如何才能为Nonvolatile我也不清楚)。
FFTtf位置这点我之前在调试时对计算FFT影响很大,因为twiddlefactor若因存储冲突肯定会造成查表值不准确,那计算FFT时肯定就不对了。
有个论坛帖子作者说一定要放在origin=0x008000开始位置,其实也不对,大家可自己去试验;后面我也会给出我的存储配置文件(即.cmd文件内容)。
无图无真相,下面给出计算实例。
假设有一信号包含两个谐波频率值,分别为413.0Hz(幅值设为1.00V)和287.0Hz(幅值设为0.400V),利用函数发生器产生这两路信号再混合,加上偏置后送入AD采样。
设采样频率1024Hz,共采样2048点。
图1的采样点均是右对齐的12位采样结果值。
图1AD采样得到的采样点图
利用TI的FFT库进行计算,查看mag数组,得到图2.
以1024Hz采样2048点,采样时间2s,对两个谐波频率可采样到整数倍周期;从另一个角度理解,此时最小频率分辨率为0.5Hz,413.0Hz与287.0Hz均是其整数倍数,故不会发生频谱展宽或混叠情况,计算得到的频谱图应该为两根尖峰线。
从图2结果也能看出这一点。
查看mag数组,可知第一根尖峰线下标574,第二根尖峰线下标826,故真实频率值分别为:
574*0.5=287Hz,826*0.5=413Hz
若要计算幅值,按照输出的Q30格式除相应系数即可。
注意最好另外定义浮点数组来做除运算,因为整型数据做除运算(或者右移位操作)会丢失小数位数据。
若要验证DSP的FFT计算结果,可将AD数据从CCS导入到MATLAB中做对比计算。
有些同学不清楚如何导入导出,下面说一下步骤。
1.点菜单栏file-data-save,选择保存类型Integer(若是其他进制还需在MATLAB中转换),点“确定”后,在“address”栏填入要保存数据的起始地址(填变量数组名或真实存储器地址皆可),在“length”栏内填入数据长度,“page”肯定选“data”页了;全部设好后点“OK”
2.找到保存的数据文件,将后缀改成.txt,再用记事本打开(也可以不改后缀直接用记事本打开),删去第一行数据;
3.打开MATLAB,点开“loaddatafile…”,选中刚才的数据文件,然后按照提示一步步往下导入即可。
最后不妨用变量temp来保存这些数据。
在MATLAB内运行如下代码:
%%%%%%%%%%%%%%%%%%
f_sample=1024;%采样频率
N=2048;%采样点N
y=temp;%temp即为导入数据的变量名
n=0:
N-1;
t=n/f_sample;
%做采样点图
plot(t,y);figure;
stem(t,y,'.');figure;
%fft变换并作图
fft_result=fft(y);
mag=abs(fft_result)/(N/2);
mag
(1)=0;%为观察谐波分量,此处特意将直流分量置为0
%求解真实频率值
f=(0:
length(fft_result)-1)*f_sample/length(fft_result);
%作频谱图
stem(f(1:
N/2),mag(1:
N/2),'b.');grid;
%%%%%%%%%%%%%%%%%%
可将MATLAB计算结果与CCS内的结果做些对比。
附:
存储器分配(.cmd文件内容),我用的是2808的板子,2812或其他的稍作改动即可。
MEMORY
{
PAGE0:
BEGIN:
origin=0x000000,length=0x000002
RAMM0:
origin=0x000002,length=0x0003FE
PRAMH0:
origin=0x3FA000,length=0x002000
RESET:
origin=0x3FFFC0,length=0x000002
BOOTROM:
origin=0x3FF000,length=0x000FC0
TESARAM:
origin=0x008000,length=0x002000
PAGE1:
BOOT_RSVD:
origin=0x000400,length=0x000080
RAMM1:
origin=0x00A000,length=0x001000
L0L1RAM:
origin=0x00B000,length=0x001000
HL0SARAM:
origin=0x3F8000,length=0x001010
DRAMH0:
origin=0x3F9010,length=0x000ff0
}
SECTIONS
{
/*Setupfor"boottoSARAM"mode:
Thecodestartsection(foundinDSP28_CodeStartBranch.asm)
re-directsexecutiontothestartofusercode.*/
codestart:
>BEGIN,PAGE=0
ramfuncs:
>RAMM0PAGE=0
.text:
>PRAMH0,PAGE=0
.cinit:
>RAMM0,PAGE=0
.pinit:
>RAMM0,PAGE=0
.switch:
>RAMM0,PAGE=0
.reset:
>RESET,PAGE=0,TYPE=DSECT/*notused,*/
FFTipcbALIGN(4096):
{}>HL0SARAM,PAGE=1
FFTtf:
>DRAMH0,PAGE=1
FFTmag:
>L0L1RAM,PAGE=1
.const:
>DRAMH0,PAGE=1
.bss:
>DRAMH0,PAGE=1
.stack:
>RAMM1,PAGE=1
.sysmem:
>RAMM1,PAGE=1
.ebss:
>DRAMH0,PAGE=1
.econst:
>PRAMH0,PAGE=0
.esysmem:
>L0L1RAM,PAGE=1
IQmath:
>PRAMH0,PAGE=0
IQmathTables:
>BOOTROM,type=NOLOAD,PAGE=0
}
另外,有些同学一开始运行FFT附带的例程,会遇上两个问题:
一是缺少函数文件,这个到TI官网下那个sprc083_SGEN包就好了,或找我我用邮箱发给大家也行;
二是提示编译不成功,找到上述配置文件的这一行
.reset:
>RESET,PAGE=0,TYPE=DSECT/*notused,*/
在后面加上TYPE=DSECT/*notused,*/。
这些是我凭印象写的,因好久没做这个FFT,可能会写出错误,欢迎大家提出来:
-)
有问题就解决问题之TI库实现2048点或更大点数FFT
2010-05-1813:
20星期二
TI提供的FFT库(sprc081)使用起来很方便,只要直接调用即可;更重要的是其中做了好些优化,而且是用汇编写成,运算速度肯定要比自己写的C代码更快好多。
这个库最大能运算1024点的实数FFT,但实际中我们可能还需要更大点数的FFT,那如何去修改呢?
其实修改的地方并不多,以前我曾花了一段时间浏览了一下那些头疼的汇编代码,发现模块化写得很好啊(佩服TI的工程师!
),这样修改起来就快多了。
1.2048点FFT
直接在工程中加入原来的库文件fft.lib,并在代码中作如下修改:
第一步:
在有关FFT代码的头文件fft.h中增添如下一段代码:
#defineRFFT32_2048P_DEFAULTS{(long*)NULL,\
(long*)NULL,\
1024,\
10,\
(long*)NULL,\
(long*)NULL,\
0,\
0,\
1,\
(void(*)(void*))CFFT32_init,\
(void(*)(void*))CFFT32_calc,\
(void(*)(void*))RFFT32_split,\
(void(*)(void*))RFFT32_mag,\
(void(*)(void*))RFFT32_win}
结构体RFFT32关键说明:
1.第一个(long*)NULL,\,是计算缓冲区(ipcb)的指针;
2.第二个(long*)NULL,\,是查表旋转因子(ft)的指针;
3.1024指的是FFT的大小,注意这里是2048点实数FFT被打包成1024点的复数,这是这个库的实数FFT运算的特点,这样节省一半运算时间和存储空间;
4.10即1024是2的10次方,因为这里采用的是基2时间抽取算法;
5.下面的两个(long*)NULL,\分别是运算结果中存放幅值(mag)以及加窗表格(win)的指针;
6.下面的两个0分别是计算结果中的峰值幅值与峰值频率,这里初始化0;
7.下面的1指的是查表的步进值,因在这个sprc081文档的库已经有Q31格式(以及Q30格式)的1024点的正弦因子表,因此对于1024点复数FFT查表步进值设为1即可。
8.下面其他的就是要调用的函数指针了,不再赘述。
明白上述结构体的定义后,改写其他点数就很容易了,比如这个库最少能做到128点的实数FFT,我们还可以相应改成64点或者更低点数的实数FFT,比如64点的结构体就是如下的一段:
#defineRFFT32_64P_DEFAULTS{(long*)NULL,\
(long*)NULL,\
32,\
5,\
(long*)NULL,\
(long*)NULL,\
0,\
0,\
32,\
(void(*)(void*))CFFT32_init,\
(void(*)(void*))CFFT32_calc,\
(void(*)(void*))RFFT32_split,\
(void(*)(void*))RFFT32_mag,\
(void(*)(void*))RFFT32_win}
第二步:
修改源程序代码里的变量定义:
#include"fft.h"
#defineN2048//定义N的点数,2048点
//分别在数据存储区为ipcb和mag分配足够的地址,具体见下面的连接器文件
#pragmaDATA_SECTION(ipcb,"FFTipcb");
#pragmaDATA_SECTION(mag,"FFTmag");
//定义变量
RFFT32fft=RFFT32_2048P_DEFAULTS;
longipcb[N+2];
longmag[N/2+1];
//其他部分函数调用都是一样的,这里不再列出
第三步:
在连接器(.cmd)文件里分配存储地址:
这里给出我在2808上的分配方法:
MEMORY
{
PAGE1:
L0L1RAM:
origin=0x00B000,length=0x001000
HL0SARAM:
origin=0x3F8000,length=0x001010
DRAMH0:
origin=0x3F9010,length=0x000ff0
}
SECTIONS
{
FFTipcbALIGN(4096):
{}>HL0SARAM,PAGE=1
FFTtf:
>DRAMH0,PAGE=1
FFTmag:
>L0L1RAM,PAGE=1
}
好,到此为止就修改好了,可以使用这个FFT库实现2048点实数FFT了。
2.4096点或更大点数的FFT
第一步:
修改汇编代码并制成库文件
其实要实现4096点的FFT需要改动的地方也不多,主要因为原来的这个库只提供了1024点的正弦值,我们需要自己算出2048点的Q31(或Q30)格式的一周期正弦值,并将其代替原来汇编代码里的1024点的表格。
汇编代码的其他地方不需要改动,不然可能产生不可预知的后果。
在汇编代码文件sel_q.asm里查看设置:
;SelectTwiddlefactorQformat
TF_QFMAT.setQ30
若是TF表设为Q30,则只需将计算出的正弦值移位成Q30格式即可。
另外需注意的是,2048点的FFT需要2048点的TF表,但这个库将其优化成只需2048×3/4=1536点,即1.5K大小的表格(见sprc081文档第10页,注意到正弦的第二部分与余弦的第一部分是相重合的,这样便节省了1/4周期的点数)。
计算好正弦值并加入到汇编代码里,然后将所有的汇编代码导进CCS里做成一个新的库,到时直接添加到工程里调用。
为了方便这里我已将其做成一个新的库命名为FFT_new.lib,大家直接用即可。
第二步:
修改头文件fft.h结构体定义
增添如下一段代码:
#defineRFFT32_4096P_DEFAULTS{(long*)NULL,\
(long*)NULL,\
2048,\
11,\
(long*)NULL,\
(long*)NULL,\
0,\
0,\
1,\
(void(*)(void*))CFFT32_init,\
(void(*)(void*))CFFT32_calc,\
(void(*)(void*))RFFT32_split,\
(void(*)(void*))RFFT32_mag,\
(void(*)(void*))RFFT32_win}
第三步:
修改变量定义和变量存储地址
这里修改的方法与上述相同,不再赘述。
必须注意的是对于32位的4096点的FFT,其缓冲区占用8K连续地址,因此分配地址需考虑周全。
这里给出我在2808上的一种分配方法。
F2808的片内18K单周期访问RAM按地址连续而言分成两块,即M0、M1的0x0000-0000~0x0000-07FF以及L0、L1、H0的0x0000-8000~0x0000-BFFF(同时映射到地址0x003F-8000~0x003F-BFFF)。
M0、M1地址用作未初始化段及堆栈段等,因此只能选用L0等大段RAM空间;另外FFT查表所需的旋转因子常数也需保存在不受冲突的RAM或Flash中。
这两段所需的地址大小就有11K×16,同时还需为程序代码和其他数据变量分配足够的存储段。
为了优化并节省空间,这里没有采用加窗方法,并且将计算缓冲区(ipcb)和结果区(mag)设为同一地址(不然的话还真没法分配):
MEMORY
{
PAGE0:
PRAMH0:
origin=0x00AC04,length=0x0013FC
FFT_TF:
origin=0x00A004,length=0x000C00
PAGE1:
FFT_IPCB:
origin=0x008000,length=0x002004
}
SECTIONS
{
.text:
>PRAMH0,PAGE=0
FFTipcbALIGN(8192):
{}>FFT_IPCB,PAGE=1
FFTtf:
>FFT_TF,PAGE=0
}
我用的2808片内只含18K的RAM,只用这些RAM对
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- TI 实现 更大 点数 FFT
![提示](https://static.bdocx.com/images/bang_tan.gif)