francexuqu 发表于 2007-12-16 13:51

[分享]用matlab编的fft变换程序

本帖最后由 牛小贱 于 2014-4-19 19:03 编辑

用蝶形算法和码位倒置法编写的fft变换程序,并验证之,有兴趣的可以看看。
function data=myfft(datat,nn,isign)
datat=mybitrevorder(datat,nn);                %这个地方也可以用matlab自带的bitrevorder,有兴趣的还是自己编一下
for i=0:length(datat)-1
    data(2*i+1)=datat(i+1); data(2*i+2)=0;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
n=nn.*2; mmax=2;
while n>mmax;
    istep=2.*mmax; theta=6.28318530717959/(isign*mmax); wtemp=sin(0.5*theta);
    wpr = -2.0*wtemp*wtemp; wpi=sin(theta); wr=1.0; wi=0.0;
    for m=1:2:mmax
      for i=m:istep:n
            j=i+mmax; tempr=wr*data(j)-wi*data(j+1); tempi=wr*data(j+1)+wi*data(j);
         data(j)=data(i)-tempr; data(j+1)=data(i+1)-tempi;
         data(i)=data(i)+tempr; data(i+1)=data(i+1)+tempi;
      end
      wtemp=wr; wr=wtemp*wpr-wi*wpi+wr; wi=wi*wpr+wtemp*wpi+wi;
    end
    mmax=istep;
end
%%%%%%%%%%%%%码位倒置程序%%%%%%%%%%%%%%%%%%%%
functionA=mybitrevorder(A,N)
j=N./2;
for i=2:1:N-1
    if i<j
      temp=A(i); A(i)=A(j+1); A(j+1)=temp;
    end
    k=N/2;
    while j>=k
      j=j-k; k=k/2;
    end
    j=j+k;
end

%%%%%%%%%%%%%进行验证(matlab自带的fft,自己编写的fft,和函数的解析fft三者的比较%%%%%%%
clear all;
isign=1; nn=256; t=linspace(0,3,nn); f=2.*exp(-3.*t);
data=myfft(f,nn,isign);
Ts=t(2)-t(1); Ws=2.*pi./Ts; W=Ws.*(0:nn./2)./nn;
for i=0:nn-1
    datar(i+1)=data(2.*i+1).*Ts; datai(i+1)=data(2.*i+2).*Ts;
end
for i=1:nn./2+1
    co(i)=datar(i); si(i)=datai(i); power(i)=(co(i).^2+si(i).^2).^0.5;
end
F=fft(f); Fp=F(1:nn./2+1)*Ts; Fpr=real(Fp); Fpi=imag(Fp);
W=Ws.*(0:nn./2)./nn; A=abs(Fp); Fa=2./(3+j.*W); B=abs(Fa);
plot(W,B,'o',W,A,'+r',W,power,'g');%%%%%%结果见附件%%%%%%%%%%%

花如月 发表于 2007-12-16 15:21

是自己写的,那可真是高手了。我都写不出来、、、、

sigma665 发表于 2007-12-16 15:39

回复 #1 francexuqu 的帖子

这里卧虎藏龙

francexuqu 发表于 2007-12-16 15:50

其实fft的蝶形算法都已经成行好久了
C++程序都有现成的;
只要搞懂蝶形算法的来龙去脉,和码位倒置的编写就很简单了;
但是因该注意的是C++中的脚标是从0开始;
而matlab是从1开始的;
ps:本人也是刚学习matlab,希望能在这里和大家共同进步;

hyl2323 发表于 2007-12-16 16:41

我用C++写了蝶形算法和另外一种算法,这是matlab版,能发不?

sigma665 发表于 2007-12-16 16:57

回复 #5 hyl2323 的帖子

发吧
哈哈

hyl2323 发表于 2007-12-17 08:41

4#举重若轻,佩服佩服!可能我智商平平,我看懂原理后,也是费了好多脑力才编写出来的。晚上来发!

我以主题帖发了我用C++写的FFT,两种算法,大伙顶一顶哈,第一次发原创!

[ 本帖最后由 ChaChing 于 2010-6-25 21:16 编辑 ]

无水1324 发表于 2007-12-17 19:15

回复 #8 hyl2323 的帖子

支持原创的东西,:lol

hyl2323 发表于 2007-12-17 19:18

无水是湖南长沙的?老乡哦,见过,握握手!

gagmeng 发表于 2007-12-27 17:22

支持原创的东西,

cwt1128 发表于 2008-6-24 17:22

真是太强了,很有用的东西,谢!

鹿鸣 发表于 2010-6-25 19:46

高手如云。

ztx106 发表于 2011-6-8 10:42

我想请问一下,如果对功率谱进行IFFT变换,再对IFFT结果进行FFT变换,那么这两种情况下的功率谱曲线是完全重合的吗?还是FFT变换结果的曲线应该有波动?

yuejunefong 发表于 2011-6-11 02:11

看不明白 还是菜鸟一个

歪猴 发表于 2014-4-19 15:37

{:{39}:}
页: [1] 2
查看完整版本: [分享]用matlab编的fft变换程序