FFT的实现不明白,希望指点下
在学FFT,有两个地方没有弄明白1 蝶形计算没看明白,希望能够详细解释下
dataR=dataR+dataR*cos_tab+dataI*sin_tab;
dataI=dataI-dataR*sin_tab+dataI*cos_tab;
dataR=TR-dataR*cos_tab-dataI*sin_tab;
dataI=TI+temp*sin_tab-dataI*cos_tab;
2 makewave函数中
INPUT=sin(PI*2*i/SAMPLENUMBER*3)*1024;
乘以3再乘1024 是什么意思?
***********************************************************************
#include "myapp.h"
#include "csedu.h"
#include "scancode.h"
#include <math.h>
#define PI 3.1415926
#define SAMPLENUMBER 128
void InitForFFT();
void MakeWave();
int INPUT,DATA;
float fWaveR,fWaveI,w;
float sin_tab,cos_tab;
main()
{
int i;
InitForFFT();
MakeWave();
for ( i=0;i<SAMPLENUMBER;i++ )
{
fWaveR=INPUT;
fWaveI=0.0f;
w=0.0f;
}
FFT(fWaveR,fWaveI);//地址传递,数组首地址传递给dataR,dataI
for ( i=0;i<SAMPLENUMBER;i++ )
{
DATA=w;
}
while ( 1 ); // break point
}
void FFT(float dataR,float dataI)
{
int x0,x1,x2,x3,x4,x5,x6,xx;
int i,j,k,b,p,L;
float TR,TI,temp;
/********** following code invert sequence ************/
/*位倒序算法:128点FFT可用7位二进制表示顺序数,这里依次为x6-x0(高到低)
x0到x6依次与0x01相与,得到顺序数i用7位二进制表示,xx则是计算倒序二进制的十进制数
dataR数组是存放顺序数的(输入),dataI数组则是存放倒序数的
注意 i/(2,4,8,16等)相当于i二进制右移1,2,3,4位等*/
for ( i=0;i<SAMPLENUMBER;i++ )
{
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;
xx=x0*64+x1*32+x2*16+x3*8+x4*4+x5*2+x6;
dataI=dataR; //dataIa 数组中存放的就是位倒序后的值
}
for ( i=0;i<SAMPLENUMBER;i++ )
{
dataR=dataI; dataI=0; //把位倒序数组中的值逐位的赋值给原输入数组,这样就少用一个数组
}
/************** following code FFT *******************/
for ( L=1;L<=7;L++ ) //第一层循环:蝶形的级数
{ /* for(1) */
b=1; i=L-1;
while ( i>0 )
{
b=b*2; i--;
} /* b= 2^(L-1) */
//第二层循环:获得旋转因子 P,第L级有2的(L-1)个不同的旋转因子,同时要计算每一个旋转因子对应的蝶形
for ( j=0;j<=b-1;j++ ) /* for (2) */
{
p=1; i=7-L;
while ( i>0 ) /* p=pow(2,7-L)*j; */
{
p=p*2; i--;
}
p=p*j;
/*第三层循环:计算蝶形 ,第L级旋转因子的个数同每一级中每个蝶形的两个输入数据相差的点数 相同*/
/*实际上这里b的值每一级中都是固定的,与循环的次数并没有关系*/
for ( k=j;k<128;k=k+2*b ) /* for (3) */
{
TR=dataR; //已经是位倒序过了
TI=dataI; // dataI中数据都赋值为0了,这里还有什么作用?
temp=dataR; // temp是蝶形中另一个输入
/****************************************************************/
dataR=dataR+dataR*cos_tab+dataI*sin_tab;
dataI=dataI-dataR*sin_tab+dataI*cos_tab;
dataR=TR-dataR*cos_tab-dataI*sin_tab;
dataI=TI+temp*sin_tab-dataI*cos_tab;
} /* END for (3) */
} /* END for (2) */
} /* END for (1) */
for ( i=0;i<SAMPLENUMBER/2;i++ )
{
w=sqrt(dataR*dataR+dataI*dataI);
}
} /* END FFT */
void InitForFFT()
{
int i;
for ( i=0;i<SAMPLENUMBER;i++ )
{
sin_tab=sin(PI*2*i/SAMPLENUMBER);
cos_tab=cos(PI*2*i/SAMPLENUMBER);
}
}
void MakeWave()
{
int i;
for ( i=0;i<SAMPLENUMBER;i++ )
{
INPUT=sin(PI*2*i/SAMPLENUMBER*3)*1024;
}
}
没深入研究 向你学习啦 本帖最后由 aresc 于 2014-5-13 22:39 编辑
1. INPUT=sin(PI*2*i/SAMPLENUMBER*3)*1024;
*1024实际没啥特别的含义,只是把一个范围为[-1.0,1.0]的正弦信号放大到了[-1024.0,1024.0]的范围.
那个3是正弦信号的频率,采样率为SAMPLENUMBER.
2. 关于蝶形运算建议你找本信号处理的书好好看看,你用的是时域抽取的算法,和频域抽取的蝶形运算是有差别的。几句话写不清楚,呵呵。关键就是旋转因子然后用欧拉恒等式展开成(cosA+jsinA)做复数乘、加。
3. DataI为0的为什么:因为你用的是复数运算的FFT,但你给的输入信号是一个实数的正弦信号,相当于复数信号只有实部,虚部为0。 楼上说的不错,想要了解找本数字信号处理好好看看,就懂了,这东西一两句说不清 aresc 发表于 2014-5-13 22:34
1. INPUT=sin(PI*2*i/SAMPLENUMBER*3)*1024;
*1024实际没啥特别的含义,只是把一个范围为[-1.0,1.0]的正弦 ...
很给力啊 NWPU_CHEN 发表于 2014-5-14 09:20 static/image/common/back.gif
楼上说的不错,想要了解找本数字信号处理好好看看,就懂了,这东西一两句说不清 ...
嗯,就是看书来理解的,但是理论和具体的代码实现还是有区别的。。。。。 FFT算法比较难理解的是最终算法实现是和算法的推导是反过来的。比如N点FFT,按奇偶分解成两个N/2点,...最后拆成N/2个2点。最后算法实现的时候是反过来的,即先做最里层的N/2个两点的FFT,最后是两个N/2点的FFT。如果按算法推导的过程写程序只能用递归函数来实现,用非递归的方法就是现在大家都在采用的程序。 aresc 发表于 2014-5-13 22:34 static/image/common/back.gif
1. INPUT=sin(PI*2*i/SAMPLENUMBER*3)*1024;
*1024实际没啥特别的含义,只是把一个范围为[-1.0,1.0]的正弦 ...
又分析了下,但还是要追问几个问题:
(1) INPUT中的3好像不是频率,应该是周期,因为按照对模拟输入信号进行采样,得到离散序列,进而进行DFT变换。Xa(t)=sin(2π f t) 令t=nT,T=1/fs,那么有Xa(t)=sin(2π f n/ fs),这样看来,3就是周期,而不是频率了。
(2) 我没有明白的地方应该是这里的 旋转因子用欧拉公式表示,后代入到蝶形运算中(之前不疑惑的地方同样是这个地方),比如
dataR=dataR+dataR*cos_tab+dataI*sin_tab; // 这个可以对应起来
dataR=TR-dataR*cos_tab-dataI*sin_tab; //这个按照欧拉展开,sin-tab 不应该是加吗?
dataI=dataI-dataR*sin_tab+dataI*cos_tab;
dataI=TI+temp*sin_tab-dataI*cos_tab; // 这两句不明白,或者说dataI 这里表示的是什么?
简单来说,这里的复数运算 没有看明白,希望前辈有针对性的解释下。 1. sin(PI*2*i/SAMPLENUMBER*3) = sin(2*PI*3*i/SAMPLENUMBER) = sin(2*PI*f*i/fs), 所以f/fs=3/SAMPLENUMBER.
2. 正向FFT时,欧拉公式展开后是cos(-A)+j*sin(-A) = cosA-jsinA.
时域抽取的一个蝶形包含运算: a' = a+b*(cosA-jsinA), b' = a - b*(cosA-jsinA). 其中a、b都是复数。展开之后就是分别求a',b'的实部、虚部,即4个运算式。
aresc 发表于 2014-5-14 22:23 static/image/common/back.gif
FFT算法比较难理解的是最终算法实现是和算法的推导是反过来的。比如N点FFT,按奇偶分解成两个N/2点,...最 ...
是的,不结合具体的FFT流图来理解,很难看懂,我是结合8点的FFT流图来分析的,和代码是一致的。 最近也在研究FFT,这个基数为2的算法真心不错,网上很多人都在用,花费了半天时间才看懂。 icekoor 发表于 2014-8-8 15:39 static/image/common/back.gif
最近也在研究FFT,这个基数为2的算法真心不错,网上很多人都在用,花费了半天时间才看懂。 ...
看懂了抽时间总结分享一下啊 你好
请问你现在弄明白了吗
页:
[1]