任意2的整数次幂点数的基-2DIT-FFT和DIF-FFT通用cc++程序

编制任意2的整数次幂点数的基-2DIT-FFT和DIF-FFT通用c/c++程序，验证其正确性，并与直接计算DFT比较点数为2^N（N=10,…,16）时运行时间的差异。要求：编制任意2的整数次幂点数的基-2DIT-FFT和DIF-FFT通用C/C++程序，验证其正确性，并与直接计算DFT比较点数2^N(其中N=10,....16)时运行时间的差异。#includeiostream#includefstream#includecomplex#includetime.h#includemath.h#definepi3.14159typedefunsignedintUINT;usingnamespacestd;UINTRevBit(UINTn,intr)//比特倒序函数{UINTtemp=n;temp=(temp&0x55555555)1|(temp&0xAAAAAAAA)1;temp=(temp&0x33333333)2|(temp&0xCCCCCCCC)2;temp=(temp&0x0F0F0F0F)4|(temp&0xF0F0F0F0)4;temp=(temp&0x00FF00FF)8|(temp&0xFF00FF00)8;temp=(temp&0x0000FFFF)16|(temp&0xFFFF0000)16;temp=temp(32-r);returntemp;}boolditfft(complexdouble*TS,complexdouble*FS,intr){//时间抽取算法,TS为时间序列，FS为频率序列，序列长度为N=2^rinti,j,k;intN=1r;//N=2^rdoublealpha;//旋转因子的相位角complexdouble*W;//旋转因子W=newcomplexdouble[N/2];for(i=0;iN/2;i++)//计算旋转因子{alpha=-i*2*pi/N;W[i]=complexdouble(cos(alpha),sin(alpha));}for(i=0;iN;i++)//时间序列倒位序FS[i]=TS[RevBit(i,r)];for(i=0;ir;i++){intix=1(i+1);//第i级每组的元素个数intiy=1(r-i-1);//第i级的总组数for(j=0;jiy;j++)//蝶形计算过程{for(k=ix/2;kix;k++)//每组的后半部分乘以旋转因子FS[j*ix+k]=FS[j*ix+k]*W[(k-ix/2)*iy];for(k=0;kix/2;k++){complexdoubleTEMP1;complexdoubleTEMP2;TEMP1=FS[j*ix+k]+FS[j*ix+k+ix/2];TEMP2=FS[j*ix+k]-FS[j*ix+k+ix/2];FS[j*ix+k]=TEMP1;FS[j*ix+k+ix/2]=TEMP2;}}}returntrue;}booldiffft(complexdouble*TS,complexdouble*FS,intr){//频率抽取算法,TS为时间序列，FS为频率序列，序列长度为N=2^rinti,j,k;intN=1r;//N=2^rdoublealpha;//旋转因子的相位角complexdouble*W;//旋转因子W=newcomplexdouble[N/2];for(i=0;iN/2;i++)//计算旋转因子{alpha=-i*2*pi/N;W[i]=complexdouble(cos(alpha),sin(alpha));}memcpy(FS,TS,sizeof(complexdouble)*N);for(i=0;ir;i++){intix=1(r-i);//第i级每组的元素个数intiy=1i;//第i级的总组数for(j=0;jiy;j++)//蝶形计算过程{for(k=0;kix/2;k++){complexdoubleTEMP1;complexdoubleTEMP2;TEMP1=FS[j*ix+k]+FS[j*ix+k+ix/2];TEMP2=FS[j*ix+k]-FS[j*ix+k+ix/2];FS[j*ix+k]=TEMP1;FS[j*ix+k+ix/2]=TEMP2;}for(k=ix/2;kix;k++)FS[j*ix+k]=FS[j*ix+k]*W[(k-ix/2)*iy];}}complexdouble*X;X=newcomplexdouble[N];memcpy(X,FS,sizeof(complexdouble)*N);for(i=0;iN;i++)FS[i]=X[RevBit(i,r)];delete[]X;returntrue;}booldft(complexdouble*TS,complexdouble*FS,intr){//正常dft算法,TS为时间序列，FS为频率序列，序列长度为N=2^rinti,j;intN=1r;//N=2^rdoublealpha;//旋转因子的相位角complexdouble*W;//旋转因子W=newcomplexdouble[N];for(i=0;iN;i++)//计算旋转因子{alpha=-i*2*pi/N;W[i]=complexdouble(cos(alpha),sin(alpha));}for(i=0;iN;i++)//根据定义计算{FS[i]=0;for(j=0;jN;j++)FS[i]=FS[i]+TS[j]*W[(i*j)%(N)];}returntrue;}voidmain(void){clock_tstart,finish;inti=0,r=0,N=0;doubletemp;ifstreamfin(10.txt);//打开测试数据文件finr;//读入所要测试的阶数N=1r;//通过一位得到Ncomplexdouble*TS;TS=newcomplexdouble[N];while(fintemp)//读入时间序列x[n]{TS[i].real(temp);TS[i].imag(0);//因为本次设置的数据都是实数，故虚部置为0i++;}complexdouble*FS;//变换结果序列FS=newcomplexdouble[N];ofstreamfout_1(output10_dit.txt);start=clock();//计时开始ditfft(TS,FS,r);//调用dit函数finish=clock();//计时结束fout_1ditfft用时：finish-startendlendl;for(i=0;iN;i++)fout_1FS[i]endl;ofstreamfout_2(output10_dif.txt);start=clock();//计时开始diffft(TS,FS,r);//调用dif函数finish=clock();//计时结束fout_2diffft用时：finish-startendlendl;for(i=0;iN;i++)fout_2FS[i]endl;ofstreamfout_3(output10_dft.txt);start=clock();//开始计时dft(TS,FS,r);//调用dft函数finish=clock();//计时结束fout_3dft用时：finish-startendlendl;for(i=0;iN;i++)fout_3FS[i]endl;}