快速傅立叶c语言函数,快速傅立叶变换

求基2、基4、基8FFT(快速傅里叶变换)的c语言程序,要能运行得出来的

1.FFT:

创新互联-专业网站定制、快速模板网站建设、高性价比银海网站开发、企业建站全套包干低至880元,成熟完善的模板库,直接使用。一站式银海网站制作公司更省心,省钱,快速模板网站建设找我们,业务覆盖银海地区。费用合理售后完善,10多年实体公司更值得信赖。

// data为输入和输出的数据,N为长度

bool CFFT::Forward(complex *const Data, const unsigned int N)

{

if (!Data || N  1 || N  (N - 1))

return false;

//   排序

Rearrange(Data, N);

//   FFT计算:const bool Inverse = false

Perform(Data, N);

return true;

}

2.IFFT:

// Scale 为是否缩放

bool CFFT::Inverse(complex *const Data, const unsigned int N,

const bool Scale /* = true */)

{

if (!Data || N  1 || N  (N - 1))

return false;

//   排序

Rearrange(Data, N);

//   FFT计算,ture表示是逆运算

Perform(Data, N, true);

//   对结果进行缩放

if (Scale)

CFFT::Scale(Data, N);

return true;

}

3.排序:

void CFFT::Rearrange(complex *const Data, const unsigned int N)

{

//   Swap position

unsigned int Target = 0;

//   Process all positions of input signal

for (unsigned int Position = 0; Position  N; ++Position)

{

//   Only for not yet swapped entries

if (Target  Position)

{

//   Swap entries

const complex Temp(Data[Target]);

Data[Target] = Data[Position];

Data[Position] = Temp;

}

//   Bit mask

unsigned int Mask = N;

//   While bit is set

while (Target  (Mask = 1))

//   Drop bit

Target = ~Mask;

//   The current bit is 0 - set it

Target |= Mask;

}

}

4.FFT计算:

void CFFT::Perform(complex *const Data, const unsigned int N,

const bool Inverse /* = false */)

{

const double pi = Inverse ? 3.14159265358979323846 : -3.14159265358979323846;

//   Iteration through dyads, quadruples, octads and so on...

for (unsigned int Step = 1; Step  N; Step = 1)

{

//   Jump to the next entry of the same transform factor

const unsigned int Jump = Step  1;

//   Angle increment

const double delta = pi / double(Step);

//   Auxiliary sin(delta / 2)

const double Sine = sin(delta * .5);

//   Multiplier for trigonometric recurrence

const complex Multiplier(-2. * Sine * Sine, sin(delta));

//   Start value for transform factor, fi = 0

complex Factor(1.);

//   Iteration through groups of different transform factor

for (unsigned int Group = 0; Group  Step; ++Group)

{

//   Iteration within group 

for (unsigned int Pair = Group; Pair  N; Pair += Jump)

{

//   Match position

const unsigned int Match = Pair + Step;

//   Second term of two-point transform

const complex Product(Factor * Data[Match]);

//   Transform for fi + pi

Data[Match] = Data[Pair] - Product;

//   Transform for fi

Data[Pair] += Product;

}

//   Successive transform factor via trigonometric recurrence

Factor = Multiplier * Factor + Factor;

}

}

}

5.缩放:

void CFFT::Scale(complex *const Data, const unsigned int N)

{

const double Factor = 1. / double(N);

//   Scale all data entries

for (unsigned int Position = 0; Position  N; ++Position)

Data[Position] *= Factor;

}

一个关于128点的快速傅立叶的C语言程序

这是我写的1024点的快速傅里叶变换程序,下面有验证,你把数组

double

A[2049]={0};

double

B[1100]={0};

double

powerA[1025]={0};

改成

A[256]={0};

B[130]={0};

power[129]={0};就行了,

void

FFT(double

data[],

int

nn,

int

isign)

的程序可以针对任何点数,只要是2的n次方

具体程序如下:

#include

iostream.h

#include

"math.h"

#includestdio.h

#includestring.h

#include

stdlib.h

#include

fstream.h

#include

afx.h

void

FFT(double

data[],

int

nn,

int

isign)

