tianyoume 发表于 2007-9-11 15:52

自己写的FFT总是得不出正确结果。

function Y=fft(x)
n=size(x);
for l=0:12,
    if 2^l>=n
      N=2^l;
      L=l;
      break
    end
end
for t=n(2)+1:N,%不足2^L的数在末尾补0
    x(t)=0;
end
x                %x为补0后的序列
nv2=N/2;         %倒位序
nm1=N-1;
I=0;
J=0;
while I<nm1
    if I<J
      t   =x(J+1);
      x(J+1)=x(I+1);
      x(I+1)=t;
    end
    K=nv2;
   while K<=J
      J=J-K;
      K=K/2;
    end
    J=J+K;
    I=I+1;
end
X=x            
M=1;
while M<=L
    LE=2^M;
    LE1=LE/2;
    U=1;
    W=exp(-j*pi/LE1);
    J=0;
    while J<=LE1-1
      I=J;
      while I<=N-1
            IP=I+LE1;
            T=x(IP+1)*U;
            x(IP+1)=x(I+1)-T;
            x(I+1)=x(I+1)+T;
            I=I+LE;
      end
      U=U*W;
      J=J+1;
    end
    M=M+1;
end
Y=x;

eight 发表于 2007-9-11 16:00

本帖最后由 wdhd 于 2016-9-8 14:18 编辑

原帖由 tianyoume 于 2007-9-11 15:52 发表
function Y=fft(x)
n=size(x);
for l=0:12,
    if 2^l>=n
      N=2^l;
      L=l;
      break
    end
end
for t=n(2)+1:N,%不足2^L的数在末尾补0
    x(t)=0;
end
x                %x ...
用 matlab 自带的 fft 函数不行吗?

zhangnan3509 发表于 2007-9-11 16:03

回复 #2 eight 的帖子

现在好像流行自己写。

tianyoume 发表于 2007-9-11 16:03

要提取循环平稳频率,需要自己写FFT

tianyoume 发表于 2007-9-13 10:23

上面那段代码是转自别人的,但我试了一下好象与MATLAB带的FFT有出入,不知道源创者是否注意到了。

graduate 发表于 2007-11-30 17:12

%***************具体的FFT的蝶形运算*************************************
for L=1:M                  
    B=2^(L-1);
    for J=0:B-1
      P=J*2^(M-L);
      for k=J:2^L:N-1
            xn(k+1)=xn(k+1)+xn(k+1+B)*(exp((-1)*i*2*pi/N*P)); % DIT:先复乘后加减
            xn(k+1+B)=xn(k+1)-xn(k+1+B)*(exp((-1)*i*2*pi/N*P));
      end
    end
end
页: [1]
查看完整版本: 自己写的FFT总是得不出正确结果。