fft.txt

一个用VC++实现傅立叶变换的程序,可以

void FFT(double fr[], double fi[], int n, int flag)
{
	int mp,arg,cntr,p1,p2;
	int i,j,a,b,k;
	double sign,pr,pi,harm,t;
	double* ca=new double[n],*sa=new double[n];
	j=0;
	if(flag!=0)
	{
		sign=1.0;
	}
	else
	{
		sign=-1.0;
	}
	
	for(i=0;i<= n-2 ;++i)
	{
		if(i<j)
		{
			t=fr[i];
			fr[i]=fr[j];
			fr[j]=t;
			////
			t=fi[i];
			fi[i]=fi[j];
			fi[j]=t;
		}
		k=n/2;
		while(k<= j)
		{
			j-=k;
			k/=2;
		}
		j+=k;
	}
	mp=0;  i=n;
	while(i!=1)
	{
		mp+=1;
		i/=2;
	}
	harm=2*PI/n;
	for(i=0;i<=n-1;++i)
	{
		sa[i]=sign*sin(harm*i);
		ca[i]=cos(harm*i);
	}
	a=2;  b=1;
	for(cntr=1;cntr<= mp; ++cntr)
	{
		p1=n/a;     p2=0;
		for(k=0  ;k<= b-1 ;++k)
		{
			i=k;
			while(i< n)
			{
				arg=i+b;
				if(k==0)
				{
					pr=fr[arg];   pi= fi[arg];
				}
				else
				{
					pr=fr[arg]*ca[p2]-fi[arg]*sa[p2];
					pi=fr[arg]*sa[p2]+fi[arg]*ca[p2];
				}
				fr[arg]=fr[i]-pr;   fi[arg]=fi[i]-pi;  fr[i]+=pr;   fi[i]+=pi;   i+=a;  
			}
			p2+=p1;
		}
		a*= 2;  ;b*= 2;
	}
	if(flag==0)
	{
		for(i=0; i< n;++i){fr[i]=fr[i]/n;  fi[i]=fr[i]/n; }
	}
	delete []ca;
	delete []sa;
}