C語(yǔ)言實(shí)現(xiàn)FFT變換

1、#include #include /* 快速福利葉變換C函數(shù)函數(shù)簡(jiǎn)介：此函數(shù)是通用的快速傅里葉變換C語(yǔ)言函數(shù)，移植性強(qiáng)，以下部分不依 賴(lài)硬件。此函數(shù)采用聯(lián)合體的形式表示一個(gè)復(fù)數(shù)，輸入為自然順序的復(fù) 數(shù)（輸入實(shí)數(shù)是可令復(fù)數(shù)虛部為0），輸出為經(jīng)過(guò)FFT變換的自然順序的 復(fù)數(shù)使用說(shuō)明：使用此函數(shù)只需更改宏定義FFT_N的值即可實(shí)現(xiàn)點(diǎn)數(shù)的改變，F(xiàn)FT_N的 應(yīng)該為2的N次方，不滿足此條件時(shí)應(yīng)在后面補(bǔ)0函數(shù)調(diào)用：FFT(s);時(shí) 間：2010-2-20版 本：Ver1.0參考文獻(xiàn)： */#include#define PI 3.1415926535897932384626433832795028841

2、971 /定義圓周率值#define FFT_N 128 /定義福利葉變換的點(diǎn)數(shù)struct compx float real,imag; /定義一個(gè)復(fù)數(shù)結(jié)構(gòu)struct compx sFFT_N; /FFT輸入和輸出：從S1開(kāi)始存放，根據(jù)大小自己定義/*函數(shù)原型：struct compx EE(struct compx b1,struct compx b2) 函數(shù)功能：對(duì)兩個(gè)復(fù)數(shù)進(jìn)行乘法運(yùn)算輸入?yún)?shù)：兩個(gè)以聯(lián)合體定義的復(fù)數(shù)a,b輸出參數(shù)：a和b的乘積，以聯(lián)合體的形式輸出*/struct compx EE(struct compx a,struct compx b) struct compx

3、c; c.real=a.real*b.real-a.imag*b.imag; c.imag=a.real*b.imag+a.imag*b.real; return(c);/*函數(shù)原型：void FFT(struct compx *xin,int N)函數(shù)功能：對(duì)輸入的復(fù)數(shù)組進(jìn)行快速傅里葉變換（FFT）輸入?yún)?shù)：*xin復(fù)數(shù)結(jié)構(gòu)體組的首地址指針，struct型*/void FFT(struct compx *xin) int f,m,nv2,nm1,i,k,l,j=0; struct compx u,w,t; nv2=FFT_N/2; /變址運(yùn)算，即把自然順序變成倒位序，采用雷德算法 nm1=F

4、FT_N-1; for(i=0;inm1;i+) if(ij) /如果ij,即進(jìn)行變址 t=xinj; xinj=xini; xini=t; k=nv2; /求j的下一個(gè)倒位序 while(k=j) /如果k=j,表示j的最高位為1 j=j-k; /把最高位變成0 k=k/2; /k/2，比較次高位，依次類(lèi)推，逐個(gè)比較，直到某個(gè)位為0 j=j+k; /把0改為1 int le,lei,ip; /FFT運(yùn)算核，使用蝶形運(yùn)算完成FFT運(yùn)算 f=FFT_N; for(l=1;(f=f/2)!=1;l+) /計(jì)算l的值，即計(jì)算蝶形級(jí)數(shù) ; for(m=1;m=l;m+) / 控制蝶形結(jié)級(jí)數(shù) /m表示第

5、m級(jí)蝶形，l為蝶形級(jí)總數(shù)l=log（2）N le=2(m-1); /le蝶形結(jié)距離，即第m級(jí)蝶形的蝶形結(jié)相距l(xiāng)e點(diǎn) lei=le/2; /同一蝶形結(jié)中參加運(yùn)算的兩點(diǎn)的距離 u.real=1.0; /u為蝶形結(jié)運(yùn)算系數(shù)，初始值為1 u.imag=0.0; w.real=cos(PI/lei); /w為系數(shù)商，即當(dāng)前系數(shù)與前一個(gè)系數(shù)的商 w.imag=-sin(PI/lei); for(j=0;j=lei-1;j+) /控制計(jì)算不同種蝶形結(jié)，即計(jì)算系數(shù)不同的蝶形結(jié) for(i=j;i=FFT_N-1;i=i+le) /控制同一蝶形結(jié)運(yùn)算，即計(jì)算系數(shù)相同蝶形結(jié) ip=i+lei; /i，ip分別表示

