fft演算法
『壹』 什麼是FFT
計算離散傅里葉變換的一種快速演算法,簡稱FFT(Fast Fourier Transform)。快速傅里葉變換是1965年由J.W.庫利和T.W.圖基提出的。採用這種演算法能使計算機計算離散傅里葉變換所需要的乘法次數大為減少,特別是被變換的抽樣點數N越多,FFT演算法計算量的節省就越顯著。
FFT 的出現,使信號分析從時域分析向頻域分析成為可能,極大地推動了信號分析在各領域的實際應用。
http://ke..com/view/1006229.htm
『貳』 FFT演算法如何加窗
沒錯。
截取一段信號相當於是加矩形窗。所以加漢明窗就不難理解了。
『叄』 怎麼用C語言實現FFT演算法 呀
float ar[1024],ai[1024];/* 原始數據實部,虛部 */
float a[2050];
void fft(int nn) /* nn數據長度 */
{
int n1,n2,i,j,k,l,m,s,l1;
float t1,t2,x,y;
float w1,w2,u1,u2,z;
float fsin[10]={0.000000,1.000000,0.707107,0.3826834,0.1950903,0.09801713,0.04906767,0.02454123,0.01227154,0.00613588,};
float fcos[10]={-1.000000,0.000000,0.7071068,0.9238796,0.9807853,0.99518472,0.99879545,0.9996988,0.9999247,0.9999812,};
switch(nn)
{
case 1024: s=10; break;
case 512: s=9; break;
case 256: s=8; break;
}
n1=nn/2; n2=nn-1;
j=1;
for(i=1;i<=nn;i++)
{
a[2*i]=ar[i-1];
a[2*i+1]=ai[i-1];
}
for(l=1;l<n2;l++)
{
if(l<j)
{
t1=a[2*j];
t2=a[2*j+1];
a[2*j]=a[2*l];
a[2*j+1]=a[2*l+1];
a[2*l]=t1;
a[2*l+1]=t2;
}
k=n1;
while (k<j)
{
j=j-k;
k=k/2;
}
j=j+k;
}
for(i=1;i<=s;i++)
{
u1=1;
u2=0;
m=(1<<i);
k=m>>1;
w1=fcos[i-1];
w2=-fsin[i-1];
for(j=1;j<=k;j++)
{
for(l=j;l<nn;l=l+m)
{
l1=l+k;
t1=a[2*l1]*u1-a[2*l1+1]*u2;
t2=a[2*l1]*u2+a[2*l1+1]*u1;
a[2*l1]=a[2*l]-t1;
a[2*l1+1]=a[2*l+1]-t2;
a[2*l]=a[2*l]+t1;
a[2*l+1]=a[2*l+1]+t2;
}
z=u1*w1-u2*w2;
u2=u1*w2+u2*w1;
u1=z;
}
}
for(i=1;i<=nn/2;i++)
{
ar[i]=4*a[2*i+2]/nn; /* 實部 */
ai[i]=-4*a[2*i+3]/nn; /* 虛部 */
a[i]=4*sqrt(ar[i]*ar[i]+ai[i]*ai[i]); /* 幅值 */
}
}
(http://..com/question/284943905.html?an=0&si=2)
『肆』 怎樣用C語言實現FFT演算法啊
1、二維FFT相當於對行和列分別進行一維FFT運算。具體的實現辦法如下:
先對各行逐一進行一維FFT,然後再對變換後的新矩陣的各列逐一進行一維FFT。相應的偽代碼如下所示:
for (int i=0; i<M; i++)
FFT_1D(ROW[i],N);
for (int j=0; j<N; j++)
FFT_1D(COL[j],M);
其中,ROW[i]表示矩陣的第i行。注意這只是一個簡單的記法,並不能完全照抄。還需要通過一些語句來生成各行的數據。同理,COL[i]是對矩陣的第i列的一種簡單表示方法。
所以,關鍵是一維FFT演算法的實現。
2、常式:
#include<stdio.h>
#include<math.h>
#include<stdlib.h>
#defineN1000
/*定義復數類型*/
typedefstruct{
doublereal;
doubleimg;
}complex;
complexx[N],*W;/*輸入序列,變換核*/
intsize_x=0;/*輸入序列的大小,在本程序中僅限2的次冪*/
doublePI;/*圓周率*/
voidfft();/*快速傅里葉變換*/
voidinitW();/*初始化變換核*/
voidchange();/*變址*/
voidadd(complex,complex,complex*);/*復數加法*/
voidmul(complex,complex,complex*);/*復數乘法*/
voidsub(complex,complex,complex*);/*復數減法*/
voidoutput();
intmain(){
inti;/*輸出結果*/
system("cls");
PI=atan(1)*4;
printf("Pleaseinputthesizeofx: ");
scanf("%d",&size_x);
printf("Pleaseinputthedatainx[N]: ");
for(i=0;i<size_x;i++)
scanf("%lf%lf",&x[i].real,&x[i].img);
initW();
fft();
output();
return0;
}
/*快速傅里葉變換*/
voidfft(){
inti=0,j=0,k=0,l=0;
complexup,down,proct;
change();
for(i=0;i<log(size_x)/log(2);i++){/*一級蝶形運算*/
l=1<<i;
for(j=0;j<size_x;j+=2*l){/*一組蝶形運算*/
for(k=0;k<l;k++){/*一個蝶形運算*/
mul(x[j+k+l],W[size_x*k/2/l],&proct);
add(x[j+k],proct,&up);
sub(x[j+k],proct,&down);
x[j+k]=up;
x[j+k+l]=down;
}
}
}
}
/*初始化變換核*/
voidinitW(){
inti;
W=(complex*)malloc(sizeof(complex)*size_x);
for(i=0;i<size_x;i++){
W[i].real=cos(2*PI/size_x*i);
W[i].img=-1*sin(2*PI/size_x*i);
}
}
/*變址計算,將x(n)碼位倒置*/
voidchange(){
complextemp;
unsignedshorti=0,j=0,k=0;
doublet;
for(i=0;i<size_x;i++){
k=i;j=0;
t=(log(size_x)/log(2));
while((t--)>0){
j=j<<1;
j|=(k&1);
k=k>>1;
}
if(j>i){
temp=x[i];
x[i]=x[j];
x[j]=temp;
}
}
}
/*輸出傅里葉變換的結果*/
voidoutput(){
inti;
printf("Theresultareasfollows ");
for(i=0;i<size_x;i++){
printf("%.4f",x[i].real);
if(x[i].img>=0.0001)printf("+%.4fj ",x[i].img);
elseif(fabs(x[i].img)<0.0001)printf(" ");
elseprintf("%.4fj ",x[i].img);
}
}
voidadd(complexa,complexb,complex*c){
c->real=a.real+b.real;
c->img=a.img+b.img;
}
voidmul(complexa,complexb,complex*c){
c->real=a.real*b.real-a.img*b.img;
c->img=a.real*b.img+a.img*b.real;
}
voidsub(complexa,complexb,complex*c){
c->real=a.real-b.real;
c->img=a.img-b.img;
}
『伍』 什麼是FFT演算法DSP是什麼
FFT是快速傅里葉變換( Fast Fourier Transform )
DSP是數字信號處理 ( Digital Signal Processing )
『陸』 fft演算法為什麼算的快
離散Fourier變換按照最平凡的方式來算需要Θ(N^2)的運算量, 但FFT可以把復雜度降到Θ(N*logN), 並且NlogN之前的系數也挺小, 所以就要快很多
『柒』 如何實現3點4點和5點的FFT演算法有知道的麻煩講一下感激不盡啊!
FFT 用於 2 的 整數次方 的 點數 的 時序 信號 做快速 傅里葉 變換。一般要求輸入 大量的 時序點 例如幾百幾千點,計算 譜。 沒聽說過 用於 3,4,5 點的。
當然,你可以用 內插法,把它插成 512 點,1024點 後再 做 FFT. 不過,有無實用價值不知道了。
『捌』 基數為2的FFT演算法
從上節所述,FFT演算法快速的關鍵在於將原來傅氏矩陣分解為每一行僅含有兩個非零項l與Wi的矩陣的乘積。下面用基數為2,即N=2n的情形討論矩陣的分解過程.並主要按時間分解的情況討論。
按時間分解的FFT演算法
設N=2n,n為正整數。考慮輸入序列x0(l)的DFT
物探數字信號分析與處理技術
將l與m用二進製表示
物探數字信號分析與處理技術
將(7-2-2)代入(7-2-1)中,得到
物探數字信號分析與處理技術
為了說明問題,我們取N=8,於是從(7-2-2)得到
物探數字信號分析與處理技術
從(7-2-4)和(7-2-3)得到
物探數字信號分析與處理技術
將(7-2-5)中的W的指數按時間l進行分解,有
物探數字信號分析與處理技術
因為
物探數字信號分析與處理技術
故從(7-2-6)得到
物探數字信號分析與處理技術
將上式代入(7-2-5)中得到
物探數字信號分析與處理技術
物探數字信號分析與處理技術
我們在公式(7-2-7)中由里往外求和,並置
物探數字信號分析與處理技術
於是得到
物探數字信號分析與處理技術
首先寫出(7-2-8)的所有式子
物探數字信號分析與處理技術
將方程組(7-2-12)寫成矩陣形式,得到
物探數字信號分析與處理技術
我們看到(7-2-13)中的方陣,正好是(7-1-13)中的方陣,也就是(7-1-12)中被分解出來的第3個矩陣,只不過這里的x1(l)與x0(l)中的l是用二進制數表示而已。
再寫出(7-2-9)的所有式子,得到
物探數字信號分析與處理技術
將方程組寫成矩陣形式,則有
物探數字信號分析與處理技術
顯然(7-2-15)中的矩陣,正好是(7-1-14)中的方陣,也就是(7-1-12)中被分解出來的第二個矩陣,只不過這里的x2(l)與x1(l)是用二進制數表示而已,最後將(7-2-10)中的全部式子寫出,得到
物探數字信號分析與處理技術
將方程組(7-2-16)寫成矩陣形式,則有
物探數字信號分析與處理技術
顯然,(7-2-17)中的方陣正是(7-1-15)中的方陣,也就是(7-1-12)中被分解出來的第1個矩陣,只不過這里的x3(l)與x2(l)中的l是用二進制數表示。
由此可見,(7-2-7)中由里往外的三個求和式(7-2-8)、(7-2-9)及(7-2-10),完全確定了(7-1-12)中三個被分解的矩陣因子。求和得到的最終結果x3(m0,m1,m2),與我們所要求的X(m2,m1,m0)正好是逆序的。
到此為止,我們就看到(7-1-11)中的方陣是怎樣被分解成三個方陣因子的。對於N=8,方程(7-2-8)~(7-2-11)就是計算DFT的FFT演算法。為了對FFT演算法有一個直觀的了解並便於編製程序,我們以N=8為例,畫出其流程圖(圖7-2-1)。
根據(7-2-13),將其中的W4用-W0代替,畫出從x0(r)到x1(r)的流程圖。這一迭代過程用符號#1表示;再根據(7-2-15),將其中的與W4、W6分別換成-W0與-W2,畫出從x1(r)到x2(r)的流程圖,這一迭代過程記為#2;最後,根據(7-2-17),將其中的W4、W6、W5、W7分別換成-W0、-W2、-W1、-W3,畫出流程圖7-2-3合並圖7-2-1~7-2-3,就得到從x0(r)到x3(r)的完整流程圖7-2-4。
圖7-2-1 第一次遞推
圖7-2-2 第二次遞推
在圖7-2-5中,畫出N=16=24的FFT演算法流程圖:
根據從x0到譜X的FFT演算法流程圖7-2-4與圖7-2-5,我們總結出如下幾點:
(1)從x0到終值xr的最大迭代次數r,由r=log2N確定。
例如,N=8,最大迭代次數r=3;N=16,最大迭代次數r=4。
(2)在第r次迭代中,要乘的W因子為
圖7-2-3 第三次遞推
例如,N=8,在第一次迭代中,要乘因子W0;在第二次迭代中要乘因子W0,W1,W2,W3。
(3)在第r次迭代中,包含2r-1個組,每個組
包含 。例如N=8,第一次迭代r=1,有
一個組,每組包含8個x(s);在第二次迭代中包含2個組,每組包含4個x(s);第三次迭代中包含4個組,每組2個x(s)。
圖7-2-4 x0(r)到x3(r)的完整流程圖
(4)在第r次迭代中,各組包含的W因子各不相同,且每一組僅包含一種類型的因子 ,此因子在組中一半數為正,另一半數為負。例如N=8,第二次迭代中,第一組包含因子W0,且在該組中一半數為正,另一半數為負;第二組包含因子W2,在該組中也是一半數為正,另一半數為負。
(5)在第r次迭代中,各組包含的W因子除正負號外類型均相同。所以只須確定每組中第一個因子,之後按半數反號,即得到所求W因子。具體做法是,在每組第一個因子WSN2r對應的xr(k)中,將k表示成n位的二進制數,n=log2N,並把這個二進制數右移n-r位,左邊空出的地方添零補足n位,之後再將此n位二進制數逆位,即得到所求W因子的指數。例如,N=8,迭代#2包含兩組,每組包含4個x2(k),第二組第一個因子W對應於x2(4)。將4表示成3位的二進制數為100,把它右移1位成10,右邊添零補成3位為010,逆位仍為010,故所求因子為W2,第2組全部W因子為W2,W2,-W2,-W2。又如,N=16,迭代#3中包含4個組,每組包含4個x3(k),第4組第1個因子W對應於x3(12)。將12表示成4位的二進制數為1100。右移1位變成110,將左邊空處添零補成4位為0110,逆位仍為0110,故所求因子為W6,從而第4組全部W因子為W6,W6,-W6與-W6。
圖7-2-5 N=16=24的FFT演算法流程圖
(6)如果已知N=2的FFT演算法,容易從它求得n=2的FFT演算法。具體作法是,在n=2n-1的FFT演算法中,將所有xr(l)的個數加倍,所有W的個數及其乘冪加倍,就得N=2n中前n-1次迭代的全部結果。之後,將新得到的第n-1次迭代中乘冪相同的W個數減半,就是第n次迭代中前2n/2個W,將這些W的乘冪依次加1,就得到後2n/2個W。例如N=16的前3次迭代,都是N=8的三次迭代中所有xr(l)的個數加倍,W的個數及其乘冪加倍的結果。再將N=16的第三次迭代中乘冪相同的W個數減半,就是第4次迭代中前8個W。
(7)在第r-1次迭代中的xr-1(i)與xr-1(j)僅用於計算r次迭代中的xr(i)與xr(j),而不會用於計算任何其它的xr(k)與xr(l)。例如N=16的第二次迭代中的x2(0)與x2(2),只用於計算第三次迭代中的x3(0)與x3(2);第三次迭代中的x3(8)與x3(9)也只用於計算第四次迭代中的x4(8)與x4(9)。因此,我們可以把第r次迭代中的xr(i)與xr(j),存放到r-1次迭代xr-1(i)與xr-1(j)所佔用的原存儲單元中去,這樣,所需要的計算機存儲容量就只限於原數據序列占據的單元數。如果是復序列的話,存儲單元要加倍。
(8)上述FFT演算法也可用於計算逆離散傅氏變換(IDFT)(圖7-2-6),只不過在計算時要將上述FFT演算法中的W因子用其共軛W*代替,並將最後迭代計算的結果統統乘以1/N.
圖7-2-6 N=8的逆離散富氏變換流程圖