蝶形运算
七百字作文-疯狂作文
2023年2月19日发(作者：瀑布图)【Java】快速傅⾥叶变换FFT的程序实现（时间抽取的基-2FFT、倒位计算、蝶
形运算）
Java——快速傅⾥叶变换（FFT）的程序实现
好久没来更新了，阿汪⼤三了。
这学期阿汪要学习两门课《数字信号处理》和《Java程序设计》，刚好前⼏天⽼师告诉我们不久后会有个实验，要求我们编写⼀个程序实现
快速傅⾥叶变换（FFT）。
所以，阿汪⽤Java写了⼀下。
//废话和原理就不多说了，直接上程序
//阿汪先⽣的博客.ws
r;
.*;
//对数函数
classLog{
staticpublicdoublelog(doublevalue,doublebase){
(value)/(base);
}
}
//乘⽅运算
classPow{
publicstaticintpowNM(intn,intm){
if(m==1)
returnn;
else
returnn*powNM(n,m-1);
//(m,n);
}
}
publicclassDFT_FFT{
publicstaticvoidmain(String[]args){
//初始化变量
finaldoublePI=3.141;
intN1=0;//序列的长度
intM=0;//序列的级数，数组的长度为2^M
intarraylength;//数组长度
inti=0;//输⼊数组的循环次数
intt=0;//输出数组的循环次数
///
//输⼊合法化
Scannersc=newScanner();
while(N1<2){
n("请确定你输⼊序列的长度！");
N1=t();
}
///
//判断长度并初始化数组
//确定N=2^M的级数M
M=((N1,2)==(int)(N1,2))?(int)(N1,2):((int)(N1,2))+1;
//确定数字长度
arraylength=(int)(2,M);
("序列长度N1为%d时，数组长度为%d，级数M为%d!n",N1,arraylength,M);
("序列长度N1为%d时，数组长度为%d，级数M为%d!n",N1,arraylength,M);
float[]Array=newfloat[arraylength];//实部
float[]IArray=newfloat[arraylength];//虚部
///
//输⼊序列
floatjc;//暂存转换后的数据
Stringxy;
("请依次输⼊序列的%d个值：n",N1);
n("");
while(i ("x[%d]实部：n",i); xy=(); if(((0)57)){//输⼊不是数字 n("阿汪先⽣(a_wang_xian_sheng)的博客！"); break;//跳出循环输⼊ } if(((0)>=48|(0)<=57)){//输⼊数字 jc=loat(xy);//String转化成float类型 Array[i]=jc;//保存String转存下来的数字 } ("x[%d]虚部：n",i);//输⼊虚部 IArray[i]=oat();//录⼊数字 i++; } n(""); n("输⼊结果如下："); while(t ("x[%d]=",t); (Array[t]); ("+j(%f)n",IArray[t]); t++; } // intw=0;//倒序处理时的外循环次数 intj=0;//倒序处理时的内循环次数 floatws=0;//倒序处理时的数组间值传递变量 intp=0;//倒序处理时下⼀倒位序的确定变量 intw1=0;//倒序处理时的外循环次数 intj1=0;//倒序处理时的内循环次数 floatws1=0;//倒序处理时的数组间值传递变量 intp1=0;//倒序处理时下⼀倒位序的确定变量 //倒序处理(实部) for(w=1,j=arraylength/2;w { if(w { ws=Array[j]; Array[j]=Array[w]; Array[w]=ws; } p=arraylength/2;//求j的下⼀个倒位序 while(p<=j)//如果k<=j,表⽰j的最⾼位为1 { j=j-p;//把最⾼位变成0 p=p/2;//k/2，⽐较次⾼位，依次类推，逐个⽐较，直到某个位为0 } j=j+p;//把0改为1 } //倒序处理(虚部) for(w1=1,j1=arraylength/2;w1 { if(w1 { { ws1=IArray[j1]; IArray[j1]=IArray[w1]; IArray[w1]=ws1; } p1=arraylength/2;//求j的下⼀个倒位序 while(p1<=j1)//如果k<=j,表⽰j的最⾼位为1 { j1=j1-p1;//把最⾼位变成0 p1=p1/2;//k/2，⽐较次⾼位，依次类推，逐个⽐较，直到某个位为0 } j1=j1+p1;//把0改为1 } // //蝶形运算变量定义并初始化 intL=1;//蝶形运算级数，⽤于循环 intN=2;//蝶形运算数据量，⽤于循环，WN(k)*X(k+N/2)中的N intdistance=1;//蝶形运算两节点间的距离，⽤于循环（distance=N/2） intgroup=arraylength/2;//蝶形运算的组数，长度为数组长度的⼀半 doubletempr1,tempr2,tempi1,tempi2;//临时变量 intk=0;//WN(k)*X(k+N/2)中的k //蝶形运算 for(;L<=M;L++){//第L级的蝶形运算 for(i=0;i for(k=0;k floattheta=(float)(-2*PI*k/N);//旋转因⼦,WN(k) tempr1=Array[N*i+k]; tempi1=IArray[N*i+k]; tempr2=(theta)*Array[N*i+k+distance]-(theta)*IArray[N*i+k+distance];//CB的实数 tempi2=(theta)*Array[N*i+k+distance]+(theta)*IArray[N*i+k+distance];//CB的复数 Array[N*i+k]=(float)(tempr1+tempr2);//(A+CB)的实数 IArray[N*i+k]=(float)(tempi1+tempi2);//(A+CB)的复数 Array[N*i+k+distance]=(float)(tempr1-tempr2);//(A-CB)的实数 IArray[N*i+k+distance]=(float)(tempi1-tempi2);//(A-CB)的复数 } } N=N*2;//下⼀级蝶形运算中循环分母N*2 distance*=2;//下⼀级蝶形运算两节点间的距离增加 group/=2;//下⼀级蝶形运算的组数减半 } // //输出运算结果 intt1=0;//数组输出的循环次数 n(""); n("FFT运算结果如下："); while(t1 ("X[%d]=%ft+j(%f)n",t1,Array[t1],IArray[t1]); t1++; } } } //阿汪先⽣的博客.ws 程序效果 1、复数输⼊ 2、输⼊实部时，输⼊任⼀⾮数字字符（串），⾃动结束录⼊过程，剩余序列位⾃动补零； 在Matlab中验算 x=[1230]; X=fft(x,4); X ⼀些说明： 1、可复数输⼊，输⼊位数不⾜2^N次，⾃动补位； 2、输⼊实部时，输⼊任⼀⾮数字字符，⾃动结束录⼊过程，剩余序列位⾃动补零； 3、倒位运算采⽤了雷德算法； 4、计算机运算时位数的限制导致上述的两段程序结果的误差，问题不⼤； 5、修改因⼦W_N^K中的N可算任意N点的FFT运算，本程序中的N为(2^M)，具有⼀定局限性。 //只要在输⼊时单独赋值，即可实现序列任意N点（可奇数点）的FFT运算。