6、參加蝶形運(yùn)算的兩個(gè)節(jié)點(diǎn) t=EE(xinip,u); /蝶形運(yùn)算，詳見(jiàn)公式 xinip.real=xini.real-t.real; xinip.imag=xini.imag-t.imag; xini.real=xini.real+t.real; xini.imag=xini.imag+t.imag; u=EE(u,w); /改變系數(shù)，進(jìn)行下一個(gè)蝶形運(yùn)算 /*函數(shù)原型：void main() 函數(shù)功能：測(cè)試FFT變換，演示函數(shù)使用方法輸入?yún)?shù)：無(wú)輸出參數(shù)：無(wú)*/void main() int i; for(i=0;iFFT_N;i+) /給結(jié)構(gòu)體賦值 si.real=sin(2*3.14159

7、2653589793*i/FFT_N); /實(shí)部為正弦波FFT_N點(diǎn)采樣，賦值為1 si.imag=0; /虛部為0 FFT(s); /進(jìn)行快速福利葉變換 for(i=0;iFFT_N;i+) /求變換后結(jié)果的模值，存入復(fù)數(shù)的實(shí)部部分 si.real=sqrt(si.real*si.real+si.imag*si.imag); while(1); #include #include /* 快速福利葉變換C程序包函數(shù)簡(jiǎn)介：此程序包是通用的快速傅里葉變換C語(yǔ)言函數(shù)，移植性強(qiáng)，以下部分不依 賴(lài)硬件。此程序包采用聯(lián)合體的形式表示一個(gè)復(fù)數(shù)，輸入為自然順序的復(fù) 數(shù)（輸入實(shí)數(shù)是可令復(fù)數(shù)虛部為0），輸出為經(jīng)過(guò)

8、FFT變換的自然順序的 復(fù)數(shù).此程序包可在初始化時(shí)調(diào)用create_sin_tab()函數(shù)創(chuàng)建正弦函數(shù)表， 以后的可采用查表法計(jì)算耗時(shí)較多的sin和cos運(yùn)算，加快可計(jì)算速度使用說(shuō)明：使用此函數(shù)只需更改宏定義FFT_N的值即可實(shí)現(xiàn)點(diǎn)數(shù)的改變，F(xiàn)FT_N的 應(yīng)該為2的N次方，不滿足此條件時(shí)應(yīng)在后面補(bǔ)0。若使用查表法計(jì)算sin值和 cos值，應(yīng)在調(diào)用FFT函數(shù)前調(diào)用create_sin_tab()函數(shù)創(chuàng)建正弦表函數(shù)調(diào)用：FFT(s);時(shí) 間：2010-2-20版 本：Ver1.1參考文獻(xiàn)： */#include#define FFT_N 128 /定義福利葉變換的點(diǎn)數(shù)#define PI 3.14

9、15926535897932384626433832795028841971 /定義圓周率值struct compx float real,imag; /定義一個(gè)復(fù)數(shù)結(jié)構(gòu)struct compx sFFT_N; /FFT輸入和輸出：從S0開(kāi)始存放，根據(jù)大小自己定義float SIN_TABFFT_N/2; /定義正弦表的存放空間/*函數(shù)原型：struct compx EE(struct compx b1,struct compx b2) 函數(shù)功能：對(duì)兩個(gè)復(fù)數(shù)進(jìn)行乘法運(yùn)算輸入?yún)?shù)：兩個(gè)以聯(lián)合體定義的復(fù)數(shù)a,b輸出參數(shù)：a和b的乘積，以聯(lián)合體的形式輸出*/struct compx EE(struc

10、t compx a,struct compx b) struct compx c; c.real=a.real*b.real-a.imag*b.imag; c.imag=a.real*b.imag+a.imag*b.real; return(c);/*函數(shù)原型：void create_sin_tab(float *sin_t)函數(shù)功能：創(chuàng)建一個(gè)正弦采樣表，采樣點(diǎn)數(shù)與福利葉變換點(diǎn)數(shù)相同輸入?yún)?shù)：*sin_t存放正弦表的數(shù)組指針輸出參數(shù)：無(wú)*/void create_sin_tab(float *sin_t) int i; for(i=0;i=0&n=FFT_N/2&n2*PI) pi2-=2*P

