fft.cpp

实现快速傅立叶变换

#include "iostream.h"


#include<fstream>//容易和#include "fstream.h"混淆
#include "stdlib.h"
#include <math.h>
#include<complex>
using namespace std;


void  FFT(complex<double> * input, complex<double> * output, int r)
{
	// 付立叶变换点数
	int	count;
	
	// 循环变量
	int		i,j,k;
	
	// 中间变量
	int		bfsize,p;
	
	// 角度
	double	angle;
	
	complex<double> *W,*X1,*X2,*X;
	
	// 计算付立叶变换点数
	count = 1 << r;
	
	// 分配运算所需存储器
	W  = new complex<double>[count / 2];
	X1 = new complex<double>[count];
	X2 = new complex<double>[count];
	
	// 计算加权系数
	for(i = 0; i < count / 2; i++)
	{
		angle = -i * 3.1415926 * 2 / count;
		W[i] = complex<double> (cos(angle), sin(angle));
	}
	
	// 将时域点写入X1
	memcpy(X1, input, sizeof(complex<double>) * count);
	
	// 采用蝶形算法进行快速付立叶变换
	for(k = 0; k < r; k++)
	{
		for(j = 0; j < 1 << k; j++)
		{
			bfsize = 1 << (r-k);
			for(i = 0; i < bfsize / 2; i++)
			{
				p = j * bfsize;
				X2[i + p] = X1[i + p] + X1[i + p + bfsize / 2];
				X2[i + p + bfsize / 2] = (X1[i + p] - X1[i + p + bfsize / 2]) * W[i * (1<<k)];
			}
		}
		X  = X1;
		X1 = X2;
		X2 = X;
	}
	
	// 重新排序
	for(j = 0; j < count; j++)
	{
		p = 0;
		for(i = 0; i < r; i++)
		{
			if (j&(1<<i))
			{
				p+=1<<(r-i-1);
			}
		}
		output[j]=X1[p];
	}
	
	// 释放内存
	delete W;
	delete X1;
	delete X2;
}

void main()
{
//	rand();
double temp;
char ch;
int size;
int	n = 1;	//可执行fft的有用输入个数
int	r = 0;//fft的层数
int j;	
cout<<"input the size of fft"<<endl;
cin>>size;
	while(n * 2 <= size)
	{
		n *= 2;
		r++;
	}
	size=n;
	cout<<"size="<<size<<endl;
	double  *mstart=new double[size];
	double  *Result=new double[size];
	complex<double> *input=new complex<double>[size];
	complex<double> *output=new complex<double>[size];	
/*	ofstream  fout1("input.text");
	  for( j=0;j<size;j++)
		{
		  cout<<"input the "<<j<<"num"<<endl;
		 cin>>temp;
		 fout1<<temp<<" ";
		}
					 
	 fout1.close();
*/	
	 ifstream fin("input.text");
	  	if(!fin)
	{
		cout<<" no  "<<endl;
	}
	  for( j =0;j<size;j++)
	 {
			fin>>mstart[j];
	    	input[j]= complex<double>(mstart[j],0);	
			cout<<mstart[j]<<endl;
	  }

	fin.close();
    FFT(input,output,r);
	cout<<"复数的幅度"<<endl;
	for(int v=0;v<size;v++)
	{
		Result[v]= sqrt(output[v].real()*output[v].real()+output[v].imag()*output[v].imag());
	
	}
    ofstream  file("output.text");
		for( j=0;j<size;j++)
	{
		cout<<Result[j]<<endl;		
	}		 
	for( j=0;j<size;j++)
	{
		file<<output[j].real()<<"+"<<output[j].imag()<<"i"<<"  ";
		cout<<output[j].real()<<"+"<<output[j].imag()<<"i"<<endl;
	}
	file.close();
			
}