codepad
[
create a new paste
]
login
|
about
Language:
C
C++
D
Haskell
Lua
OCaml
PHP
Perl
Plain Text
Python
Ruby
Scheme
Tcl
#include<stdio.h> #include<math.h> #define N 256 /* FFT関数 */ void fft(double xRe[],double xIm[],int p) { int i,ip,j,k,m,me,me1,n,nv2; double uRe,uIm,vRe,vIm,wRe,wIm,tRe,tIm; double sin(),cos(); n=1<<p; nv2=n/2; /* ビットリバース */ j=0; for(i=0;i<n-1;i++){ if(j>i){ tRe=xRe[j]; tIm=xIm[j]; xRe[j]=xRe[i]; xIm[j]=xIm[i]; xRe[i]=tRe; xIm[i]=tIm; } k=nv2; while(j>=k){ j-=k; k/=2; } j+=k; } /* 主計算 */ for(m=1;m<=p;m++){ me=1<<m; me1=me/2; uRe=1.0; uIm=0.0; wRe=cos(M_PI/me1); wIm=-sin(M_PI/me1); for(j=0;j<me1;j++){ for(i=j;i<n;i+=me){ ip=i+me1; /* バタフライ演算 */ tRe=xRe[ip]*uRe-xIm[ip]*uIm; tIm=xRe[ip]*uIm+xIm[ip]*uRe; xRe[ip]=xRe[i]-tRe; xIm[ip]=xIm[i]-tIm; xRe[i]+=tRe; xIm[i]+=tIm; } vRe=uRe*wRe-uIm*wIm; vIm=uRe*wIm+uIm*wRe; uRe=vRe; uIm=vIm; } } } /* main関数 */ int main(){ int d,i,p; double xRe[N],xIm[N]; FILE *fpr,*fpw; p=log10(N) / log10(2); /* ファイル入力 */ fpr=fopen("./CYMBAL.txt","r"); for(i=0;i<N;i++){ fscanf(fpr,"%d",d); xRe[i]=(double)d; xIm[i]=0.0; //printf("%lf\n",xRe[i]); } fclose(fpr); /* FFT関数呼び出し */ fft(xRe,xIm,p); for(i=0;i<N;i++){ printf("%d,%lf,%lf\n",i,xRe[i],xIm[i]); } /* ファイル出力 */ fpw=fopen("./test.txt","w"); for(i=0;i<N;i++){ fprintf(fpw,"%lf,%lf\n",xRe[i],xIm[i]); } fclose(fpw); }
Private
[
?
]
Run code
Submit