void FftExe(IN_TYPE *pIn, OUT_TYPE *pRe, OUT_TYPE *pIm)函数中有实现FFT的一段，能讲解下这三级循环吗？？？
我看的书是西安电子科技大学出版社·丁美玉·《数字信号处理》，上面的程序段和书P103页说的程序流程图不一样····
想搞清楚FFT，所以请大概讲解一下这段代码··Thanks````
---------------------------------------------------------------------
//先计算2点的
for(j=0;j<LENGTH;j+=2)
{
tr=pIn[j+1];
pRe[j+1]=(pIn[j]-tr);
pIm[j+1]=0;
pRe[j]=(pIn[j]+tr);
pIm[j]=0;
}
for(BlockSize=4;BlockSize<=LENGTH;BlockSize<<=1) //再一层层计算
{
for(j=0;j<LENGTH;j+=BlockSize)
{
for(i=0;i<BlockSize/2;i++)
{
#ifndef USE_TABLE
c=(long)(1024*cos(2*PI*i/BlockSize));
s=(long)(1024*sin(2*PI*i/BlockSize));
#else
OffSet0=LENGTH/BlockSize*i;//
c=COS_TABLE[OffSet0];//
s=SIN_TABLE[OffSet0];//
#endif
OffSet1=i+j;
OffSet2=OffSet1+BlockSize/2；
tr=(OUT_TYPE)((c*pRe[OffSet2]+s*pIm[OffSet2])/1024);
ti=(OUT_TYPE)((c*pIm[OffSet2]-s*pRe[OffSet2])/1024);
#ifdef UNITARY //如果要对结果归一化处理，则每次运算要除以2
pRe[OffSet2]=(pRe[OffSet1]-tr)/2;
pIm[OffSet2]=(pIm[OffSet1]-ti)/2;
pRe[OffSet1]=(pRe[OffSet1]+tr)/2;
pIm[OffSet1]=(pIm[OffSet1]+ti)/2;
#else
pRe[OffSet2]=(pRe[OffSet1]-tr);
pIm[OffSet2]=(pIm[OffSet1]-ti);
pRe[OffSet1]=(pRe[OffSet1]+tr);
pIm[OffSet1]=(pIm[OffSet1]+ti);
#endif
}
}
}
---------------------------------------------------------------------
先自己理解的点点是：
for(j=0;j<LENGTH;j+=2)
{
tr=pIn[j+1];
pRe[j+1]=(pIn[j]-tr);
pIm[j+1]=0;
pRe[j]=(pIn[j]+tr);
pIm[j]=0;
} 这段是完成第一级2点的FFT，采用实数运算，所以虚部取0.这段代码的表达式可以用书P102上实数的碟形运算规律解释。
|
楼主
标题置顶
标题高亮
