如何将AD采集到的数据,用FFT进行变换
发布网友
发布时间:2022-04-26 10:13
我来回答
共3个回答
热心网友
时间:2022-06-04 13:18
展开1全部#include<math.h>
kkfft(double pr[], double pi[], int n, int k, double fr[], double fi[], int l, int il)
{
int it,m,is,i,j,nv,l0;
double p,q,s,vr,vi,poddr,poddi;
for (it=0; it<=n-1; it++)
{
m = it;
is = 0;
for (i=0; i<=k-1; i++)
{
j = m/2;
is = 2*is+(m-2*j);
m = j;
}
fr[it] = pr[is];
fi[it] = pi[is];
}
pr[0] = 1.0;
pi[0] = 0.0;
p = 6.283185306/(1.0*n);
pr[1] = cos(p);
pi[1] = -sin(p);
if (l!=0)
pi[1]=-pi[1];
for (i=2; i<=n-1; i++)
{
p = pr[i-1]*pr[1];
q = pi[i-1]*pi[1];
s = (pr[i-1]+pi[i-1])*(pr[1]+pi[1]);
pr[i] = p-q;
pi[i] = s-p-q;
}
for (it=0; it<=n-2; it=it+2)
{
vr = fr[it];
vi = fi[it];
fr[it] = vr+fr[it+1];
fi[it] = vi+fi[it+1];
fr[it+1] = vr-fr[it+1];
fi[it+1] = vi-fi[it+1];
}
m = n/2;
nv = 2;
for (l0=k-2; l0>=0; l0--)
{
m = m/2;
nv = 2*nv;
for (it=0; it<=(m-1)*nv; it=it+nv)
for (j=0; j<=(nv/2)-1; j++)
{
p = pr[m*j]*fr[it+j+nv/2];
q = pi[m*j]*fi[it+j+nv/2];
s = pr[m*j]+pi[m*j];
s = s*(fr[it+j+nv/2]+fi[it+j+nv/2]);
poddr = p-q;
poddi = s-p-q;
fr[it+j+nv/2] = fr[it+j]-poddr;
fi[it+j+nv/2] = fi[it+j]-poddi;
fr[it+j] = fr[it+j]+poddr;
fi[it+j] = fi[it+j]+poddi;
}
}
if (l!=0)
{
for (i=0; i<=n-1; i++)
{
fr[i] = fr[i]/(1.0*n);
fi[i] = fi[i]/(1.0*n);
}
}
if (il!=0)
{
for (i=0; i<=n-1; i++)
{
pr[i] = sqrt(fr[i]*fr[i]+fi[i]*fi[i]);
if (fabs(fr[i])<0.000001*fabs(fi[i]))
{
if ((fi[i]*fr[i])>0)
pi[i] = 90.0;
else
pi[i] = -90.0;
}
else
pi[i] = atan(fi[i]/fr[i])*360.0/6.283185306;
}
}
return;
}
void main()
{
double s1[128],s2[128],s3[128],s4[128],TT,TM;
int i,f,ff;
f=128;
ff=4;
for(i=0;i<128;i++)
{
TT=i;
TM=TT/f;
s1[i]=2*sin(2*3.14*ff*TM)+10;
}
kkfft(s1, s2, 128, 7, s3, s4, 0, 1);
printf("s1=%f\n",s1[0]/128) ;
getch();
for(i=1;i<12;i++)
{
printf("s1=%f,i=%d\n",s1[i]/64,i) ;
}
getch();
}
这是一段用wintc下编译的fft,double pr[], double pi[],double fr[], double fi[]是入口数组,在主函数中定义好,点数自定,你把ad采集的数存入
pr[]对应的数组中,记住如程序中的
for(i=0;i<128;i++)
{
TT=i;
TM=TT/f;
s1[i]=2*sin(2*3.14*ff*TM)+10;
}
s1对应pr
运行后可以看到直流分量(i=0点),交流分量(i=4点)
说明: kkfft(s1, s2, 128, 7, s3, s4, 0, 1);
s1对应ad才来的值,s2是ad采样的虚部一般为0,128采样点数一定要是2的
n次方,7就是那个n了,s3,s4背管了是副角返回值,0,1是程序内定的标志
fft结果的意义:她是一个数组如S1[n],那么s1[0]/采样点数是直流分量值,
而实际信号的频率在数组中的反应是f=n*采样频率/采样点数;
kkfft(double pr[], double pi[], int n, int k, double fr[], double fi[], int l, int il)
在51(020,120,060)单片机上运行通过;
热心网友
时间:2022-06-04 13:18
FFT 对离散数据作时域到频域的变换。
FFT 利用函数对称性使计算加速,所以输入数据必须是 2 的整数次方,例如1024,2048,4096 。。。
所以你要对 采到的数据要做长度处理。
一种是 保留时间长度,改变 dt, 重新采样(内插),变成 N=2^k 个数。
一种是 保留dt,加长时间,用尾部加0,补充到 N=2^k
一种是 保留dt,假定信号重复,把头部信号接到尾部,补充到 N=2^k
对时域信号要做光滑(修匀)处理,去掉毛刺
FFT 后的频域 结果 要对 高频和低频做 截断处理。
目的是消除 离散 和 时间截断 引起的 泄漏和旁瓣效应。
函数参数随编程人员而定。通常有 点数,采样时间步长,时系值数组,低频截断值,高频截断值。也许有输出数组,也许用时系值数组存放变换后的值。
热心网友
时间:2022-06-04 13:19
fft变换是方法不是目的, 问题就你到底想要取得什么样的东西。
对无限时间的序列来说,随便取一段都可以做fft, 而这个fft只能表示取的这段时间里的频率分布, 如果频率分布是稳定的, 你可以随意取足够长的一段来fft, 如果频率分布是不稳定的, 这就成了时频分析, 方法太多了, 各有特点, 其中比较容易想到的就是对每连续的n个点(比如1024)都做一次fft。