(完整word版)快速傅里葉變換(FFT)原理及源程序(20210413200648)

1、測(cè)試信號(hào)分析及處理課程作業(yè) 快速傅里葉變換 一、程序設(shè)計(jì)思路 快速傅里葉變換的目的是減少運(yùn)算量，其用到的方法是分級(jí)進(jìn)行運(yùn)算。全部 計(jì)算分解為M級(jí)，其中M =|og2N ；在輸入序列Xi中是按碼位倒序排列的， 輸出序列X k是按順序排列；每級(jí)包含 N個(gè)蝶形單元，第i級(jí)有個(gè)群，每個(gè) 2 2 群有21 個(gè)蝶形單元；每個(gè)蝶形單元都包含乘 WN和-wN系數(shù)的運(yùn)算，每個(gè)蝶形 單元數(shù)據(jù)的間隔為2iJ ,i為第i級(jí)；同一級(jí)中各個(gè)群的系數(shù) W分布規(guī)律完全相 同。 將輸入序列Xi按碼位倒序排列時(shí)，用到的是倒序算法一一雷德算法。自然 序排列的二進(jìn)制數(shù)，其下面一個(gè)數(shù)總比上面的數(shù)大 1,而倒序二進(jìn)制數(shù)的下面一 個(gè)數(shù)是上

2、面一個(gè)數(shù)在最高位加1并由高位向低位僅為而得到的。 若已知某數(shù)的倒序數(shù)是J ,求下一個(gè)倒序數(shù)，應(yīng)先判斷J的最高位是否為0， 與k = N進(jìn)行比較即可得到結(jié)果。如果k J，說(shuō)明最高位為0,應(yīng)把其變成1, 2 即J -，這樣就得到倒序數(shù)了。如果J，即J的最高位為1，將最高位化為 2 0,即J - N ,再判斷次高位；與k =進(jìn)行比較，若為0,將其變位1,即J , 244 即得到倒序數(shù)，如果次高位為1,將其化為0，再判斷下一位即從高位到低 位依次判斷其是否為1,為1將其變位0,若這一位為0,將其變位1,即可得到 倒序數(shù)。若倒序數(shù)小于順序數(shù)，進(jìn)行換位，否則不變，防治重復(fù)交換，變回原數(shù)。 注：因?yàn)?的倒序