{

//复数的快速傅里叶变换

int

n,j,i,m,mmax,istep;

double

tempr,tempi,theta,wpr,wpi,wr,wi,wtemp;

n

=

2

*

nn;

j

=

1;

for

(i

=

1;

i=n

;

i=i+2)

//这个循环进行的是码位倒置。

{

if(

j

i)

{

tempr

=

data[j];

tempi

=

data[j

+

1];

data[j]

=

data[i];

data[j

+

1]

=

data[i

+

1];

data[i]

=

tempr;

data[i

+

1]

=

tempi;

}

m

=

n

/

2;

while

(m

=

2

j

m)

{

j

=

j

-

m;

m

=

m

/

2;

}

j

=

j

+

m;

}

mmax

=

2;

while(

n

mmax

)

{

istep

=

2

*

mmax;

//这里表示一次的数字的变化。也体现了级数,若第一级时,也就是书是的第0级,其为两个虚数,所以对应数组应该增加4,这样就可以进入下一组运算

theta

=

-6.28318530717959

/

(isign

*

mmax);

wpr

=

-2.0

*

sin(0.5

*

theta)*sin(0.5

*

theta);

wpi

=

sin(theta);

wr

=

1.0;

wi

=

0.0;

for(

m

=

1;

m=mmax;

m=m+2)

{

for

(i

=

m;

i=n;

i=i+istep)

{

j

=

i

+

mmax;

tempr=double(wr)*data[j]-double(wi)*data[j+1];//这两句表示蝶形因子的下一个数乘以W因子所得的实部和虚部。

tempi=double(wr)*data[j+1]+double(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;

}

wtemp

=

wr;

wr

=

wr

*

wpr

-

wi

*

wpi

+

wr;

wi

=

wi

*

wpr

+

wtemp

*

wpi

+

wi;

}

mmax

=

istep;

}

}

void

main()

{

//本程序已经和MATLAB运算结果对比,准确无误,需要注意的的是,计算中数组都是从1开始取得,丢弃了A[0]等数据

double

A[2049]={0};

double

B[1100]={0};

double

powerA[1025]={0};

char

line[50];

char

dataA[20],

dataB[20];

int

ij;

char

ch1[3]="\t";

char

ch2[3]="\n";

int

strl1,strl2;

CString

str1,str2;

ij=1;

//********************************读入文件data1024.txt中的数据,

其中的数据格式见该文件

FILE

*fp

=

fopen("data1024.txt","r");

if(!fp)

{

cout"Open

file

is

failing!"endl;

return;

}

while(!feof(fp))

//feof(fp)有两个返回值:如果遇到文件结束,函数feof(fp)的值为1,否则为0。

{

memset(line,0,50);

//清空为0

memset(dataA,0,20);

memset(dataB,0,20);

fgets(line,50,fp);

//函数的功能是从fp所指文件中读入n-1个字符放入line为起始地址的空间内

sscanf(line,

"%s%s",

dataA,

dataB);

//我同时读入了两列值,但你要求1024个,那么我就只用了第一列的1024个值

//dataA读入第一列,dataB读入第二列

B[ij]=atof(dataA);

//将字符型的dataA值转化为float型

ij++;

}

for

(int

mm=1;mm1025;mm++)//A[2*mm-1]是实部,A[2*mm]是虚部,当只要输入实数时,那么保证虚部A[mm*2]为零即可

{

A[2*mm-1]=B[mm];

A[2*mm]=0;

}

//*******************************************正式计算FFT

FFT(A,1024,1);

//********************************************写入数据到workout.txt文件中

for

(int

k=1;k2049;k=k+2)

{

powerA[(k+1)/2]=sqrt(pow(A[k],2.0)+pow(A[k+1],2.0));//求功率谱

FILE

*pFile=fopen("workout.txt","a+");

//?a+只能在文件最后补充,光标在结尾。没有则创建

memset(ch1,0,15);

str1.Format("%.4f",powerA[(k+1)/2]);

if

(A[k+1]=0)

str2.Format("%d\t%6.4f%s%6.4f

%s",(k+1)/2,A[k],"+",A[k+1],"i");//保存fft计算的频谱,是复数频谱

else

str2.Format("%d\t%6.4f%6.4f

%s",(k+1)/2,A[k],A[k+1],"i");

strl1=strlen(str1);

strl2=strlen(str2);

//

法:fwrite(buffer,size,count,fp);

//

buffer:是一个指针,对fwrite来说,是要输出数据的地址。

//

size:要写入的字节数;

//

count:要进行写入size字节的数据项的个数;

//

fp:目标文件指针。

fwrite(str2,1,strl2,pFile);

fwrite(ch1,1,3,pFile);

fwrite(ch1,1,3,pFile);

fwrite(str1,1,strl1,pFile);

fwrite(ch2,1,3,pFile);

fclose(pFile);

}

cout"计算完毕,到fft_test\workout.txt查看结果"endl;

}

快速傅立叶变换在c中怎么实现?

快速傅氏变换,是离散傅氏变换的快速算法,它是根据离散傅氏变换的奇、偶、虚、实等特性,对离散傅立叶变换的算法进行改进获得的。它对傅氏变换的理论并没有新的发现,但是对于在计算机系统或者说数字系统中应用离散傅立叶变换,可以说是进了一大步。

设x(n)为N项的复数序列,由DFT变换,任一X(m)的计算都需要N次复数乘法和N-1次复数加法,而一次复数乘法等于四次实数乘法和两次实数加法,一次复数加法等于两次实数加法,即使把一次复数乘法和一次复数加法定义成一次“运算”(四次实数乘法和四次实数加法),那么求出N项复数序列的X(m), 即N点DFT变换大约就需要N2次运算。当N=1024点甚至更多的时候,需要N2=1048576次运算,在FFT中,利用WN的周期性和对称性,把一个N项序列(设N=2k,k为正整数),分为两个N/2项的子序列,每个N/2点DFT变换需要(N/2)2次运算,再用N次运算把两个N/2点的DFT 变换组合成一个N点的DFT变换。这样变换以后,总的运算次数就变成N+2(N/2)2=N+N2/2。继续上面的例子,N=1024时,总的运算次数就变成了525312次,节省了大约50%的运算量。而如果我们将这种“一分为二”的思想不断进行下去,直到分成两两一组的DFT运算单元,那么N点的 DFT变换就只需要Nlog2N次的运算,N在1024点时,运算量仅有10240次,是先前的直接算法的1%,点数越多,运算量的节约就越大,这就是 FFT的优越性.

傅里叶变换(Transformée de Fourier)是一种积分变换。因其基本思想首先由法国学者傅里叶系统地提出,所以以其名字来命名以示纪念。

应用

傅里叶变换在物理学、数论、组合数学、信号处理、概率论、统计学、密码学、声学、光学、海洋学、结构动力学等领域都有着广泛的应用(例如在信号处理中,傅里叶变换的典型用途是将信号分解成幅值分量和频率分量)。

概要介绍

傅里叶变换能将满足一定条件的某个函数表示成三角函数(正弦和/或余弦函数)或者它们的积分的线性组合。在不同的研究领域,傅里叶变换具有多种不同的变体形式,如连续傅里叶变换和离散傅里叶变换。最初傅里叶分析是作为热过程的解析分析的工具被提出的(参见:林家翘、西格尔著《自然科学中确定性问题的应用数学》,科学出版社,北京。原版书名为 C. C. Lin L. A. Segel, Mathematics Applied to Deterministic Problems in the Natural Sciences, Macmillan Inc., New York, 1974)。

傅里叶变换属于谐波分析。

傅里叶变换的逆变换容易求出,而且形式与正变换非常类似;

正弦基函数是微分运算的本征函数,从而使得线性微分方程的求解可以转化为常系数的代数方程的求解.在线性时不变的物理系统内,频率是个不变的性质,从而系统对于复杂激励的响应可以通过组合其对不同频率正弦信号的响应来获取;

卷积定理指出:傅里叶变换可以化复杂的卷积运算为简单的乘积运算,从而提供了计算卷积的一种简单手段;

离散形式的傅里叶变换可以利用数字计算机快速的算出(其算法称为快速傅里叶变换算法(FFT)).


文章名称:快速傅立叶c语言函数,快速傅立叶变换
链接分享:http://scyanting.com/article/hcdisp.html