11、I; a=sin_tab(pi2); return a;/*函數(shù)原型：void FFT(struct compx *xin,int N)函數(shù)功能：對(duì)輸入的復(fù)數(shù)組進(jìn)行快速傅里葉變換（FFT）輸入?yún)?shù)：*xin復(fù)數(shù)結(jié)構(gòu)體組的首地址指針，struct型輸出參數(shù)：無(wú)*/void FFT(struct compx *xin) int f,m,nv2,nm1,i,k,l,j=0; struct compx u,w,t; nv2=FFT_N/2; /變址運(yùn)算，即把自然順序變成倒位序，采用雷德算法 nm1=FFT_N-1; for(i=0;inm1;i+) if(ij) /如果ij,即進(jìn)行變址 t=xinj;

12、 xinj=xini; xini=t; k=nv2; /求j的下一個(gè)倒位序 while(k=j) /如果k=j,表示j的最高位為1 j=j-k; /把最高位變成0 k=k/2; /k/2，比較次高位，依次類(lèi)推，逐個(gè)比較，直到某個(gè)位為0 j=j+k; /把0改為1 int le,lei,ip; /FFT運(yùn)算核，使用蝶形運(yùn)算完成FFT運(yùn)算 f=FFT_N; for(l=1;(f=f/2)!=1;l+) /計(jì)算l的值，即計(jì)算蝶形級(jí)數(shù) ; for(m=1;m=l;m+) / 控制蝶形結(jié)級(jí)數(shù) /m表示第m級(jí)蝶形，l為蝶形級(jí)總數(shù)l=log（2）N le=2(m-1); /le蝶形結(jié)距離，即第m級(jí)蝶形的蝶形

13、結(jié)相距l(xiāng)e點(diǎn) lei=le/2; /同一蝶形結(jié)中參加運(yùn)算的兩點(diǎn)的距離 u.real=1.0; /u為蝶形結(jié)運(yùn)算系數(shù)，初始值為1 u.imag=0.0; /w.real=cos(PI/lei); /不適用查表法計(jì)算sin值和cos值 / w.imag=-sin(PI/lei); w.real=cos_tab(PI/lei); /w為系數(shù)商，即當(dāng)前系數(shù)與前一個(gè)系數(shù)的商 w.imag=-sin_tab(PI/lei); for(j=0;j=lei-1;j+) /控制計(jì)算不同種蝶形結(jié)，即計(jì)算系數(shù)不同的蝶形結(jié) for(i=j;i=FFT_N-1;i=i+le) /控制同一蝶形結(jié)運(yùn)算，即計(jì)算系數(shù)相同蝶形結(jié)

14、 ip=i+lei; /i，ip分別表示參加蝶形運(yùn)算的兩個(gè)節(jié)點(diǎn) t=EE(xinip,u); /蝶形運(yùn)算，詳見(jiàn)公式 xinip.real=xini.real-t.real; xinip.imag=xini.imag-t.imag; xini.real=xini.real+t.real; xini.imag=xini.imag+t.imag; u=EE(u,w); /改變系數(shù)，進(jìn)行下一個(gè)蝶形運(yùn)算 /*函數(shù)原型：void main() 函數(shù)功能：測(cè)試FFT變換，演示函數(shù)使用方法輸入?yún)?shù)：無(wú)輸出參數(shù)：無(wú)*/void main() int i; create_sin_tab(SIN_TAB); for

15、(i=0;iFFT_N;i+) /給結(jié)構(gòu)體賦值 si.real=sin(2*3.141592653589793*i/FFT_N); /實(shí)部為正弦波FFT_N點(diǎn)采樣，賦值為1 si.imag=0; /虛部為0 FFT(s); /進(jìn)行快速福利葉變換 for(i=0;iFFT_N;i+) /求變換后結(jié)果的模值，存入復(fù)數(shù)的實(shí)部部分 si.real=sqrt(si.real*si.real+si.imag*si.imag); while(1); #include #include /* 快速福利葉變換C程序包函數(shù)簡(jiǎn)介：此程序包是通用的快速傅里葉變換C語(yǔ)言函數(shù)，移植性強(qiáng)，以下部分不依 賴(lài)硬件。此程序包采用

