快速傅立叶算法的C语言实现0001

1、快速傅立叶算法的语言实现#include math.hvoid kfft(pr,pi,n,k,fr,fi,l,il)int n,k,l,il;double pr,pi,fr,fi; 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; frit=pris; fiit=piis;pr0=1.0; pi0=0.0; p=6.283185306/(1.0*n); pr1=cos(p

2、); pi1=-sin(p);if (l!=0) pi1=-pi1;for (i=2; i=n-1; i+) p=pri-1*pr1; q=pii-1*pi1; s=(pri-1+pii-1)*(pr1+pi1); pri=p-q; pii=s-p-q;for (it=0; it=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=prm*j*frit+j+nv/2; q=pim*j*fiit+j+nv/2; s=prm*j+pim*j; s=s*(frit+j+nv/2+fi

3、it+j+nv/2); poddr=p-q; poddi=s-p-q; frit+j+nv/2=frit+j-poddr; fiit+j+nv/2=fiit+j-poddi; frit+j=frit+j+poddr; fiit+j=fiit+j+poddi;if (l!=0)for (i=0; i=n-1; i+) fri=fri/(1.0*n);fii=fii/(1.0*n);if (il!=0)for (i=0; i=n-1; i+) pri=sqrt(fri*fri+fii*fii);if (fabs(fri)0) pii=90.0;else pii=-90.0;elsepii=atan

4、(fii/fri)*360.0/6.283185306;return;这是在 PIC18F458 里用的 FFT 程序/* * 函 数 名:void ADCHD()* 功能描述 :单通道高频采样 ,数据处理 ,显示 .*/void ADCHD()long D,SHU; / 数据 .int n_x,k_x,i; / 循环参数 .float Ur,Ui,Urn,Uin; / 数据处理中间变量 . ADWORDS=ADCAN+7; / 确定 MAX197 控制字 . for(ADD=0;ADD64;ADD+) / 采样 .OUTADC(ADWORDS); / 送控制字 .DELAY(0); / 延时

5、 .DA TA0ADD=INADC(); / 读取数据 . for(ADD=0;ADD64;ADD+) / 显示波形 .D=(int)(40.0*DATA0ADD/4096.0); / 取数据 .LCDPIEX(ADD+64,40-D); / 显示 .for(n_x=0;n_x5;n_x+) / 计算 .Urn=0.0; / 实部 .Uin=0.0; / 虚部 .for(k_x=0;k_x32;k_x+) /n_x 次谐波 .D=DATA0k_x; / 取数据计算 .Urn=Urn+D/409.6*cos(2*n_x+1)*(k_x+1)*0.196);Uin=Uin+D/409.6*sin(

6、2*n_x+1)*(k_x+1)*0.196); Ur=Urn/16.0; /Ui=Uin/16.0; /SHU=(long)(100*sqrt(Ur*Ur+Ui*Ui);UI0n_x=SHU; / 第 n_x 次谐波幅值 .UI05=SHU*SHU+UI05; /UI05=(long)sqrt(UI05); / 总幅值 . for(i=0;i5;i+) /SHU=UI0i*1000; /SHU=SHU/(UI05); /UP0i=SHU; / 第 i 次谐波占有率 . SHU=1000*(UI05-UI00); /UP05=SHU/UI05;/畸变率 .LCDCLEAN(12,2,126,7

7、); / 清除数据显示区 . for(i=0;i6;i+)OUTNUM(UI0i,1,28,2+i); / 显示幅值 .OUTNUM(UP0i,3,56,2+i); / 显示占有率 . ADWORDS=ADCAN; / 还原控制字 .* 函 数 名:void XBFX()* 功能描述 : 全周波傅立叶积分计算各次谐波的幅值,并返回结果 void XBFX()long D,SHU; / 数据 . int n_x,k_x,i,n; / 循环参数 .float Ur,Ui,Urn,Uin,UIL2; / 数据处理中间变量 for(n=0;n2;n+) / 两路信号 .for(n_x=0;n_x5;n

8、_x+) / 计算 .Urn=0.0; / 实部 .Uin=0.0; / 虚部 . for(k_x=0;k_x32;k_x+) /n_x 次谐波 .D=DATACHAN*2+nk_x; / 取数据计算 . Urn=Urn+D/409.6*cos(2*n_x+1)*(k_x+1)*0.196);Uin=Uin+D/409.6*sin(2*n_x+1)*(k_x+1)*0.196); Ur=Urn/16.0; /Ui=Uin/16.0; SHU=(long)(760*sqrt(Ur*Ur+Ui*Ui); / UInn_x=SHU; /UIn5=SHU*SHU+UIn5; / 第 n_x 次谐波幅值 if(n_x=0)UILn=atan(Ui/Ur);/相位 . UIn5=(long)sqrt(UIn5); / 总幅值 . for(i=0;i5;i+) /SHU=UIni*1000; /SHU=SHU/(UIn5); /UPni=SHU; / 第 i 次谐波占有率 . SHU=1000*(UIn5-UIn0); / UPn5=SHU/UIn5;/畸变率 .功率因数角 LCHAN=abs(long)