Physics 版 (精华区)

发信人: zjliu (秋天的萝卜), 信区: Physics
标  题: FFT Made by Visaul C++
发信站: 哈工大紫丁香 (2003年06月07日18:30:09 星期六), 站内信件


找到这个东西和大家SHARE：
FFT Made by Visaul c++

//small fft program,using complex.


#include "complex"
#include <iostream>
#include <conio.h>
using namespace std; //very important
#define PI 3.1415926
typedef complex<double> complexd; //use template
void myfft(double y[][2], int N) //N is the power of 2 (as 8,64, 1024)
{
  complexd t1,t2,t;
  int u,mu,i,j,k,n,l,s;
  complexd *w =  new complexd[N/2];
  complexd *x =  new complexd[N];
  mu=0; u=N-1;
  while (u>0) {mu++; u/=2;};光), 信区: MathTools             讨论区 
[MathTools]
  for (i=0; i<N; i++) x[i =complexd(y[i][0],y[i][1]);
  for (i=0; i<N; i++) { //sort input number array
    k=i; l=0;
    for (j=2; j<=N; j=j*2) {
      l+=(N/j)*(k%2); k=k/2;
    }
    if (l>i) {
      t=x[i];
      x[i =x[l];
      x[l =t;
    }
  }
  for (i=0; i<N/2; i++) w[i =complexd(cos(2*PI/N*i),sin(2*PI/N*i));
  s=1;
  for (i=0; i<mu; i++) { //fft by binary dividing
    for (j=0; j<N; j+=s*2) {
      for (k=0; k<s; k++) {
  t=x[j+k+s]*(w[k*N/s/2]);
  t1=x[j+k]+t;
  t2=x[j+k]-t;
  x[j+k =t1;
  while (u>0) {mu++; u/=2;};光), 信区: MathTools             讨论区 
[MathTools]
  for (i=0; i<N; i++) x[i =complexd(y[i][0],y[i][1]);
  for (i=0; i<N; i++) { //sort input number array
    k=i; l=0;
    for (j=2; j<=N; j=j*2) {
      l+=(N/j)*(k%2); k=k/2;
    }
    if (l>i) {
      t=x[i];
      x[i =x[l];
      x[l =t;
    }
  }
  for (i=0; i<N/2; i++) w[i =complexd(cos(2*PI/N*i),sin(2*PI/N*i));
  s=1;
  for (i=0; i<mu; i++) { //fft by binary dividing
    for (j=0; j<N; j+=s*2) {
      for (k=0; k<s; k++) {
  t=x[j+k+s]*(w[k*N/s/2]);
  t1=x[j+k]+t;
  t2=x[j+k]-t;
  x[j+k =t1;
  x[j+k =t1;
  x[j+k+s =t2;
      }
    }
    s*=2;
  }
  for (i=0; i<N; i++) {y[i][0 =x[i].real(); y[i][1 =x[i].imag();}
  delete w; //very important. If not, the system will lost memory
  delete x;
}
void main()
{
 int i;
 double u[8][2];
 for (i=0; i<8; i++) u[i][0 =i+1; //the input
 myfft(u,8);
 for (i=0; i<8; i++) cout<<'('<<u[i][0]<<','<<u[i][1]<<')'<<endl;
 _getch();
}
--
╔═══════════════════╗
║★★★★★友谊第一  比赛第二★★★★★║
╚═══════════════════╝

※ 来源:·哈工大紫丁香 bbs.hit.edu.cn·[FROM: 202.118.229.86]