16、聯(lián)合體的形式表示一個(gè)復(fù)數(shù)，輸入為自然順序的復(fù) 數(shù)（輸入實(shí)數(shù)是可令復(fù)數(shù)虛部為0），輸出為經(jīng)過(guò)FFT變換的自然順序的 復(fù)數(shù).此程序包可在初始化時(shí)調(diào)用create_sin_tab()函數(shù)創(chuàng)建正弦函數(shù)表， 以后的可采用查表法計(jì)算耗時(shí)較多的sin和cos運(yùn)算，加快可計(jì)算速度.與 Ver1.1版相比較，Ver1.2版在創(chuàng)建正弦表時(shí)只建立了1/4個(gè)正弦波的采樣值， 相比之下節(jié)省了FFT_N/4個(gè)存儲(chǔ)空間使用說(shuō)明：使用此函數(shù)只需更改宏定義FFT_N的值即可實(shí)現(xiàn)點(diǎn)數(shù)的改變，F(xiàn)FT_N的 應(yīng)該為2的N次方，不滿足此條件時(shí)應(yīng)在后面補(bǔ)0。若使用查表法計(jì)算sin值和 cos值，應(yīng)在調(diào)用FFT函數(shù)前調(diào)用create_si

17、n_tab()函數(shù)創(chuàng)建正弦表函數(shù)調(diào)用：FFT(s);時(shí) 間：2010-2-20版 本：Ver1.2參考文獻(xiàn)： */#include#define FFT_N 128 /定義福利葉變換的點(diǎn)數(shù)#define PI 3.1415926535897932384626433832795028841971 /定義圓周率值struct compx float real,imag; /定義一個(gè)復(fù)數(shù)結(jié)構(gòu)struct compx sFFT_N; /FFT輸入和輸出：從S0開(kāi)始存放，根據(jù)大小自己定義float SIN_TABFFT_N/4+1; /定義正弦表的存放空間/*函數(shù)原型：struct compx EE(s

18、truct compx b1,struct compx b2) 函數(shù)功能：對(duì)兩個(gè)復(fù)數(shù)進(jìn)行乘法運(yùn)算輸入?yún)?shù)：兩個(gè)以聯(lián)合體定義的復(fù)數(shù)a,b輸出參數(shù)：a和b的乘積，以聯(lián)合體的形式輸出*/struct compx EE(struct compx a,struct compx b) struct compx c; c.real=a.real*b.real-a.imag*b.imag; c.imag=a.real*b.imag+a.imag*b.real; return(c);/*函數(shù)原型：void create_sin_tab(float *sin_t)函數(shù)功能：創(chuàng)建一個(gè)正弦采樣表，采樣點(diǎn)數(shù)與福利葉變換

19、點(diǎn)數(shù)相同輸入?yún)?shù)：*sin_t存放正弦表的數(shù)組指針輸出參數(shù)：無(wú)*/void create_sin_tab(float *sin_t) int i; for(i=0;i=0&nFFT_N/4&n=FFT_N/2&n=3*FFT_N/4&n2*PI) pi2-=2*PI; a=sin_tab(pi2); return a;/*函數(shù)原型：void FFT(struct compx *xin,int N)函數(shù)功能：對(duì)輸入的復(fù)數(shù)組進(jìn)行快速傅里葉變換（FFT）輸入?yún)?shù)：*xin復(fù)數(shù)結(jié)構(gòu)體組的首地址指針，struct型輸出參數(shù)：無(wú)*/void FFT(struct compx *xin) int f,m,n

20、v2,nm1,i,k,l,j=0; struct compx u,w,t; nv2=FFT_N/2; /變址運(yùn)算，即把自然順序變成倒位序，采用雷德算法 nm1=FFT_N-1; for(i=0;inm1;i+) if(ij) /如果ij,即進(jìn)行變址 t=xinj; xinj=xini; xini=t; k=nv2; /求j的下一個(gè)倒位序 while(k=j) /如果k=j,表示j的最高位為1 j=j-k; /把最高位變成0 k=k/2; /k/2，比較次高位，依次類(lèi)推，逐個(gè)比較，直到某個(gè)位為0 j=j+k; /把0改為1 int le,lei,ip; /FFT運(yùn)算核，使用蝶形運(yùn)算完成FFT運(yùn)算 f=FFT_N; for(l=1;(f=f/2)!=1;l+) /計(jì)算l的值，即計(jì)算蝶形級(jí)數(shù) ; for(m=1;m=l;m+) / 控制蝶形結(jié)級(jí)數(shù) /m表示第m級(jí)蝶形，l為蝶形級(jí)總數(shù)l=log（2）N le=2(m-1); /le蝶形結(jié)距離，即第m級(jí)蝶形的蝶形結(jié)相距l(xiāng)e點(diǎn) lei=le/2; /同一蝶形結(jié)中參加運(yùn)算的兩點(diǎn)的距離 u.real=1.0; /u為蝶形結(jié)運(yùn)算系數(shù)，初始值為1 u.imag=0.0; /w.real=cos(PI/lei); /不適用查表法計(jì)算sin值和cos值 / w.imag