#include<stdio.h>
#include<math.h>
#include<stdlib.h>
#define N 1024
typedef struct{
double real;
double imag;
}complex;
complex x[N], *W;
int size=1024;
double l[5000],h[5000];
double PI=4.0*atan(1);
void output()
{
int i;
for(i=0;i<size;i++)
{
printf("%lf %lf %.4f",h[i],l[i],x[i].real);
if(x[i].imag>=0.0001)
{
printf("+%.4fj\n",x[i].imag);
}
else if(fabs(x[i].imag)<0.0001)
{
printf("\n");
}
else
{
printf("%.4fj\n",x[i].imag);
}
}
}
void change()
{
complex temp;
unsigned short i=0,j=0,k=0;
double t;
for(i=0;i<size;i++)
{
k=i;
j=0;
t=(log(size)/log(2));
while( (t--)>0 )
{
j=j<<1;
j|=(k & 1);
k=k>>1;
}
if(j>i)
{
temp=x[i];
x[i]=x[j];
x[j]=temp;
}
}
}
void transform()
{
int i;
W=(complex *)malloc(sizeof(complex) * size);
for(i=0;i<size;i++)
{
W[i].real=cos(2*PI/size*i);
W[i].imag=-1*sin(2*PI/size*i);
}
}
void add(complex a,complex b,complex *c)
{
c->real=a.real+b.real;
c->imag=a.imag+b.imag;
}
void sub(complex a,complex b,complex *c)
{
c->real=a.real-b.real;
c->imag=a.imag-b.imag;
}
void mul(complex a,complex b,complex *c)
{
c->real=a.real*b.real - a.imag*b.imag;
c->imag=a.real*b.imag + a.imag*b.real;
}
void fft()
{
int i=0,j=0,k=0,m=0;
complex q,y,z;
change();
for(i=0;i<log(size)/log(2) ;i++)
{
m=1<<i;
for(j=0;j<size;j+=2*m)
{
for(k=0;k<m;k++)
{
mul(x[k+j+m],W[size*k/2/m],&q);
add(x[j+k],q,&y);
sub(x[j+k],q,&z);
x[j+k]=y;
x[j+k+m]=z;
}
}
}
}
int main()
{
int i;
double fs=1000;
double t[5000];
for(i=0;i<size;i++)
{
t[i]=i/fs;
x[i].real=3+2*cos(2*PI*850*t[i])+3.6*cos(2*PI*300*t[i])+2.6*sin(2*PI*400*t[i])+5.5*cos(2*PI*200*t[i]+PI/2)+3.8*sin(2*PI*80*t[i]);//+4*cos(2*PI*150*t[i]+PI/2)+16*sin(2*PI*50*t[i]);//输入连续时间信号,并对其进行采样
x[i].imag=0;
h[i]=i*(fs/(double)size);
}
transform();
fft();
for(i=0;i<size;i++)
{
l[i]=sqrt(x[i].imag*x[i].imag+x[i].real*x[i].real) ;
}
printf("输出FFT后的结果\n");
output();
return 0;
}
改一下变量名 让其能在STM32上运行
刚开始学 怎么都弄不出来,麻烦了 salute