作者weast (kuan)
看板C_and_CPP
标题[问题] FFT
时间Sun Mar 29 19:14:47 2009
我这边有一段FFT程式 可是他规定是只能用256点去跑
请问各位 有没有办法再这段程式的基础上改成可以给定任意点呢?
void FFT(double *segment, double **X)
{
int M = 8 , iN = Dim;
int i, j, k, limit_1, limit_2, limit_3, ij, il;
double farg, fpi = 3.1415926, fCos, fSin, fbr, fbi, fcr, fci, fx1, fx2;
double *fxr = new double[window + 1];
double *fxi = new double[window + 1];
//****************** 将segment[256] => fxr[window] *************//
for (i = 1; i <= window; i++)
fxr[i] = fxi[i] = 0; //初始化fxr[window], fxi[window]
for (i = 1; i <= window; i++)
fxr[i] = segment[i - 1]; //将原始信号的256点个值 => fxr[window]
// ***************** 傅立叶主程式开始 *************//
limit_1 = iN - 1;
limit_2 = iN / 2;
ij = 1;
// **************** 重排输入值 ************//
for (i = 1; i <= limit_1; i++)
{
if ((float) i <= ((float) ij - 0.01))
{
fx1 = fxr[ij];
fx2 = fxi[ij];
fxr[ij] = fxr[i];
fxi[ij] = fxi[i];
fxr[i] = fx1;
fxi[i] = fx2;
}
il = limit_2;
while(1)
{
if ((float)il > ((float) (ij) - 0.01))
break;
ij -= il;
il /=2;
}
ij += il;
}
// **************** 傅立叶转换 *************** //
for(i = 1; i <= M; i++) // (limit_1, limit_2) =>(1, 256) (2, 128) (4, 64) (8, 32)..(256,1)
{
limit_1 = (int) pow(2, i - 1);
limit_2 = (int) pow(2, (M - i));
for(j = 1; j <= limit_2; j++) //(1:256) (1:128) (1:64)..........(1:1)
{
for (k = 1; k <= limit_1; k++)
{
limit_3 = (k - 1) + (j - 1) * 2 * limit_1 + 1; //Odd Points
//*******************
fbr = fxr[limit_3];
fbi = fxi[limit_3];
fcr = fxr[limit_3 + limit_1];
fci = fxi[limit_3 + limit_1];
//*******************
farg = 2.0 * fpi * (k - 1) * limit_2 / iN;
fCos = cosl (farg);
fSin = sinl (farg);
//*******************
fx1 = fcr * fCos + fci * fSin ;
fx2 = fci * fCos + (0 - fcr) * fSin ; //
fxr[limit_3] = fbr + fx1;
fxi[limit_3] = fbi + fx2;
fxr[limit_3 + limit_1] = fbr - fx1;
fxi[limit_3 + limit_1] = fbi - fx2;
}
}
}
//
for (i = 1; i <= window; i++)
{
X[0][i - 1] = fxr[i];
X[1][i - 1] = fxi[i];
}
//cout << "X[] = " << X[0][0];
delete[] fxr, fxi;
}
谢谢各位了!!
--
※ 发信站: 批踢踢实业坊(ptt.cc)
◆ From: 140.120.90.193
1F：推 yoco315:我可以先问一下，你希望的数字是多少吗？ 03/29 20:14
2F：推 zhi5566:上网随边抓都有现成的CODE 就256 512 1024 最普遍 03/30 00:22
3F：→ weast:我是自己可以找到512 1024的 可是我是要用128~256中间的 03/30 00:23
4F：→ weast:数字阿 ORZ 03/30 00:23
5F：→ weast:目前需要用到192~~ 03/30 22:45
6F：→ dcgkii:FFT没有192的= =也不可能有 把192後面的值补零 补到256 03/31 13:05
7F：→ weast:可是为什麽matlab的FFT可以用任意数去转换呢? 03/31 13:28
8F：→ zhi5566:去看FFT的演算法 你就知道为什麽一定是2^n点 03/31 22:18
9F：推 ohiammatt:想要任意点数就只能跑DFT了吧 04/02 14:23
10F：→ weast:恩 DFT我会写 只是速度太慢了 ORZ 04/02 16:32