3、數(shù)為0,所以可從1開(kāi)始進(jìn)行求解。 二、程序設(shè)計(jì)框圖 (1)倒序算法一一雷德算法流程圖 三、FFT源程序 void fft(x, n) int n; double x; int i,j,k,l,m,n1,n2; double c,c1,e,s,s1,t,tr; for(j=1,i=1;i n/2;i+) m=i; /得到流程圖的共幾級(jí) /如果ij,即進(jìn)行變址 j=2*j; if(j=n )break; n仁n-1; for(j=0,i=0;i n1;i+) if(ivj) tr=xj; xj=xi; xi=tr; k=n/2; while(k(j+1) /求j的下一個(gè)倒位序 /如果k(j+1),

4、表示j的最高 位為1 j=j-k;/把最高位變成0 k=k/2;/k/2，比較次高位，依次類(lèi)推，逐個(gè)比較，直到某個(gè)位為0 j=j+k;/ 把 0 改為 1 for(i=0;i n ;i+=2) tr=xi; xi=t 葉xi+1; xi+1=tr-xi+1; n2=1; for(l=1;l=m;l+)/ 控制蝶形結(jié)級(jí)數(shù) n 4=n2; n 2=2* n4; n1=2* n2; e=6.28318530718/n1; for(i=0;in;i+=n1)/控制同一蝶形結(jié)運(yùn)算，即計(jì)算系數(shù)相同蝶形結(jié) tr=xi; xi=t 葉xi+n2; xi+n 2=tr-xi+n2; xi+n2+n4=-xi+n

5、2+n4; a=e; for(j=2;j=(n4-1);j+)/控制計(jì)算不同種蝶形結(jié)，即計(jì)算系數(shù)不同 的蝶形結(jié) i仁 i+j; i2=i-j+n2; i3=i+j+n2; i4=i-j+n1; cc=cos(a); ss=s in( a); a=a+e; t1=cc*xi3+ss*xi4; t2=ss*xi3-cc*xi4; xi4=xi2-t2; xi3=-xi2-t2; xi2=xi1-t1; xi1=xi1+t1; 四、計(jì)算實(shí)例及運(yùn)行結(jié)果 設(shè)輸入序列x(i)為 x(i) =sin200二i：t,(i =0,1,2,，n1) 其離散傅里葉變換為 N 二 X(k)八 x(i)wNk,(k =

6、0,1,2,n-1) i ：5 .j 這里W =e n。選n=512,計(jì)算離散傅里葉變換X(k) 所用軟件為T(mén)urbo c 2.0,操作界面如圖1所示 圖1 Turbo c 2.0操作界面 程序運(yùn)行結(jié)束后的界面如圖2所示 回|2 國(guó)3 C:Wi ndc32c - . exe 5.1911638+J 3-6933914+J 2.B497513+J 2.3099713+J 1-5360773+J 1 .6627248+J 1-4549735+J 1.2924530+J L.1 double x; int i,j,k,l,i1,i2,i3,i4,n4,m,n1,n2; double a,e,cc,s

7、s,tr,t1,t2; for(j=1,i=1;i n/2;i+) m=i; j=2*j; if(j=n )break; n仁n-1; for(j=0,i=0;i n1;i+) if(ivj) tr=xj; xj=xi; xi=tr; k=n/2; while(k(j+1) j=j-k; k=k/2; j=j+k; for(i=0;i n ;i+=2) tr=xi; xi=t 葉xi+1; xi+1=tr-xi+1; n2=1; for(l=1;l=m;l+) n 4=n2; n 2=2* n4; n1=2* n2; e=6.28318530718/n1; for(i=0;i n ;i+=n1

8、) tr=xi; xi=t 葉xi+n2; xi+n 2=tr-xi+n2; xi+n2+n4=-xi+n2+n4; a=e; for(j=2;j=( n4-1);j+) i仁 i+j; i2=i-j+n2; i3=i+j+n2; i4=i-j+n1; cc=cos(a); ss=s in( a); a=a+e; t1=cc*xi3+ss*xi4; t2=ss*xi3-cc*xi4; xi4=xi2-t2; xi3=-xi2-t2; xi2=xi1-t1; xi1=xi1+t1; mai n() FILE *p; int i,j, n; double dt=0.001; double x512

9、; p=fope n(d:123.c,w); n=512; for(i=0;i n;i+) xi=si n(200*pi*i*dt); for(i=0;i n;i+) fprin tf(p,%10.7f,xi); fprin tf(p,n); prin tf(%10.7f,xi); prin tf(n); fft(x, n); fprin tf(p,n DISCRETE FOURIER TRANSFORMS); prin tf(n DISCRETE FOURIER TRANSFORMS); fprin tf(p,%10.7f,x0); prin tf(%10.7f,x0); fprin tf(

10、p,%10.7f+J%10.7fn,x1,x n-1); prin tf(%10.7f+J%10.7fn,x1,x n-1); for(i=2;i n/2;i+=2) fprin tf(p,%10.7f+J%10.7f,xi,x n-i); fprin tf(p,%10.7f+J%10.7f,xi+1,x n-i-1); fprin tf(p,n); prin tf(%10.7f+J%10.7f,xi,x n-i); prin tf(%10.7f+J%10.7f,xi+1,x n-i-1); prin tf(n); fprin tf(p,%10.7f,x n/2); prin tf(%10.7f,x n/2); fprin tf(p,%10.7f+J%10.7fn,x n/2-1,-x n/2+1); for(i=2;i n/2;i+=2) fprin tf(p,%10.7f+J%10.7f,x n/2-i,-x n/2+i); fprin tf(p,%10.7f+J%10.7f,x n/2-i-1,-x n/2+i+1); fprin tf(p,n); prin tf(%10.7